Can hm_source orient a palindromic SNV?

GH #371. The short answer is no, for every value of hm_source, and the reason is structural rather than statistical: the harmonization stage that writes hm_source never looks at an allele.

---

0. Why this was asked

Gate C defers palindromic SNVs — A/T and C/G pairs, whose strand complement is the pair itself — because nothing in the alleles says which strand they were written on, and guessing wrong inverts the dosage silently. It is the one failure mode where being wrong is both undetectable downstream and maximally harmful. At the Gate E2 scope that deferral costs 64,810,321 rows, 2.99 % of the fleet, concentrated in the small curated pre-2015 scores: 1,687 production-eligible scores carry at least one palindromic row, and among those the median score loses 10.41 % of its rows.

The assumed way out is a frequency comparison — the scoring file’s declared allelefrequency_effect against a reference panel’s frequency at that locus. Concordant means same strand, complementary means flipped. But only 167 of the 1,687 affected scores declare an effect-allele frequency at all, capping that route at 9.5 % of the deferred rows.

GH #371 proposed a second route. Harmonized files carry an hm_source column, and a 150-file sample found 44.2 % of palindromic rows carrying hm_source = ENSEMBL. If ENSEMBL meant “this variant was matched to an Ensembl/dbSNP record”, and dbSNP records are strand-defined, then those rows would already be forward-strand and would need no frequency at all — lifting the recoverable share from ~9.5 % to ~45 %.

That hypothesis is false. This document says why, on four independent lines of evidence, and prices what is left.

---

1. The answer

hm_source	What it actually means	Can it orient a palindrome?
ENSEMBL	The variant’s position was looked up from its rsID in the Catalog’s Ensembl coordinate map.	No.
liftover	The position was converted from the author’s build with pyliftover.	No.
Author-reported	The author’s own coordinates were kept, because the file was already on the target build.	No.
Unknown	No position could be resolved at all. Emitted by the harmonizer but absent from the published column documentation.	No — such a row is unplaceable anyway.

The verdict is the same for every value for one reason: hm_source is the provenance of the coordinate, not of the alleles. The PGS Catalog’s own column definition says so in as many words — “Data source of the variant position” — and the harmonization stage that writes it does not read, write, or check effect_allele or other_allele.

There is therefore no hm_source-based route, and no partial one. H1, H2 and H3 from the issue are all answered “no”, and H1 is answered “no” for a different reason than the issue supposed — not because ENSEMBL matching is unreliable, but because ENSEMBL never claimed anything about alleles in the first place.

§2 is the mechanism. §4.1 is the measurement that would have settled it without reading a line of the harmonizer: in one scoring file, hm_source = ENSEMBL labels 26,901 forward-strand rows and 26,801 reverse-complemented ones.

---

2. The evidence

2.1 The harmonizer has two stages, and only the second one touches strand

PGScatalog/pgs-harmonizer runs in two passes, and the distinction between them is the whole answer:

Stage	What it does	What it writes
HmPOS	Resolves each variant’s coordinate on the target build, by rsID lookup, by liftover, or by keeping the author’s.	hm_source, hm_rsID, hm_chr, hm_pos, hm_inferOtherAllele
HmVCF	Looks the variant up in a reference VCF and checks the alleles against it: matches, palindromic, reverse-complement.	hm_code, hm_vid, and on a correction hm_fixedStrandFlip + hm_reported_effect_allele

variant_HmPOS() returns exactly five values — hm_source, hm_rsID, hm_chr, hm_pos, hm_inferOtherAllele — and effect_allele is not among them. It is passed through untouched. hm_source = 'ENSEMBL' is set on the single condition that the row’s rsID was found in the coordinate map; the alleles are never compared to the Ensembl record’s.

Strand lives entirely in HmVCF. DetermineHarmonizationCode() maps the triple (matches-VCF, is-palindromic, is-flipped) to an hm_code, and FixStrandFlips() rewrites alleles only where hm_code == -4, preserving the original in hm_reported_effect_allele. The codes that matter here:

hm_code	Meaning
5	Variant matches the reference VCF
3	Matched, but with ambiguous orientation (A/T, C/G)
-4	Strand flip suspected; corrected by default
-5	Not present in the reference VCF

