QuickloadGenome-class {rtracklayer} | R Documentation |
Quickload Genome Access
Description
A Quickload data source is a collection of tracks and sequences,
separated by genome. This class, QuickloadGenome
provides
direct access to the data for one particular genome.
Constructor
-
QuickloadGenome(quickload, genome, create = FALSE, seqinfo = seqinfo(genome), title = toString(genome))
: Constructs a newQuickloadGenome
object, representinggenome
in the repositoryquickload
(a URI string or aQuickload
object).The
genome
argument can be an ID corresponding to a genome (potentially) inquickload
or an installedBSgenome
package. It can also be any instance of a class which has methods fororganism
andreleaseDate
. A good example isBSgenome
or any other derivative ofGenomeDescription
. Those items are necessary for constructing the canonical Quickload genome string (G_Species_Month_Year).If
create
isTRUE
, and the genome does not already exist, the genome will be created, usingseqinfo
for the sequence lengths andtitle
for the display name of the genome in a UI. Creation only works if the repository is local and writeable. Reasonable defaults are used forseqinfo
andtitle
when the necessary methods are available (and they are forBSgenome
).
Accessor Methods
In the code snippets below, x
and object
represent a
Quickload
object.
-
seqinfo(x)
,seqinfo(x) <- value
: Gets or sets theSeqinfo
object indicating the lengths of the sequences in the genome. No circularity information or genome identifier is stored. -
quickload(x)
: Get the Quickload object that contains this genome. -
uri(x)
: Get the uri pointing to the genome directory in the Quickload repository -
genome(x)
: Get the name of the genome, e.g. “H_sapiens_Feb_2009”. -
releaseDate(x)
: Get the release portion of the genome name, e.g., “Feb_2009”. -
organism(object)
: Get the organism portion of the genome name, e.g., “H sapiens”.
Data Access
-
length(x)
: number of datasets -
names(x), trackNames(x)
: names of the datasets -
mcols(x)
: merged metadata on the datasets -
track(x, name), x$name
: get the track calledname
-
track(x, name, format = bestFileFormat(value), ...) <- value, x$name <- value
: store the trackvalue
undername
. Note that track storing is only supported for local repositories, i.e., those with afile://
URI scheme.Currently, supported
value
types include aGenomicRanges
,GRangesList
, or a file resource (copied to the repository). The file resource may be represented as a path, URL,BiocFile
orRsamtoolsFile
. If not a file name,value
is written informat
. For generic interval data, this means a BigWig file (if there is a numeric “score” column) or a BED file otherwise. AnRleList
(e.g., coverage) is output as BigWig. ForUCSCData
values, the format is chosen according to the type of track line. ForRsamtoolsFile
objects, the file and its index are copied.The arguments in
...
become attributes in the XML metadata. The “description” attribute is standard and is a blurb for describing the track in a UI. For the rest, the interpretation is up to the client. IGB supports an ever-growing list; please see its documentation. -
referenceSequence(x)
: Get the reference sequence, as aDNAStringSet
. -
referenceSequence(x) <- value
: Set the reference sequence, as aDNAStringSet
. It is written as a 2bit file. This only works on local repositories.
Author(s)
Michael Lawrence
Examples
tests_dir <- system.file("tests", package = "rtracklayer")
ql <- Quickload(file.path(tests_dir, "quickload"))
qlg <- QuickloadGenome(ql, "T_species_Oct_2011")
seqinfo(qlg)
organism(qlg)
releaseDate(qlg)
names(qlg)
mcols(qlg)
if (.Platform$OS.type != "windows") { # temporary
qlg$bedData
}
## Not run:
## populating the test repository
ql <- Quickload(file.path(tests_dir, "quickload"), create = TRUE)
reference_seq <- import(file.path(tests_dir, "test.2bit"))
names(reference_seq) <- "test"
qlg <- QuickloadGenome(ql, "T_species_Oct_2011", create = TRUE,
seqinfo = seqinfo(reference_seq))
referenceSequence(qlg) <- reference_seq
test_bed <- import(file.path(tests_dir, "test.bed"))
names(test_bed) <- "test"
qlg$bedData <- test_bed
test_bedGraph <- import(file.path(tests_dir, "test.bedGraph"))
names(test_bedGraph) <- "test"
start(test_bedGraph) <- seq(1, 90, 10)
width(test_bedGraph) <- 10
track(qlg, "bedGraphData", format = "bw") <- test_bedGraph
## End(Not run)