VcfFile {VariantAnnotation} | R Documentation |
Use VcfFile()
to create a reference to a Vcf file (and its
index). Once opened, the reference remains open across calls to
methods, avoiding costly index re-loading.
VcfFileList()
provides a convenient way of managing a list of
VcfFile
instances.
## Constructors
VcfFile(file, index = paste(file, "tbi", sep="."), ..., yieldSize=NA_integer_) VcfFileList(..., yieldSize=NA_integer_)
## Accessors index(object) path(object, ...) isOpen(con, rw="") yieldSize(object, ...) yieldSize(object, ...) <- value show(object)
## Opening / closing open(con, ...) close(con, ...)
conAn instance of VcfFile
.
fileA character(1) vector to the Vcf file path; can be remote (http://, ftp://).
indexA character(1) vector of the Vcf file index (.tbi file).
yieldSizeNumber of records to yield each time the file is read
from using scanVcf
or readVcf
.
...Additional arguments. For VcfFileList
, this can
either be a single character vector of paths to Vcf files, or
several instances of VcfFile
objects.
rwcharacter() indicating mode of file.
Objects are created by calls of the form VcfFile()
.
VcfFile
and VcfFileList
classes inherit fields from the
TabixFile
and TabixFileList
classes.
VcfFile
and VcfFileList
classes inherit methods from the
TabixFile
and TabixFileList
classes.
Opening / closing:
Opens the (local or remote) path
and
index
. Returns a VcfFile
instance.
yieldSize
determines the number of records parsed during
each call to scanVcf
or readVcf
; NA
indicates
that all records are to be parsed.
Closes the VcfFile
con
; returning
(invisibly) the updated VcfFile
. The instance may be
re-opened with open.VcfFile
.
Accessors:
Returns a character(1) vector of the Vcf path name.
Returns a character(1) vector of Vcf index (tabix file) name.
Return or set an integer(1) vector indicating yield size.
Valerie Obenchain
fl <- system.file("extdata", "chr7-sub.vcf.gz", package="VariantAnnotation", mustWork=TRUE) vcffile <- VcfFile(fl) vcffile param <- GRanges("7", IRanges(c(55000000, 55900000), width=10000)) vcf <- readVcf(vcffile, "hg19", param=param) dim(vcf)