GH #396 asked for the two stages that feed Gate D’s normalizer. D1 reads the reference panel’s allele frequencies over the locus union. D2 scores every panel sample under the patient’s own scoring policy, so her raw polygenic sums can be placed in a distribution of the same quantity. D3, the regression and percentile, had landed the previous evening. The issue was specific about the contract, what was out of scope, and the exit criteria. It ended up as PR #426, merged.

The arithmetic turned out to be the easy part. What I want on record is the order in which things went wrong, because that order is the lesson.

The requirement the issue did not state

The issue said D2 must use the patient’s ScoringPolicy, and a comment from the #404 spike added that missing panel genotypes must be dropped, not imputed, to match Gate C. Both make sense. But the primer had a third line, in §5.7: “the patient and the panel must be scored over the same locus set, so the deflation cancels by construction rather than by hope.” The issue did not say how to satisfy it.

I spent longer here than anywhere else, and most of that time went on models that did not survive. Dropping missing calls on both sides looks symmetric. With random-sign weights it is not: her particular missing loci remove whole terms from her sum, while the panel’s missingness is smeared thinly across samples. My back-of-envelope put the resulting offset near a third of a standard deviation for genome-wide scores. For curated, positive-weight scores it nearly cancels. So “same policy” alone does not make the percentile like-for-like.

The version that settled it was cruder. A locus the panel does not carry is dropped from the panel for everyone, because absence is never zero. Her Gate C sum still counts it. At a row whose effect allele is the reference, a homozygous-reference patient carries two copies, so her side gains 2w that no panel sample can gain. That is the REF-effect bug’s shape, one level up.

The fix was to score her in the same pass, from her Gate A callset, as one more column under identical rules. run_gate_d.py then places that panel-aligned score and says so in the payload, and her stored raw score is left alone. This went beyond the issue, so the PR flags it as a decision for review rather than burying it. The residual asymmetry, her specific missing loci against the panel’s scattered ones, is now measured per score rather than argued about.

Making “same arithmetic” falsifiable

“The panel was scored under the patient’s policy” is exactly the kind of claim this project has learned to distrust: a check that shares the scorer’s assumptions cannot fail. So the centrepiece test builds a Gate A-shaped callset for every fixture panel sample from its own genotypes. It runs each through Gate C’s own score_variants and requires D2’s vectorized number to match to 1e-12. It covers 140 unrelated fixture samples, split multiallelic sites, haploid males on chrX and chrY, two deliberately inconsistent split calls, and five policies, including both per-file rules. It passed on its first real run, which made me suspicious. The suspicion was misplaced; what failed was elsewhere.

What went wrong, in order

• The first batch design did not fit in memory. Each matched row carried a site array as well as its slot, and scipy’s COO→CSR conversion took several more copies. That was 24–36 bytes a row, against a fleet of ~2.2 billion rows. I rewrote it before any test ran: 12 bytes a row held, one stable sort to build the CSR, and the dropped-call accounting riding on W’s own sparsity instead of a second matrix.
• BufferError on the first run. np.frombuffer over an array('q') holds a view, and clearing the array under it is refused. A copy fixed it. Trivial, but it cost a round trip.
• My test was wrong, not the code. I asserted 3 deferred palindromes in the kitchen-sink scoring file. The code said 5. I had written two C/G pairs into rows I meant as ordinary multiallelic and never-observed cases. The golden test had already agreed with Gate C on those rows. My hand accounting had not.
• My test helper skipped a normalization production does. The patient-column test expected her panel-absent contribution to be 6.5 and got 5.9. The missing 0.6 was an alt-contig row: I had built her callset with raw chromosome names, not Gate C’s parser, so 1_KI270706v1_random never matched its uppercase lookup. The fix was to use the production parser in the test.
• A merge that would always refuse. Rereading my own code after everything was green, I found the matrix header’s patient block carried a timing, elapsed_s. Shard headers must be identical to merge, so every real multi-shard D2 run would have been refused. The tests passed because both fixture shards rounded to the same tenth of a second.

The real data had its own opinions

The issue deferred running against the real panel, rightly. But one step is cheap and read-only: building D1’s index on the host from the real .pvar.zst and the v4 union, with no .pgen reads. It returned three things no fixture could have told me.

• It peaked at 3.9 GB. I had held each kept record’s REF and ALT as Python strings until assembly, ~16M records’ worth. Interning to ids while streaming dropped the same run to ~1.4 GB. There are only 87,495 distinct allele strings among 16,055,496 records.
• The panel carries 90.7 % of the v4 union, not 98.93 %. That is 15,200,494 of 16,755,919 positions. The 98.93 % in the primer, and in two vector memories, was measured on the older, smaller union. The v4 promotion added ~5.4M positions, and only ~73 % of them are in the panel. This does not change the design, but it makes the patient-aligned column matter more than the old number implied.
• The format held. All 84,325,213 records streamed in 221 s with no multi-ALT records and no REF conflicts, and the position-adjacency guard passed. The KING selection reproduced 3,330 unrelated samples to the individual super-population counts.

Two estimates of mine were also wrong without any run needed to show it. I had put the panel index at “~0.5 GB”; computing it from its own dtypes gives ~1 GB, so a D2 task at the default batch peaks around 8–9 GB. And my rule for when a block is read as one contiguous range was far too timid: 2× span, where #404’s measured 16.6k records/s in order against ~1k/s scattered says anything within 16× wins.

Where it stands

D1 and D2 are built and tested; neither has run. The real pass waits on three things. The #339 worker-image rebuild now has to add zstandard, because the .pvar.zst is unreadable in every container. D2 waits on the E3 re-score: ScoringPolicy.from_dict refuses to rebuild the E2 rows’ policy, because they predate two fields, and I would rather it refuse at the first moment the mismatch is knowable than let run_gate_d refuse it row by row later. And nothing yet computes the ancestry PC table D3 needs. run_gate_d now at least joins that table to the matrix by sample id and refuses a missing row, where before it would have fitted pairs by position.

What I take from it

The golden test was the right investment: it made “same arithmetic” a measurement, and it held. But every defect that mattered sat where a fixture’s size or shape hid it. Two shards’ timings happened to match. A 23-record panel cannot show a 3.9 GB footprint. And a coverage figure was inherited from an older artifact. The cheap read-only run on real data should come before quoting any memory figure for a worker task, not after. And a header that shards must agree on should never hold a timing.