Can a file’s own unambiguous SNVs orient its palindromes?
GH #404, the successor to
#371. The short answer is: for most files, yes, with a
measured error under 1 %. The rule cannot tell which files are the exceptions, though, and
the exceptions are real. So it is licensed per class of file, and the class that holds most
of the rows cannot be measured with the data we have.
0. Why this was asked
Gate C defers palindromic SNVs (A/T and C/G pairs, whose strand complement is the pair itself). A wrong strand guess inverts the dosage silently. At E2 scope that deferral costs 64,834,179 rows across 1,687 production-eligible scores.
#371 closed the hm_source route and proved something
more useful along the way: the PGS Catalog’s harmonizer never touches an allele. So
every scoring file’s alleles are its author’s, verbatim, and its unambiguous rows (which
Gate C checks against the pinned GRCh38 reference) are a direct readout of the author’s
strand convention. #404 proposed to lean on that corollary: if a file’s unambiguous rows
all match forward, read its palindromes forward too.
The issue named the risk precisely. A pipeline that orients the variants it can check and leaves the palindromes alone produces a file that looks forward everywhere we can look and is unoriented exactly where we cannot. The complement-only rate is blind to that by construction. Only ground truth can price it, and the 167 scores that declare an effect-allele frequency (EAF) are the only ground truth there is.
1. The answer
| Class of file | What the frequency test found | Measured flipped share (95 % bound) |
|---|---|---|
| Genome-wide (≥ 100k rows), complement-only ≤ 0.1 % | 11 files, 6 publications, 5,353,342 decidable palindromes | ≤ 0.12 %, and it tracks each file’s own complement-only rate |
| Small curated (< 1k rows), zero complement-only rows | 61 files, 49 publications, 953 decidable palindromes | ≤ 0.96 % |
| Mid-size (1k–100k rows) | one file, PGS000020: unambiguous rows forward, ~40 % of palindromes reversed | refuted |
| Strand-mixed, A/B-coded (the control) | PGS000714: 5 of 5 palindromes reversed | as expected |
So the forward rule is right almost everywhere it was measured. Three things stop that from becoming a fleet-wide licence:
- The counterexample exists. PGS000020 has a complement-only rate of 0.175 %, the same order as the clean genome-wide files (up to 0.098 %). No threshold separates them except a post hoc one.
- The measured genome-wide files are not the ones the rule would mostly touch. They hold 9.7 % of the genome-wide class’s palindromic rows. 87 % sit in three publications that declare no frequency: Ma 2022, Truong 2024 (PRSmix) and Bakker 2023 (metaGRS). Nothing we hold can test them.
- The EAF-carrying set under-samples mid-size files. It has one, and it fails; the fleet has 535.
The recommendation (§7) is to license the small-curated class, which the measurement covers well, and to put the genome-wide class to Richard as an explicit extrapolation. That class is 98.8 % of the rows, so its answer decides the arithmetic.
2. The instrument, and why agreement alone is not the answer
The test. For an A/T pair, reading the effect allele on the other strand turns it into
the other allele, so the reverse reading’s panel frequency is exactly 1 − p. Comparing
the author’s declared EAF with the panel’s
allele frequency for the effect base read
forward is therefore a sign test. Same side of 0.5 is forward, opposite sides is
flipped, decisive only when both frequencies are clear of 0.5 by a margin (0.08 here, the
usual MAF < 0.42 convention). This is exactly the rule
mods_prs_normalization.resolve_palindrome_strand implements for Gate D, so this spike also
validates the code D1 will run.
The panel. HGDP+1kGP (pgsc_HGDP+1kGP_v1), read from its .pgen with pgenlib: allele
counts over the 3,330 unrelated samples, overall and per super-population. Each file’s
frequencies are compared with the population that sits closest to its EAF on its own
unambiguous rows (European for 114 of 167 files).
Why agreement is not the answer. A cohort frequency and a panel frequency can straddle
0.5 for reasons that have nothing to do with strand: ancestry mismatch, sampling noise, an
EAF column that means something else. So raw agreement sits below 100 % even if every
palindrome is forward. The same file’s unambiguous SNVs calibrate the instrument: their
strand is known from the panel’s own alleles, so its error rate on them is measured. Under
“every palindrome is forward” the palindromic discordance should equal that error rate, cell
for cell in (|EAF − 0.5|, |p − 0.5|). Each flipped palindrome raises its cell from e to
1 − e, so
f_hat = (observed discordant palindromes − expected) / (n − 2·expected)
estimates the flipped share. On synthetic data with known flip rates of 0, 1 %, 5 % and 50 % it recovers each within its interval.
The bound that needs no calibration. Observed discordance is r = e + f(1 − 2e). For
any instrument error e below 0.5, f ≤ r. So the raw discordance
upper-bounds the
flipped share without trusting the calibration at all. §5 leads with these bounds wherever
the calibration is shaky.
3. Three ways the ground truth could lie, and the check for each
A population mismatch. An East Asian cohort’s EAF compared with a European panel frequency straddles 0.5 far more often, which would inflate discordance. Check: each file is compared with its closest population (above). The calibration then absorbs whatever mismatch remains, cell by cell.
An EAF column that is not an effect-allele frequency. Some files report a minor-allele frequency, a constant, or the other allele’s frequency. Such a column would read as a wholesale “flip”. Check: a file’s instrument must call ≥ 80 % of its own known-strand SNVs correctly before its palindromes count. 8 files fail, set aside rather than pooled. PGS001775 calls 0 % of them right: its column is the other allele’s frequency, which without this check would have counted as 100 % flipped palindromes.
A frequency that was looked up, not measured. This is the circularity threat, and the reason it matters is structural. An EAF copied from a forward-strand reference by allele letter reports the forward frequency of whatever letter the file wrote. It therefore calls every palindrome forward whether it is or not: a check that shares the assumption it is checking, the lesson of the REF-effect bug. Check: the commonest such reference is 1000 Genomes, whose samples are inside the panel. So a copied EAF sits almost exactly on those subsets’ frequencies, far closer than an independent cohort can (every publication judged independent here sits at a median of 0.0055 or more from the nearest 1000 Genomes subset). One publication fails: PGP000324 (PGS002294) sits a median 0.0015 from 1000 Genomes EAS against 0.0095 from the full EAS group. Every genome-wide file is measured (medians 0.0065–0.0125, no closer to 1000 Genomes than to the full group).
4. The ground-truth set
167 EAF-carrying scores from 103 publications, streamed through the connector’s own row parser and classified by the same branch order Gate B uses. The palindromic counts reconcile file by file against the manifest. 6,180,181 palindromic and 35,688,115 unambiguous rows carry a frequency and a position, and 99.31 % of their positions are present in the panel.
The set is two populations wearing one label, and every result is reported in that shape:
| Files | Publications | Decidable palindromes (margin 0.08) | |
|---|---|---|---|
| Genome-wide (≥ 10k palindromic EAF rows) | 11 | 6 | 5,353,342 |
| Curated | 156 (118 with a verified instrument) | 97 | 1,831 |
A pooled row rate over all 167 would be a statement about 6 pipelines. This is the row-aggregate-versus-per-score trap of `prs_primer.md` §5.4, arriving again.
The reconciliation that anchors everything: the complement-only bands recomputed here match the issue’s table exactly (784 / 298 / 325 / 81 / 9 / 190 scores; 64,834,179 palindromic rows).
5. Results
5.1 The instrument, on rows whose strand is known
| Rows | Genome-wide | Curated |
|---|---|---|
| Forward unambiguous SNVs, called forward | 99.982 % of 30,997,612 | 97.735 % of 17,261 |
| Complement-only SNVs, called reversed | 99.905 % of 23,103 | 85.3 % of 68 |
The second row is the positive control. Complement-only rows are genuinely written on the reverse strand, and the frequency test sees them as reversed. So it has the power to see a real flip, on real data, at the rate that matters.
5.2 Genome-wide files
| Score | Publication | Complement-only rate | Decidable palindromes | Agreement | Ceiling | Flipped share |
|---|---|---|---|---|---|---|
| PGS003576–85 (×4) | Dahl 2022 | 0.098 % | 756,112 each | 99.817 % | 99.958 % | 0.141 % |
| PGS005349 | Tanha 2026 | ~0 | 721,600 | 99.998 % | 99.998 % | 0.001 % |
| PGS005241 | Tanha 2025 | ~0 | 528,018 | 99.999 % | 99.999 % | 0.000 % |
| PGS000909 | Campos 2021 | 0.055 % | 367,893 | 99.924 % | 99.998 % | 0.075 % |
| PGS000908 | Campos 2021 | 0.056 % | 350,186 | 99.925 % | 99.998 % | 0.073 % |
| PGS000907 | Campos 2021 | 0.009 % | 218,452 | 99.962 % | 99.998 % | 0.036 % |
| PGS000882 | Revez 2020 | 0.057 % | 71,597 | 99.950 % | 99.999 % | 0.050 % |
| PGS004926 | Tsepilov 2023 | 0.057 % | 71,148 | 99.949 % | 99.999 % | 0.050 % |
Pooled: agreement 99.883 % against a ceiling of 99.976 %. The calibration-free bound is 6,258 / 5,353,342 = 0.117 %, ≤ 0.119 %.
The column worth reading is the complement-only rate beside the flipped share: they move together, from ~0 with ~0 up to 0.098 % with 0.141 %. That is what uniform residual strand error looks like: whatever put 0.1 % of the unambiguous rows on the wrong strand put about the same share of the palindromes there, where no allele check can see them. It is also flat across the margin. The flipped share is 0.090 % at a margin of 0.20 and 0.104 % at 0.02, with no sign of a pipeline that aligned the easy palindromes by frequency and left the near-0.5 ones unaligned.
5.3 Curated files
| Set | Files | Discordant / decidable | Flipped share, 95 % bound |
|---|---|---|---|
| Rate = 0 | 62 (61 with a decidable palindrome) | 4 / 953 = 0.420 % | ≤ 0.958 % |
| Rate ≤ 1 %, excluding the two outliers of §5.4 | 66 | 4 / 1,068 = 0.375 % | ≤ 0.855 % |
| Fewer than 20 orientable rows | 24 | 2 / 36 | ≤ 16.5 % (too small to judge) |
For curated files the calibration runs slightly pessimistic: the rate-0 band’s f_hat is
−1.4 %, meaning palindromes agree better than the files’ unambiguous rows predict. So the
calibration-free bounds are the numbers to use here. Of the curated files with any discordant
palindrome, 9 were observed against 21.5 expected from instrument noise alone.
5.4 The counterexample: PGS000020
PGS000020 is a type 2 diabetes score (Läll et al., Genet Med 2016; pruning and
thresholding on external GWAS summary statistics). It has 7,502 rows. Its unambiguous rows
match GRCh38 forward: 12 complement-only of 6,867, 0.175 %, with a balanced allele
profile (top_fraction 0.483, so not an A/B convention). Its instrument calls its own
unambiguous SNVs almost perfectly (expected palindromic discordance 0.10). Its palindromes:
169 of 425 decidable are reversed, across both pair types (A/T 75 of 169, C/G 94 of 256).
The examples are unambiguous: 10:112989181 A declared at 2 % where the panel has A at 97 %.
That is the failure mode the issue described, observed. The forward rule would have silently inverted roughly 40 % of this score’s palindromes. Nothing in the file’s unambiguous behaviour warns of it, and its complement-only rate is within a factor of two of the cleanest genome-wide files’. A clean complement-only rate is necessary and not sufficient.
It matters beyond its own rows because of what it is. It is the only mid-size (1k–100k rows) file in the ground truth, and the fleet has 535 of them.
5.5 The control
The issue asked for the comparison to be run on files that fail the complement-rate gate,
expecting no better than chance. The one genuinely strand-mixed EAF file, PGS000714
(Sipeky 2020, all-A/x allele convention, 47.8 % complement-only), has 5 of 5 decidable
palindromes reversed (binomial p ≈ 4 × 10⁻⁸ against its calibration). The gate is measuring
what we think. The band-level “gate-failing” pool (302 decidable, 6 discordant) is mostly
small files in the 1–5 % band, where one or two complement-only rows out of fifty is not
strand-mixing. Their palindromes are fine, which is itself a finding about how coarse the
rate is for small files.
6. Where the rule would land
The fleet’s deferred palindromic rows by size of file and complement-only band, beside the ground truth covering each cell:
| File size | Band | Scores | Palindromic rows | Share | Ground truth |
|---|---|---|---|---|---|
| < 1k | rate = 0 | 527 | 14,250 | 0.022 % | 62 files → ≤ 0.96 % |
| < 1k | 0.1–1 % | 156 | 7,739 | 0.012 % | 4 files → ≤ 2.6 % |
| < 1k | 1–5 % | 81 | 2,279 | 0.004 % | 24 files → ≤ 1.7 % |
| < 1k | > 5 % | 9 | 46 | — | the control |
| < 1k | < 20 orientable | 190 | 394 | 0.001 % | too small to judge |
| 1k–100k | rate = 0 | 182 | 277,032 | 0.43 % | none |
| 1k–100k | 0–0.1 % | 184 | 399,563 | 0.62 % | none |
| 1k–100k | 0.1–1 % | 169 | 55,092 | 0.09 % | PGS000020 → ~40 % reversed |
| ≥ 100k | rate = 0 | 75 | 12,412,327 | 19.1 % | none |
| ≥ 100k | 0–0.1 % | 114 | 51,665,457 | 79.7 % | 11 files → ≤ 0.12 % |
And the genome-wide class (≥ 100k rows, ≤ 0.1 %) by publication, since a pipeline’s strand handling belongs to its publication:
| Publication | Files | Palindromic rows | Share | Method (PGS Catalog) | EAF? |
|---|---|---|---|---|---|
| Ma et al. 2022, AJHG | 74 | 37,998,118 | 59.3 % | DBSLMM, lassosum, P+T and others | no |
| Truong et al. 2024, Cell Genom | 64 | 12,057,352 | 18.8 % | PRSmix | no |
| Bakker et al. 2023, Stroke | 6 | 5,696,454 | 8.9 % | elastic-net metaGRS | no |
| the 6 measured publications | 11 | 6,184,026 | 9.7 % | PRSice-2, SBayesR(C), SCT | yes |
| 10 others | 34 | 2,141,834 | 3.3 % | mixed | no |
Two of the three dominant pipelines deserve a sentence each, labelled as inference from the method rather than measurement. PRSmix and metaGRS files are weighted mixtures of existing PGS Catalog scores, so their palindromes carry whatever orientation their components had. If any component was a PGS000020, the mixture inherited it. Ma 2022 is a multi-method benchmark built on external summary statistics. The measured set does include external consortium pipelines (Tanha on BCAC; Campos), which is some comfort, but it is not these.
7. Recommendation
| Option | Scope | Scores | Palindromic rows | Measured flipped share |
|---|---|---|---|---|
| A. Refuse | nothing | 0 | 0 | n/a |
| B. Small curated (recommended) | < 1k rows, zero complement-only over ≥ 20 orientable | 527 | 14,250 | ≤ 0.96 %, 61 files from 49 publications |
| C. B + genome-wide | adds ≥ 100k rows at ≤ 0.1 % complement-only | 716 | 64,092,034 | ≤ 0.12 % on 9.7 % of the added rows; unmeasured on the rest |
| never | 1k–100k rows | 535 | 731,687 | counterexample |
Recommended rule (B). Read a file’s palindromic SNVs on the forward strand when it has fewer than 1,000 rows, zero complement-only rows over at least 20 orientable rows, and is not strand-mixed by the two-signal test. Measured error: ≤ 0.96 % of palindromes reversed (95 %), point estimate 0.42 % raw, indistinguishable from zero after calibration.
It is small in rows and large where it counts. The primer’s §1.4 point is that these scores are ML features, consumed one score at a time. For the 527 scores in B, palindromes are a median 11.3 % of each score’s rows (IQR 8.0–14.3 %), the curated, interpretable scores where the deferral hurts most per score.
C is Richard’s call, and it is a real one. For the pipelines measured, the genome-wide rule is excellent: ≤ 0.12 %, error tracking each file’s own complement-only rate. The case for C is that it recovers ~10 % of the loci (median 10.1 % per score) in 189 genome-wide scores, at a probable error that small. The case against is that 87 % of those rows are in pipelines built differently from anything measured. Two of them mix existing scores and inherit their components’ orientation, and PGS000020 shows the unseen failure is not hypothetical. My lean is to hold C until one of the two routes in §8 can measure the dominant publications. For a genome-wide score, the deferral mostly costs variance that Gate D’s normalization carries, while a systematic inversion is exactly what §1.4 calls “a confident, learnable lie”.
What landed, inert until someone chooses:
ScoringPolicy.palindrome_forward_scoreshas the shape ofstrand_flip_recovery_scores: a set of ids, empty by default, named in the run’s provenance. Each score recordspalindrome_forward_ruleandn_palindromic_read_forwardin itsgate_cpayload, so Gate F can gate on them and a reader can see how much of a score rests on the rule.scripts/prs/palindrome_forward_candidates.pyderives the set per class. B is the default; C only with--include-genome-wide, which prints the caveat. The mid-size class is not offered at all, andtests/test_palindrome_forward_candidates.pykeeps it that way. Like the strand-recovery script, it exits non-zero while any manifest row is on an older classifier.- Populating the set changes the policy fingerprint, so it rides E3 with the patient and D2 with the panel, never a pass of its own (GH #360).
8. What this does not settle
The genome-wide class, and two ways to measure it. Both are feasible and neither was in scope.
- The source GWAS as the second frequency. #371 §5.1
argued that a bigger panel cannot help, because orientation needs the author’s
frequency and a file without EAF has none. The source summary statistics do have one. The
Catalog records the GWAS accession (Ma 2022’s scores cite
GCSTsources), and most harmonized summary statistics carry effect-allele frequencies and β signs. Joining a scoring file to its source by position recovers the author-side frequency, and comparing weight sign with β sign on unambiguous rows calibrates the join. That would measure Ma 2022 directly. - Mixtures through their components. A PRSmix or metaGRS file’s palindromes come from named component scores. Where a component is EAF-carrying or already measured, its verdict transfers.
The weight-sign mode. A pipeline can relabel a palindrome’s letters to its reference and attach the other allele’s β. The letters then look perfect and the contribution is inverted anyway. The frequency test catches this only when the EAF travels with the β, i.e. when it was measured in the author’s cohort, which is what §3’s provenance check establishes. Letter checks, including the REF/ALT convention of a file’s unambiguous rows (strict in Truong, Bakker, Tanha and Revez), cannot catch it. So they are not offered as evidence here.
The last 0.02. Palindromes within 0.02 of 0.5 are invisible to any frequency test. §5.2’s flat margin curve says the genome-wide flip rate does not climb toward 0.5 through 0.02. Below that it is unmeasurable, and the rule reads those palindromes forward on the strength of the neighbouring bins.
The EAF route and the forward rule, together. 73 of the 527 scores in B also carry an EAF, so D1’s frequency route could orient their palindromes row by row. The frequency route is the better instrument where it is decisive. When both are implemented, it should take precedence per row and the forward rule should cover what it leaves undecided. That wiring belongs to D1.
A trap in the fingerprint, noted rather than fixed. Today’s default policy fingerprint
(b96ec1786f6f19b4) already differs from E2’s (c6a513393403a5cb), because #376 added a
field after E2 ran. This change adds another. That is by design: a new knob must move the
fingerprint. But it means a D2 panel scored with current code will refuse every E2 row,
which is right, since E3 is planned before normalization anyway, and it would be confusing
to discover mid-run.
9. Reproducing it
backend/scripts/prs/spike_palindrome_forward_rule.py. Four stages, each resumable. The
module docstring carries the exact commands, including the throwaway-container step that
reads the .pgen (pgenlib has no Windows wheel and lives only in the live worker’s writable
layer, not in the maple-worker image). The stages were measured at 102 s, 175 s, 7.3 min
over 8,041,370 panel records, and ~4 min. The report prints the reconciliation first: if the
recomputed bands do not match the issue’s table, every number under it is void. Outputs land
in backend/data/prs_palindrome_404/ (host-only): report.txt, results.json, and the
per-file target rows.
Further reading: `palindrome_orientation_spike.md` (#371, the harmonizer and the EAF ceiling) · `prs_primer.md` §5.4 and §5.7b (why Gate C defers, and the road to release) · `compatibility_classification.md` (how a palindromic row is counted in the first place).