# K19 – Independent check of the Greenland–Longnecker slopes in update2026

Controller: K19 (AI agent), 30.10.2026. Files: `GL_ANTAKELSER.md` (inputs and assumptions per set), `gl_check.py` (own code; imports no pipeline module), `gl_check_out.json` (all numbers). Nothing outside `kontroll/K19/` was changed.

## Verdict in brief

1. **All GL slopes reproduce exactly once the pipeline's rules are written down.** 15 of 15 category sets match to 4 dp (and 2 dp). 10 of 10 update2026 rows built from GL (7 main pool, 2 multiverse alternatives, Vulcan) match to 4 dp. The linear two-stage slope matches to 4 dp: 1.2363 (1.0603–1.4415), printed as 1.24 (1.06–1.44).
2. **Why K16a did not reproduce them: an undocumented covariance variant, not missing-data assumptions.** The pipeline does not use the textbook Greenland–Longnecker/Orsini covariance. It takes the *correlation* implied by the GL pseudo-counts and rescales it to the published standard errors (S_ij = r_ij·sd_i·sd_j), where textbook GL keeps the pseudo-count covariance itself (S_ij = 1/A0*, …). The module docstring says "kovarians 1/A0*", so the code and its own description disagree. K16a used textbook GL. My textbook implementation reproduces K16a's numbers (Pietinen 0.7955, Tiemersma 1.9205, Ollberding 0.9656, Knuppel 1.4976, Yiannakou 1.6817 vs K16a 1.6815, linear slope 1.2262 vs K16a 1.2262). So both codes are right for what they each implement.
3. **Classification: a real ambiguity in method, wrongly labelled, not a coding bug.** The rescaled-correlation variant is coded on purpose (inline comment "korrelasjon fra pseudo, skalert til publisert varians"). It always gives a positive-definite matrix, which textbook GL does not guarantee. But it is not Greenland–Longnecker as the protocol (R4) and the docstring promise, and it is not Hamling 2008 either. It has one unwanted property (§3): with equal denominators the case counts drop out completely.
4. **Two small real bugs, both without visible effect:** (a) case-control margins: n_controls is used as the category total, so non-cases B = controls − cases where it should be controls (Tiemersma 1.9195 → 1.9191); (b) total-case lookup is not sex-specific, so each Al Rajabi sex stratum gets the whole-cohort total of 265 cases (0 effect in the pipeline variant, small effect under textbook GL). Plus one mislabel: n_participants goes through the odds model (1/A0 + 1/B0), not the risk model the docstring calls "kumulativ insidens" (1/A0 − 1/N0). The effect is ≤ 0.0001 because the outcome is rare.
5. **Impact on results:** main pool 1.0967 (1.0447–1.1513) with the pipeline's GL rows, 1.0955 (1.0433–1.1503) with textbook GL rows. Linear two-stage slope 1.24 (1.06–1.44) vs 1.23 (1.05–1.43). No conclusion changes. The variant matters most for single studies with equal or assumed denominators: Larsson 1.23 vs 1.17, Pietinen 0.83 vs 0.80, Al Rajabi 1.25 vs 1.22. It also matters for Yiannakou (full data) 1.63 vs 1.68.

## 1. Reproduction with the documented rules (variant PIPELINE, own code)