So the Catalog does identify palindromes — it just labels them and moves on. Its own source carries the unimplemented intent beside that code: a ToDo to use allele frequency to resolve the hm_code=3 variants. The frequency route is an open item in the Catalog’s harmonizer too, which is a useful calibration on how tractable it is.

One corollary worth stating, because it closes a tempting side door: hm_inferOtherAllele carries no strand information either. infer_OtherAllele() requires the author’s effect allele to appear verbatim in the Ensembl record’s allele list, with no complement check. For a palindromic A/T variant the author’s A matches whichever strand it came from, and T comes back either way.

2.2 The files we hold never went through HmVCF

That is the mechanism; here is our data. Scanning every cached scoring file’s header and column names (spike_palindrome_hm_source.py headers):

	files
cached GRCh38 scoring files scanned	5,377
carrying hm_source	5,377
carrying HmPOS_build / HmPOS_date header keys	5,377
carrying any HmVCF artifact — hm_code, hm_vid, hm_info, hm_fixedStrandFlip, hm_reported_effect_allele, hm_reported_other_allele, or an HmVCF_* header key	0

Not one file in the catalog we hold has been through the stage that checks alleles. The hm_code column that would flag a palindrome, and the strand correction that would fix a flip, are both absent from every row we have ever scored.

This is not a gap in our mirror. The Catalog’s FTP publishes exactly two harmonized flavours per score — PGS######_hmPOS_GRCh37.txt.gz and PGS######_hmPOS_GRCh38.txt.gz — and no hmVCF flavour at all. hm_code is not available to any consumer of the PGS Catalog, so there is nothing to fetch and no version of this route that works by downloading a better file.

2.3 Harmonization does not change alleles — measured, not inferred

The two arguments above are a code reading. This one is our data. For nine scores chosen to span all three hm_source values, the author-submitted scoring file was downloaded from the Catalog’s FTP and joined to our harmonized copy row by row:

scores	rows joined	effect_allele differs	other_allele differs
9 (3 ENSEMBL, 3 liftover, 3 Author-reported)	401	0	0

Harmonization moved coordinates — in the liftover files by up to megabases — and left every allele exactly as the author wrote it. Whatever strand a scoring file arrived on is the strand we score it on.

2.4 What the witness does

pgsc_calc — the Catalog’s own calculator, and our reconciliation witness (see `prs_primer.md` §5.4) — labels A/T and C/G SNPs ambiguous in its match step and excludes them by default (keep_ambiguous = false). It has no frequency-based resolution path anywhere; the only alternative it offers is to keep them as written, and its own parameter documentation describes that as assuming the scoring file and the target data are on the same strand.

Its match code is blunter still: even on the same build, it notes, data can arrive improperly strand-oriented.

So Gate C’s deferral is not a conservative deviation from the field. It is the field’s default, arrived at independently, and the one opt-out that exists is documented as an assumption rather than a resolution.

---

3. The fleet-wide distribution

The issue’s table was a 150-file sample. This is the full measurement over all 1,687 production-eligible scores carrying a palindromic row — every row of every one of them, 645,528,195 rows in total, of which 64,834,179 are palindromic.

hm_source on a palindromic row	rows	share	dominant in	of those rows, carrying an EAF
liftover	38,674,476	59.65 %	616 files	1,076
ENSEMBL	14,072,760	21.71 %	984 files	6,179,102
Author-reported	12,078,169	18.63 %	87 files	20
Unknown	8,774	0.01 %	0 files	6,266

and over the same files’ 577,298,945 non-palindromic biallelic SNV rows, for contrast: Author-reported 43.88 %, liftover 38.58 %, ENSEMBL 17.54 %, Unknown 0.00 %.

The measured palindromic total reconciles exactly against the manifest’s sum(n_palindromic) over the same scores — 64,834,179 both ways — which is what establishes that this scan’s allele-pair classification is the same one Gate B and Gate C use rather than a re-implementation that drifted.

Two things worth reading off that table.

