解釋頁 03 · ISM 內部機制 · 給教授的細節層

ISM 怎麼「讀」甲基化、又怎麼「篩」？

從 BAM 檔的 MM/ML 標籤，到一張 read × CpG 的可分析矩陣，中間經過一連串解析與篩選。這頁把每一步拆開，每個門檻都標出源碼行號，並澄清「篩選」在 ISM 裡其實是兩件完全不同的事。

L0 一眼看懂

ISM 把每條 read 上 ONT basecaller 標記的每個 CpG 甲基化機率抽出來，對齊到參考基因組座標，組成一張「哪條 read 在哪個 CpG 是否甲基化」的矩陣。讀＝解析 MM/ML 標籤＋CIGAR 對位；篩＝沿路丟掉不可信／不合格的 read、位點、與模糊訊號。

▸ 給教授一句話：甲基「讀取」是 BAM 標籤解析（5mC-only、CpG-context、機率＝ML/255）；甲基「篩選」要分兩層 — 處理層（QC＋覆蓋＋二值化灰帶丟棄，決定什麼進矩陣）與應用層（已測的甲基 variant-filter formulation 為 TESTED NEGATIVE；新假說仍可經 pre-decision audit 重開）。兩者常被混為一談。

● 本頁最重要的一條：兩種「篩選」不要混

處理層篩選（data cleaning）＝ read QC、SNV 覆蓋、去重、CpG-context、二值化灰帶丟棄、C_min 覆蓋門檻。目的：讓矩陣乾淨可算。這是每次跑都在做的正常前處理。
應用層篩選（用甲基當變異 filter）＝拿甲基特徵去判斷某個 somatic 呼叫是 TP 還是 FP。目的：提升 caller 精度。已測 formulation 有多個 convergent negative tests；這不是不可推翻的 universal negative，重開須 materially new hypothesis + pre-decision audit。
教授若問「甲基化能不能幫忙篩變異？」——目前答案是「能讀、能 characterize；已測 formulation 未形成可用 filter」。這不是 universal negative；materially new hypothesis 可經 pre-decision audit 重開。§⑦ 給完整證據。

① 先建立語言（7 個關鍵詞）

完整定義在名詞地基；這裡放本頁會反覆用到的 7 個，點開看一句話定義。

MM tag 修飾位置（delta 編碼）

直覺：不存「第幾個鹼基被甲基化」，而存「跳過幾個未修飾的 C 再標一個」。例 C+m?,2,1,0＝跳 2 個 C → 標、再跳 1 個 → 標、再跳 0 個 → 標。C+m? 代表 5mC，C+h? 代表 5hmC。

ML tag 修飾機率（uint8 陣列）

直覺：一串 0–255 的整數，每個對應 MM 裡一個被標記的位置。甲基化機率 = ML[i] / 255。255→1.0（幾乎確定甲基化）、0→0.0。

5mC-only 只用 5-甲基胞嘧啶

直覺：BAM 可能同時帶 5mC（C+m?）與 5hmC（C+h?）。ISM 只取 5mC 當訊號；5hmC 的機率值只用來「對齊 ML 索引」，不進矩陣。

CpG-context 只認 CpG 雙核苷酸

直覺：一個 C 必須後面接 G（CpG）才算數。解析時會拿參考序列核對 ref[i]=='C' && ref[i+1]=='G'，不在 CpG 脈絡的修飾呼叫直接丟棄。

binary 0 / 1 / −1 二值化三態

直覺：把連續機率壓成三態：1＝甲基化（≥0.8）、0＝未甲基化（≤0.2）、−1＝缺失或灰帶模糊（0.2–0.8）。−1 在二值距離裡被當缺失丟掉。

C_min 最少共同 CpG 數

直覺：兩條 read 要算距離，至少要有 C_min 個「兩者都覆蓋到」的 CpG，否則這對距離視為無效，effective current CLI 預設用 SKIP 排除。生產預設 C_min = 3。

reliability gate cluster × label 列聯表可靠度閘

