Skip to content

How to Convert Genomic Coordinates Between Assemblies

Oak
An iridescent winged reptile specimen suspended on a chain of glowing glass links between two brass cradles against a black backdrop.

By the end of this you will have a reproducible command-line path for moving coordinates between GRCh37/hg19, GRCh38/hg38, and T2T-CHM13v2.0, for the four file shapes that matter: VCF, BED/interval, per-position lookups, and transcript (HGVS) coordinates. You will also have a rejection log you can read, because the useful output of any lift is the list of things that did not lift and the reason. You need: bcftools ≥1.20 with the liftover plugin, UCSC liftOver, CrossMap (Python), samtools, roughly 30 GB of disk for reference FASTAs and their indices, and the chain files. Everything below assumes a POSIX shell.

1. Establish which assembly and which convention your file is already in

Most bad lifts start here. Three things vary independently:

  • Assembly. GRCh37 vs GRCh38 vs CHM13v2.0.
  • Contig naming. UCSC style (chr1, chrM) vs Ensembl/NCBI style (1, MT) vs RefSeq accessions (NC_000001.11).
  • Coordinate convention. BED is 0-based half-open. VCF, GFF/GTF, SAM, and HGVS g. are 1-based inclusive. A BED line chr1 100 200 and a VCF record at chr1:101 describe the same first base.

The fastest assembly check is contig length. chr1 is 249,250,621 in GRCh37, 248,956,422 in GRCh38, and 248,387,328 in CHM13v2.0.

bcftools view -h in.vcf.gz | grep '^##contig' | head -3
samtools faidx ref.fa && head -3 ref.fa.fai

If the VCF header has no ##contig lines and no ##reference, do not guess from the filename. Take a handful of common variants and check the REF alleles against each candidate FASTA (step 9 shows how). hg19 and GRCh37 are the same nuclear sequence with different names, but the mitochondrion differs: hg19 chrM is the older 16,571 bp NC_001807 sequence, GRCh37 MT is the 16,569 bp rCRS. Any tool that lifts chrM between those two without special handling will produce garbage offsets past position ~3,100.

2. Fetch chain files and references

mkdir -p chains refs && cd chains
wget https://hgdownload.soe.ucsc.edu/goldenPath/hg19/liftOver/hg19ToHg38.over.chain.gz
wget https://hgdownload.soe.ucsc.edu/goldenPath/hg38/liftOver/hg38ToHg19.over.chain.gz
wget https://hgdownload.soe.ucsc.edu/goldenPath/hg38/liftOver/hg38ToHs1.over.chain.gz

hs1 is UCSC’s name for T2T-CHM13v2.0. Chains are directional: hg19ToHg38 cannot be run backwards. The Ensembl equivalents (GRCh37_to_GRCh38.chain.gz) use unprefixed contig names, so pick the chain that matches your file’s naming rather than renaming contigs twice.

Get the target reference and index it. You need the actual FASTA for VCF lifting, because the lift has to re-derive REF alleles.

cd ../refs
samtools faidx GRCh38_full_analysis_set_plus_decoy_hla.fa

3. Lift a VCF with bcftools +liftover

This is the tool we use. It is the only common one that handles REF/ALT swaps, strand flips, and the resulting need to update AF, AC, and genotype phase in a single pass.

export BCFTOOLS_PLUGINS=/path/to/bcftools/plugins

bcftools +liftover --no-version -Ou in.GRCh37.vcf.gz -- \
  -s refs/human_g1k_v37.fasta \
  -f refs/GRCh38_full_analysis_set_plus_decoy_hla.fa \
  -c chains/hg19ToHg38.over.chain.gz \
  --reject rejected.GRCh37.vcf.gz --reject-type z \
  --write-src \
| bcftools sort -Oz -o out.GRCh38.vcf.gz -T /tmp/bcfsort
bcftools index -t out.GRCh38.vcf.gz

Points that matter:

  • -s (source FASTA) is not optional in practice. Without it the plugin cannot distinguish a genuine REF/ALT swap from a mismatch, and you lose records.
  • --write-src copies the original CHROM/POS/REF/ALT into INFO fields. Keep it. When a downstream result looks wrong, you want the pre-lift coordinate on the same line.
  • Output is unsorted by construction, because the chain reorders segments. Piping into bcftools sort is mandatory, not stylistic.
  • Records where REF and ALT trade places get a SWAP tag, reverse-strand lifts get FLIP. Filter on these when you compare allele frequencies against a population reference.

Picard’s LiftoverVcf is the alternative:

java -jar picard.jar LiftoverVcf I=in.vcf.gz O=out.vcf.gz \
  CHAIN=chains/hg19ToHg38.over.chain.gz \
  REJECT=rejected.vcf.gz R=refs/GRCh38.fa \
  RECOVER_SWAPPED_REF_ALT=true WARN_ON_MISSING_CONTIG=true

It works and it is well tested, but it is slower on whole-genome VCFs and it rejects more indels near chain boundaries. Use it if you are already inside a GATK pipeline and want matching sequence dictionaries.

