Skip to content

Getting started with MTBC

This end-to-end tutorial takes you from a fresh install to a lineage call for the Mycobacterium tuberculosis complex (MTBC), using the published pathotypr marker panel and pre-trained model. You will classify an assembled genome, genotype raw reads, run the machine-learning predictor, and read every output file.

What you'll do

  1. Install pathotypr
  2. Get the MTBC panel, model, and reference
  3. Classify an assembly with pathotypr classify
  4. Genotype raw reads with pathotypr split-fastq
  5. Predict a lineage with pathotypr predict
  6. Read the output

The MTBC panel resolves lineages L1–L10 and animal-adapted A1–A4. Two complementary routes are shown: marker-based calling (classify / split-fastq, exact known SNPs) and model-based prediction (predict, a Random Forest over k-mer features). Run either or both.

Your data Command Route
Assembled genome (FASTA) classify Marker-based
Raw reads (FASTQ) split-fastq Marker-based
Assembled genome (FASTA) predict Model-based

1. Install pathotypr

Pick any install route from the Installation guide — the desktop app, Bioconda, or a build from source. The fastest for the CLI is Bioconda:

conda create -n pathotypr -c bioconda pathotypr
conda activate pathotypr

Confirm the binary is on your PATH:

pathotypr --version
pathotypr --help

Verbose logging

Every subcommand accepts -v (debug) or -vv (trace). Add it to any command below if you need to see what pathotypr is doing under the hood.

2. Get the MTBC panel, model, and reference

The marker panel and pre-trained model live on Zenodo under the concept DOI 10.5281/zenodo.19210043, which always resolves to the newest version. Open it and download the files, or fetch them from the command line as below; the desktop app resolves the same deposit on its own. Filenames carry the catalogue version and therefore change between releases, so the commands read them from the deposit instead of assuming one. The marker positions are numbered in H37Rv coordinates (NC_000962.3), but the reference sequence you pass to -r must be the MTBC ancestor genome the panel was built against — the same file the desktop app downloads.

# Newest marker deposit: lineage panel, drug-resistance panel and pre-trained model
curl -sL "https://zenodo.org/api/records/19210043/versions/latest" \
  | jq -r '.files[] | "\(.links.self) \(.key)"' \
  | while read -r url name; do curl -L -o "$name" "$url"; done

# Reference the panel was built against (single-record, separate deposit)
curl -L -o MTB_ancestor_reference.fasta \
  "https://zenodo.org/records/3497110/files/MTB_ancestor_reference.fasta?download=1"

# Newest marker deposit: lineage panel, drug-resistance panel and pre-trained model
wget -qO- "https://zenodo.org/api/records/19210043/versions/latest" \
  | jq -r '.files[] | "\(.links.self) \(.key)"' \
  | while read -r url name; do wget -O "$name" "$url"; done

# Reference the panel was built against (single-record, separate deposit)
wget -O MTB_ancestor_reference.fasta \
  "https://zenodo.org/records/3497110/files/MTB_ancestor_reference.fasta?download=1"

Use the ancestor reference, not H37Rv

classify and split-fastq splice each marker's allele into flanking bases taken from whatever you pass to -r, so the reference must be the one the panel was built against. That is MTB_ancestor_reference.fasta — the file pathotypr itself downloads (see defaults.rs). It is single-record and shares H37Rv's coordinates, but carries ancestral alleles, so substituting the H37Rv sequence silently produces different diagnostic k-mers and markers may stop matching. The pre-trained model is self-contained and does not need a reference.

After this step your working directory holds the deposit's files plus MTB_ancestor_reference.fasta. You'll also need your own sample(s): an assembled genome (sample.fasta) and, for the reads route, a FASTQ pair (sample_R1.fastq.gz, sample_R2.fastq.gz). See Input formats for the exact specs.

3. Classify an assembly

Call the known lineage SNPs in an assembled genome. --nested-classification turns on hierarchical calling so the multi-level marker columns resolve sub-lineages, not just the top level.

pathotypr classify \
  -m pathotypr_lineage_markers_*.tsv \
  -r MTB_ancestor_reference.fasta \
  -i sample.fasta \
  -o sample_lineage \
  --nested-classification \
  --excel

# Same command, swap in the DR panel to scan for resistance mutations
pathotypr classify \
  -m pathotypr_dr_markers_ancestor_*.tsv \
  -r MTB_ancestor_reference.fasta \
  -i sample.fasta \
  -o sample_dr \
  --excel

This writes sample_lineage.tsv (one row per marker hit), sample_lineage_summary.tsv (columns genome, lineage:count, major_lineage) — one row per genome; with -i each FASTA record counts as its own genome, and — because of --excel — the matching .xlsx files.