直覺：計算 HP label × algorithmic cluster 的 Cramér's V 時，另檢查列聯表夠不夠「飽」（Cochran：≥80% 格期望值 ≥5）。不夠時 reliability 旗標為 false；raw V 仍保留，只有 summary projection 寫 0，legacy gate 仍直接讀 raw V。三個口徑不可混用，且都不判定 cluster、label 或 cellular group 真值。詳見 ISM 方法本體。

basecaller 鹼基判讀（甲基標籤的來源）

直覺：把 ONT 測序的電訊號轉成 ACGT 序列、並順便標出每個 CpG 的甲基化機率的軟體（如 Dorado）。它把這些機率寫成 BAM 的 MM/ML 標籤 — ISM 讀的就是它的輸出，不自己重判。

CIGAR 比對格式碼（座標換算用）

直覺：描述一條 read 怎麼貼到參考的一串碼：M＝比對（讀與參考都前進）、I＝插入（讀有、參考無）、D＝缺失（參考有、讀無）、S＝soft-clip、N＝skip。ISM 用它把「read 上第幾個鹼基」換算成「基因組第幾號座標」。

somatic vs germline 後天突變 vs 遺傳變異

直覺：somatic＝腫瘤後天才出現（只在腫瘤細胞）；germline＝生殖系遺傳（全身都有）。ISM 以 somatic SNV 當錨點（要驗的目標）、germline 單倍型（HP）當骨幹（可信的分群依據）。

PERMANOVA 指定 label 的置換式 read-distance association 檢定

直覺：拿 read-read 距離矩陣，在指定分組（如 HP label）下比較 between-label 與 within-label variation，再隨機洗 label建立置換 null。顯著表示 read-distance variation／centroid separation 與 label 在該 null 下有 association；須與 PERMDISP 合讀，不能證明 unsupervised cluster existence/truth、cellular group 或因果。ISM 生產跑 999 次置換、seed=42；現行 p 解析度下限為 0.001。

LOSO 留一樣本交叉驗證

直覺：leave-one-sample-out — 每次留一個樣本當測試、其餘當訓練。若模型只在「訓練過自己」時才有效、留一就崩，代表它記住答案（循環）而非學到真規律。甲基 filter 在 LOSO 下顯示 100% 循環。

Cohen ribbon 效應量的「實質意義」門檻

直覺：一個改善要大到一定程度才算「有意義」而非雜訊。本專案對 ΔF1 設約 +0.005 的 ribbon；即使 ΔF1>0，只要低於 ribbon 就視為實質無改善。甲基 filter 的 ΔF1=+0.00242 正落在 ribbon 之下。

② 甲基讀取：MM/ML 標籤 → read×CpG 矩陣

核心檔案 MethylationParser::parse_read()。輸入一條 BAM read，輸出它在每個 CpG 的甲基化機率。下圖是五個步驟。

圖 2‑1　一條 read 的甲基訊號怎麼變成矩陣裡的一格（示意）

關鍵：甲基讀取是確定性的 — 同一條 read 同一參考永遠得到同一結果，沒有隨機性、沒有訓練。所有「智慧」都在後面的距離/聚類/檢定，不在讀取。5mC-only 與 CpG-context 是兩個會「丟資料」的設計選擇（見 §④ 與名詞地基的 ism-evolution 稽核）。
圖中 CIGAR 碼：M＝比對、I＝插入、D＝缺失、S＝soft-clip（完整定義見 §① CIGAR 卡）。

逐步拆解 MM delta 編碼（具體例子）

Worked Example · MM/ML 解碼（輸入為示意）假設一條 forward read，序列裡的 C 依序編號 C0 C1 C2 C3 C4 …，標籤為 MM: C+m?,2,1,0、ML: [230, 51, 200]：

delta=2 → 跳過 C0,C1 → 標記 C2，機率 = 230/255 = 0.902 delta=1 → 跳過 C3 → 標記 C4，機率 = 51/255 = 0.200 delta=0 → 跳過 0 個 → 標記 C5，機率 = 200/255 = 0.784