Do not lift a gVCF. Reference blocks (END=) span chain breaks and the result is incoherent. Re-call from the alignment instead.

4. Intervals, BED, GTF, BAM, bigWig

For plain intervals, UCSC liftOver is fine and fast:

liftOver -minMatch=0.95 -bedPlus=6 -tab \
  regions.hg19.bed chains/hg19ToHg38.over.chain.gz \
  regions.hg38.bed regions.unmapped.bed

-minMatch defaults to 0.95: at least 95% of bases in the interval must map. For long intervals (whole genes, CNV calls, 100 kb+ windows) that default silently drops a lot. Drop it to 0.8 and then inspect what changed, rather than accepting either number blindly. -multiple returns all mapping locations for regions that map more than once, which is what you want in segmental duplications. -minBlocks controls how fragmented a mapping may be.

For anything with structure inside the record, use CrossMap, which understands file semantics:

CrossMap bed   chains/hg19ToHg38.over.chain.gz in.bed out.bed
CrossMap gff   chains/hg19ToHg38.over.chain.gz in.gtf out.gtf
CrossMap bam   chains/hg19ToHg38.over.chain.gz in.bam out.bam
CrossMap bigwig chains/hg19ToHg38.over.chain.gz in.bw out

CrossMap bam rewrites alignment coordinates without re-aligning the reads. It is useful for a quick browse. It is not a substitute for re-alignment if you plan to call variants, because the CIGAR strings still reflect alignment to the old reference and any base that changed between assemblies will produce a spurious mismatch.

A GTF lift needs one extra check afterwards: exons of a transcript can land on different contigs or in a different order. Group by transcript_id and assert monotonic coordinates and a single chromosome before you use the file.

5. Single positions, and why an rsID beats a coordinate

If you are carrying a handful of variants between assemblies, do not lift them. Look them up by rsID. The dbSNP RefSNP identifier is assembly-independent by design, and dbSNP remaps every RefSNP onto each new assembly using the same alignments it curates for the whole map1. That gives you the assembly’s own answer rather than a chain-file approximation.

curl -s "https://api.ncbi.nlm.nih.gov/variation/v0/refsnp/334" \
| jq '.primary_snapshot_data.placements_with_allele[]
      | select(.seq_id | startswith("NC_"))
      | {seq_id, pos: .alleles[0].allele.spdi.position}'

Or via Ensembl, which will give you both assemblies:

curl -s "https://rest.ensembl.org/variation/human/rs334?content-type=application/json" | jq '.mappings[]'
curl -s "https://grch37.rest.ensembl.org/variation/human/rs334?content-type=application/json" | jq '.mappings[]'

Caveats worth holding in mind: RefSNP IDs merge (an old ID redirects to a newer one), some map to multiple locations, and a small number are withdrawn. Always record both the rsID and the assembly-qualified coordinate in your notes, so a future merge does not silently repoint your record.

6. Transcript coordinates (HGVS c. and p.)

c. coordinates are defined relative to a transcript, not a genome. NM_000546.6:c.215C>G means the same base whichever assembly you are in, as long as the transcript version is pinned. The conversion you need is projection, not liftover: genome → transcript → genome.

import hgvs.parser, hgvs.dataproviders.uta, hgvs.assemblymapper

hp  = hgvs.parser.Parser()
hdp = hgvs.dataproviders.uta.connect()
am37 = hgvs.assemblymapper.AssemblyMapper(hdp, assembly_name="GRCh37", alt_aln_method="splign")
am38 = hgvs.assemblymapper.AssemblyMapper(hdp, assembly_name="GRCh38", alt_aln_method="splign")

v37 = hp.parse_hgvs_variant("NC_000017.10:g.7578475C>G")
c   = am37.g_to_c(v37, "NM_000546.6")
v38 = am38.c_to_g(c)
print(str(c), str(v38))

This routes through curated transcript alignments, so it fails loudly when a transcript has no alignment on the target assembly instead of producing a plausible wrong number. VariantValidator and Mutalyzer expose the same idea over HTTP. Mutalyzer’s Position Converter does not normalize or validate the description it is given, so feed it already-normalized input.

Version matters. NM_000546.5 and NM_000546.6 can differ in UTR length, which shifts every c. coordinate in the 5’ UTR. Pin the version in your records.

7. Lifting into T2T-CHM13

hg38ToHs1.over.chain.gz works for the ~95% of the genome where both assemblies agree. It cannot help you in the regions where CHM13 adds sequence that GRCh38 does not have: centromeric higher-order repeats, the acrocentric short arms, several segmental duplication families, and the rDNA arrays. There is no source coordinate to lift from.

If your interest is in those regions, lift nothing and realign. The complete diploid benchmark built on the CHM13 and HG002 assemblies covers substantially more of the genome than GRCh38-based benchmarks did, including regions that were previously unassayable by short-read calling, which is exactly the set you cannot reach by chain-file projection2. Treat a CHM13 analysis as a separate alignment run with its own variant calls, and use the lift only to cross-reference annotations.

8. Array and probe manifests