| Category set | Pipeline (derived_per100g.csv) | Own, same variant | 2 dp | 4 dp | Own, textbook GL |
|---|---|---|---|---|---|
| pietinen1999 men | 0.8338 (0.4313–1.6119) | 0.8338 (0.4313–1.6119) | yes | yes | 0.7955 (0.4110–1.5398) |
| tiemersma2002 | 1.9195 (0.8541–4.3141) | 1.9195 (0.8541–4.3141) | yes | yes | 1.9205 (0.8524–4.3267) |
| larsson2005 women | 1.2252 (0.8978–1.6721) | 1.2252 (0.8978–1.6721) | yes | yes | 1.1652 (0.8631–1.5730) |
| ollberding2012 | 0.9665 (0.8292–1.1264) | 0.9665 (0.8292–1.1264) | yes | yes | 0.9656 (0.8271–1.1273) |
| knuppel2020 m2 | 1.4956 (1.0964–2.0402) | 1.4956 (1.0964–2.0402) | yes | yes | 1.4976 (1.0970–2.0445) |
| alrajabi2021 men m2 | 1.4102 (0.6885–2.8880) | 1.4102 (0.6885–2.8880) | yes | yes | 1.3756 (0.7032–2.6912) |
| alrajabi2021 women m2 | 0.9363 (0.3077–2.8489) | 0.9363 (0.3077–2.8489) | yes | yes | 0.9091 (0.3251–2.5425) |
| yiannakou2022 | 1.6262 (0.8412–3.1436) | 1.6262 (0.8412–3.1436) | yes | yes | 1.6817 (0.8689–3.2547) |
| wada2017 men | 2.1759 (1.2272–3.8581) | 2.1759 (1.2272–3.8581) | yes | yes | 2.1714 (1.2291–3.8359) |
| wada2017 women | 1.0905 (0.4884–2.4351) | 1.0905 (0.4884–2.4351) | yes | yes | 1.0624 (0.4773–2.3644) |
| bradbury2020 | 1.3871 (0.9939–1.9360) | 1.3871 (0.9939–1.9360) | yes | yes | 1.3912 (0.9867–1.9614) |
| vulcan2017 men | 1.0043 (0.6554–1.5390) | 1.0043 (0.6554–1.5390) | yes | yes | 1.0135 (0.6590–1.5587) |
| vulcan2017 women | 1.0262 (0.6659–1.5815) | 1.0262 (0.6659–1.5815) | yes | yes | 1.0243 (0.6658–1.5757) |
| norat2005 m2 (spline set only) | (not printed per set) | 1.2279 (1.0036–1.5023) | – | – | 1.2280 (0.9896–1.5239) |
| mehta2020 (spline set only) | (not printed per set) | 1.0838 (0.4523–2.5969) | – | – | 1.0723 (0.4406–2.6099) |

Norat 2005 and Mehta 2020 per-set slopes are not written anywhere by the pipeline. They are verified indirectly, because the 12-set linear pool they enter reproduces to 4 dp.

Rows in update2026 (sex strata combined with fixed effects):

| update2026 row | Role | Published | Own (strata rounded to 4 dp, as the pipeline does) | 4 dp |
|---|---|---|---|---|
| Al Rajabi 2021 | main pool | 1.2506 (0.6845–2.2848) | 1.2506 (0.6845–2.2848) | yes |
| Wada 2017 (Takayama only) | multiverse | 1.7240 (1.0815–2.7483) | 1.7240 (1.0815–2.7483) | yes |
| Vulcan 2017 | sensitivity | 1.0150 (0.7491–1.3754) | 1.0150 (0.7491–1.3754) | yes |

Without the rounding these rows differ in the 4th decimal (2.2847, 2.7484, 1.0151). The pipeline combines the strata from `derived_per100g.csv`, which is written with 4 dp. This is harmless but worth knowing for future 4-dp controls. The other 7 single-set rows equal their set values.

## 2. Linear two-stage slope (12 spline sets, DL)

| Variant | RR per 100 g (95 % CI) | τ² | Q |
|---|---|---|---|
| Pipeline (json) | 1.2363 (1.0603–1.4415) | – | – |
| Own, PIPELINE variant | **1.2363 (1.0603–1.4415)** – 4 dp match | 0.0204 | 17.01 |
| Own, PIPELINE + case-control margin fix | 1.2364 (1.0604–1.4417) | 0.0204 | 17.02 |
| Own, + risk model for n_participants | 1.2364 (1.0604–1.4417) | 0.0204 | 17.02 |
| Own, textbook GL | 1.2262 (1.0488–1.4335) | 0.0220 | 17.35 |

The CI here uses z = 1.96. Everywhere else the pipeline uses 1.959964. This is a cosmetic inconsistency.

## 3. What exactly differs, and why it matters where it does

Textbook GL (Greenland & Longnecker 1992; Orsini et al. 2006, `dosresmeta`/`glst`) works like this. Diagonal = published var(log RR_i). Off-diagonal = c0, the pseudo-count covariance of the shared reference: 1/A0* for person-time, 1/A0* − 1/N0 for cumulative incidence, 1/A0* + 1/B0* for case-control.

The pipeline works differently. Off-diagonal = c0 · √(v_i v_j) / √(pv_i pv_j), where pv_i = 1/A_i* (+1/B_i*) + c0 is the pseudo-count variance. It keeps the pseudo correlation and inflates or deflates the covariance by the ratio of published to pseudo SEs. When adjustment makes the published variances larger than the crude pseudo variances (the usual case), the pipeline's covariances are larger than textbook ones. That changes the GLS weights.