每個被標記的 C 還要通過 ④ CpG-context（後面必須是 G）才會留下；通過後以 參考座標當欄位寫進矩陣。公式 prob=ML/255：MethylationParser.cpp:140,180；delta 推進 next += deltas[i]+1：:150,188。

輸入 / 輸出一覽

項目	內容	來源 / 行號
輸入	一條 `bam1_t`（含 MM/ML/CIGAR/HP）＋該區參考序列＋ ref 起點座標	`parse_read()` MethylationParser.cpp:15
目標修飾	只取 `C+m?`（5mC）；`C+h?`（5hmC）跳過但計入 ML 偏移	MethylationParser.cpp:60, 71‑82
機率	`prob = ML[ml_offset+i] / 255.0`，型別 float ∈ [0,1]	MethylationParser.cpp:140, 180
座標對位	CIGAR 走訪建 seq→ref 映射；插入(I)/soft-clip(S) 標 −1（無參考座標）	build_seq_to_ref_map :232‑281
CpG 驗證	forward 找 C 且 ref 下一位是 G；reverse 反向找 G 且前一位是 C	:176‑178（fwd）, :138‑139（rev）
輸出	`vector<MethylCall>`，每個 = (ref_pos 1‑based, probability)	MethylationParser.hpp:14‑22

🧠 高認知負荷：為什麼反向股要「反著掃」？

MM 標籤是按原始 read 的 5'→3' 順序列出 delta 的。但比對到負股（reverse strand）的 read，BAM 裡存的 SEQ 是反向互補（3'→5' 的原 read）。所以要從 seq_len−1 往 0 掃，才能對上 MM 的順序；而且原本的 C 在反向互補後變成 G，故負股要找的目標鹼基是 'G'（CpG 的 G 端）。

原始註解與邏輯：MethylationParser.cpp:104‑157。這是甲基讀取最容易出錯的一段 — 弄反會讓負股所有甲基呼叫錯位。

③ 二值化：把連續機率壓成 0 / 1 / −1（以及「灰帶」去哪了）

矩陣先以原始機率存（raw），要算二值距離時再用兩個門檻壓成三態。這一步會悄悄丟掉「模糊」的位點，是個常被忽略的篩選。

圖 3‑1　二值化的三個區間（門檻 0.2 / 0.8）

關鍵設計：灰帶（0.2–0.8）在二值路徑被當缺失（−1）丟掉 — basecaller「不確定」的位點不參與二值距離。有效 CLI 在未指定 metric 時預設走 NHD，因此灰帶會被排除；只有明示請求 BERNOULLI 才會保留原始機率，並用信心權重把 p≈0.5 的貢獻壓低（見 §⑤）。

Worked Example · 同一條 read 的二值化沿用 §② 解碼出的三個機率：

0.902 → ≥0.8 → binary 1（甲基化） 0.200 → ≤0.2 → binary 0（未甲基化） 0.784 → 0.2<0.784<0.8 → binary −1（灰帶，二值距離中丟棄）

第三個位點明明有資料（basecaller 給了 0.784），卻因為「不夠確定」在二值化被丟。指派邏輯 RegionProcessor.cpp:1417-1422。

原始機率 raw	raw 矩陣	binary 矩陣	含義
無覆蓋（−1.0）	`NaN`	`−1`	這條 read 根本沒讀到這個 CpG
≥ 0.8	原值	`1`	甲基化
≤ 0.2	原值	`0`	未甲基化
0.2 – 0.8（灰帶）	原值（保留）	`−1`	模糊：二值距離丟、Bernoulli 軟降權（詳見 §⑤）

④ 處理層篩選：一條 read 要過幾關才進矩陣

在甲基讀取「之前與之後」，ISM 有一串守門。下圖是漏斗，下表是每一關的門檻與源碼。

圖 4‑1　read 進入甲基矩陣的漏斗（每關都可能淘汰）

