AssemblyMapper
In one line: the high-level projector — moves a variant between genomic, transcript and protein coordinates for a chosen assembly, and normalizes as it goes.
Construction
import hgvs.assemblymapper, hgvs.dataproviders.uta
hdp = hgvs.dataproviders.uta.connect(pooling=True)
am = hgvs.assemblymapper.AssemblyMapper(
hdp,
assembly_name="GRCh38", # or "GRCh37"
alt_aln_method="splign", # for RefSeq transcripts
normalize=True, # normalize output (default True)
replace_reference=True, # fill in ref bases from the sequence
prevalidation_level="INTRINSIC", # None | INTRINSIC | EXTRINSIC
in_par_assume="X", # pseudoautosomal region → assume X
)
Build it once. Each construction opens database connections and caches transcript data. Constructing it per variant is the most common performance mistake in hgvs code.
The parameters worth understanding:
replace_reference=Truefills in the reference bases from the actual sequence. This turnsc.76_78delintoc.76_78delCTTand — more importantly — surfaces mismatches between what the variant claims and what the sequence says.normalize=Trueapplies 3'-shifting to the output. Almost always what you want. See hgvs Normalizer.in_par_assumeresolves the pseudoautosomal regions, where a locus exists on both X and Y with identical sequence. Without it, PAR variants raise.
The projection methods
am.c_to_g(var_c) # transcript coding → genomic
am.g_to_c(var_g, tx_ac) # genomic → a SPECIFIC transcript
am.c_to_p(var_c) # coding → protein (a PREDICTION)
am.c_to_n(var_c); am.n_to_c(var_n)
am.g_to_n(var_g, tx_ac); am.n_to_g(var_n)
am.g_to_t(var_g, tx_ac) # dispatches to g_to_c or g_to_n
am.t_to_g(var_t)
am.t_to_p(var_t)
am.relevant_transcripts(var_g) # which transcripts overlap this genomic variant?
am.assembly_name; am.data_version()
g_to_c requires a transcript accession — that is the point. A genomic variant has one position; its coding consequence depends entirely on which transcript you ask about, and a gene can have twenty. There is no "the" transcript.
The transcript selection problem
var_g = hp.parse_hgvs_variant("NC_000017.11:g.7676154G>C")
am.relevant_transcripts(var_g)
# ['NM_000546.6', 'NM_001126112.3', 'NM_001126113.3', ...]
Which one do you report? This is a policy decision, not a technical one, and it must be made explicitly:
| Strategy | Rationale |
|---|---|
| MANE Select | current best practice; one agreed transcript per gene, RefSeq+Ensembl matched |
| RefSeq Select | NCBI's representative transcript |
| Longest CDS | crude but deterministic |
| Clinically established | some genes have a legacy transcript everyone uses (e.g. specific BRCA1 and CFTR transcripts) |
| All of them | most complete; report every consequence |
MANE Select is the right default in 2026. Where a clinical convention exists for a gene, it overrides MANE — because a variant reported against a different transcript than the literature uses will not match anything.
MANE = pd.read_csv("MANE.GRCh38.v1.4.summary.txt.gz", sep="\t")
MANE_BY_GENE = dict(zip(MANE["symbol"], MANE["RefSeq_nuc"]))
def preferred_tx(gene, available):
mane = MANE_BY_GENE.get(gene)
if mane in available:
return mane
return sorted(available)[0] if available else None
Record which transcript you used in every output row. A c. variant without its accession is not interpretable.
GRCh37 vs GRCh38
am37 = AssemblyMapper(hdp, assembly_name="GRCh37", alt_aln_method="splign")
am38 = AssemblyMapper(hdp, assembly_name="GRCh38", alt_aln_method="splign")
am37.c_to_g(var_c) # NC_000017.10:g.7579472G>C
am38.c_to_g(var_c) # NC_000017.11:g.7676154G>C
Same variant, ~97,000 bp apart. A genomic coordinate without an assembly is meaningless, and mixing assemblies within one analysis is a silent, serious error.
Note the accession itself encodes the assembly: NC_000017.10 is GRCh37 chr17, NC_000017.11 is GRCh38 chr17. That versioning is a genuine safety feature — unlike chr17, which tells you nothing.
Projecting through a transcript is not liftover. Going g.(GRCh37) → c. → g.(GRCh38) works only where the transcript aligns to both assemblies consistently. For genuine coordinate conversion, use CrossMap or UCSC liftOver with the official chain files.
Protein projection is a prediction
var_p = am.c_to_p(var_c)
str(var_p) # 'NP_000537.3:p.(Pro72Arg)'
The parentheses are HGVS's marker for a predicted consequence. c_to_p applies the reference codon table to the reference sequence; it does not know about:
- nonsense-mediated decay
- alternative splicing that skips the exon
- the actual translated product
- readthrough or reinitiation
For a frameshift it will predict a new stop position, which is a mechanical calculation, not an observation. Treat these as annotations, not facts, and preserve the parentheses in your reports.
Error handling
from hgvs.exceptions import (
HGVSError, HGVSDataNotAvailableError, HGVSInvalidVariantError,
HGVSInvalidIntervalError, HGVSUsageError, HGVSUnsupportedOperationError,
)
def project(am, var_c):
try:
return am.c_to_g(var_c), None
except HGVSDataNotAvailableError as e:
return None, f"no transcript data: {e}" # UTA lacks this accession
except HGVSInvalidIntervalError as e:
return None, f"position out of range: {e}" # position beyond the transcript
except HGVSInvalidVariantError as e:
return None, f"invalid variant: {e}" # ref base mismatch
except HGVSError as e:
return None, f"hgvs error: {e}"
Catch the specific exception types — each means a different data problem, and lumping them together loses the diagnostic information you need.
HGVSUnsupportedOperationError from c_to_p typically means a variant type the protein predictor cannot handle (some complex indels, variants spanning the start codon).
Batch pattern
def project_all(df, am, hp):
rows = []
for r in df.itertuples():
row = {"input": r.hgvs, "g": None, "p": None, "error": None}
try:
v = hp.parse_hgvs_variant(r.hgvs)
row["g"] = str(am.c_to_g(v))
row["p"] = str(am.c_to_p(v))
except HGVSError as e:
row["error"] = f"{type(e).__name__}: {e}"
rows.append(row)
out = pd.DataFrame(rows)
print(out["error"].notna().value_counts())
return out
Never let one bad variant kill a batch of 50,000. Collect the errors, report the counts, and inspect the patterns — they usually cluster around one systematic cause.
AssemblyMapper vs VariantMapper
import hgvs.variantmapper
vm = hgvs.variantmapper.VariantMapper(hdp)
vm.c_to_g(var_c, "NC_000017.11") # you supply the genomic accession
VariantMapper is the lower-level engine; AssemblyMapper wraps it with assembly awareness, transcript filtering and normalization. Use AssemblyMapper unless you need the extra control.
Common mistakes
- Constructing
AssemblyMapperper variant. - Not recording which transcript was used.
- Mixing GRCh37 and GRCh38.
- Treating
p.(...)as observed. - Using transcript projection as liftover.
- Catching bare
Exception, losing the diagnostic type. - Assuming one transcript per gene.
- Forgetting
replace_reference=True, so ref-base mismatches go undetected. - Ignoring PAR regions on X/Y.
See also
Projection c to g to p · UTA and Data Providers · hgvs Normalizer · hgvs Validator · hgvs Pitfalls · Applied - Annotating a Variant Table