A consequence I verified numerically: with equal denominators (Pietinen, Mehta, Larsson, Al Rajabi), the fitted pseudo-cases are A_i* = RR_i·A0*, so r_ij = 1/√((1+1/RR_i)(1+1/RR_j)). **The case counts cancel out completely.** Pietinen's slope is 0.8338 whether the cases are the printed 55/35/50/45, an absurd 10/500/3/90, or ×7. Under textbook GL the same three inputs give 0.7955, 0.6743 and 0.6510. Two things follow:
- in the pipeline, the "cases assumed" rule (total × RR/ΣRR) and the non-sex-specific total for Al Rajabi have **zero** influence on the results;
- the printed case counts for Pietinen and Mehta are silently ignored. That is hard to defend as "GL", since GL's whole information about the covariance comes from the counts.

The variant is closer in spirit to the Hamling 2008 idea (covariance consistent with the published variances) than to GL, but it is not Hamling. Hamling solves for effective counts that reproduce both RRs and published variances. The pipeline does not use Hamling, so I did not implement it.

Checks that came out clean: the dose is natural (g/day), not log; log RR is natural log; x = dose − reference dose (GLS through the origin, no intercept); CI→SE with z = 1.959964; the reference is always rank 1 with RR = 1 (sets with another reference are skipped, not re-referenced); z is used, not t (no df correction); open categories are handled only through the dose rules (CUP rule: open top = cut-point + half the neighbour width; lowest category from 0); zero cells are not handled, but none of the 15 sets has one and the skipped log has no GL failures.

## 4. Bugs and inaccuracies found (pipeline)

| # | Finding | Type | Effect |
|---|---|---|---|
| B1 | `gl_doseresponse.py` docstring and `mvmeta.gl_cov` docstring describe GL/Orsini covariance (1/A0*); the code uses pseudo correlation rescaled to published SEs | Method mislabelled / undocumented choice | Explains every K16a mismatch. Main pool 1.0967 vs 1.0955 (textbook); linear 1.24 vs 1.23; single studies up to ±0.06 |
| B2 | n_controls is used as the category total M_i in the odds model (B_i* = controls − A_i*); it should be M_i = cases + controls | Coding bug | Tiemersma 1.9195 → 1.9191; linear slope +0.0001 |
| B3 | n_participants is fitted with the odds model and c0 = 1/A0 + 1/B0; the docstring promises a cumulative-incidence (risk) model, c0 = 1/A0 − 1/N0 | Mislabel / minor model choice | ≤ 0.0001 (rare outcome) |
| B4 | The total-case lookup (`total_cases`) is keyed only under sex = all, so sex strata with missing cases get the whole-cohort total (Al Rajabi men and women: 265 each) | Coding bug | 0 in the pipeline variant (counts cancel); would matter under textbook GL |
| B5 | Sex-stratum combination uses rows rounded to 4 dp | Rounding | ≤ 0.0001 |
| B6 | z = 1.96 in the linear two-stage CI vs 1.959964 elsewhere | Cosmetic | none at 4 dp |

No bug was found in my own code. My three variants agree with each other where they should (rate model sets are identical across PIPELINE/FIXED_CC/FIXED_ALL), and the textbook variant independently reproduces K16a.

## 5. Recommendations

- Decide on purpose which covariance to use. If textbook GL is wanted (as the protocol R4 says), change the off-diagonal to c0 and rerun: main pool 1.10 (1.04–1.15) either way, linear 1.23 (1.05–1.43). If the rescaled correlation is kept, document it in REANALYSE_PROTOKOLL.md (change log) and in the docstrings, and say that with equal denominators the case counts do not enter.
- Fix B2 (case-control margins) and B4 (sex-specific totals), or document them; neither changes any reported number at 2 dp.
- Report the textbook-GL alternative as one more dimension of the multiverse/sensitivity table. It is a real analyst choice and moves single cohorts by up to 0.06.
- The spline was not re-fitted here. It uses the same `gl_cov`, so it carries the same variant. A textbook-GL spline would shift slightly in the same direction as the linear slope.
- Side note (not checked): the Vulcan 2017 men and women strata have identical per-category medians, probably the sex-combined quintile medians.