注意：G5 是 tumor 專屬 — normal read 即使 UNKNOWN（沒覆蓋到 SNV）也保留，因為它的甲基資訊對 normal-anchored 基線有用（ReadAggregator.cpp:46-55）。這就是為什麼 normal BAM 是 ISM 的重要輸入。

七道守門（含位點/距離層）

FLAG 過濾：丟 secondary / supplementary / duplicate / unmapped。ReadParser.cpp:31‑42
MAPQ ≥ 20：比對品質太低不要。ReadParser.cpp:45；門檻 ReadParser.hpp:16
read 長度 ≥ 1000 bp：太短的不要（長讀才有結構資訊）。ReadParser.cpp:50‑51；門檻 :17
需有 MM + ML 標籤：沒甲基標籤的 read 進不來。ReadParser.cpp:56‑64；旗標 :19
tumor 須覆蓋 somatic SNV：alt_support==UNKNOWN（沒蓋到 SNV 位）的 tumor read 丟；normal 不丟。ReadAggregator.cpp:49‑55
去重：同 read_name 只留第一條。ReadAggregator.cpp:57‑60
位點/距離層：二值化灰帶丟（§③）＋ C_min=3 共同覆蓋門檻（§⑤）。RegionProcessor.cpp:1421；Config.hpp:37

Takeaway 處理層篩選的目的只有一個：讓矩陣每一格都可信、每一對 read 都可比。它不判斷對錯，只清掉雜訊與不可比的資料。

⚠️ 誠實標註：兩個「宣告了但沒接線/可選」的參數

① min_site_coverage = 5（「每個 CpG 至少 5 條 read 覆蓋才保留該欄」）：在 Config.hpp:36 宣告，但 grep src/ include/ 全域只有這一處宣告、無任何使用點（2026‑06‑13 驗）。屬預留 / 殘留參數 — 從程式碼無法判定是「刻意停用」或「未接線殘留」（標 undetermined，非指控疏失）。實務含義：現行矩陣建構不做 per‑CpG 覆蓋過濾，對外不可宣稱有此步。per-site cohort coverage 與 pairwise C_min 是不同 gate，對現有結論的影響為 UNVERIFIED，需接線後做 sensitivity rerun 才能判定。（若日後接線需更新本頁）

② PMD gating（pmd_gating=true, pmd_bed_path；PMD＝partially methylated domain 部分甲基化區，非 post-mortem damage）：同樣只在 Config.hpp 與 config receipt serializer 出現，核心路徑沒有 CLI／filter wiring。屬declared-but-unwired，不是可用 optional feature；預設旗標雖為 true，現行也不會實際套用 PMD BED。

這兩條是 §13.0「先有驗證過的事實才寫」的範例：與其寫「合理推測有做」，不如誠實標「宣告未接線」。

⑤ 距離計算：C_min 門檻、灰帶處理、與有效預設 NHD

矩陣建好後，算每對 read 的距離。這裡有第三種「篩選」— 共同覆蓋不足的 read 對被判無效。ISM 實作了 6 種距離；有效 CLI 無參數預設是 NHD，BERNOULLI 只在明示指定時使用。

距離	用哪個矩陣	怎麼處理灰帶	一句話
NHD 有效無參數預設	binary	灰帶＝−1 → 丟	不一致位點 / 共同位點（最直覺）
BERNOULLI 需明示指定	raw 機率	保留，用 2\|p−0.5\| 降權	期望不一致率，按信心加權
L1 / L2	raw	保留	機率差的平均 / 均方根
CORR	raw	保留（需 ≥3 點）	1 − 皮爾森相關
JACCARD	binary	灰帶＝−1 → 丟	甲基位點集合重疊度

為什麼會明示選 BERNOULLI？ 依源碼設計說明（DistanceMatrix.cpp:237‑252），Bernoulli「用原始機率 + 信心權重，降低 p≈0.5 低信心位點的影響」。好處：保留 basecaller 的機率資訊；代價是比 NHD 多算權重、較不直覺。這是可選路徑，不是目前 CLI 預設。

NHD vs BERNOULLI：灰帶位點的命運（具體例子）