The issue’s sample was right about files and wrong about rows, and the discrepancy is instructive. It reported ENSEMBL at 44.2 % of palindromic rows against the fleet’s 21.71 %. But it capped each file at 300 palindromic rows, which makes it a file-weighted measure wearing a row-weighted label — and on the file measure it was nearly exact: 84 of 150 files (56 %) ENSEMBL-dominant against the fleet’s 984 of 1,687 (58.3 %). Row-weighted, liftover dominates because the multi-million-row genome-wide scores were harmonized that way; file-weighted, ENSEMBL dominates because most files are small. This is the same row-aggregate-versus-per-score trap the primer flags at §5.4, arriving from a third direction.

Essentially every EAF-carrying palindromic row is ENSEMBL-sourced — 6,179,102 of the 6,186,464 in §5. So the one route that does work and the hm_source value the hypothesis singled out are almost perfectly correlated, which is exactly the coincidence that makes a spurious hypothesis look promising. The correlation is real; the causal claim behind it is not.

---

4. The counterexample: strand-discordant files are real

Everything above says hm_source cannot certify forward strand. A separate question is whether the catalog’s palindromes happen to be forward-strand anyway, which would make “assume forward” empirically safe even though it is formally unjustified.

It is not safe, and the demonstration is cheap.

The detector. On the forward strand, A/G and C/T are the same transition seen from the two strands and occur about equally often; likewise the A/C / G/T transversions. So for any file, top_fraction — the share of its unambiguous SNVs written as A/C or A/G — should sit at ~0.5. On three normal files measured here it lands at 0.4991, 0.5017, 0.5018, with the paired counts agreeing to within 1 %.

A value near 1.0 means something else entirely: every row was written in Illumina’s TOP allele convention, where the reported pair is always A/x. TOP is assigned per variant from the surrounding sequence, so roughly half of those rows are reverse-complemented relative to GRCh38 — and nothing in the scoring file records which.

It happens. Auditing the allele pairs of the outliers directly against the pinned GRCh38 FASTA:

score	trait	single-base rows	allele pairs present	forward-match	complement-only	palindromic
PGS005393	post-traumatic stress disorder	53,705	A/G, A/C only	50.1 %	49.9 %	0
PGS000741	prostate cancer	46	A/G, A/C only	41.3 %	45.7 %	6
PGS004237	coronary artery disease	200,000	all four	99.9 %	0.1 %	0

PGS005393 does not contain a single C/T or G/T row out of 53,705. That is not a statistical tendency, it is a coding convention, and half its rows are on the opposite strand from the genome. PGS000741 shows the same signature — and it carries six palindromic rows (three A/T, three C/G) that a forward-strand assumption would resolve by coin flip.

top_fraction alone is not enough, and checking that mattered. Three files — PGS005087, PGS005088, PGS005089 — are top_fraction = 1.000 and yet 0 % of their orientable rows are complement-only. Audited against the FASTA, other_allele == REF at 5,326 of 5,326 checkable loci in PGS005087: its alleles genuinely are the forward-strand pair, so the all-A/x shape describes its variant set rather than a relabelling of it, and its palindromes are no worse off than anyone’s. The detector therefore requires both signals — an A/B-shaped allele convention and rows that only match complemented. The first names the convention; only the second proves the file is strand-mixed. Shipping top_fraction alone would have refused three innocent files.

How much is actually at stake — the honest number. Under the two-signal rule, of the 1,687 palindrome-carrying eligible scores, 3 are strand-mixed (PGS000714, PGS000741, PGS000733), carrying 21 deferred palindromic rows between them. Relaxing to “any complement-only rate above 5 %, whatever the convention” gives 18 scores and 62 rows; above 1 %, 99 scores and 2,341 rows — 0.0036 % of the fleet’s 64.8 M. Fleet-wide the complement-only rate is 490,602 of 2,037,784,515 orientable rows = 0.0241 % (the primer’s 0.0226 % is the same numerator over all 2.166 B scoring rows).

So this finding changes the arithmetic of the palindrome question by essentially nothing. What it changes is the warrant. “Assume forward strand” cannot be defended as empirically safe, because the exceptions exist, are identifiable in advance, and fail in the one mode — silent dosage inversion — that nothing downstream catches. A policy that is right 99.996 % of the time and unfalsifiable on the remainder is not the same object as a policy that is right.

