主题
R9 第二、三轮核查：数字可追溯性与独立复算
承接 第一轮定向排查。
- 第二轮：把 PPT / docs 里断言的每个数字，反查到产生它的产物文件，逐条判 PASS/FAIL
- 第三轮：脱离流水线独立复算关键结论，验证"产物本身算得对"
结果：第二轮 79 项断言零不一致；第三轮三项复算零差异。
一、第二轮：数字可追溯性
方法
不做肉眼比对——写核对脚本把断言与产物逐条比对，输出 PASS/FAIL/UNKNOWN。 查不到来源的标 UNKNOWN 而非默认通过。脚本可重复执行， 以后任何一轮重跑后再跑一遍，就能立刻知道 PPT/docs 里哪个数字过时了。
覆盖范围与结果
| 批次 | PASS | FAIL | UNKNOWN |
|---|---|---|---|
| 主核对 | 63 | 0 | 0 |
| 补充核对 | 15 | 1（脚本断言写错，非产物问题） | 0 |
核对内容：
- 结局表：7 个结局的病例/对照/工具数/显著数、工具强度 minF/medianF、最小 F ≥ 40.6
- 共定位：评估蛋白 95、判定四分布 13/14/54/14、13 个支持蛋白的簇 PP 逐个核对
- 三法：两法 19/13/0、SuSiE 32 对 / 12 收敛 / 4 支持、合表 142 对、Undetermined 24
- 糖尿病对照：57 条与四分类 30/11/10/6
- 外部复制：三组的覆盖 / 同向 / 同向 p<0.05
- LDSC：rg(T1D, T2D) = 0.087
- 靶点分层：31 蛋白、A/B/C/D = 3/4/11/13、A 级名单
- 脱靶负担：7 个关键蛋白的
n_traits_bonf、PheWAS 关联 15,789、覆盖 92 rsid - complonly：8 个支持蛋白、两版共有 6 个
- 数据完整性：
_parts/*.csv= 95、.skip= 0、X 染色体蛋白 39、补下 tar 49、本地 tar 146
唯一的 FAIL 是我自己写错了断言
补充核对里断言 REPORT/table_S_* 应为 14 张，实测 13 张。 核实后确认：实际就是 13 张 CSV + 1 个 Excel，是断言把 Excel 混进了 CSV 计数。 且"14 张附表"这个说法只存在于对话叙述中，未落进 PPT / docs / 记忆任何交付物。
第二轮的边界
它验证的是「PPT/docs 的数字与产物一致」，不验证产物本身算得对。 后者属于第三轮。
二、第三轮：脱离流水线独立复算
原则：不调用流水线任何函数，用独立实现从上游产物重新算，再与流水线结果比对。
A. 糖尿病对照四分类
从 screen_R9_dm/*_significant.csv 独立取 β/se/p，独立实现四分类判据 （糖尿病端 p≥0.05 → 并发症特异；方向相反 → 方向背离；z 检验显著且量级更大 → 并发症增强； 否则 → 糖尿病驱动）。
结果：57/57 一致，0 差异。
B. 两法 verdict
从 hyprcoloc_shared.csv 与 coloc_abf_pairs.csv 独立重建判定 （HyPrColoc：蛋白在簇内且簇 PP>0.7；coloc.abf：PP.H4≥0.8）。
结果：142/142 一致，0 差异。
C. coloc.abf 的 PP.H4（独立实现 ABF 公式）
不调用 R 的 coloc 包：自己用 tarfile 解 tar 取蛋白区域、自己 zcat 取疾病区域、 自己按等位对齐、自己实现 Giambartolomei (2014) 的近似贝叶斯因子与后验合并。
| 位点 | 交集 | 等位同向 | 等位翻转 | 等位不一致 | 独立用 | 流水线用 | PP.H4 差 |
|---|---|---|---|---|---|---|---|
| MINDY1 × T2D | 1666 | 987 | 679 | 0 | 1666 | 1666 | 5.55e-16 |
| ERI1 × T2D | 4482 | 2410 | 2066 | 6 | 4476 | 4482 | 1.26e-03 |
| KIT × T1D | 5119 | 5094 | 0 | 25 | 5094 | 5119 | 1.61e-05 |
MINDY1 因为零等位不一致，两者精确到机器精度相同（5.55e-16）—— 这证明独立实现与 coloc 包算法等价，从而 ERI1 / KIT 的小差确实只来自 SNP 取舍， 且数字完全对得上：4482 − 6 = 4476、5119 − 25 = 5094。
D. 等位对齐审计（顺带查出的重要一点）
| 脚本 | 是否对齐等位 | 是否正确 |
|---|---|---|
05_hyprcoloc.R:195-200 | ✅ 一致保留 / 翻转变号 / 模糊丢弃 | 正确 |
26_coloc_susie.R | ✅ 同上 | 正确 |
24_coloc_abf.R:96 | ❌ 按 chr:pos 合并，不校验等位 | 对 coloc.abf 无偏（见下） |
为什么 HyPrColoc 对齐与否是生死攸关的
实测 MINDY1 有 679 个、ERI1 有 2066 个 SNP 的等位在两个数据源间是翻转的，占比很高。 HyPrColoc 按带符号的 β 做跨性状聚类——若不对齐，聚类结果会完全错误。 所幸它对齐了。这是本轮检查中少数几处「设计本来就对」的地方。
24_coloc_abf.R 不对齐不影响 PP.H4：coloc.abf 的 ABF 只依赖 z² = (β/se)²，与 β 符号无关。 但它会把同坐标的不同变异（ERI1 6 个、KIT 25 个）一并纳入，实测影响 ≤1.3e-3，不改变任何判定。
这一项不计入缺陷计数——它没有产生错误结论，与「读错列导致假分歧」那类并列会夸大问题严重性。 建议作为代码卫生项在下次重构时处理。
三、三轮核查的总结论
| 轮次 | 目的 | 结果 |
|---|---|---|
| 第一轮 | 拿已知缺陷模式定向排查其余脚本 | 查出 4 个确认缺陷，均已修；其中 PheWAS 缓存问题使关联数从 10,976 增至 15,789 |
| 第二轮 | 数字可追溯性 | 79 项断言零不一致 |
| 第三轮 | 独立复算 | 三项复算零差异；算法层与 coloc 包等价 |
结论层面：A 级候选 ERMAP / IFNAR1 / APOL1 自始至终未变。
三轮下来最值得记的一条方法论：本轮全部缺陷无一是靠通读代码发现的—— 全部来自「换数据实跑」「结果对不上时追根因」「拿已知缺陷模式定向排查」三条路径。