P3 validated: normalization bit-identical on 2141 vendor GTCs

Normalization 1.1.2 reproduces all 109,191 bins of the 1000 Genomes
Omni2.5 GTCs (AutoCall 1.6.2.2) with 0 ULP difference, reserved fields
included; 1.2.0 and 1.1.2 are bit-identical to the no-FMA gtc2vcf oracle
on 36 GSA samples. AWS run: 33 min, ~$0.15, no resources left.

README.md

@@ -9,9 +9,10 @@ Personal research toolkit for Illumina Infinium genotyping BeadChip data
 
 ## Status
 
-Phases **P0 (setup)**, **P1 (file parsers)** and **P2 (canonical model, raw X/Y)** are complete — see
-`docs/validation/`. Validated on 2141 vendor GTCs of the 1000 Genomes Omni2.5 set.
-Available commands: `simurg version`, `simurg inspect <files…> [--json]`, `simurg ingest`.
+Phases **P0 (setup)**, **P1 (file parsers)**, **P2 (canonical model, raw X/Y)** and **P3 (normalization)** are
+complete — see `docs/validation/`. Each is validated on 2141 vendor GTCs of the 1000 Genomes Omni2.5 set
+(normalization transforms bit-identical).
+Available commands: `simurg version`, `simurg inspect <files…> [--json]`, `simurg ingest`, `simurg normalize`.
 
 ## Quick start
 

docs/known_unknowns.md

@@ -4,6 +4,8 @@ Plan §20 lists the open questions; this file records what the implementation ha
 
 | ID | Date | Topic | Finding | Status |
 |---|---|---|---|---|
+| K-013 | 2026-09-27 | Normalization version of public vendor GTCs | The 1000 Genomes Omni2.5 GTCs (AutoCall 1.6.2.2) are reproduced bit-for-bit by normalization 1.1.2 (109 191/109 191 bins, reserved fields included). | resolved |
+| K-012 | 2026-09-27 | FMA and bit-exactness | gcc on aarch64 contracts a*b+c into FMA by default; the gtc2vcf oracle then differs from the vendor by up to 4 ULP in some bins. Oracles are built with `-ffp-contract=off`. Any future C/Rust port (P9) must also forbid FMA contraction in exact mode. | resolved; constraint for P9 |
 | K-011 | 2026-09-27 | Plan U1 (channel order) | Infinium II X = Red, Y = Grn confirmed on 2141 vendor GTCs (the GTC PDF's "Raw X = green" wording is wrong). | resolved |
 | K-010 | 2026-09-27 | Ingest throughput | On AWS, `ingest` of 2141 Omni samples was I/O-bound (gp3 125 MB/s baseline): 18.9 min, 4 of 16 cores busy. For P9 benchmarks provision gp3 throughput or use instance-store NVMe. | open (P9) |
 | K-009 | 2026-09-27 | Missing-intensity rule | 3 loci (in 3 of 2141 vendor GTCs, AutoCall 1.6.2.2) keep one channel's intensity where we zero both. Hypothesis: old AutoCall applies the NBeads ≥ 2 rule per channel. Confirm with the bead counts of 5767788121_R01C01 (rs178512), 5826354152_R02C01 (SNP12-29082770), 6268201099_R02C01 (rs17133206); if confirmed, add a `--missing-rule per-channel` compatibility option. | open |

docs/validation/P3.md

@@ -0,0 +1,41 @@
+# Phase P3 validation report — normalization
+
+Date: 2026-09-27. Development and oracle tests on the Pi; vendor comparison on a self-destructing AWS spot instance
+(`run-20260927T220605Z`, 33 min, ≈ $0.15). Summary: `P3_omni2141_summary.json`.
+
+## Delivered
+
+| Item | Module |
+|---|---|
+| Vendor numerics: percentile, trimmean, IQR, MAD sigma (old/new), bisquare robust fits, sequential sums | `signal/normalization/numerics.py` |
+| Six-parameter fit (outliers, offset, rotation, shear, scale, nn12/rr12) and the two-parameter single-channel model; normalization 1.2.0 and 1.1.2 | `signal/normalization/fit.py` |
+| Exact and binned (vendor) nearest-neighbour sweeps, vectorised | `signal/normalization/fit.py` |
+| Transform application and θ/R | `signal/normalization/apply.py` |
+| `simurg normalize` (per-sample transforms + normalized X/Y/θ/R Parquet) | `pipeline/normalize.py` |
+
+## Acceptance (plan §18 P3: transforms within 2 ULP for ≥ 99.9 % of bins)
+
+| Check | Data | Result |
+|---|---|---|
+| Transforms, normalization 1.2.0 | 36 GSA v4 samples × 48 bins, oracle = gtc2vcf built without FMA | **1728/1728 bins, 0 ULP** (all 10 fields incl. nn12, rr12, taa, tbb) |
+| Transforms, normalization 1.1.2 | same samples, 1.1.2 oracle | **1728/1728 bins, 0 ULP** |
+| Transforms vs **vendor GTCs** (AutoCall 1.6.2.2) | 2141 Omni2.5 samples × 51 bins | **1.1.2: 109 191/109 191 bins, 0 ULP; all 6 reserved fields identical in 2141/2141 samples** |
+| Same vendor GTCs with 1.2.0 | | 1/109 191 bins (expected: that software predates 1.2.0) |
+| `ingest` → `normalize` end to end | 18 GSA samples | bit-identical transforms |
+
+The gate is exceeded: every bin is bit-identical (0 ULP), not just within 2 ULP.
+
+## Findings
+
+- The public 1000 Genomes GTCs were produced with normalization **1.1.2**; `simurg normalize --norm-version 1.1.2`
+  reproduces them exactly.
+- Compilers that fuse `a*b+c` into FMA instructions (gcc on aarch64 by default) shift the last bits of the reference
+  C implementation by up to 4 ULP; the vendor's .NET code does not fuse. Oracle GTCs are therefore generated with
+  `-ffp-contract=off` (`scripts/make_oracle_gtcs.sh`).
+
+## Performance
+
+| Machine | Normalization fit per sample |
+|---|---|
+| Raspberry Pi 5, 1 core, GSA (683 k loci) | 1.4 s |
+| AWS m7g.4xlarge, 16 workers, Omni2.5 (2.45 M loci), both versions | 2141 samples in 8.5 min (≈ 3.8 s per sample-core for two fits) |