There is also a reason the exposure is small, and it is not luck: the large A/B-coded files carry no palindromes at all. PGS005393, the cleanest case in the catalog at 53,705 rows and 50.8 % complement-only, has zero. A submitter who knows their data is in A/B coding generally also knows the ambiguous SNPs are unrecoverable and drops them. The files that kept theirs are the ones that did not know — which is precisely why their palindromes should not be trusted.

4.1 And these files are hm_source = ENSEMBL

The same audit, split by hm_source, is the direct empirical refutation of H1 — not an argument from the harmonizer’s source, but a measurement on rows we hold:

score	ENSEMBL rows that are forward	ENSEMBL rows that are reverse-complemented	palindromic rows, all ENSEMBL-sourced
PGS005393	26,901	26,801	0
PGS000741	19	21	6
PGS000714	24	22	9
PGS000047	26	21	3

In PGS005393, hm_source = ENSEMBL labels 26,901 forward-strand rows and 26,801 reverse-strand rows in the same file. The label takes the same value on both halves of a file split down the middle by strand, which is what “carries no strand information” looks like when you measure it rather than reason about it.

And the rows H1 would have licensed us to orient — the ENSEMBL-sourced palindromes of PGS000741, PGS000714 and PGS000047 — sit in files that are 45–48 % reverse-strand. They are not merely un-certified by hm_source; they are among the most dangerous rows in the fleet, and H1 would have singled them out as safe.

(PGS000047 sits at top_fraction 0.875, just under the A/B gate, with a 44.7 % complement-only rate — a file that is strand-mixed without being cleanly A/B-shaped. It is in this table as evidence about hm_source, not as one of the three the two-signal rule flags.)

---

5. What is left of the frequency route

The 9.5 % ceiling is a count of rows that have a declared frequency. It is not a count of rows a panel could orient, because a frequency near 0.5 is as consistent with one strand as the other — that is the same ambiguity one level up.

Measured over all 167 EAF-carrying scores (spike_palindrome_hm_source.py eaf):

	rows
palindromic rows in those 167 scores	6,186,687
…declaring a parseable effect-allele frequency	6,186,464 (223 rows carry an unusable value)

and of those, the share a panel comparison could actually separate, by the margin from 0.5 required to call it:

required |EAF − 0.5|	rows orientable	of EAF-carrying palindromic rows	of the fleet’s 64.8M deferred
> 0.02	6,014,185	97.22 %	9.28 %
> 0.05	5,753,395	93.00 %	8.88 %
> 0.08	5,495,356	88.83 %	8.48 %
> 0.10	5,318,598	85.97 %	8.21 %
> 0.15	4,861,478	78.58 %	7.50 %
> 0.20	4,389,812	70.96 %	6.77 %

The curve is gentle, which is the good news: the EAF-carrying scores are not clustered at 0.5, so the margin can be set conservatively without losing much. At a margin of 0.08 — equivalent to requiring a minor allele frequency below 0.42, a common convention in GWAS harmonization — the route delivers 8.48 % of the fleet’s deferred palindromic rows rather than 9.5 %.

This is an upper bound in one further respect that only D1 can close: it counts rows whose declared frequency is far enough from 0.5 to be decisive, and assumes the panel agrees with the author’s cohort closely enough to call the side. Where a score was trained in a cohort whose allele frequencies differ materially from HGDP+1kGP, a row can clear the margin and still not be callable. That measurement needs the panel and belongs to D1.

Note what the margin does not interact with: a file being Illumina A/B-coded (§4) does not spoil the frequency route. A/B coding is a per-variant strand assignment, and a per-row frequency comparison resolves it per row, which is precisely the case it is for. What A/B coding rules out is assuming forward strand — §6 item 4.

5.1 Why a bigger frequency source would not raise the ceiling

The obvious objection is that 1,520 scores lack a frequency only because the file lacks one, and frequencies are freely available — dbSNP, gnomAD, ALFA, and the HGDP+1kGP panel already staged for Gate D all carry them. Why not look the locus up?

Because orientation needs two frequencies, not one. The test is whether the value the author attached to their effect allele matches the panel’s frequency for the same allele read forward, or its complement. A panel frequency on its own says how common each base is at that locus; it says nothing about which of them the author meant by “effect”. Swapping in a larger panel replaces one side of a comparison whose other side does not exist.

