pythonforbio.
[WASM idle]
HGVS Nomenclature07/11

hgvs Normalizer

Starting sandbox…
Beginnerlesson

hgvs Normalizer

In one line: rewrites a variant into its canonical form so that equivalent variants become identical strings — the step that makes variant comparison possible.

If you take one thing from the hgvs notes, take this one.

The problem

Consider a deletion of one A from a run of four:

reference:  ...GC AAAA TG...
                 ↑↑↑↑
observed:   ...GC AAA  TG...

Which A was deleted? All four descriptions produce the identical result:

c.100del    c.101del    c.102del    c.103del

They are the same biological event with four different names. Without a canonical form, string comparison says they are four different variants — so your ClinVar lookup misses, your deduplication fails, and your cohort frequency is wrong by a factor of four.

The HGVS rule: 3'-most

HGVS requires the most 3' position on the reference sequence. So all four become c.103del.

import hgvs.normalizer

hn = hgvs.normalizer.Normalizer(hdp)
hn.normalize(var)

VCF requires the opposite — leftmost. The same event is POS=100 in a VCF and c.103del in HGVS. Both are correct under their own standard. Converting between them without normalizing at each step loses variants silently. See hgvs Pitfalls.

The direction depends on the reference sequence

"3'-most" is relative to the reference sequence being used, not the genome.

For a gene on the minus strand:

  • On the transcript (c.), 3' means increasing c. coordinates.
  • On the genome (g.), that is decreasing genomic coordinates.

So normalizing the same event as c. and as g. shifts it in opposite genomic directions. This is correct per the standard and endlessly confusing in practice. It is also why you must normalize after projecting, in the target coordinate system — not once at the start.

Configuration

hn = hgvs.normalizer.Normalizer(
    hdp,
    shuffle_direction=3,        # 3 = 3'-most (HGVS); 5 = 5'-most (VCF-like)
    cross_boundaries=False,     # allow shifting across exon boundaries?
    replace_reference=True,     # fill in ref bases from the sequence
    validate=True,
)

hn5 = hgvs.normalizer.Normalizer(hdp, shuffle_direction=5)   # for VCF interop

shuffle_direction=5 gives you left-aligned output, which is what you want when producing a VCF from HGVS input.

cross_boundaries=False is the safer default: it prevents a shift from moving a variant out of an exon and into an intron, which would change its predicted consequence.

What normalization does

  1. Shifts ambiguously-placed indels to the 3'-most equivalent position.
  2. Rewrites ins as dup where the inserted sequence duplicates the preceding bases — HGVS requires dup.
  3. Fills in reference bases (c.76_78delc.76_78delCTT) with replace_reference=True.
  4. Simplifies delins to a simpler type where possible (a delins that is really a substitution becomes one).
  5. Trims common prefixes and suffixes from delins alleles.
hn.normalize(hp.parse_hgvs_variant("NM_000546.6:c.76_77insG"))
# → NM_000546.6:c.76dup      (if base 76 is G)

Where it matters

Deduplication. Two labs reporting the same variant with different descriptions:

a = hn.normalize(hp.parse_hgvs_variant("NM_004006.2:c.100del"))
b = hn.normalize(hp.parse_hgvs_variant("NM_004006.2:c.103del"))
str(a) == str(b)          # True after normalization; False before

Database lookup. ClinVar, gnomAD and HGMD store normalized representations. An unnormalized query silently returns nothing, which looks like "this variant is novel" rather than "I asked the wrong question". This failure mode has real clinical consequences.

Cohort counting. Unnormalized variants split one allele across several rows, so every frequency is understated.

VCF ↔ HGVS. Round-tripping without renormalizing at each step drops variants.

In the projection pipeline

am = AssemblyMapper(hdp, assembly_name="GRCh38",
                    alt_aln_method="splign",
                    normalize=True,             # ← projections normalize their output
                    replace_reference=True)

AssemblyMapper with normalize=True (the default) normalizes after each projection, in the target coordinate system. That is the correct behaviour and the reason you should generally use AssemblyMapper rather than hand-rolling projections.

The recommended order for a raw input variant:

var = hp.parse_hgvs_variant(s)      # 1. parse
hv.validate(var)                     # 2. validate against the sequence
var = hn.normalize(var)              # 3. normalize in its own coordinate system
var_g = am.c_to_g(var)               # 4. project (normalizes the output too)

Error handling

from hgvs.exceptions import HGVSInvalidVariantError, HGVSDataNotAvailableError

try:
    var = hn.normalize(var)
except HGVSInvalidVariantError as e:
    ...        # ref base mismatch — wrong accession version or wrong assembly
except HGVSDataNotAvailableError as e:
    ...        # sequence not available; check SeqRepo

With validate=True, normalization checks the reference bases as it goes, so HGVSInvalidVariantError here usually means the variant does not match the sequence at all — a genuine data problem worth investigating rather than suppressing.

Normalization depends on your data version

The 3'-shift is computed from the actual sequence. Different SeqRepo snapshots can, in principle, produce different normalized output for a variant near a corrected region. That is why hgvs asks you to pin the library to a minor version and why you must record your SeqRepo and UTA snapshot dates.

print(f"hgvs {hgvs.__version__}, UTA {hdp.data_version()}, seqrepo {os.environ['HGVS_SEQREPO_DIR']}")

Put that line in every output file. See Applied - Reproducible Environment.

Batch normalization

def normalize_column(series, hp, hn):
    out, errs = [], []
    for s in series:
        try:
            out.append(str(hn.normalize(hp.parse_hgvs_variant(s))))
            errs.append(None)
        except HGVSError as e:
            out.append(None)
            errs.append(f"{type(e).__name__}: {e}")
    return pd.DataFrame({"normalized": out, "norm_error": errs}, index=series.index)

df = df.join(normalize_column(df["hgvs"], hp, hn))
changed = (df["hgvs"] != df["normalized"]) & df["normalized"].notna()
print(f"{changed.sum()} of {len(df)} variants changed on normalization")

That last count is diagnostic. A high proportion changing means your upstream source was not normalizing — worth knowing, and worth mentioning to whoever produced the file.

Common mistakes

  • Not normalizing at all. The root cause of most variant-matching failures.
  • Normalizing once, before projection, instead of in each target coordinate system.
  • Assuming HGVS and VCF shift the same way. They are opposite.
  • Forgetting that "3'" is reference-relative, so minus-strand genes shift the other way genomically.
  • Comparing variant strings without normalizing both.
  • Not recording the SeqRepo/UTA versions.
  • Suppressing HGVSInvalidVariantError — it usually indicates a real data problem.
  • cross_boundaries=True without understanding that a shift can move a variant out of an exon.

See also

hgvs Validator · AssemblyMapper · SeqRepo · HGVS Nomenclature · hgvs Pitfalls · Applied - Variant Normalization Pipeline

scratch

No output yet — run the code to populate this drawer.