Maps between gene identifier systems (Ensembl, Entrez, HGNC symbol, UniProt, RefSeq, MANE) using AnnotationDbi, biomaRt, mygene, pyensembl, and Ensembl REST. Encodes Ensembl version stripping with GENCODE _PAR_Y preservation, the Ziemann 2016 Excel autocorrect debacle and Bruford 2020 HGNC renames (SEPT*->SEPTIN*, MARCH*->MARCHF*, MARC*->MTARC*, DEC1->DELEC1), OCT4/POU5F1 alias resolution, biomaRt archive endpoints for release pinning, the `filters` (plural) gotcha, MANE Select for clinical reporting, cross-species orthology via Ensembl Compara / OMA / OrthoDB, and tx2gene construction for tximport. Use when converting gene IDs across systems, handling renamed symbols, building tx2gene, pinning to a specific Ensembl release for reproducibility, or mapping cross-species orthologs.
Maps between gene identifier systems (Ensembl, Entrez, HGNC symbol, UniProt, RefSeq, MANE) using AnnotationDbi, biomaRt, mygene, pyensembl, and Ensembl REST. Encodes Ensembl version stripping with GENCODE _PAR_Y preservation, the Ziemann 2016 Excel autocorrect debacle and Bruford 2020 HGNC renames (SEPT*->SEPTIN*, MARCH*->MARCHF*, MARC*->MTARC*, DEC1->DELEC1), OCT4/POU5F1 alias resolution, biomaRt archive endpoints for release pinning, the `filters` (plural) gotcha, MANE Select for clinical reporting, cross-species orthology via Ensembl Compara / OMA / OrthoDB, and tx2gene construction for tximport. Use when converting gene IDs across systems, handling renamed symbols, building tx2gene, pinning to a specific Ensembl release for reproducibility, or mapping cross-species orthologs.
Before using code patterns, verify installed versions match. If versions differ:
R: packageVersion('<pkg>') then ?function_name to verify parameters
Python: pip show <package> then help(module.function) to check signatures
If code throws ImportError, AttributeError, or TypeError, introspect the installed
package and adapt the example to match the actual API rather than retrying.
Gene ID Mapping
"Convert gene IDs from X to Y" -> Query the appropriate annotation source (local org.db for speed, biomaRt for Ensembl-specific attributes, mygene for cross-database aliases, Ensembl REST for low-level access), with version pinning for reproducibility and explicit handling of one-to-many mappings, withdrawn symbols, and species-specific naming.
The Single Most Important Modern Insight -- Excel autocorrect renamed real genes
Ziemann, Eren, El-Osta 2016 Genome Biol 17:177 scanned 18 leading genomics journals and found ~20% of papers with Excel-attached supplementary gene lists had silently mangled symbols (SEPT2 -> 2-Sep, MARCH1 -> 1-Mar, ...). Five years later the problem persisted. HGNC's response (Bruford, Braschi, Denny, Jones, Seal, Tweedie 2020 Nat Genet 52:754) was to rename the affected genes:
Old
New
Affected
SEPT#
SEPTIN#
SEPT1 - SEPT14 -> SEPTIN1 - SEPTIN14
MARCH#
MARCHF#
MARCH1 - MARCH11 -> MARCHF1 - MARCHF11
MARC#
MTARC#
MARC1, MARC2 -> MTARC1, MTARC2
DEC1
DELEC1
DEC1 -> DELEC1
Code that hard-codes old symbols silently drops these genes when joined against post-2020 annotations. Detection on import: if a gene column contains ^\d{1,2}-(Jan|Feb|Mar|...|Dec)$ patterns, the file was Excel-corrupted. Always read.csv(colClasses=c(gene='character')) (R) or pd.read_csv(dtype={'gene': str}) (Python) -- but the damage is at Excel-save time, not import time.
Two related insights that determine half the practical work:
Ensembl version suffixes matter sometimes and not others.ENSG00000123456.7 is release-specific; the unversioned ENSG00000123456 is the stable cross-release ID. STRIP for cross-release joins, MSigDB lookups, gene-set databases. KEEP for intra-release reproducibility and clinical reports. CRITICAL: the naive sub('\\..*', '', x) regex ALSO strips the GENCODE _PAR_Y suffix in releases 25-43, collapsing chrY PAR duplicates onto their chrX counterparts. Use sub('\\.[0-9]+(_PAR_Y)?$', '\\1', x).
Never use HGNC symbols as the primary computational key. Symbols change. Use Ensembl or Entrez as keys; carry symbols only as display labels in the final results table.
Algorithmic Taxonomy
Tool
Source
Speed
Strength
Use for
AnnotationDbi + org.Hs.eg.db / org.Mm.eg.db
NCBI Gene snapshot, pinned at Bioc install
Fast, local
Stable, version-pinned
Default for Ensembl <-> Entrez <-> Symbol within Bioconductor
multiVals options: 'first' (silent), 'asNA' (NA for ambiguous), 'list' (preserve all). For DE results tables, 'first' is typical but the mapping rate should be reported.
For mouse: org.Mm.eg.db. For other organisms: check Bioconductor AnnotationData -> OrgDb list.
biomaRt with Version Pinning
Goal: Query Ensembl BioMart with the EXACT release version, for reproducibility.
Approach:useEnsembl(version=N) pins; listEnsemblArchives() lists available archives.
The filters= argument is PLURAL. The singular filter= may work via R's partial matching but breaks unpredictably if another argument starts with f. Always spell filters= and values= fully.
Without version=, biomaRt floats to the current release -- a script written in 2023 against Ensembl 109 produces different mappings in 2026 against Ensembl 113. ALWAYS pin for any published analysis. Cache the mapping table alongside the analysis for reproducibility.
listEnsemblArchives() shows the available historical releases.
mygene (Python)
Goal: Map between any identifier systems using the curated MyGene.info meta-database with alias fallback.
Approach:MyGeneInfo().querymany(ids, scopes, fields, species); auto-batches at 1000 IDs server-side.
import mygene
mg = mygene.MyGeneInfo()
results = mg.querymany(['ENSG00000141510', 'ENSG00000012048', 'ENSG00000141736'],
scopes='ensembl.gene', fields='symbol,entrezgene,uniprot',
species='human')
mapping = {r['query']: r.get('symbol', None) for r in results}
results = mg.querymany(['SEPT1', 'MARCH1', 'OCT4'],
scopes='symbol,alias,prev_symbol',
fields='symbol,entrezgene,ensembl.gene',
species='human')
For paper-derived gene lists where symbols may be old or aliases (OCT4 vs POU5F1, MARCH1 vs MARCHF1, SEPT2 vs SEPTIN2), scopes='symbol,alias,prev_symbol' handles the resolution. The MyGene database aggregates HGNC's prev/alias columns.
OCT4 is the common usage; POU5F1 is the official HGNC symbol; in MSigDB the gene is POU5F1; in a Western blot legend it's "Oct4". For mapping a stem-cell paper to an Ensembl-quantified matrix, scope to aliases.
pyensembl downloads and caches the release database on first use; thereafter offline and version-locked.
Apply Mapping to Count Matrix -- Handling One-to-Many
Goal: Convert the gene index of a count matrix to a different ID type, summing reads from multiple source IDs that map to the same target.
Approach: Look up mapping, replace index, aggregate duplicates by SUM (not mean -- counts add).
import pandas as pd
import mygene
defmap_count_matrix_ids(counts, from_type='ensembl.gene', to_type='symbol',
species='human'):
'''Map gene IDs in count matrix index, summing reads when multiple source map to one target.'''
mg = mygene.MyGeneInfo()
clean = [g.split('.')[0] for g in counts.index]
results = mg.querymany(clean, scopes=from_type, fields=to_type, species=species)
mapping = {r['query']: r[to_type] for r in results if to_type in r}
new_index = [mapping.get(g.split('.')[0], g) for g in counts.index]
counts_mapped = counts.copy()
counts_mapped.index = new_index
counts_mapped = counts_mapped.groupby(counts_mapped.index).sum()
return counts_mapped
mapped = map_count_matrix_ids(counts, 'ensembl.gene', 'symbol')
Counts ADD when collapsing multiple source genes to one target. Means or medians would be wrong (they understate library size for the merged target).
defrobust_id_mapping(gene_ids, from_type, to_type, species='human'):
import mygene
mg = mygene.MyGeneInfo()
clean = [g.split('.')[0] for g in gene_ids]
results = mg.querymany(clean, scopes=from_type, fields=to_type, species=species)
mapping, unmapped = {}, []
for r in results:
original = gene_ids[clean.index(r['query'])]
if to_type in r:
mapping[original] = r[to_type]
else:
mapping[original] = original
unmapped.append(original)
print(f'Mapped: {len(gene_ids) - len(unmapped)}/{len(gene_ids)}')
return mapping, unmapped
Unmapped fraction is a QC signal:
<5% unmapped: normal (rare genes, recent deprecations)
5-20% unmapped: check Ensembl release alignment; check for HGNC renames
20% unmapped: wrong annotation release, wrong species, or wrong source ID type
MANE Select for Clinical Reporting
Goal: Use the single representative transcript per gene with identical exon/CDS in RefSeq AND Ensembl for clinical variant reporting.
Approach: Download the MANE TSV; join on Ensembl_Gene -> Ensembl_nuc (transcript) and RefSeq_nuc.
Morales J, Pujar S, Loveland JE et al. 2022 Nature 604:310-315 established MANE Select. ~19,000+ protein-coding genes have a single agreed transcript with matched coordinates across RefSeq (NM_xxxxxx) and Ensembl/GENCODE (ENST00000xxxxxxx). MANE Plus Clinical adds extra transcripts at loci where Select misses clinical variants.
For clinical reports with HGVS notation like NM_000546.6:c.215C>G, use the MANE Select RefSeq accession. The MANE TSV (downloadable from NCBI) provides the Ensembl crosswalk.
Cross-Species Orthologs
Goal: Map mouse <-> human (or any pair) for cross-species integration or pathway transfer.
Approach: Ensembl Compara via biomaRt getLDS; filter to orthology type appropriate to use.
library(biomaRt)
human <- useEnsembl(biomart ='genes', dataset ='hsapiens_gene_ensembl', version =110)
mouse <- useEnsembl(biomart ='genes', dataset ='mmusculus_gene_ensembl', version =110)
orthologs <- getLDS(attributes=c('hgnc_symbol','ensembl_gene_id'),
filters ='ensembl_gene_id',
values = human_gene_ids,
mart = human,
attributesL =c('mgi_symbol','ensembl_gene_id','mmusculus_homolog_orthology_type'),
martL = mouse
)
Strategy
When
Trade-off
one2one orthologs only
Cross-species scRNA-seq integration; conservative DE comparison
Loses genes with paralog expansions; lower coverage
Include one2many
Broader gene coverage needed
Must select within group (highest confidence; highest expression)
Include many2many
Maximum inclusivity
Introduces ambiguity; use with caution
The "homology threshold" problem: no automatic threshold reliably separates true orthologs from paralogs across all gene families. For pathway transfer (mouse signature -> human), filter to one2one and accept the coverage loss.
Alternative sources: OMA (Hierarchical Orthologous Groups, cleaner one2one when present, smaller coverage); OrthoDB (hierarchical at multiple taxonomic levels). OrthoFinder for custom genomes.
PAR Gene Complications
Pseudo-autosomal region (PAR) genes exist on both X and Y with identical sequences. In GENCODE 25-43, the chrY copy has a _PAR_Y suffix. In GENCODE 44+ (Ensembl 110+), chrY PAR genes get their own ENSG accessions.
Reads from PAR regions cannot be unambiguously assigned to X or Y. Some references mask the Y-chromosome PAR to avoid double-counting; verify what the alignment reference does before building the matrix.
Build tx2gene for tximport
Goal: Create the transcript-to-gene mapping needed by tximport for gene-level summarization.
Approach: Build from the SAME GTF used to construct the Salmon/kallisto index, OR pull from biomaRt with version pinning.
library(biomaRt)
mart <- useEnsembl(biomart ='genes', dataset ='hsapiens_gene_ensembl', version =110)
tx2gene <- getBM(attributes=c('ensembl_transcript_id_version','ensembl_gene_id_version'),
mart = mart
)
colnames(tx2gene)<-c('TXNAME','GENEID')
import pandas as pd
deftx2gene_from_gtf(gtf_path):
records = []
withopen(gtf_path) as f:
for line in f:
if line.startswith('#') or'\ttranscript\t'notin line:
continue
attrs = line.strip().split('\t')[8]
gene_id = [a.split('"')[1] for a in attrs.split(';') if'gene_id'in a][0]
tx_id = [a.split('"')[1] for a in attrs.split(';') if'transcript_id'in a][0]
records.append({'TXNAME': tx_id, 'GENEID': gene_id})
return pd.DataFrame(records).drop_duplicates()
CRITICAL: the tx2gene MUST use the same versioning convention as the Salmon/kallisto index. If the index used ENST00000269305.9 and tx2gene has ENST00000269305 (unversioned), tximport drops the transcripts. Mismatched versions silently lose data.
ID Type Reference
Type
Example
Stability
Use case
Ensembl Gene
ENSG00000141510
Stable across releases; versioned
RNA-seq, GTFs, primary computational key
Ensembl Transcript
ENST00000269305
Stable; versioned
Transcript-level analysis
Entrez Gene
7157
Stable; never reused
NCBI databases, KEGG pathways
HGNC Symbol
TP53
Changes (see SEPT/MARCH renames)
Display labels only
UniProt
P04637
Stable; versioned releases
Protein databases
RefSeq mRNA
NM_000546
Stable; versioned
Clinical reports, HGVS notation
MANE Select
NM_000546.6 / ENST00000269305.9
Stable consensus
Clinical variant reporting
Per-Method Failure Modes
_PAR_Y stripped, chrY duplicates collapsed
Trigger: GENCODE v40 count matrix; rownames(counts) <- sub('\\..*', '', rownames(counts)); duplicate row indices and inflated chrY PAR gene counts.
Mechanism: Default regex strips _PAR_Y along with the version suffix. Two distinct rows (chrX and chrY copies) become the same ENSG ID; aggregate sums them.
Symptom: Counts for PAR genes double; sex check shows females expressing chrY genes; downstream rowGroupBy returns warnings.
Fix: Use the preserving regex: sub('\\.[0-9]+(_PAR_Y)?$', '\\1', x). Or upgrade quantification to GENCODE 44+ where _PAR_Y is retired.
Mechanism: R's partial matching usually resolves filter -> filters, but in some package versions or with conflicting argument names, the call silently passes nothing.
Symptom: Empty result data frame; no error.
Fix: Always spell filters= and values= fully.
HGNC SEPT/MARCH symbols silently dropped
Trigger: Code copies a pre-2020 list of septin genes (SEPT1, SEPT2, ...); current org.db / biomaRt returns no matches.
Mechanism: HGNC renamed all SEPT# to SEPTIN# in 2020.
Fix: Use mygene scopes='symbol,alias,prev_symbol'; or update the input list to current symbols.
biomaRt drift between runs
Trigger: A 2023 analysis used useEnsembl() without version=; rerun in 2026 produces 200 fewer significant genes.
Mechanism: Without version=, biomaRt floats to the current release. Symbols, biotypes, and gene boundaries change between releases.
Symptom: Non-reproducible results across runs of the same script.
Fix: Pin useEnsembl(version=N) where N is the release used in the original analysis. Cache the mapping table.
tx2gene version mismatch with Salmon index
Trigger:tximport(files, type='salmon', tx2gene) runs but the gene-level counts have far fewer genes than expected.
Mechanism: Salmon index built with versioned transcript IDs (ENST00000269305.9) but tx2gene has unversioned IDs (ENST00000269305). Transcripts silently drop during the mapping step.
Symptom: Lower-than-expected gene count; warning from tximport about missing transcript IDs.
Fix: Match versioning convention: rebuild tx2gene with the same versioning as the index. GenomicFeatures::makeTxDbFromGFF on the same GTF as the index is the safest path.
Cross-species mapping reports many2many, user picks one arbitrarily
Trigger: Mouse-to-human mapping returns 1.3 mouse genes per human gene on average; user takes the first row of each duplicate.
Mechanism: Many2many orthology is genuinely ambiguous; "first row" is unprincipled and irreproducible across biomaRt API versions.
Symptom: Different mappings on rerun; conflicting downstream gene sets.
Fix: Either filter to mmusculus_homolog_orthology_type == 'ortholog_one2one' (conservative) or aggregate via highest homology confidence score (mmusculus_homolog_perc_id_r1).
Common errors
Error / symptom
Cause
Fix
filters returns empty
Singular filter= partial-matched against another argument
Spell filters= fully
1-Mar in gene column
Excel autocorrected MARCH1
Re-import with explicit string type; map back to MARCHF1
pyensembl ValueError: gene not found
ID not in pinned release; or unversioned ID against versioned database
Strip version before lookup; verify release
Duplicate rownames after aggregate
Collapsed multiple source IDs to one target; OR _PAR_Y stripped
Sum-collapse expected; for PAR_Y use preserving regex
biomaRt timeout for >5k IDs
Query too large
Chunk into batches of 1000
Wrong species mapping
Default species='human' in mygene; mouse query returns nothing
Pass species='mouse' explicitly
ENSEMBL keytype not available
Older org.db package or non-human/mouse
keytypes(orgdb) to verify
References
Ziemann M, Eren Y, El-Osta A. 2016. Gene name errors are widespread in the scientific literature. Genome Biol 17:177. doi:10.1186/s13059-016-1044-7
Bruford EA, Braschi B, Denny P, Jones TEM, Seal RL, Tweedie S. 2020. Guidelines for human gene nomenclature. Nat Genet 52:754-758. doi:10.1038/s41588-020-0669-3
Morales J, Pujar S, Loveland JE, et al. 2022. A joint NCBI and EMBL-EBI transcript set for clinical genomics and research. Nature 604:310-315. doi:10.1038/s41586-022-04558-8
Durinck S, Spellman PT, Birney E, Huber W. 2009. Mapping identifiers for the integration of genomic datasets with the R/Bioconductor package biomaRt. Nat Protoc 4(8):1184-1191. doi:10.1038/nprot.2009.97
Carlson M, Falcon S, Pages H, Li N. 2019. org.Hs.eg.db: Genome wide annotation for Human. R package. Bioconductor.
Wu C et al. 2013. BioGPS and MyGene.info: organizing online, gene-centric information. Nucleic Acids Res 41(D1):D561-D565. doi:10.1093/nar/gks1114
Frankish A et al. 2021. GENCODE 2021. Nucleic Acids Res 49(D1):D916-D923. doi:10.1093/nar/gkaa1087
Howe KL et al. 2021. Ensembl 2021. Nucleic Acids Res 49(D1):D884-D891. doi:10.1093/nar/gkaa942
Soneson C, Love MI, Robinson MD. 2015. Differential analyses for RNA-seq: transcript-level estimates improve gene-level inferences. F1000Res 4:1521. doi:10.12688/f1000research.7563.2
Related Skills
counts-ingest - Building count matrices and tx2gene
metadata-joins - Joining annotation with sample tables
normalization - Biotype filtering before normalization
sparse-handling - Single-cell row metadata in AnnData
differential-expression/de-results - Annotating DE results; gene-symbol display
