Take a user-defined gene vector and dynamically map the input to either the object rownames or the gene names (symbols). These functions are useful for writing code that needs to handle either gene identifier or gene name input dynamically (e.g. for single-cell RNA-seq marker analysis).

mapGenesToRownames(object, ...)

mapGenesToIDs(object, ...)

mapGenesToSymbols(object, ...)

# S4 method for Gene2Symbol
mapGenesToRownames(object, genes, strict = TRUE)

# S4 method for SummarizedExperiment
mapGenesToRownames(object, genes,
  strict = TRUE)

# S4 method for Gene2Symbol
mapGenesToIDs(object, genes, strict = TRUE)

# S4 method for SummarizedExperiment
mapGenesToIDs(object, genes,
  strict = TRUE)

# S4 method for Gene2Symbol
mapGenesToSymbols(object, genes, strict = TRUE)

# S4 method for SummarizedExperiment
mapGenesToSymbols(object, genes,
  strict = TRUE)

Arguments

object

Object.

genes

character. Gene identifiers.

strict

logical(1). Require all genes to match. Recommended by default. If set FALSE, instead will return a warning to the user, and subset the genes vector to only include matches.

...

Additional arguments.

Value

character.

Note

Updated 2019-08-21.

Ambiguous gene names

Some genomes (e.g. Homo sapiens, Mus musculus) contain duplicated gene names for multiple gene identifiers. Normally we handle these ambiguous gene names by sanitizing them with make.names. If a user requests a gene name that is duplicated, these functions will return a warning.

Examples

data(RangedSummarizedExperiment, package = "acidtest") rse <- RangedSummarizedExperiment print(rse)
#> class: RangedSummarizedExperiment #> dim: 500 12 #> metadata(3): version date interestingGroups #> assays(1): counts #> rownames(500): gene001 gene002 ... gene499 gene500 #> rowData names(8): broadClass description ... geneName seqCoordSystem #> colnames(12): sample01 sample02 ... sample11 sample12 #> colData names(1): condition
rownames <- head(rownames(rse)) print(rownames)
#> [1] "gene001" "gene002" "gene003" "gene004" "gene005" "gene006"
g2s <- Gene2Symbol(rse) geneIDs <- head(g2s[["geneID"]]) print(geneIDs)
#> [1] "ENSG00000000003" "ENSG00000000005" "ENSG00000000419" "ENSG00000000457" #> [5] "ENSG00000000460" "ENSG00000000938"
geneNames <- head(g2s[["geneName"]]) print(geneNames)
#> [1] "TSPAN6" "TNMD" "DPM1" "SCYL3" "C1orf112" "FGR"
## Row names. mapGenesToRownames(rse, genes = rownames)
#> gene001 gene002 gene003 gene004 gene005 gene006 #> "gene001" "gene002" "gene003" "gene004" "gene005" "gene006"
mapGenesToRownames(rse, genes = geneIDs)
#> ENSG00000000003 ENSG00000000005 ENSG00000000419 ENSG00000000457 ENSG00000000460 #> "gene001" "gene002" "gene003" "gene004" "gene005" #> ENSG00000000938 #> "gene006"
mapGenesToRownames(rse, genes = geneNames)
#> TSPAN6 TNMD DPM1 SCYL3 C1orf112 FGR #> "gene001" "gene002" "gene003" "gene004" "gene005" "gene006"
## Gene identifiers. mapGenesToIDs(rse, genes = rownames)
#> gene001 gene002 gene003 gene004 #> "ENSG00000000003" "ENSG00000000005" "ENSG00000000419" "ENSG00000000457" #> gene005 gene006 #> "ENSG00000000460" "ENSG00000000938"
mapGenesToIDs(rse, genes = geneIDs)
#> ENSG00000000003 ENSG00000000005 ENSG00000000419 ENSG00000000457 #> "ENSG00000000003" "ENSG00000000005" "ENSG00000000419" "ENSG00000000457" #> ENSG00000000460 ENSG00000000938 #> "ENSG00000000460" "ENSG00000000938"
mapGenesToIDs(rse, genes = geneNames)
#> TSPAN6 TNMD DPM1 SCYL3 #> "ENSG00000000003" "ENSG00000000005" "ENSG00000000419" "ENSG00000000457" #> C1orf112 FGR #> "ENSG00000000460" "ENSG00000000938"
## Gene names (symbols). mapGenesToSymbols(rse, genes = rownames)
#> gene001 gene002 gene003 gene004 gene005 gene006 #> "TSPAN6" "TNMD" "DPM1" "SCYL3" "C1orf112" "FGR"
mapGenesToSymbols(rse, genes = geneIDs)
#> ENSG00000000003 ENSG00000000005 ENSG00000000419 ENSG00000000457 ENSG00000000460 #> "TSPAN6" "TNMD" "DPM1" "SCYL3" "C1orf112" #> ENSG00000000938 #> "FGR"
mapGenesToSymbols(rse, genes = geneNames)
#> TSPAN6 TNMD DPM1 SCYL3 C1orf112 FGR #> "TSPAN6" "TNMD" "DPM1" "SCYL3" "C1orf112" "FGR"