Worked Example · 一個灰帶位點，兩種距離怎麼算（輸入為示意）兩條 read 在某 CpG 的原始機率：read i = 0.6、read j = 0.7（都落在灰帶 0.2–0.8）。

NHD（二值）：0.6→−1、0.7→−1 → 兩者皆缺 → 此位點不計入（共同位點少一個） BERNOULLI（原始）： 期望不一致 delta = 0.6·(1−0.7) + (1−0.6)·0.7 = 0.18 + 0.28 = 0.46 信心權重 w = 2|0.6−0.5| · 2|0.7−0.5| = 0.2 · 0.4 = 0.08（很低） → 此位點仍貢獻，但只佔權重 0.08（幾乎被忽略）

對比一個高信心位點 read i=0.95、read j=0.05：delta=0.95·0.95+0.05·0.05=0.905，權重 w=2·0.45·2·0.45=0.81 → 貢獻是灰帶位點的 ~10 倍。公式 delta／weight：DistanceMatrix.cpp:262-301（calculate_bernoulli）。

C_min 與無效對

每對 read 必須有 ≥ C_min 個「兩者都有效」的共同 CpG 才算得出距離；不足則該對標記為無效。effective current CLI 預設使用 SKIP 排除無效對；MAX_DIST=1.0 只是 legacy/experimental 選項，不是無參數生產行為。

參數	生產值	來源	備註
C_min（min_common_coverage）	3	Config.hpp:37 → RegionProcessor.cpp:790	`DistanceConfig` internal initializer 是 5，但被 current Config 的 3 覆寫
無效對策略	SKIP	Config.hpp:39 → RegionProcessor.cpp:791	不足 C_min → 標記 invalid 並排除
legacy/experimental max_distance_value	1.0	Config.hpp	只在明示選 MAX_DIST 時視為「最遠」；非 effective default
主距離 metric	NHD（無 CLI 參數）	ArgParser.hpp:181‑194	Config.hpp:40 的 BERNOULLI initializer 會被清空；BERNOULLI 需明示指定

給教授可講的一句：「無參數時我們用 NHD，灰帶當缺失；若要依 basecaller 信心加權，必須明示選 BERNOULLI；兩者都受 C_min 共同覆蓋門檻限制。」

⑥ 矩陣之後：可靠度閘 + 四道顯著性檢定（一句話帶過）

距離矩陣 → UPGMA 聚類 → 顯著性檢定。本頁專注「讀與篩」，這裡只點出與「篩選/可靠度」相關的關卡；完整在 ISM 方法本體。

關卡	做什麼	篩選/可靠度角色	來源
reliability flag／summary projection	HP label × cluster 列聯表做 Cochran 檢查（≥80% 格期望 ≥5）	不夠飽 → reliability=false、summary V=0；raw V 與 legacy gate 另保留，只約束 cluster × label association strength	`MathUtils.cpp::cramers_v` + `RegionProcessor.cpp` summary projection
PERMANOVA	置換檢定指定 HP label 與 read-distance variation／centroid separation 的 association	n_permutations=999、seed=42；現行 p-floor=0.001；須合讀 PERMDISP，不證明 cluster truth	`src/core/RegionProcessor.cpp`（`structure_config.n_permutations`）
GlobalTest / LocalTest / PerCpgAsm	cluster × label 或 per-CpG label × methylation association	Fisher／Cramér's V 只量 association，非 truth、cellular group、causality 或單一 filter	SignificanceAnalyzer

註：這些檢定量的是指定 label 與 read-distance variation、algorithmic cluster 或 per-CpG methylation 的 association，不是 unsupervised cluster existence/truth、cellular group、因果，也不是「這個 somatic 呼叫對不對」。PERMANOVA 還須與 PERMDISP 合讀；把 association 升格為 truth，正是 §⑦ 要避免的誤讀。

⑦ 應用層「篩選」：已測甲基 variant-filter formulation＝TESTED NEGATIVE

