Skip to content

jlanej/nonhuman-screen

Repository files navigation

nonhuman-screen

PyPI Python versions CI License: MIT

Classify sequencing reads by non-human taxonomic content using kraken2, and reduce the per-read verdicts to a non-human fraction (NHF) plus per-domain breakdowns (bacterial / archaeal / fungal / protist / viral / UniVec-Core).

The headline capability is the allele-based NHF: given a BAM/CRAM and a variant, compute the non-human fraction of the reads that support that variant's ALT allele — in one call, no boilerplate:

from nonhuman_screen import classify_variant_alt_reads

v = classify_variant_alt_reads(
    "sample.bam", "kraken2_db", "chr1", 12345, "A", "T", ref_fasta="ref.fa",
)
v.nonhuman_fraction   # 0.0–1.0 over ALT-supporting reads
v.fractions.bacterial # per-domain breakdown
v.supporting_reads    # denominator

nonhuman-screen is sample-agnostic: it classifies whatever reads you hand it and never asks whose they are. The same call works on a proband, a parent, a tumour, or a plain FASTQ — so it can back higher-level analyses (e.g. flagging inherited variant calls that are actually parental contamination) without knowing anything about pedigrees.

Extracted from kmer-denovo-filter, where this methodology was originally developed, so it can be reused independently.

How it works (in one paragraph)

Reads are written to a temporary FASTQ and classified by kraken2's LCA algorithm. Each read's assigned taxid is mapped to a domain using the database's nodes.dmp taxonomy. A read counts as non-human only if it is classified and lies outside the human lineage, outside the human clade, and outside UniVec-Core. Two conservative guards reduce false positives: a human-homology guard (any k-mer voting for human, taxid 9606, drops the read from every non-human numerator — important for integrating viruses such as HBV/HPV/ERVs) and UniVec-Core exclusion (synthetic vector/adapter sequences, taxid 81077, are tracked separately and never counted as contamination). Both mechanisms are defined relative to human, so the kraken2 database must contain the human genome — PrackenDB does. See docs/methodology.md.

A failed kraken2 run raises Kraken2Error (CLI: exit 3) instead of returning nonhuman_fraction = 0.0: for a contamination screen, a broken run that looks clean is the one result you can't afford. See methodology.md §7.

Install

pip install nonhuman-screen           # core engine (standard library only)
pip install 'nonhuman-screen[bam]'    # + BAM/CRAM & allele-based NHF (pysam)

The bam extra is required for the BAM/allele helpers and for the nonhuman-screen CLI (both classify modes import pysam). The core install covers only Kraken2Runner.classify_sequences on in-memory reads.

You also need the kraken2 binary on PATH and a kraken2 database — see docs/database.md.

Usage

Allele-based NHF for many variants (batched — one kraken2 call)

from nonhuman_screen import classify_variants_alt_reads

variants = [("chr1", 12345, "A", "T"), ("chr2", 999, "G", "GAC")]  # pos 0-based
for v in classify_variants_alt_reads("sample.bam", "kraken2_db", variants,
                                      ref_fasta="ref.fa"):
    print(v.variant_key, v.nonhuman_fraction, v.supporting_reads)

The reads supporting all variants' ALT alleles are gathered and classified in a single kraken2 invocation (its database load dominates runtime), then split back per variant. You never re-implement fetch/allele-filter/classify yourself.

Classify an arbitrary set of reads

from nonhuman_screen import Kraken2Runner

result = Kraken2Runner("kraken2_db").classify_sequences(
    {"read1": "ACGT...", "read2": "TTGCA..."}
)
result.nonhuman_fraction        # 0.0–1.0
result.fractions()              # TaxonomicFractions(nonhuman=…, bacterial=…, …)
result.taxonomy_available       # False if nodes.dmp was missing (fractions unreliable)
result.nonhuman_read_names      # set of read names

Command line

# Per-variant allele NHF table (headline)
nonhuman-screen classify \
    --bam sample.bam --kraken2-db kraken2_db --ref-fasta ref.fa \
    --variants candidates.vcf.gz --out-prefix sample_contam
# -> sample_contam.variant_nhf.tsv, sample_contam.summary.json

# Whole-BAM summary
nonhuman-screen classify --bam sample.bam --kraken2-db kraken2_db

See docs/cli.md.

Coordinate convention

All pos arguments in the Python API are 0-based (pysam reference coordinates); variant keys are "{chrom}:{pos}:{ref}:{alt}". The CLI reads a VCF and handles the conversion for you.

API summary

Symbol Needs [bam] Purpose
Kraken2Runner.classify_sequences(seqs) no Classify {name: seq}ClassificationResult
ClassificationResult no Per-domain read-name sets, counts, fractions(), nonhuman_fraction, taxonomy_available, classification_failed
Kraken2Error no Raised when a kraken2 run did not produce a trustworthy tally (opt out with strict=False)
TaxonomicFractions no Per-domain fractions; from_result / over_reads
read_supports_alt(read, variant_pos, ref, alt) no Does an aligned read carry the ALT allele?
parse_kmer_votes(kmer_string) no Parse kraken2 per-read k-mer detail into taxid votes
VariantNHF no Result of the allele-based helpers: variant_key, nonhuman_fraction, supporting_reads, fractions, taxonomy_available, to_dict() (pos 0-based)
classify_variant_alt_reads(...) yes Allele-based NHF for one variant → VariantNHF
classify_variants_alt_reads(...) yes Batched allele-based NHF for many variants
classify_reads_from_bam(...) yes Classify reads by name and/or locus
reads_supporting_alt(...) yes Read names supporting an ALT allele

Docker

The included Dockerfile installs the pinned kraken2 (v2.17.1) and the package with the [bam] extra:

docker build -t nonhuman-screen .

docker run --rm -v "$PWD:/data" nonhuman-screen \
    classify --bam /data/sample.bam --kraken2-db /data/kraken2_db \
             --variants /data/calls.vcf.gz --out-prefix /data/contam

The image entrypoint is nonhuman-screen; mount your BAM and kraken2 database in. You still supply the kraken2 database yourself (see docs/database.md).

Testing

pip install -e '.[bam,dev]'
pytest

Unit tests need nothing extra. The integration and functional tests additionally need the kraken2 binary and self-skip without it. The functional tests run the CLI end-to-end against committed Streptococcus controls — a strep-injected BAM (positive), clean human BAMs (negative), and a broken database (which must exit 3, not report 0.0). Their rationale, how the fixtures are built, what they assert, and what they do not establish are documented in docs/testing.md.

Releasing

Releases are published to PyPI automatically by GitHub Actions using Trusted Publishing (OIDC — no API tokens), via .github/workflows/release.yml:

  • Continuous builds — every push to main publishes <major>.<minor>.<run-number> (the <major>.<minor> base is read from pyproject.toml), so the latest code is always installable from PyPI.
  • Tagged releases — pushing a vX.Y.Z tag publishes exactly that version:
    git tag v0.2.0
    git push origin v0.2.0
  • Bumping the minor/major — edit version in pyproject.toml on main; continuous builds then track the new base.

To switch to deliberate, semver-only releases, remove main from the workflow's push.branches trigger (leaving only tags).

License

MIT — see LICENSE.

About

Classify sequencing reads by non-human taxonomic content using kraken2, and reduce the per-read verdicts to a non-human fraction (NHF) plus per-domain breakdowns (bacterial / archaeal / fungal / protist / viral / UniVec-Core).

Resources

Stars

Watchers

Forks

Releases

Packages

Contributors

Languages