stayed exported from <code>obikindex</code> (<code>new_capped</code> too) since it's
general-purpose, but nothing in <code>obikphylo</code> calls it.</li>
<li>Splitting resolution chunks even finer (<code>/(n_workers*8)</code>, cap 1024,
instead of <code>/n_workers</code>, cap 4096) to smooth the residual sawtooth.
Measured ~10% <em>slower</em>, wider dips, not narrower. Reverted to the
original chunk sizing.</li>
</ul>
<p><strong>Known remaining limitation, not yet worth fixing:</strong> within one layer, the
four stages (sequential <code>unitigs.bin</code> read → parallel generation →
parallel resolution → sequential annex write) never overlap — confirmed by
1s-interval sampling: generation alone occupies ~17 threads evenly, but the
next layer's read/generation never starts until the current layer's
resolution and write are both done. This produces a real, periodic (~layer
duration) alternation between "many cores" and "few cores" that neither of
the fixes above touches, since both operate <em>within</em> one layer's resolution
step. The only remaining lever is overlapping consecutive layers (e.g. a
depth-2 pipeline: start layer N+1's read/generation while layer N's
resolution/write is still running) — a real restructuring, not a parameter
tweak, and explicitly <em>not</em> to be combined with the reverted
budget-capping idea above (let each phase use however many cores it
naturally wants; only the <em>scheduling</em> needs to overlap). Deferred, not
started.</p>
<h2id="cross-partition-batch-resolution-current-state-vs-the-batched-accumulator-design-discussion-2026-08-14">Cross-partition batch resolution — current state vs. the batched-accumulator design (discussion, 2026-08-14)</h2>
<p><code>family_scan.rs::scan_layer_families</code> (shared by <code>snp_pseudo_alignment</code>,
implements most of a dispatch/accumulate/resolve pipeline: generation
(cheap, CPU-only — builds <code>outgoing[dest_partition]</code> from <code>FamilyMask</code> and
buckets cross-partition queries) runs on an <code>obipipeline::throttle</code> +
<code>make_pipe!</code> stage, decoupled from resolution (I/O-bound, <code>rayon::par_iter</code>
<em>across partitions</em>, one generated batch resolved at a time, never several
concurrently — this ordering is deliberate, see the module's own docs on a
reverted concurrent-batch-resolution attempt that scattered mmap access).
The fast/slow mode gate (<code>PartitionCache::fast_mode</code>, <code>cache.rs:162-163</code>)
already exists: <code>n_layers <= 7</code> (checked once from the first non-empty
partition's <code>PartitionMeta::n_layers</code>, documented as identical across every
partition of an index — a structural, build-time property, never a
per-partition state) decides whether <code>FamilyMask</code>'s recorded <code>layer_value</code>
can be trusted to skip straight to the right layer
(<code>find_presence_batch_fast</code>) or must fall back to scanning every layer of
the destination partition (<code>find_presence_batch</code>).</p>
<p><strong>Real gap, confirmed not implemented</strong>: resolution is triggered by the
<em>source</em> batch finishing (<code>FAMILY_BATCH = 65536</code> minorants read from the
scanned layer), not by an <em>output</em> accumulator filling up. Since most
central-base variants of a family route back to the same partition being
scanned (~48% per the k=31/m=11 measurement above), a <code>FAMILY_BATCH</code>'s
<code>outgoing[dest]</code> is large for the local/self partition and thin for the
other ~255 (or however many) destination partitions — each of those gets
resolved at low query density every batch instead of being accumulated
across several source batches until resolving it is worthwhile. This is
distinct from, and not fixed by, the fast/slow layer gate above.</p>
<p>Redesign sketched (not built): per-destination accumulators decoupled from
<code>FAMILY_BATCH</code>, flushed on reaching a size threshold instead of on source-batch
completion — a "hot" accumulator for the partition being scanned (sharded
one-per-generation-worker, no lock, since all <code>n_workers</code> pipeline workers
write to it concurrently — this differs from an earlier, simpler mental
model of "one thread owns one layer's local collector," which doesn't hold
here since <code>n_workers</code> threads cooperate on scanning <em>one</em> layer at a time,
not one thread per layer) and "cold" mutex-per-partition accumulators for
the rest, low contention expected since traffic to any single cold
destination is a small fraction of total.</p>
<p>This breaks the current strict-iteration-order delivery of <code>on_family</code>
(today: a reorder buffer keyed by batch, since a whole <code>FAMILY_BATCH</code>
resolves atomically). With cross-batch accumulation, a family only becomes
complete once <em>every</em> accumulator holding one of its outgoing queries has
flushed, at unpredictable, independent times — no longer streamable
strictly in order without a large, unbounded pending buffer. Resolution
sketched: replace order-dependent consumers with coordinate-addressed
writes instead of order-dependent appends (see <code>PseudoAlignment</code> idea
below) wherever possible, since <code>sibling_annex_stats</code>'s reduction (plain
counts) is already order-independent and needs nothing here.</p>
<p><strong>Superseded 2026-08-15</strong> by the <code>--subsample</code>/<code>--shannon</code> design below,
which sidesteps the accumulator redesign for now: bounding the number of
families actually resolved per layer (via sampling) keeps per-layer
resolution volume small enough that the batch-density problem above stops
mattering in practice for these two consumers. The accumulator redesign
remains relevant for a future <em>unsampled, full-index</em> run, but is not
required to ship <code>--subsample</code>/<code>--shannon</code>.</p>
<h2id="pseudo-alignment-at-scale-pruning-is-unavoidable-discussion-2026-08-1415">Pseudo-alignment at scale — pruning is unavoidable (discussion, 2026-08-14/15)</h2>
<p>The reference run (<code>phyloskims_sal_vac</code>-scale bacterial test set,
<code>iqtree.fasta</code>) produced a dense alignment for 13 genomes × 383,965 sites
(4.8 MB) — trivially small. The in-progress plant index is expected to
carry on the order of 9 billion minorant families; a dense byte-per-cell
alignment at that column count is unbuildable regardless of genome count
(hundreds of GB even at a handful of genomes). Long-term ambition is 6,000–
8,000 genomes on a large machine, which makes the per-cell cost dominant in
the other dimension too. Pruning the retained family set before
materializing anything is mandatory, not an optimization.</p>
<p><strong>Already free</strong>: <code>family_size() < 2</code> (no sibling variant registered at
all) is a zero-cost structural filter, read directly off <code>FamilyMask</code> bits,
already applied in <code>snp_pseudo_alignment</code>. Insufficient alone — per the
~80-85% mono-family estimate from earlier discussion, this only brings 9
billion down to roughly 1.3-1.8 billion, still unusable.</p>
<p><strong>Entropy definition — settled 2026-08-15, correcting an earlier wrong
turn.</strong> The project does <strong>not</strong> encode families as IUPAC ambiguity codes
interpreted the classical way (Fitch-parsimony subset-compatibility, or
ML's "one true state, uncertain which"); see
<code>docmd/theory/evolutionary_distances.md</code> ("Why the IUPAC/DNA encoding used
for the first <code>--snp</code> test was wrong") and the Sankoff resolution that
followed it. The real model is a genuine 16-state alphabet (the powerset of
<code>{A,C,G,T}</code>, <code>∅</code> included as a real state) scored with a <em>calibrated
<p><strong>Why not implemented</strong>: <code>coverage[slot]</code> is a per-member marginal —
summing members' coverages to approximate a family's entropy silently
assumes no genome carries more than one member at once. It cannot
represent or detect joint co-occurrence (a genome carrying both <code>A</code> and
<code>C</code> at once, i.e. a combined 15-symbol state) at all, which is exactly the
phenomenon <code>family_entropy</code>'s 15-state definition exists to capture (see
<code>CardinalityTally</code>/<code>cardinality_transition_probs</code>, the project's own
existing machinery for this same co-occurrence structure, built for the
Sankoff matrix calibration). A family that is in truth uniformly <code>AC</code>
across every carrying genome would look like a well-balanced 2-state split
under the marginal approximation (entropy ≈ 1) while its true 15-state
entropy is 0 — i.e. the marginal proxy's failure mode lands on exactly the
"saturated, uninformative" tail this pre-filter would need to catch,
undermining the point. It stays plausible as a coarse filter for the
<em>low</em> tail only (a dominant single member's marginal share reliably
predicts low true entropy too), but not as a stand-in for the high tail —
not pursued further for now.</p>
<h2id="-free-loss-tnt-pipeline-four-independent-scans-three-of-them-unsampled-found-2026-08-15-fixed-2026-08-15-see-implemented-below"><code>--free-loss</code>/<code>--tnt</code> pipeline: four independent scans, three of them unsampled (found 2026-08-15, fixed 2026-08-15 — see "Implemented" below)</h2>
<p>Measured on <code>phyloskims_sal_vac</code> (91 genomes): <code>obikmer phylo --subsample 500000
--free-loss --tnt</code> logs four sequential stages —
So today, <code>--subsample 500000</code> only bounds the pseudo-alignment step —
the SNP-distance matrix, the Sankoff base-pair calibration, and the
cardinality histogram are always computed over the <strong>full, unsampled</strong>
index regardless of the flag. This is not merely "different subsamples
per stage" (which would already be a problem worth fixing) — it's that
three of the four stages never subsample, which explains most of the
~10x runtime gap against <code>snp_pseudo_alignment</code> on its own.</li>
</ol>
<p><strong>Decided requirement</strong>: all four stages consume <strong>one shared selection</strong>,
computed once, not each stage either scanning everything or drawing its
own independent sample. Per-genome-pair SNP counts, the Sankoff base-pair
calibration, the cardinality histogram, and the pseudo-alignment all
describe the same set of families — the calibration and the alignment it
calibrates are now guaranteed to agree on which sites exist.</p>
<p><strong>Correction to the "single pass" framing above</strong>: <code>base_pair_tally</code>/
<code>cardinality_tally</code> both need <code>raw_snp_distance</code>'s <em>complete</em> aggregate
SNP/shared counts before they can derive <code>included[i,j]</code> (the
<code>ratio_ceiling</code> filter) — a genuine sequential dependency (<code>included</code>
can't be known until every pair's aggregate count is final), not an
artifact of the old code's structure. So the fix is <strong>two</strong> passes over
the shared selection, not one: pass A computes the aggregate counts (and
derives <code>included</code>); pass B fuses <code>base_pair_tally</code> + <code>cardinality_tally</code>
+ the pseudo-alignment (mutually independent once <code>included</code> is known)
into a single scan. Still a 4→2 reduction, and — per the clarification
that settled this — pass A itself now runs over the <em>same shared
selection</em> pass B uses (not the full unsampled index): "les ratios, on
les fait sur les sites sélectionnés, c'est tout, les autres sites
n'existent pas" — once a selection is chosen, both passes are bounded by
it, so on a real <code>--subsample</code>/<code>--entropy</code> run pass A is cheap too, not
just pass B.</p>
<p><strong>Entropy-biased selection's own resolution to the "forward-looking
complication"</strong> (entropy must be known before selection, but selection
happens during the same scan that would resolve it): see "Entropy-biased
selection" below — resolved via a persisted per-layer entropy annex, not
by restructuring the scan into an inline pre-pass.</p>
<h2id="entropy-biased-selection-soft-gaussian-weighting-not-a-hard-cutoff-decided-and-implemented-2026-08-15">Entropy-biased selection: soft Gaussian weighting, not a hard cutoff (decided and implemented 2026-08-15)</h2>
<p>Refines the "forward-looking complication" above with a concrete
mechanism. Instead of a hard <code>[low, high]</code> entropy band (or any other
exact cutoff) deciding which non-monomorphic minorant families are
eligible, selection is weighted by an <strong>unnormalized Gaussian kernel</strong>
centered on a target entropy: <code>w(entropy) = exp(-(entropy - μ)² / (2σ²))</code>
— deliberately not the normalized Gaussian density (which would peak
below 1 and complicate the "probability" reading) — this kernel form
equals 1 exactly at <code>entropy = μ</code> and decays smoothly to 0 away from it,
so it reads directly as an acceptance weight: the further a family's
entropy from the target, the less likely it is picked, with no hard
in/out boundary — a few "bad" sites can still get in, by design. <code>μ</code>
(default ~1.0) and <code>σ</code> (default ~0.5) are meant to be user-tunable.</p>
<p><strong>Mechanism: joint probability, not a weighted reservoir</strong>. Not
Efraimidis–Spirakis weighted reservoir sampling (an earlier, more complex
proposal, superseded before implementation) — instead, a single
independent accept/reject draw per qualifying candidate: draw
<code>u ~ Uniform(0,1)</code>, accept iff <code>u < p₀ · w(entropy)</code>. <code>p₀</code> is a single
<strong>index-wide</strong> rate, <code>N / total_count</code> (<code>total_count</code> = the sum of
<code>non_monomorphic_counts</code> across every layer, <code>N</code> = <code>--subsample</code>'s
target), applied identically at every layer — this alone already gives
each layer its proportional share (the same effect the old uniform
reservoir's explicit per-layer <code>n_layer = N · count_layer / total_count</code>
computation achieved, but without needing to compute it: applying one
rate uniformly is mathematically the same as proportioning per layer).
<code>p₀ = 1.0</code> when there is no <code>--subsample</code> at all — <code>--entropy</code> alone is a
pure soft entropy filter over the whole index, no size target. Properties:
(1) <strong>strictly generalizes the existing uniform sampler</strong> — with <code>σ</code> large
enough that <code>w ≈ 1</code> everywhere, this reduces to the old uniform <code>N/total_count</code>
draw; (2) yields "approximately N", not exactly N — expected accepted
count is <code>N · mean(w)</code>, always <code>≤ N</code> — an intentional relaxation, matching
"sous-échantillonnage à environ n" rather than the old reservoir's exact-N
guarantee; (3) one streaming pass, one random draw per candidate, no
reservoir state.</p>
<p><strong>Resolving "entropy must be known before selection, but selection
happens during the same resolving scan"</strong>: solved with a <strong>persisted,
<p><strong>Resolved</strong>: the existing hard "non-monomorphic minorant" eligibility
filter stays a hard gate upstream of the Gaussian weighting — only
qualifying families ever get a stored entropy value or a weighted draw.</p>
<p><strong>CLI, implemented</strong>: two <code>phylo</code> options, <code>--entropy <μ></code> and
<code>--entropy-sd <σ></code> (<code>obikmer/src/cmd/phylo/args.rs</code>). The entropic filter
activates as soon as <em>either</em> is given (<code>mod.rs</code>, computed once into an
<code>Option<EntropyBias></code> threaded through <code>--snp</code>/<code>--family-overlap</code>/
<code>--shannon</code>/the fused sankoff pipeline below). If active but one or both
are unset, defaults are <code>μ = 1.0</code>, <code>σ = 0.5</code>.</p>
<h2id="two-entropy-definitions-kept-side-by-side-for-comparison-2026-08-15">Two entropy definitions kept side by side, for comparison (2026-08-15)</h2>
<p><code>--shannon</code>'s CSV carries both <code>entropy15</code> (<code>family_entropy</code> — the settled
15-non-empty-state definition, see above) and <code>entropy4</code> (<code>family_entropy_4</code>
— plain nucleotide reduction), computed from the same already-resolved
<code>genome_mask</code>, not from the marginal approximation above. A genome carrying
several bases at once contributes to <em>each</em> base's count (counted once per
base present, not fractionally split, not folded into one combined state)
— a genome polymorphic for the family is present at more than one base by
construction, so it is expected to count more than once; the denominator is
the total base-occurrence count, not the genome count (the two coincide
only when no genome carries more than one base). Kept side by side
specifically to measure, on real data, how much the two diverge — not yet
analyzed.</p>
<h2id="persistentsparsebitmatrix-implemented-and-measured-2026-08-15"><code>PersistentSparseBitMatrix</code> — implemented and measured (2026-08-15)</h2>
<p>A row-major (k-mer-major), deduplicated sparse alternative to
<code>obicompactvec::PersistentBitMatrix</code>, motivated by the same sparsity that
drove <code>--subsample</code>/<code>--shannon</code> above, but pursued as a foundational
storage-layer change rather than an index-level workaround. Full design
history, rationale, and rejected alternatives (external Elias-Fano crates,
<code>cacheline-ef</code>, a single unsplit <code>dict_id</code> array) are in the dedicated
implementation plan (<code>vivid-mapping-tiger.md</code> at the time of writing — the
content below is the durable summary, not a pointer to a session-scoped
file). Also directly informed by Alanko, Bille, Gørtz, Navarro, Puglisi,
"Compact Data Structures for Collections of Sets" (2025,
<code>biblio/Alanko et al. - Compact Data Structures for Collections of
Sets.pdf</code>) — this design implements only their exact-duplicate special
case (a plain dedup dictionary), not their full subset-containment
hierarchy.</p>
<p><strong>Design</strong>: four on-disk components, each mmap-backed, built once per
layer (matching how the rest of the build pipeline already works — never
the whole multi-billion-row index at once): an <code>is_multi</code> rank-capable
flag per row (singleton vs. multi-genome), a fixed-bit-width array for
singleton rows (genome index directly, <code>ceil(log2(n_cols))</code> bits), a
separate fixed-bit-width array for multi-genome rows (<code>dict_id</code>,
<code>ceil(log2(n_distinct_multi_sets))</code> bits — kept apart from the singleton
array specifically because <code>n_distinct_multi_sets</code> can be large in
absolute terms even when multi-genome rows are a small <em>fraction</em> of all
rows, and a single shared array would force every row, singletons
included, to pay the wider width), and a deduplicated dictionary of
distinct multi-genome sets (Elias-Fano-encoded byte offsets + a
varint-encoded values blob). New low-level primitives added to
<code>obicompactvec</code> to build this: <code>PersistentFixedIntVec</code> (arbitrary,
runtime-parameterized bit width, width 0 included — needed once a real
bug surfaced, see below), <code>PersistentRankSelectBitVec</code> (rank1/rank0/select1
on top of the crate's existing <code>count_ones</code>, using
<code>common_traits::SelectInWord</code>), <code>EliasFano</code> (composes the two). A new
<scriptid="__config"type="application/json">{"annotate":null,"base":"../..","features":[],"search":"../../assets/javascripts/workers/search.2c215733.min.js","tags":null,"translations":{"clipboard.copied":"Copied to clipboard","clipboard.copy":"Copy to clipboard","search.result.more.one":"1 more on this page","search.result.more.other":"# more on this page","search.result.none":"No matching documents","search.result.one":"1 matching document","search.result.other":"# matching documents","search.result.placeholder":"Type to start searching","search.result.term.missing":"Missing","select.version":"Select version"},"version":null}</script>