The same answer disposes of “look the rsID up in dbSNP, which is strand-defined.” dbSNP reports an A/T variant as A/T — the allele set is anchored to the forward strand, but a palindromic set is its own complement, so the file’s A still matches under either reading. That is the definition of the problem, not a way around it. It is also exactly what our own reference check already does against the pinned FASTA, and why the scorer defers instead of resolving.

---

6. Recommended policy

1. Keep deferring palindromic SNVs, and stop treating hm_source as a possible reprieve. There is no partial credit here: the column is positional provenance and promoting any of its values would be exactly the “decide an ambiguous mapping that isn’t crystal clear” that CLAUDE.md forbids — with dosage inversion as the failure.

2. Do not run the Catalog’s HmVCF stage ourselves to get hm_code. It would tell us which rows are palindromic — which Gate B already computes — and its answer for those rows is 3, meaning unresolved. The one thing it adds, hm_code = -4 strand-flip detection, we already measure ourselves as n_complement_only_match, against a pinned reference rather than against Ensembl’s VCF.

3. Resolve palindromes in D1 only for the EAF-carrying subset, and only above a stated margin from 0.5. §5 prices it: 8.48 % of the fleet’s deferred rows at a margin of 0.08, against the 9.5 % headline. Everything else becomes a Gate F coverage gate, not a pending decision — there is no instrument that would settle it later.

