Skip to content

BrucellaID 2.0.0

Choose a tag to compare

@duceppemo duceppemo released this 08 Nov 02:26
· 24 commits to main since this release

The extended reference set is now the default, and that is the breaking change. A run with no --reference-set identifies against 432 curated genomes covering all 30 named species of the Brucella clade — the former Ochrobactrum included — instead of vSNP's 16 genomes and 8 species. Over every usable public Brucella read set in ENA, 4,042 runs, it matches the submitter on 98.7% against the 16-genome set's 96.2%; the two disagree on 146 runs, where the extended set is right on 91 of them and the smaller set on none. Across the 2,018 B. melitensis and 1,118 B. abortus runs that carry most real traffic the two agree, so the change adds coverage without disturbing the calls the smaller set already made well. The numbers are in the wiki.

What this costs a 1.x pipeline. No table changed shape — not one column was added or removed anywhere — but three things a parser might read did change:

  • ClosestReference is a GCA file name now (GCA_000742195.1_ASM74219v1_genomic.fna.gz), and MashDistance and SharedHashes are measured against a different, usually closer, genome.
  • The code column of --xlsx and species_code in the JSON hold a species epithet (abortus, abortus5) rather than a vSNP code. The break is partial, which is worse than total: 7 of vSNP's 16 codes — canis, microti, neotomae, ovis, suis1, suis2, suis3 — spell the same thing in both tables and keep matching, while ab1, ab3, mel1, mel1b, mel2, mel3, ceti1, ceti2 and suis4 silently stop. A code-to-species lookup will half work, per species. --xlsx now warns when it is writing epithets.
  • Notes is populated on nearly every sample, where in 1.x it was usually empty, because the dense set almost always has another genome of the same species inside the close-match margin. On the four tutorial read sets it goes from 1 of 4 to 4 of 4. A rule of "a non-empty Notes means look at this one" now fires on everything.

brucellaid_mash.tsv holds every distance, so it grows with the set: 5 kB to 170 kB for four samples, and brucellaid.json 24 kB to 132 kB.

--reference-set vsnp gives back 1.x's reference data — the same 16 genomes, the same sketch, the same threshold, and therefore the same species, closest reference, distance, shared hashes and sequence type on every sample we have checked. It is not a byte-for-byte 1.x run: the 2.0 features that do not depend on the reference set still apply, so a library past --loose-distance gains a note, inserted first, and Loose match beyond appears in the run parameters. --reference-set default, which only ever existed in development builds, is rejected with exit 2 rather than quietly accepted.

The default moved at the library level too, not only on the command line: reference_sketch() and load_references() resolve to the extended files when given none, so a script that calls analyse_samples() directly measures against the 432 genomes without knowing the flag exists.

Two other things are new. A closest reference that is not close at all now says so: past --loose-distance, 0.005 by default, the call carries a note with the distance and the shared hashes, without changing the species. The threshold comes from the ENA sweep, where nineteen in twenty calls that agreed with the submitter sat within 0.005. And a nearly-as-close reference of the same species is counted rather than listed, so the note reads "4 other Brucella melitensis genome(s)" instead of four file names that change nothing.

One of the 30 species took some finding: nothing in INSDC is labelled B. papionis, but its genome is there under "Brucella sp.", and how we identified it is in the repository, along with every mislabelled public record the validation turned up, all of them reported to the archives.

Fixed while reviewing this release, none of it a consequence of the new default: the legacy --xlsx report wrote a value beginning with = as a live Excel formula, so a sample name could display something other than the sample; a named pipe or a link to a device in an input folder could hang a run or end it before any report was written, losing every other sample; and a species resting on the sequence type alone, with no reference behind it, now counts as a warning rather than "ok" in the quality assurance report.

One caveat worth knowing, which predates 2.0 and is unchanged by it: the derived k-mer copy threshold is still low for some libraries, and on two of the four tutorial read sets it leaves MashDistance inflated roughly 25-fold — --min-kmer-copies 12 takes SRR8550484 from 775/1000 hashes at 0.0065 to 990/1000 at 0.00024. The species and the sequence type are right either way, and --reference-set vsnp gives the same figures, so read a distance from a default run as a floor on the match rather than a measurement of it. The evidence is on the Validation page.

See an example report, and the full list is in the changelog.