transcripts {GenomicFeatures} | R Documentation |
Generic functions to extract genomic features from an object. This page documents the methods for TranscriptDb objects only.
transcripts(x, ...) ## S4 method for signature 'TranscriptDb' transcripts(x, vals=NULL, columns=c("tx_id", "tx_name")) exons(x, ...) ## S4 method for signature 'TranscriptDb' exons(x, vals=NULL, columns="exon_id") cds(x, ...) ## S4 method for signature 'TranscriptDb' cds(x, vals=NULL, columns="cds_id") genes(x, ...) ## S4 method for signature 'TranscriptDb' genes(x, vals=NULL, columns="gene_id", single.strand.genes.only=TRUE) #promoters(x, upstream=2000, downstream=200, ...) ## S4 method for signature 'TranscriptDb' promoters(x, upstream=2000, downstream=200, ...) disjointExons(x, ...) ## S4 method for signature 'TranscriptDb' disjointExons(x, aggregateGenes=FALSE, includeTranscripts=TRUE, ...) microRNAs(x) ## S4 method for signature 'TranscriptDb' microRNAs(x) tRNAs(x) ## S4 method for signature 'TranscriptDb' tRNAs(x)
x |
A TranscriptDb object. |
... |
Arguments to be passed to or from methods. |
vals |
Either |
columns |
Columns to include in the output.
Must be
If the vector is named, those names are used for the corresponding column in the element metadata of the returned object. |
single.strand.genes.only |
|
upstream |
For |
downstream |
For |
aggregateGenes |
For |
includeTranscripts |
For |
These are the main functions for extracting transcript information
from a TranscriptDb object. With the exception of
microRNAs
, these methods can restrict the output based on
categorical information. To restrict the output based on interval
information, use the transcriptsByOverlaps
,
exonsByOverlaps
, and cdsByOverlaps
functions.
The promoters
function computes user-defined promoter regions
for the transcripts in a TranscriptDb
object. The return object is a
GRanges
of promoter regions around the transcription start
site the span of which is defined by upstream
and downstream
.
For additional details on how the promoter range is computed and the
handling of +
and -
strands see
?'promoters,GRanges-method'
.
disjointExons
creates a GRanges
of non-overlapping
exon parts with metadata columns of gene_id and exonic_part.
Exon parts that overlap more than 1 gene can be dropped with
aggregateGenes=FALSE
. When includeTranscripts=TRUE
a tx_name
metadata column is included that lists all
transcript names that overlap the exon fragment. This function
replaces prepareAnnotationForDEXSeq
in the DEXSeq
package.
A GRanges object. The only exception being
when genes
is used with single.strand.genes.only=FALSE
,
in which case a GRangesList object is returned.
M. Carlson, P. Aboyoun and H. Pages. disjointExons
was contributed
by Mike Love and Alejandro Reyes.
transcriptsBy
and transcriptsByOverlaps
for more ways to extract genomic features
from a TranscriptDb object.
select-methods for how to use the simple "select" interface to extract information from a TranscriptDb object.
id2name
for mapping TranscriptDb internal ids
to external names for a given feature type.
The TranscriptDb class.
## transcripts(), exons(), genes(): txdb <- loadDb(system.file("extdata", "UCSC_knownGene_sample.sqlite", package="GenomicFeatures")) vals <- list(tx_chrom = c("chr3", "chr5"), tx_strand = "+") transcripts(txdb, vals) exons(txdb, vals=list(exon_id=1), columns=c("EXONID", "TXNAME")) exons(txdb, vals=list(tx_name="uc009vip.1"), columns=c("EXONID", "TXNAME")) genes(txdb) # a GRanges object cols <- c("tx_id", "tx_chrom", "tx_strand", "exon_id", "exon_chrom", "exon_strand") single_strand_genes <- genes(txdb, columns=cols) ## Because we've returned single strand genes only, the "tx_chrom" ## and "exon_chrom" metadata columns are guaranteed to match ## 'seqnames(single_strand_genes)': stopifnot(identical(as.character(seqnames(single_strand_genes)), as.character(mcols(single_strand_genes)$tx_chrom))) stopifnot(identical(as.character(seqnames(single_strand_genes)), as.character(mcols(single_strand_genes)$exon_chrom))) ## and also the "tx_strand" and "exon_strand" metadata columns are ## guaranteed to match 'strand(single_strand_genes)': stopifnot(identical(as.character(strand(single_strand_genes)), as.character(mcols(single_strand_genes)$tx_strand))) stopifnot(identical(as.character(strand(single_strand_genes)), as.character(mcols(single_strand_genes)$exon_strand))) all_genes <- genes(txdb, columns=cols, single.strand.genes.only=FALSE) all_genes # a GRangesList object multiple_strand_genes <- all_genes[elementLengths(all_genes) >= 2] multiple_strand_genes mcols(multiple_strand_genes) ## microRNAs() : ## Not run: library(TxDb.Hsapiens.UCSC.hg19.knownGene) library(mirbase.db) microRNAs(TxDb.Hsapiens.UCSC.hg19.knownGene) ## End(Not run) ## promoters() : head(promoters(txdb, 100, 50))