assignGeneID {CAGEfightR} | R Documentation |
Annotate a set of ranges in a GRanges object with gene IDs (i.e. Entrez Gene Identifiers) based on their genic context. Features overlapping multiple genes are resolved by distance to the nearest TSS. Genes are obtained from a TxDb object, or can manually supplied as a GRanges.
assignGeneID(object, geneModels, ...) ## S4 method for signature 'GenomicRanges,GenomicRanges' assignGeneID( object, geneModels, outputColumn = "geneID", swap = NULL, upstream = 1000, downstream = 100 ) ## S4 method for signature 'RangedSummarizedExperiment,GenomicRanges' assignGeneID(object, geneModels, ...) ## S4 method for signature 'GenomicRanges,TxDb' assignGeneID( object, geneModels, outputColumn = "geneID", swap = NULL, upstream = 1000, downstream = 100 ) ## S4 method for signature 'RangedSummarizedExperiment,TxDb' assignGeneID(object, geneModels, ...)
object |
GRanges or RangedSummarizedExperiment: Ranges to be annotated. |
geneModels |
TxDb or GRanges: Gene models via a TxDb, or manually specified as a GRangesList. |
... |
additional arguments passed to methods. |
outputColumn |
character: Name of column to hold geneID values. |
swap |
character or NULL: If not NULL, use another set of ranges contained in object to calculate overlaps, for example peaks in the thick column. |
upstream |
integer: Distance to extend annotated promoter upstream. |
downstream |
integer: Distance to extend annotated promoter downstream. |
object with geneID added as a column in rowData (or mcols).
Other Annotation functions:
assignMissingID()
,
assignTxID()
,
assignTxType()
data(exampleUnidirectional) # Obtain gene models from a TxDb-object: library(TxDb.Mmusculus.UCSC.mm9.knownGene) txdb <- TxDb.Mmusculus.UCSC.mm9.knownGene # Assign geneIDs assignGeneID(exampleUnidirectional, geneModels=txdb, outputColumn='geneID') # Assign geneIDs using only TC peaks: assignGeneID(exampleUnidirectional, geneModels=txdb, outputColumn='geneID', swap='thick')