Markov State ModelMarkov state models
用一堆短軌跡拼出狀態之間的轉換速率。
01一句話只有一句
用一堆短軌跡估出狀態之間的轉換機率，再從中解出速率。
02Why should I care?我為什麼要讀這頁
因為它承諾了一件你非常想要的事：不挑 CV、不加偏壓，也能拿到速率——而且它的算力形狀剛好長得像你手上的資源。
自由能地景那條路要你先押注一個 collective variable，押錯就前功盡棄；而且它要一條又一條的長 job。MSM 走的是完全相反的路：跑很多條短軌跡，事後用統計把它們縫起來。共用叢集上同時跑得動的 job 只有個位數，長 job 排隊久、還會撞時限被砍——「很多條短的」在排隊系統裡反而好活。
先講結論免得你抱錯期待：MSM 便宜的是「job 的形狀」，不是「總取樣量」。 這一頁大半篇幅在講這句話的後果。
用你熟悉的語言講，MSM 的三個核心步驟你其實都做過：
- featurization 加 TICA 就是選基底。 你在 DFT 裡挑 plane-wave cutoff、挑投影到哪組局域軌域；這裡你挑哪些距離／二面角當特徵，再用 TICA 投影到最慢的幾個方向。基底選不好，後面全錯。
- implied timescale 對 lag time 收斂，就是 k-mesh 收斂測試的時間版本。 你不是在找「對的答案」，你是在找答案不再隨參數變的那個 plateau。看到 plateau 才有資格報數字，跟你看到能量對 k-point 密度變平才敢報吸附能是同一件事。
- 轉移矩陣對角化就是解本徵問題。 特徵值接近 1 的模式對應最慢的過程，就像低能激發態；PCCA+ 把大量 microstate 併成少數 macrostate，就像把細密的能帶併成幾條有物理意義的 band。
03What這是什麼
MSM 的假設只有一句話： 在某個 lag time τ 之下，系統從一個離散狀態跳到另一個狀態的機率，只跟現在在哪個狀態有關，跟怎麼來的無關。
把這個矩陣估出來，你就同時拿到熱力學與動力學兩樣東西：
- 穩態分布 π（特徵值 λ1 = 1 對應的左特徵向量）給你各狀態的平衡族群，因此也給你自由能：Fi = −kBT ln πi。
- 其餘特徵值給你慢過程的時間尺度：
- MFPT（mean first passage time） 給你狀態之間的平均轉換時間，倒數就是速率。
它是你的「狀態劃分 + lag time」這一組選擇的性質。記憶效應藏在狀態內部：把 τ 拉長，或把狀態切細，都能讓模型更接近 Markov。所以「這個系統是不是 Markovian」是個假問題，真問題是「在我這組劃分下，多長的 τ 才夠」。
04Why為什麼重要
因為它是少數幾個能同時給你 π（熱力學）與 MFPT（動力學）的框架，而這正是 funnel 缺的那一對。
以平衡結合 ΔG 排序的 screening funnel，排的是 π 那一側。你缺的是速率那一側，而 MSM 給的東西不只是速率，是整張狀態圖：有幾個構形狀態、族群各佔多少、彼此之間怎麼跑、哪一步最慢。Dagdas 與 Sternberg 兩組用單分子 FRET 看到的正是這種東西——Cas9 在 DNA 結合與切割之間有離散的構形狀態，HNH 要 docking 到催化構形才放行。Dagdas 那篇直接解出三個態（EFRET 分別是 0.19 ± 0.02、0.34 ± 0.03、0.97 ± 0.01），smFRET 的 dwell-time 分析其實就是一個實驗版的三態 MSM，這意味著你的計算輸出跟實驗輸出可以直接對話——這是本站少數幾個「算出來的東西有對照答案」的地方。
另一個吸引力是：它不需要你事先知道反應座標。慢過程是從資料裡浮出來的，不是你押注的。對 Cas9 檢查點這種「大家都知道 HNH 在動、但沒人確定該用哪個座標描述」的問題，這是很實在的優勢。
還有一個更務實的用法，甚至可能是這個專案最該做的：把 MSM 當成 CV 探勘工具。 就算你建不出可信的 MSM，TICA 的第一根慢分量本身就是資料告訴你的「最慢的方向」。把它拿去餵給 umbrella sampling 當 CV，比你憑直覺挑一個距離有說服力得多。
05How怎麼做
完整流程有七步，每一步都有一個標準的失敗模式：
- Featurization。 選一組平移旋轉不變的特徵：成對 Cα 距離的子集（不要全部，維度會爆）、骨架二面角的 sin/cos、接觸指示函數，或領域專屬的距離（HNH 催化殘基到 scissile phosphate，或 L1／L2 這兩段把 HNH 掛在蛋白主體上的 linker 的長度）。要選就要比：用 VAMP-2 score 客觀比較幾組候選特徵，不要憑感覺。
- 降維：TICA。 找出自相關衰減最慢的線性組合。這裡有一個必須內化的區別：PCA 找的是變異最大的方向，TICA 找的是最慢的方向。對動力學而言，變異大不等於慢——一段大幅甩動的末端 loop 變異很大，但它快得要命。用 PCA 當 MSM 的降維幾乎必然失敗，見 PCA。
- 分群。 在 TICA 空間用 k-means 或 regular-space clustering 切出數百到上千個 microstate。狀態要夠細，細到每個狀態內部的鬆弛比 τ 快。
- 估計轉移矩陣。 在 lag τ 下數轉換次數，通常施加可逆性（detailed balance）約束。注意估計器只會用最大連通集——沒連上的狀態被靜默丟掉，這是一個你必須主動去查的陷阱。
- 選 τ：implied timescale 圖。 把 ti 對 τ 作圖，選在 plateau 上的 τ。圖上一定要畫出 t = τ 這條線，落在線下方的時間尺度沒有意義。硬限制：τ 不能超過單條軌跡的長度——這是短軌跡策略最直接的代價。
- 驗證：Chapman-Kolmogorov test。 用 τ 建的模型去預測 T(kτ)，再跟直接在 kτ 估出來的比。兩者要在信賴區間內吻合。這是 MSM 唯一像樣的自我檢驗，沒有 CK test 的 MSM，reviewer 有充分理由不信。
- 粗粒化與解讀。 用 PCCA+ 把 microstate 併成 3 到 10 個 macrostate，算 macrostate 之間的 MFPT。誤差用 Bayesian MSM 或對軌跡（不是對幀）做 bootstrap。
軟體： PyEMMA 與 MSMBuilder 是把 MSM 推廣開來的經典套件，近年維護重心移到 deeptime。三者的流程概念一致，換套件不會改變上面七步的邏輯。
Adaptive sampling（決定你花不花得起的關鍵）： 用當前的 MSM 判斷哪些狀態取樣不足，下一輪就從那裡重新播種。這件事跟「只能少量並行、但可以跑很多短軌跡」的排隊環境非常相容——一輪跑幾條短的，分析完再決定下一輪從哪開始。它能大幅降低所需的總取樣量，但代價是穩態分布 π 的估計變得依賴重權，而且你多少得先知道轉換區大概在哪。
06When什麼時候用
- 當你要速率、但押不出一個可辯護的 CV 時。 這是 MSM 相對於 umbrella／metadynamics 最大的優勢。
- 當你的算力是「很多小塊」而不是「一大塊」時。 短 job 好排隊、被砍了損失也小、天然平行。
- 當目標過程比單條軌跡慢 10 到 1000 倍，而不是慢 106 倍時。 差距太大就不可行，見 When NOT。
- 在截短系統上。 只留 HNH 與它的直接環境，總取樣量的門檻才有機會落在可負擔的範圍。
- 當你想重複利用已經跑過的軌跡時。 featurization 加 TICA 加 implied timescale 這一段幾乎不用額外算力，跑一次就有資訊。
- 當你要當 CV 探勘工具時。 就算 MSM 本身建不起來，TICA 找出來的慢方向可以直接送去做 umbrella sampling。
07When NOT什麼時候別用
- 當總取樣量遠小於目標過程的時間尺度時，不要宣稱建出了 MSM。 為什麼會壞：MSM 確實可以估出比單條軌跡更長的時間尺度，因為它把跨越屏障的轉換事件縫起來——但你必須真的觀察到那些事件。從平衡態自然等待的話，所需總取樣量大約是「你想觀察的轉換次數 × MFPT」。給個定位：你目前是 4 replica × 約 87 ns，總量約 0.35 μs。而 East 等人用溶液 NMR 的 CPMG relaxation dispersion，在 HNH 上直接量到 ms 尺度的構形交換；整體切割更慢，pre-steady-state 動力學給的 intrinsic 切割速率是 4 到 6 s−1（stopped-flow 另外量到 RuvC 4.2 s−1），observed 速率約 1 s−1。
把數字接起來（注意這是把 NMR 量到的交換尺度當成檢查點轉換 MFPT 的代理，不是有人直接量過 docking 的 MFPT）：取 ms 當 MFPT，光是看到一次轉換就要約 103 倍於現有總量（這跟 HNH 距離 那頁講的三到四個數量級是同一件事）；而 MFPT 要有統計意義得看到上百次，於是差距變成 105 到 106 倍。就算樂觀地假設檢查點交換其實只在 μs，也還差 102 到 103 倍。怎麼看出它壞了：永遠先報原始轉換次數。 如果兩個 macrostate 之間只數到個位數次轉換，那個 MFPT 沒有統計意義。另一個特徵是 implied timescale 永遠不 plateau，只是隨 τ 線性往上爬；再加軌跡，它還是往上爬。
-
不要用 PCA 當 MSM 的降維。 為什麼會壞：PCA 最大化變異，而變異會被最柔軟的區域主宰。你自己的系統就有現成的地雷——建模補出來的 5′ 單股末端（res 1 到 11）RMSF 有 6 到 8 Å，晶體本體只有 3.9 Å；那段東西會直接吃掉前幾個主分量，而它跟 HNH 檢查點一點關係都沒有。怎麼看出它壞了：檢查主分量的載荷是不是集中在末端或補建區段；並且比較用 PCA 與用 TICA 在同一個 τ 下估出的時間尺度——估出來比較短的那個就是比較差的劃分，因為變分原理保證所有估計都是真值的下限。
-
不要跳過 Chapman-Kolmogorov test。 為什麼會壞：implied timescale 的 plateau 可以被一個壞劃分假造出來（切太粗時，好幾個過程被混成一個看似穩定的數字）。CK test 檢查的是模型自己的預測能力，比 plateau 嚴格。怎麼看出它壞了：CK 圖上預測的存活機率系統性地高於或低於直接估計值，而且超出信賴帶——這代表你的 τ 還不夠長，或狀態切得不夠細。
-
不要忽略連通性。 為什麼會壞：估計器只在最大連通集上工作，沒連上的狀態連同它們的幀一起被丟掉，你的 π 與 MFPT 描述的其實是一個子系統。怎麼看出它壞了：報出 active set 佔總幀數的比例。掉了三成就必須說明，而不是安靜地把它算掉。
-
不要拿 MSM 當「便宜的替代品」去掃 100 個變體。 為什麼會壞：MSM 不會在變體之間攤提成本，每個變體都要自己的一整套總取樣量。配上只能少量並行的排隊配額，這是死路。怎麼看出它壞了：你會先撞到磁碟再撞到科學。535,265 原子的系統要建 MSM 就得存高頻率的幀，而磁碟配額很緊，這個限制會在你發現統計不足之前先發作。
-
當你其實已經有一個好 CV 時，不要硬用 MSM。 為什麼會壞：無偏 MSM 要靠自然發生的稀有事件，而增強取樣可以直接把系統推過去。同一道屏障，umbrella sampling 通常便宜得多。怎麼看出它壞了：你辛苦建完 MSM，發現最慢的 TICA 分量跟你一開始就想手選的那個距離高度相關——那就是預算白花了的訊號。
-
不要對 adaptive sampling 出來的資料直接讀 π。 為什麼會壞：重新播種破壞了平衡取樣，族群比例被你的播種策略決定。怎麼看出它壞了：π 隨播種輪次系統性漂移。正解是明確使用重權估計器，並且把 π 與速率的可信度分開陳述——速率通常還撐得住，族群比例往往撐不住。
08Project Lens跟我的 project 多相關
為什麼不是更低（不是一星或二星）：
- 它的形狀跟你的資源吻合。 共用叢集上同時跑得動的 job 只有個位數；一堆短軌跡在這種環境下比一條超長軌跡好活得多，被砍了也只損失一小塊。
- 分析本身幾乎免費。 你已經有 4 replica × 約 87 ns 的 5F9R 完整三元軌跡。在上面跑一次 featurization 加 TICA 加 implied timescale，成本是 CPU 分鐘等級。
- 就算建不出 MSM，那張 ITS 圖也有價值。 一條永遠不 plateau 的 implied timescale 曲線，是一張誠實而且量化的「我們的取樣離目標還差多遠」證據。這種圖在 methodology paper 裡是加分的，因為它把「我們沒做到」變成「我們量出了差距有多大」。
- 它是最好的 CV 探勘工具。 TICA 的慢分量可以直接送去當 自由能地景的 CV，這個接力比兩者各自單獨做都有價值。
為什麼不是四星或五星：總取樣量。
要對 HNH 檢查點建出一個能通過 CK test 的 MSM，需要的 aggregate 大約是你現有的 105 到 106 倍——這是最粗的量級估計（NMR 量到的 HNH 構形交換在 ms 尺度，而你需要的是上百次轉換，現有總量只有 0.35 μs），可能因為 adaptive sampling 而降低，但降不到同一個量級。在 535,265 原子的完整三元上，這不可行；磁碟配額甚至會比 GPU 時數更早成為瓶頸。可行的路是先把系統截小，但截小之後動力學還算不算數，本身就是一整篇 paper 的論證，不是一個 gate 裡的一步。
具體怎麼用：
當成 diagnostic 與 open question，不要當成 gate。三個做得起的動作，按 CP 值排序：
- 在既有軌跡上跑 TICA，把最慢的分量畫出來，看它跟 HNH 距離相不相關。這一步幾乎不用錢，而且直接回答「我手選的 CV 有沒有道理」。
- 畫 implied timescale 圖並附上跨 replica 的 spread，當作取樣充分性的量化證據。這跟你既有的 CV 可靠性分析是同一套精神——你已經知道聚合量 CV 在 2 到 10%、逐殘基與 DNA RMSF 在 24 到 146%，ITS 圖是這個故事在時間軸上的延伸。
- 把「完整三元的 MSM 目前不可行、需要多少量級的取樣才可行」寫進 open questions。這是誠實的 negative result，不是失敗。
要老實承認的一點： 這個專案目前的 replica 數（4 條）對 MSM 而言太少，而且你已經吃過小 n 的虧——2 replica 時逐殘基一致性看起來 6/8 很漂亮，跑到 4 replica 掉到 3/8。MSM 的統計量需求比逐殘基分析高得多，所以在把 replica 數推上去之前，任何 MSM 的定量結論都不該寫進 paper。
09Reviewer Thinkingreviewer 會問什麼
你的兩個 macrostate 之間，實際觀察到幾次轉換？點開看參考答案
這是 MSM 最致命的一問，而且答案是一個赤裸裸的整數，藏不住。個位數次轉換的 MFPT 沒有統計意義，相對誤差大致隨 1/√N轉換 走。主動把 count matrix 的關鍵元素放進補充材料，比等人問出來好。如果數字難看，正確的做法是把結論降級成定性描述，並用 ITS 圖說明還差多少取樣。
implied timescale 收斂了嗎？CK test 過了嗎？誤差怎麼估？點開看參考答案
三樣缺一不可。ITS 圖要畫出 t = τ 的禁區線與跨 replica 的 spread；CK test 要附信賴帶；誤差要用 Bayesian MSM 或對軌跡做 bootstrap——對幀 bootstrap 是假的，因為幀之間高度相關。只給一條 ITS 曲線就宣稱收斂，是這個領域最常見的敷衍。
你怎麼選 feature？有沒有客觀比較過？點開看參考答案
說「我們用了成對距離」不夠。要用 VAMP-2 score 比較至少兩三組候選特徵，並說明為什麼最後選這一組。同時要說清楚有沒有把建模補出來的區段（例如 5′ 端 res 1 到 11）排除——那些區段沒有實驗座標約束，它們的運動是建模產物，一旦進了 feature 就會汙染整個 TICA 空間。
你的 MFPT 跟實驗速率差幾個數量級？點開看參考答案
必須正面回答，不能迴避。實驗上 Cas9 的 intrinsic 切割速率是 4 到 6 s−1、observed 約 1 s−1，而 HNH 本身的構形交換由 NMR CPMG 量在 ms 尺度；你的總取樣量則是 μs 以下等級。誠實的框架是：「本模型的 MFPT 是在現有取樣下的下限估計，與實驗值差 N 個數量級，因此我們只用它做相對比較與慢過程的辨識，不宣稱絕對速率。」把限制講成可量化的距離，比含糊帶過可信得多。
MSM 找到的最慢過程，跟你宣稱的 HNH docking 是同一件事嗎？點開看參考答案
不能假設。要把最慢的 TICA 分量或最慢的本徵向量投影回結構空間，明確指出是哪些原子在動、對應哪個 domain 的什麼運動，並附上代表性構形。如果最慢的過程其實是某個表面 loop 的翻轉或末端的甩動，那你的 MSM 是對的，只是它在回答另一個問題——而這件事必須自己說出來。
你用了 adaptive sampling，那你的平衡族群還算數嗎？點開看參考答案
不直接算數。重新播種破壞了平衡取樣，必須用重權估計器（例如 reversible / Bayesian 估計搭配軌跡權重）才能談 π。誠實的處理是把速率與族群分開講：轉換速率通常還撐得住，平衡族群比例則要標明依賴重權，並附上不同播種輪次下的穩定性檢查。
10Common Mistakes最常犯的錯
- 只報 MFPT 不報轉換次數。 後果：一個由三次事件撐起來的速率被當成定量結果。正解：一律附 count matrix 的關鍵元素。
- 用 PCA 取代 TICA。 後果：慢過程被高變異的柔軟區蓋掉，估出來的時間尺度系統性偏短。正解：用 TICA，並用 VAMP-2 比較特徵組。
- 憑一張 ITS 圖就宣稱 Markovian。 後果：粗劃分假造出來的 plateau 被誤認為收斂。正解：加做 CK test，並附信賴帶。
- 對幀做 bootstrap 估誤差。 後果：假精度，誤差棒小一個數量級。正解：對獨立軌跡 bootstrap，或用 Bayesian MSM。
- 不查連通性。 後果：三成的幀被靜默丟掉，π 描述的是子系統。正解：報 active set 佔比。
- 把補建區段留在 feature 裡。 後果：無實驗約束的人工運動主宰 TICA 空間。正解：排除，並在方法段寫明排除了哪些殘基。
- 用 τ 大於單條軌跡長度。 後果：根本估不出來，或估出無意義的結果。正解：軌跡長度要先滿足 τ 的需求，這是規劃階段就要決定的事。
- 把 MSM 當成便宜的 workflow 主幹。 後果：磁碟先爆，接著是沒有一個變體有足夠統計量。正解：定位成 diagnostic 與 CV 探勘工具，或只在截短系統上做。
11Further Reading讀哪幾篇
12Summary三句話
Markov state model 把大量短軌跡縫成一個轉移矩陣，從中同時解出平衡族群與狀態之間的轉換速率，因此是少數不必押注 collective variable 就能碰到動力學的框架，而且它「很多小 job」的計算形狀正好適合只能少量並行的排隊環境。它的雙面刃在於便宜的是 job 的形狀而不是總取樣量，要對 HNH 檢查點建出能通過 Chapman-Kolmogorov test 的模型，所需 aggregate 比目前的 4 replica × 約 87 ns 高出好幾個數量級，在 535,265 原子的完整三元上並不可行。所以它在這個專案的正確定位是診斷工具與 CV 探勘工具，用 TICA 找出資料自己認定的慢方向再送去做增強取樣，並把「還差多少取樣」量化成一張 implied timescale 圖寫進 open questions。