Gene and amino-acid change columns

Pass a GFF3 annotation with --gff (for -i) to populate the Gene, Gene_Start, Gene_End, AA_Pos, and AA_Change columns — especially useful with the DR panel to name the affected gene and mutation. Without a GFF, Gene_Start, Gene_End and AA_Pos stay empty; Gene and AA_Change fall back to the marker TSV's own gene and mutation columns when the panel provides them. See classify for GFF handling.

Many genomes at once

Swap -i sample.fasta for -l genomes.tsv (a sample_name<TAB>fasta_path[<TAB>gff_path] list) or --input-files a.fasta b.fasta … to run a whole batch under one output prefix.

4. Genotype raw reads

No assembly? Genotype straight from FASTQ. split-fastq counts REF vs ALT diagnostic k-mers at each marker — no alignment step. Paired files named _R1/_R2 (or _1/_2) are auto-detected; --paired makes the pairing explicit.

pathotypr split-fastq \
  -m pathotypr_lineage_markers_*.tsv \
  -r MTB_ancestor_reference.fasta \
  -i sample_R1.fastq.gz -i sample_R2.fastq.gz \
  --paired \
  -o sample_genotype \
  --nested-classification \
  --excel

Outputs are sample_genotype_<sample>_mutations.tsv (per-marker detail: position, alleles, counts, ALT fraction, lineage) and sample_genotype_summary.tsv (lineage counts plus the major call).

Tuning calls for coverage

A variant is called when depth ≥ --min-depth (default 10) and the ALT fraction ≥ --min-alt-percent (default 95). Lower --min-depth for shallow sequencing, or lower --min-alt-percent if you expect mixed / heteroresistant samples.

Indels are skipped from reads

Short reads give unreliable k-mer matches across MTBC repeats (PE/PPE, IS elements), so split-fastq calls SNPs and MNVs only. For indel-bearing markers, use the assembly route in Step 3.

5. Predict a lineage with the model

predict is the model-based alternative: it vectorizes the assembly with the same feature hasher used at training time and lets the Random Forest vote. It needs only your FASTA and the downloaded model — no markers, no reference.

pathotypr predict \
  -i sample.fasta \
  -m pathotypr_rf_model_*.pathotypr \
  -o sample_prediction.tsv \
  --excel

This writes sample_prediction.tsv (plus sample_prediction.xlsx) with one row per input sequence.

classify or predict?

classify reports exactly which known markers were hit — auditable, and it also flags resistance mutations with the DR panel. predict needs no reference or markers and returns a lineage with a calibrated confidence. Running both is a good cross-check.

6. Read the output

classify — detailed TSV

sample_lineage.tsv has one row per marker match, with these columns:

Column Meaning
genome Sample / contig the hit came from
k-mer Diagnostic k-mer that matched
k-merPOS Position of the match within the genome
SNPgenome 1-based position of the variant allele in the query
SNPreference 1-based marker position in the reference genome
REF / ALT Marker reference and alternate alleles
lineage Lineage label carried by the marker
Gene, Gene_Start, Gene_End Gene context (populated only with --gff)
AA_Pos, AA_Change Amino-acid position and change (only with --gff)

The companion sample_lineage_summary.tsv collapses this to a single row: the genome and its major lineage call (with sub-lineage when --nested-classification is on).

split-fastq — mutations + summary

sample_genotype_<sample>_mutations.tsv lists each called marker with its position, REF/ALT alleles, REF and ALT read counts, ALT fraction, and lineage. sample_genotype_summary.tsv gives one row per sample: the per-lineage marker counts and the final major-lineage call.

predict — predictions TSV

sample_prediction.tsv has one row per input sequence — except sequences shorter than the model's k-mer size, which are skipped with a warning:

Column Meaning
Header Sequence header from the input FASTA
Predicted_Lineage Most-voted lineage label
Confidence Fraction of trees voting for the winner (0–1)
Confidence_Margin Gap between winner and runner-up (0–1)
Other_Votes Top 3 alternatives with their vote fractions

Reading a confidence score

A Confidence of 0.95 means 95% of the trees agreed on the call; a Confidence_Margin of 0.40 means the winner drew 40% of the whole ensemble more than the runner-up (e.g. 60% of trees vs 20%). High values on both signal a clean, unambiguous prediction.

Excel everywhere

Every command above used --excel, so each .tsv has a side-by-side .xlsx. Drop the flag for TSV-only output.

Next steps

  •   classify

    Full options, GFF annotation, and the lineage-calling logic.

  •   split-fastq

    Read-based genotyping, depth/frequency thresholds, and pairing modes.

  •   predict

    Model-based prediction, confidence metrics, and output details.

  •   train

    Build your own Random Forest model from labeled genomes.