Skip to content

Add refilter_profile.py: apply stricter SNV thresholds without re-profiling - #230

Merged
MrOlm merged 1 commit into
masterfrom
claude/clever-dirac-xruhka
Oct 20, 2026
Merged

MrOlm merged 1 commit into
masterfrom
claude/clever-dirac-xruhka

Conversation

@MrOlm

@MrOlm MrOlm commented Oct 20, 2026

Copy link
Copy Markdown
Owner

Adds auxiliary_scripts/refilter_profile.py for #229. It re-applies stricter --min_cov / --min_freq to an existing inStrain profile without re-running inStrain profile, which for large studies can take months.

python auxiliary_scripts/refilter_profile.py -i sample.IS -o sample.c10_f10.IS --min_cov 10 --min_freq 0.1 [-g genes.fna]

How it works

Raising the thresholds can only remove calls, and a profile already stores everything needed to recalculate them:

  • raw_snp_table has every site called at the original thresholds, with its A/C/T/G counts at each mismatch level. Each site is re-called with the new thresholds, following the same per-mm logic as update_snp_table, including "cryptic" SNVs.
  • covT has per-position coverage at each mismatch level. It's used for breadth_minCov, and to tell "not called because of low coverage" apart from "not called because only the reference base was present".
  • clonT has per-position nucleotide diversity for positions at ≥ the original min_cov. Nucleotide diversity doesn't depend on min_freq, so positions below the new min_cov are simply masked out.

The scaffold, genome and (with -g) gene tables are then rebuilt with inStrain's own functions (make_coverage_table, genomeLevel_from_IS, profile_genes). The original profile isn't modified; a new one is written to -o.

Limitations:

  • Linkage needs the reads, so it's left out.
  • Only stricter thresholds work: lowering them raises an error, because sites that would be gained aren't stored.
  • Profiles made before v1.12 are re-called with the old SNV thresholds by default, so results stay consistent with the user's original run. --current_snv_thresholds overrides this.

Testing

test_refilter_profile profiles the test data at the defaults and at -c 10 -f 0.1, then re-filters the default profile to 10×/10%. The re-filtered SNVs.tsv, scaffold_info, genome_info and gene_info tables are identical to the fresh 10×/10% run. The only exclusions are linkage and rarefied nucleotide diversity, which is a random subsample in every run and doesn't depend on these thresholds.

This doesn't change the inStrain package itself, so merging won't publish a new version.


Generated by Claude Code

Re-applies higher --min_cov / --min_freq to an existing inStrain profile
using the data it already stores: sites in raw_snp_table are re-called
from their A/C/T/G counts at each mm level (stricter thresholds can only
remove calls), clonT is masked to positions with the new min_cov using
covT, and the scaffold / genome / gene tables are rebuilt with inStrain's
own functions. Linkage needs the reads, so it's left out.

The test profiles the test data at the defaults and at -c 10 -f 0.1, and
checks that re-filtering the default profile gives identical SNVs,
scaffold_info, genome_info and gene_info tables (except rarefied
diversity, which is a random subsample in every run).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01X4KRcccF9fQPJkzpN9YaXw
@MrOlm
MrOlm merged commit e8de867 into master Oct 20, 2026
1 check passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants