get_MNV

Codon-level variant annotation

A VCF lists bases.
get_MNV reads codons.

Two substitutions in one codon make one amino acid, not two. get_MNV translates every codon whole, places each indel against its gene and, given a BAM or CRAM, counts the reads that carry the changes together.

In
a VCF, BCF or iVar table, a FASTA, a GFF or gene table, and optionally a BAM or CRAM
Out
a TSV, a VCF with ANN, and a report that opens offline
Figure 1

4 of 4 · The whole codon

example/G35894.var.snp.vcf

104941CT

104942CG

104944AG

Rv0095c, as the gene reads minus strand: positions count down

codon 91 · Val1 change

GCC Ala · one SNP row

codon 92 · Gly2 changes

104942 aloneCGGArg

104941 aloneGAGGlu

both, one moleculeCAGGln

codon 93 · Hisno change

  1. The VCF lists three substitutions, one record each, written on the reference's forward strand.
  2. Grouped by codon. The gene runs on the minus strand, so C>T at 104941 is G>A in the gene. 104941 and 104942 share codon 92; 104944 sits in codon 91 on its own.
  3. One base at a time, codon 92 reads Arg and Glu. That is what a base-by-base annotator reports: the input's own snpEff ANN says p.Gly92Arg and p.Gly92Glu.
  4. The whole codon reads Gln: the amino acid the protein carries when both changes sit on one molecule, which a BAM can confirm read by read. get_MNV writes it as one SNP/MNV row, p.Gly92Gln, with the two single readings beside it.
Rv0095c in the M. tuberculosis sample that ships with get_MNV. Every position, base and amino acid above comes from running it on example/G35894.var.snp.vcf.

Try it: one codon, two readings

Change any base. Each change is read on its own, the way a base-by-base annotator reads it, and then all together, the way get_MNV reads the codon.

The script that makes this interactive did not run. Here is the first example, worked out.

Read asCodonAmino acid
ReferenceCGTArg
Position 1 alone, C>TTGTCys missense
Position 2 alone, G>CCCTPro missense
Both, as get_MNV reads itTCTSer missense

One base at a time this codon gives Cys and Pro. Read whole, it gives Ser, which neither of them names.

Translated with the bacterial code (table 11), get_MNV's default; nine NCBI tables are supported. This is a codon inside a gene: the first codon of a gene is read as an initiator instead.

What it does with a VCF

Steps 1 to 3 need no reads. A BAM or CRAM only adds step 4, which turns "these two changes are in one codon" into "these two changes are on one molecule".

  1. 1

    Decompose

    Each record's REF and ALT into the changes it holds: SNVs, MNVs, insertions, deletions, delins.

  2. 2

    Place

    Each change in its gene and codon, from a GFF or a gene table, on either strand and across exons.

  3. 3

    Translate

    The whole codon, in one of nine NCBI genetic codes, and each indel against the protein.

  4. 4

    Count

    With reads: how many carry each change, on which strand, and whether the changes travel together.

  5. 5

    Write

    A TSV, a VCF with ANN, BCF on request, and a report that opens offline.

How get_MNV tells one molecule from two · Indels and local haplotypes

Install

One binary, no C dependencies, on macOS, Linux and Windows. A desktop app does the same with drag and drop.

Bioconda

conda install -c bioconda get_mnv

From source

git clone https://github.com/PathoGenOmics-Lab/get_MNV.git
cd get_MNV
cargo install --path .

Rust 1.91 or newer.

Desktop app

macOS, Linux and Windows builds on the Releases page, with a guided tutorial.

Not signed with an Apple certificate: on the first launch, right-click the app and choose Open.

Check what you downloaded

Every release carries a SHA256SUMS file covering all of its assets. Download it next to the file you want and check them together; each file it finds prints OK.

sha256sum --ignore-missing -c SHA256SUMS

On macOS, shasum -a 256 --ignore-missing -c SHA256SUMS does the same.

Your first run

The repository ships a M. tuberculosis sample: a VCF, its reference and a gene table. From a clone:

get_mnv \
  --vcf example/G35894.var.snp.vcf \
  --fasta example/MTB_ancestor.fas \
  --genes example/anot_genes.txt \
  --both --output-dir out \
  --report report.html
  • 950 records in, 941 rows out: 797 SNP, 10 SNP/MNV codons and 134 intergenic
  • out/G35894.var.snp.MNV.tsv and out/G35894.var.snp.MNV.vcf
  • report.html, the example report

The same run, line by line · Common recipes · Every option

What you get

Codon 92 of Rv0095c from Figure 1, in each of the three outputs.

TSV one row per codon or indel

Positions
104941, 104942
AA Changes
Gly92Gln
SNP AA Changes
Gly92Glu, Gly92Arg
Variant Type
SNP/MNV
HGVS c.
c.[274G>C;275G>A]
DBS Class
CC>TG

VCF each base, and the codon

104941  C   T   …AA=Gly92Glu…
104942  C   G   …AA=Gly92Arg…
104941  CC  TG  …AA=Gly92Gln;…
  ANN=TG|missense_variant|MODERATE|
  Rv0095c_Rv0095c|…|p.Gly92Gln|
  …|INFO_COMPOUND_ANNOTATION

ANN is in snpEff's format, so bcftools +split-vep reads it.

Report one HTML file, no server

The example report: summary counts, consequences, mutational spectra and the variant table The example report: summary counts, consequences, mutational spectra and the variant table Open the example report

Every column and INFO key