VariantExperiment 1.2.0
This package includes 2 vignettes. This One is about the class
definitions, basic methods and operations of VariantExperiment that
are defined or directly inherited from SummarizedExperiment. For
users who want to apply VariantExperiment directly to their analysis
high-throughput genotyping or DNA-seq data sets, please refer to the
other vignette.
As the sequencing data (DNA-seq, RNA-seq, single cell RNA-seq…) gets increasingly larger-profile, the memory space in R has been an obstable for fast and efficient data processing, because most available R or Bioconductor packages are developed based on in-memory data manipulation. SingleCellExperiment has achieved efficient on-disk saving/reading of the large-scale count data as HDF5Array objects. However, there was still no such light-weight containers available for high-throughput DNA-seq / genotyping data.
We have developed VariantExperiment, a Bioconductor package for
saving the VCF/GDS format data into RangedSummarizedExperiment
object. The high-throughput genetic/genomic data are saved in
GDSArray objects, a direct extension of DelayedArray with GDS
back-end. In addition to the light-weight Assay data, We have also
realized the on-disk saving of annotation data for both
features/samples (corresponding to rowData/colData) by implementing
the DelayedDataFrame data structure. The on-disk representation
of both assay data and annotation data realizes on-disk reading and
processing and saves R memory space significantly. The interface of
RangedSummarizedExperiment data format enables easy and common
manipulations for high-throughput genetic/genomic data with common
SummarizedExperiment metaphor in R and Bioconductor.
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("VariantExperiment")Or install the development version of the package from Github.
BiocManager::install(“Bioconductor/VariantExperiment”) library(VariantExperiment)GDSArray is a Bioconductor package that represents GDS files as
objects derived from the DelayedArray package and DelayedArray
class. It converts GDS nodes into a DelayedArray-derived data
structure. The rich common methods and data operations defined on
GDSArray makes it more R-user-friendly than working with the GDS
file directly. The array data from GDS files are always returned with
the first dimension being features (genes/variants/snps) and the
second dimension being samples. This feature is consistent with the
assay data saved in SummarizedExperiment, and makes the GDSArray
package more interoperable with other established Bioconductor data
infrastructure.
The GDSArray() constructor takes 2 arguments: the file path and the
GDS node name inside the GDS file.
file <- SeqArray::seqExampleFileName("gds")
#> GDSArray(file, "genotype") #> deactivate temporarily 3/5/20
GDSArray(file, "sample.id")## <90> array of class GDSArray and type "character":
## NA06984 NA06985 NA06986 . NA12891 NA12892
## "NA06984" "NA06985" "NA06986" . "NA12891" "NA12892"More details about GDS or GDSArray format can be found in the
vignettes of the gdsfmt, SNPRelate, SeqArray, GDSArray
and DelayedArray packages.
DelayedDataFrame is a Bioconductor package that implements
delayed operations on DataFrame objects using standard DataFrame
metaphor. Each column of data inside DelayedDataFrame is represented
as 1-dimensional GDSArray with on-disk GDS file. Methods like
show,validity check, [, [[ subsetting, rbind, cbind are implemented for
DelayedDataFrame. The DelayedDataFrame stays lazy until an
explicit realization call like DataFrame() constructor or
as.list() incurred. More details about DelayedDataFrame data
structure could be found in the vignette of DelayedDataFrame
package.
VariantExperiment classVariantExperiment classVariantExperiment class is defined to extend
RangedSummarizedExperiment. The difference would be that the assay
data are saved as GDSArray, and the annotation data are saved by
default as DelayedDataFrame (with option to save as ordinary
DataFrame), both of which are representing the data on-disk with GDS
back-end. There are coercion methods defined for VCF and GDS files into
VariantExperiment objects (Check the method vignette for more details).
Here we take an example for illustration.
ve <- makeVariantExperimentFromGDS(gds)
veassay data are in GDSArray format, and could be retrieve by the
assays()/assay() function.
assays(ve)
assay(ve, 1)The rowData() of the VariantExperiment is saved in
DelayedDataFrame format. We can use rowRanges() / rowData() to
retrieve the feature-related annotation file, with/without a
GenomicRange format.
rowRanges(ve)
rowData(ve)sample-related annotation is in DelayedDataFrame format, and could
be retrieved by colData().
colData(ve)The gdsfile() will retrieve the gds file path associated with the
VariantExperiment object.
gdsfile(ve)Some other getter function like metadata() will return any metadata
that we have saved inside the VariantExperiment object.
metadata(ve)The VariantExperiment object supports subsetting methods with [,
$ and subsetByOverlap etc. as in SummarizedExperiment.
ve[1:10, 1:5]$ subsettingThe $ subsetting could be operated directly on colData() columns,
for easy sample extraction. Note that the colData/rowData are in the
DelayedDataFrame format, with each column saved as GDSArray. So
when doing subsetting, we need to use as.logical() to convert the
1-dimensional GDSArray into ordinary vector.
colData(ve)
ve[, as.logical(ve$family == "1328")]subsetting by rowData() columns.
rowData(ve)
ve[as.logical(rowData(ve)$REF == "T"),]VariantExperiment objects support all of the findOverlaps()
methods and associated functions. This includes subsetByOverlaps(),
which makes it easy to subset a VariantExperiment object by an
interval.
ve1 <- subsetByOverlaps(ve, GRanges("22:1-48958933"))
ve1VariantExperiment objectNote that after the subsetting by [, $ or Ranged-based operations, and you feel satisfied with the data for downstream
analysis, you need to save that VariantExperiment object to
synchronize the gds file (on-disk) associated with the subset of data
(in-memory representation) before any statistical analysis. Otherwise,
an error will be returned.
For example, after we subset the ve by GRanges("22:1-48958933"),
and we want to calculate the hwe based on the 23 variants, an error
will be generated indicating that we need to sync the on-disk and
in-memory representations.
hwe(ve1)
## Error in .saveGDSMaybe(gdsfile) : use
## 'saveVariantExperiment()' to synchronize on-disk and
## in-memory representationsVariantExperiment objectUse the function saveVariantExperiment to synchronize the on-disk
and in-memory representation. This function writes the processed data
as ve.gds, and save the R object (which lazily represent the
backend data set) as ve.rds under the specified directory. It
finally returns a new VariantExperiment object into current R
session generated from the newly saved data.
a <- tempfile()
ve2 <- saveVariantExperiment(ve1, dir=a, replace=TRUE)VariantExperiment objectYou can alternatively use loadVariantExperiment to load the
synchronized data into R session, by providing only the file
directory. It reads the VariantExperiment object saved as ve.rds, as lazy
representation of the backend ve.gds file under the specific
directory.
ve3 <- loadVariantExperiment(dir=a)
gdsfile(ve3)
all.equal(ve2, ve3)Now we are all set for any downstream analysis as needed.
head(hwe(ve2))sessionInfo()## R version 4.0.0 (2020-04-24)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.4 LTS
##
## Matrix products: default
## BLAS: /home/biocbuild/bbs-3.11-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.11-bioc/R/lib/libRlapack.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=C
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] parallel stats4 stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] VariantExperiment_1.2.0 DelayedDataFrame_1.4.0
## [3] GDSArray_1.8.0 gdsfmt_1.24.0
## [5] SummarizedExperiment_1.18.0 DelayedArray_0.14.0
## [7] matrixStats_0.56.0 Biobase_2.48.0
## [9] GenomicRanges_1.40.0 GenomeInfoDb_1.24.0
## [11] IRanges_2.22.0 S4Vectors_0.26.0
## [13] BiocGenerics_0.34.0 BiocStyle_2.16.0
##
## loaded via a namespace (and not attached):
## [1] tidyselect_1.0.0 xfun_0.13 purrr_0.3.4
## [4] splines_4.0.0 lattice_0.20-41 generics_0.0.2
## [7] vctrs_0.2.4 GWASExactHW_1.01 htmltools_0.4.0
## [10] yaml_2.2.1 mgcv_1.8-31 rlang_0.4.5
## [13] pillar_1.4.3 glue_1.4.0 GenomeInfoDbData_1.2.3
## [16] lifecycle_0.2.0 stringr_1.4.0 zlibbioc_1.34.0
## [19] Biostrings_2.56.0 evaluate_0.14 knitr_1.28
## [22] SeqArray_1.28.0 broom_0.5.6 Rcpp_1.0.4.6
## [25] backports_1.1.6 BiocManager_1.30.10 XVector_0.28.0
## [28] digest_0.6.25 stringi_1.4.6 bookdown_0.18
## [31] dplyr_0.8.5 grid_4.0.0 tools_4.0.0
## [34] bitops_1.0-6 magrittr_1.5 RCurl_1.98-1.2
## [37] tibble_3.0.1 mice_3.8.0 tidyr_1.0.2
## [40] crayon_1.3.4 SNPRelate_1.22.0 pkgconfig_2.0.3
## [43] ellipsis_0.3.0 Matrix_1.2-18 assertthat_0.2.1
## [46] rmarkdown_2.1 logistf_1.23 R6_2.4.1
## [49] SeqVarTools_1.26.0 nlme_3.1-147 compiler_4.0.0