Pretrained deep-learning ITS classifiers read the flanking regions

Harness, intermediate tables and figure code for the audit manuscript.

Two pretrained MycoAI classifiers are benchmarked against best-hit alignment, SINTAX and HiTaC on 5,222 identical UNITE queries, one per genus, at two loci — full-length ITS and the ITS2 subregion of the same records — with a flank-grafting ablation to establish why the two differ.

Layout

code/                    the harness, in pipeline order below
data/truth/              truth tables, donor assignments, ITSxRust extraction QC
data/predictions/        raw classifier and search output, one file per method and arm
data/scored/             per-query scored tables; every number in the manuscript
figures/                 figA_methods.pdf
taxa/                    per-rank label sets recovered from the distributed models

Not included, and why

Everything omitted is either large or rebuilt by a script that is here. The scored tables in data/scored/ are the ones every number in the manuscript comes from, and they can be re-derived from data/predictions/ without rebuilding anything upstream.

• The manuscript source. Added at submission.
• The MycoAI models (246 MB and 1.1 GB). Zenodo record 10904344.
• UNITE+INSD. The release distributed with dnabarcoder; see the manuscript.
• Query sets and the five ablation arms (~19 MB of FASTA). Rebuilt by code/build_pilot_queries.py and code/build_flank_arms.py. Their truth tables are kept, since the scorers read those directly.
• Reference FASTAs and taxonomy maps (~8 MB). Rebuilt by code/build_align_reference.py and code/build_sintax_reference.py at each per-genus depth.
• Fitted HiTaC classifiers (*.pkl). Regenerate with code/hitac_refit.sh; a fit takes about ten minutes on eight threads.

Rebuilding from UNITE reproduces the query set only if the same release is used; the release is named in the manuscript.

Pipeline

# 1. recover each model's label space from its per-rank encoders
python code/extract_labels.py --model MycoAI-CNN.pt --outdir taxa

# 2. one query per genus, stratified by whether that genus is in the label space
python code/build_pilot_queries.py --classification unite.classification \
    --fasta unite.fasta --genus-file taxa/genus.txt \
    --family-file taxa/family.txt --max-per-genus 1

# 3. references mirroring the label space, at three depths
python code/build_align_reference.py --max-per-genus 20 ...
python code/build_sintax_reference.py ...

# 4. flank ablation arms (i)-(iv)
python code/build_flank_arms.py ...

# 5. searches and classification
bash code/its2_arm.sh          # ITS2 arm, alignment + SINTAX + both classifiers
bash code/sweep_depth.sh       # reference-depth sensitivity
bash code/hitac_refit.sh align_ref_sintax.fasta fullits_queries.fasta pilot_queries.fasta

# 6. scoring
python code/score_matched.py --hits ... --ref-map ... --truth ... --classified ...
python code/score_sintax.py ...
python code/score_hitac.py --fullits ... --its2 ...
python code/score_donor_attribution.py ...
python code/mcnemar_tests.py ... --out mcnemar_results_v2.tsv
python code/auroc_ci.py --table matched_matched.tsv --score identity

# 7. figure
python code/patch_figA.py && python code/make_figA.py

code/make_figA.py reads every plotted value from data/scored/results_summary.tsv, so the figure and the tables cannot diverge.

Notes on two things a reader may trip over

HiTaC reference headers need a terminal semicolon. Its parser expects one and silently truncates the final rank label without it, so a SINTAX-format reference that ends at genus loses a character from every genus name. code/hitac_refit.sh adds it. code/hitac_semicolon_check.sh verifies the fix on a subsample before a full fit.

Two arms are constructed inputs. The padded and inter-phylum-chimera arms are diagnostic and no sequencing workflow would produce them; they identify what the models attend to and are not estimates of field performance.

Environment

vsearch v2.31.0, ITSxRust v0.2.2, mycoai-its v0.0.5 with torch v2.2.2, scikit-learn v1.3.2, numpy v1.26.4, biopython v1.87 on Python 3.10.20. Several are pinned by MycoAI's dependency constraints rather than being current releases.

License

MIT, see LICENSE.