Genotyping and methylation manifests ship with coordinates baked in, usually hg19 for older products. Lift the probe coordinate with liftOver on a BED built from the manifest, then re-derive strand and allele orientation from the probe sequence rather than trusting the manifest’s Strand column across the lift.

For methylation arrays there is a second, unrelated problem that the lift will not fix: a meaningful fraction of EPIC probes cross-hybridize to multiple genomic locations, and others overlap common SNPs, so their signal does not report the CpG their coordinate names3. Filter those probes using a published exclusion list before any positional analysis, whatever assembly you end up in. A correct coordinate on a bad probe is still a bad measurement.

9. Validate the lift before you trust it

Four checks, in order of how often they catch something:

# 1. REF allele agreement with the new reference
bcftools norm -f refs/GRCh38.fa -c w out.GRCh38.vcf.gz -Ou 2> ref_check.log >/dev/null
grep -c 'REF_MISMATCH' ref_check.log

# 2. How much was rejected, and why
bcftools query -f '%INFO/FAILED_LIFTOVER\n' rejected.GRCh37.vcf.gz | sort | uniq -c

# 3. Counts in and out
bcftools index -n in.GRCh37.vcf.gz; bcftools index -n out.GRCh38.vcf.gz

# 4. Sanity on a known site
bcftools view -r chr11:5227002 out.GRCh38.vcf.gz | grep -v '^##'

A whole-genome SNV+indel VCF lifted GRCh37→GRCh38 typically loses on the order of 0.1–0.5% of records, concentrated in segmental duplications, the MHC, and near assembly patches. If you are losing several percent, something upstream is wrong: chain direction, contig naming, or a missing source FASTA.

One more check that catches subtle errors: compare your post-lift allele frequencies against gnomAD on the target assembly. Systematic REF/ALT confusion shows up as a cloud of variants with AF near 1 − expected.

10. When to skip lifting entirely

If the data will underpin analysis you revisit for years, realign rather than lift. From CRAM or FASTQ:

samtools fastq -@8 -1 r1.fq.gz -2 r2.fq.gz -0 /dev/null -s /dev/null -n in.cram
bwa-mem2 mem -t 32 -R '@RG\tID:1\tSM:sample\tPL:ILLUMINA' refs/GRCh38.fa r1.fq.gz r2.fq.gz \
| samtools sort -@8 -o out.GRCh38.bam

A 30x human genome realigns in a few hours on 32 cores. That is cheap next to the cost of tracking down a result that turned out to hinge on a chain-file artifact in a segmental duplication. Lift when you need a quick cross-reference or when you only have a summary-level file. Realign when the coordinates are load-bearing.

Common problems

Chain and file contig names disagree. liftOver reports everything as unmapped and exits 0. Check the first column of your input against the chain’s chr fields before blaming the data.

Output VCF is unsorted. Symptoms: bcftools index fails, or tabix queries return partial results. Always pipe through bcftools sort.

Reversed chain. hg38ToHg19 applied to a GRCh37 file will map some fraction of records to plausible-looking wrong positions, because the assemblies share long stretches of identical sequence. The REF-mismatch count in step 9 is what catches this.

REF/ALT swap not handled. Roughly 0.5% of sites differ in which allele GRCh38 calls reference. If your tool rejects rather than swaps, you lose real variants. Use --write-src and grep the SWAP tag to see how many were affected.

Ambiguous A/T and G/C variants on arrays. Strand flips cannot be resolved from the alleles alone. Resolve using flanking sequence or allele frequency against a reference panel, not by assumption.

chrM. Lift it separately or not at all, and check whether you are in rCRS or NC_001807 coordinates first.

PAR regions on chrX/chrY. GRCh37 and GRCh38 differ in how PAR sequence is represented and masked. Variants there may lift to chrX, chrY, both, or neither, depending on the reference build you used. Decide on a convention and apply it explicitly.

Interpretation. Coordinates and variant calls are measurements. Any clinical meaning attached to a variant, including whether a result warrants follow-up, belongs with a clinician or a certified diagnostic laboratory.

Oak builds longitudinal molecular profiles of individuals: whole-genome sequencing, RNA sequencing, proteomics, blood biomarkers, and continuous glucose data, integrated into one model of you. Build your profile.

Footnotes

  1. The International SNP Map Working Group, Cold Spring Harbor Laboratories:, Ravi Sachidanandam, et al. A map of human genome sequence variation containing 1.42 million single nucleotide polymorphisms. Nature, 2001. https://doi.org/10.1038/35057149 ↩

  2. Nancy F. Hansen, Nathan Dwarshuis, Hyun Joo Ji, et al. A complete diploid human genome benchmark for personalized genomics. Cell, 2026. https://doi.org/10.1016/j.cell.2026.06.016 ↩

  3. Ruth Pidsley, Elena Zotenko, Timothy J. Peters, et al. Critical evaluation of the Illumina MethylationEPIC BeadChip microarray for whole-genome DNA methylation profiling. Genome biology, 2016. https://doi.org/10.1186/s13059-016-1066-1 ↩