前面六節都是「處理層」。這節是完全不同的問題：已測的甲基特徵與資料切分，能否形成可用的 somatic TP/FP filter？ 多重驗證在該 scope 下為 negative。這限制的是已測 formulation，不排除 materially new hypothesis。

Verdict 甲基化「能讀、能 characterize ASM、能輔助 phasing」；但在已測 formulation 與資料 scope 下，未形成可用的 somatic variant filter。多道驗證收斂為 negative；新假說可另案重開。對外引用守「存在 ≠ 已證明能判別」。

+0.00242

ΔF1（甲基增強 filter pilot）
低於 +0.005 Cohen ribbon → 不算改善

第 5 順位

多軸 LR 的 deviance 貢獻（非機率，不求和=1）
caller_af 0.393 ≫ CN 0.211 ≫ HP 0.063 ≫ LOH 0.038 ≫ 甲基

100%

已測 LOSO 留一設計顯示的循環性
不外推到 materially new hypothesis

多個 convergent negative tests（限已測 formulation、scope與日期）

#	驗證	結果	機制 / 為何失敗	來源
1	甲基增強 FP filter pilot（HCC1395 單樣本）	ΔF1=+0.00242 < ribbon	在此模型中甲基訊號與 caller_af 高度重疊，未證明獨立增量軸	報告 20260518_V6_Methyl_Filter_Pilot
2	TP 救援（用甲基救回被丟的 TP）	TESTED NEGATIVE	95.2% 流失 TP 是低 AF subclone（caller_af<0.3）；已測救援規則同時重新引入 FP	同上 Step 5c
3	多軸 logistic regression 消融	甲基排第 5	caller_af deviance 0.393 主導；甲基在該模型的增量貢獻最低	同上 H 機制
4	LOSO 留一 + 全位點存在性掃描（6 樣本 28 萬位點）	100% 循環 / 非 usable	TP 顯著率 3.95%(11655/294723) > FP 1.07%(34/3177) 但絕對靈敏度 ~4%、COLO829 TP≈FP	報告 20260604_ISM_complete_TPFPFN

為什麼「scope-bound association」不等於「能判別」？（教授常問）

在具名資料與既定比較軸中觀察到 scope-bound allele-associated methylation — BRCA2 的 observed HP-axis Δβ=−0.122、6/6 biological IDs 的 excess-over-null >0 跨 3 癌種（見 INDEX 的 ASM 卡）。這些 association 不證明 causal cis 或 tumor-acquired；而「觀察到 association」與「能告訴你這個 somatic 呼叫是真是假」也是兩個正交問題：

association：具名 scope 的 label × methylation 統計訊號（顯著率 TP 略 > FP），不是 biological group truth 或因果。
判別：要在單一位點給可操作的 TP/FP 判斷 — 但靈敏度 ~4%、且 strong‑ASM 反而在 FP 富集（OR=0.194），方向甚至相反。

這個區分是本站「4 條鐵則」第 1 條，貫穿所有 ASM/甲基結論：association 不可升格為 truth 或 causality。

⑧ 收束：兩種「篩選」對照表

一張表把整頁收斂。教授若只看一張圖，看這張。

處理層 data cleaning

決定「什麼資料進矩陣 / 可比」。每次跑都在做。

read QC（FLAG/MAPQ/長度/MM-ML）
tumor 須覆蓋 SNV、去重
5mC-only、CpG-context
二值化灰帶丟 / Bernoulli 降權
C_min=3 共同覆蓋

狀態：正常前處理，確定性，無爭議。

應用層 variant filtering

想用甲基「判 somatic 呼叫真假」。提升 caller 精度。

甲基特徵 → TP/FP 分類器
pilot ΔF1 低於 ribbon
TP 救援 NEGATIVE
LR 消融：甲基第 5 順位
LOSO 循環、跨樣本不一致

狀態：已測 formulation concluded negative — 不升格成 universal claim；重開須新假說與 pre-decision audit。

One Line 「讀」是標籤解析；「篩」要分清是清資料（處理層）還是判變異（應用層；已測 formulation 為 negative，可用新假說重開）。