Annotates RLRanges with entrez ids for overlapping genes. See details.

geneAnnotation(object, txdb = NULL)

Arguments

object

An RLRanges object.

txdb

The TxDb or EnsDb object containing gene annotations. If not supplied, annotations will be automatically downloaded from AnnotationHub. See also GenomicFeatures::TxDb.

Value

An RLRanges object with gene overlaps included. The results are available via rlresult(object, "geneAnnoRes"). The result object is a tbl with a mapping of peak_name (peak names from names(object)) to gene_id (entrez gene IDs).

Details

The geneAnnotation function provides a simple procedure for annotating RLRanges with gene IDs by overlap.

Annotations

First. gene annotations are automatically downloaded using AnnotationHub::query with the following pattern:

AnnotationHub::query(
    x = ah,
    pattern = c("TxDb", "UCSC", "knownGene", genome)
)

Where genome is the UCSC genome id for the RLRanges object. If these annotations are unavailable, they should be provded using the txdb parameter. See also GenomicFeatures::TxDb.

Overlaps

The annotations are subsequently overlapped with the ranges in the supplied RLRanges object using valr::bed_intersect and saved in the RLResults object as a tbl with a mapping of peak names to gene_id (entrez gene IDs).

Examples


# Example RLRanges data
rlr <- readRDS(system.file("extdata", "rlrsmall.rds", package = "RLSeq"))

# Perform gene annotation
rlr <- geneAnnotation(rlr)
#> snapshotDate(): 2022-08-04
#> Loading required package: GenomicFeatures
#> Loading required package: GenomicRanges
#> Loading required package: AnnotationDbi
#> Loading required package: Biobase
#> Welcome to Bioconductor
#> 
#>     Vignettes contain introductory material; view with
#>     'browseVignettes()'. To cite Bioconductor, see
#>     'citation("Biobase")', and for packages 'citation("pkgname")'.
#>   1662 genes were dropped because they have exons located on both strands
#>   of the same reference sequence or on more than one reference sequence,
#>   so cannot be represented by a single genomic range.
#>   Use 'single.strand.genes.only=FALSE' to get all the genes in a
#>   GRangesList object, or use suppressMessages() to suppress this message.

# Supply custom TxDb if needed
if (GenomeInfoDb::genome(rlr)[1] == "hg19") {
    library(TxDb.Hsapiens.UCSC.hg19.knownGene)
    rlr <- geneAnnotation(rlr, txdb = TxDb.Hsapiens.UCSC.hg19.knownGene)
}