4. Gate any assumption-based orientation on the per-score complement-only rate, which we already store. Gate C writes n_complement_only_match on every score; a rate above ~5 % of orientable rows means the file is strand-mixed and no forward-strand assumption may be applied to it — including the TOP/BOT reconstruction in §7. No new manifest column is needed for the guard, because the scorer computes the signal anyway. top_fraction is worth adding to Gate B only to explain which convention caused it, and §4 is the reason it cannot stand alone as the guard.

   Landed (GH #377): classifier v5 counts n_snv_pair_top / n_snv_pair_bot (migration 108), strand_profile() joins them to Gate C’s measured complement-only rate and returns strand_mixed — null, not false, when either signal is unmeasured — and the Manifest tab renders both signals side by side rather than collapsing them into a verdict.

   This does not gate the frequency route. A frequency comparison resolves strand per row and is precisely the right tool for an A/B-coded file.

5. Separately, the strand-mixed files are losing rows. In an A/B-coded file about half the non-palindromic rows fail the reference check and are dropped as ref_allele_mismatch — which reads as “the score and the genome describe different variants” when the real cause is a known, correctable convention. A per-file, evidence-gated flip would recover 26,885 rows across 5 scores, without weakening the global check that allow_strand_flip=True would destroy.

   (An earlier revision said 56,411 rows in 37 scores. That was the complement-rate cut alone, before the two-signal rule existed. Requiring the allele convention too removes PGS005285 — 29,389 rows at a 5.8 % rate but a balanced allele profile, so its mismatches are not an A/B convention and a flip there would be a guess, not a repair. The rule halving its own headline number is the rule working.)

   Landed (GH #376): ScoringPolicy.strand_flip_recovery_scores is a set of pgs_ids, not a second boolean — the global flag is an assumption about the catalog, this is a finding about named files, and the run’s provenance records the ids so the difference is auditable. scripts/prs/strand_recovery_candidates.py derives the set from the two signals and exits non-zero while any manifest row is still at an older classifier, so a provisional answer cannot be piped into a scoring policy by accident. Empty by default: turning it on obliges a re-score of the affected scores, which is Gate F’s call and can ride E3. It does not block Gate D, it cannot touch a palindromic row, and it is coverage rather than correctness.

---

7. What this does not settle

Added 2026-09-27 — this section originally read as though the frequency route were the only one left. It is not, and this document’s own §2 is what makes the other one legitimate. Because HmPOS never touches alleles, a scoring file’s alleles are the author’s, verbatim and uniformly — so its non-palindromic rows are a direct readout of that author’s strand convention, measured against the pinned reference on thousands of rows per file. §6 item 4 uses that signal as a guard; its positive side is a licence, and I wrote the first without drawing the second.

Banded over all 1,687 palindrome-carrying eligible scores: 99.899 % of the 64,834,179 deferred palindromic rows sit in files whose non-palindromic rows match the reference forward at ≥ 99.9 %. Forty-six rows are in strand-suspect files. Against the declared-EAF ceiling of 8.48 % in §5, that is roughly 12× the coverage, and it needs no scoring file to have volunteered a frequency.

It is a hypothesis with one named failure mode — a pipeline that orients non-palindromic variants and leaves palindromes alone, which the GWAS Catalog and OpenTargets sumstats harmonisers genuinely do — and the 167 EAF-carrying scores are the instrument that decides it. GH #404 specifies the experiment.

Measured the same day (`palindrome_forward_rule.md`): the licence holds per class of file, not for the 99.899 %. The failure mode is real — PGS000020 has ~40 % of its palindromes reversed at a 0.175 % complement-only rate — and 87 % of the genome-wide rows sit in pipelines with no frequency to test.

Whether a TOP-coded file’s palindromes could be oriented by reconstructing Illumina’s TOP/BOT rule. The designation for an ambiguous SNP is decided by walking outward through the flanking sequence — which we have, in the pinned FASTA. So in principle, a file known to be TOP-coded could have its A/T and C/G rows oriented by computing what TOP would have meant at that locus.

This is not recommended, for two reasons. The first is that any convention reporting the A-bearing strand produces A/C and A/G pairs on unambiguous SNVs, so top_fraction ≈ 1 identifies the family of conventions, not Illumina’s specifically — and only Illumina’s full rule defines a designation for the ambiguous pairs. Assuming the submitter implemented the flanking-sequence walk rather than a naive “A-first” relabelling is exactly the judgement call that must be tested before it ships, not shipped because it is plausible.

The second is that it cannot be tested on what we hold. Validating it needs a file that is both A/B-shaped and EAF-carrying, so the reconstruction can be scored against a frequency comparison. Of the 15 A/B-shaped scores in the palindrome-carrying set, two declare a frequency, and one of those has three unambiguous SNVs — its top_fraction is noise. The single usable candidate, PGS000714, carries 9 palindromic rows. That is not a ground-truth set; it is an anecdote. The hypothesis stays unproven and unprovable here, which is the strongest available reason to leave those rows deferred rather than the weakest.

The EAF route’s actual yield against a real panel. This spike measures how many rows could be oriented if a panel frequency were available. Whether the declared frequency and the panel’s agree — and how often a score’s declared cohort frequency is simply too far from HGDP+1kGP to be usable — is a D1 measurement, not a spike one.

Measured since, by GH #404 (`palindrome_forward_rule.md`), against the panel’s per-super-population frequencies. On unambiguous rows, where strand is known, the declared frequency and the panel agree in 99.982 % of decisive cases in genome-wide files and 97.7 % in curated ones. 8 of the 167 files carry a frequency column that fails that test, one because it reports the other allele’s frequency, and one publication’s frequencies were copied from 1000 Genomes rather than measured. #404 also tested §6 item 4’s inversion, reading palindromes forward where a file’s unambiguous rows match. That is licensed for small curated files and refuted by one mid-size file.

---

8. Reproducing it

backend/scripts/prs/spike_palindrome_hm_source.py, which needs DATABASE_URL (a git worktree has no .env; export the canonical checkout’s value):

python backend/scripts/prs/spike_palindrome_hm_source.py headers --out headers.json
python backend/scripts/prs/spike_palindrome_hm_source.py rows --shards 6 --shard N --out-dir out/
python backend/scripts/prs/spike_palindrome_hm_source.py eaf --out eaf.json
python backend/scripts/prs/spike_palindrome_hm_source.py report --in-dir out/ --with-e2

rows reads 10.7 GB of gzip; six shards, longest-processing-time-first on n_variants, balance to 1.00×. report prints the reconciliation first — if the measured palindromic total does not equal the manifest’s, every number under it is void.

Further reading: `prs_primer.md` §5.4 (why Gate C defers rather than guesses) and §5.7b (the staged road to release) · `compatibility_classification.md` (how a palindromic row is counted in the first place).