鹽橋分析Salt bridge analysis
帶電殘基與 DNA 磷酸骨架的靜電握手，蛋白與核酸介面的主力。
01一句話只有一句
帶正電殘基的氮與 DNA 磷酸的氧靠得夠近，就算一座鹽橋。
02Why should I care?我為什麼要讀這頁
因為 eSpCas9(1.1) 這個變體，本質上就是一個「把三座鹽橋拆掉」的實驗。
Slaymaker 等人的做法是：找出 Cas9 上那條容納 non-target 股（R-loop 形成時被頂開、沒有跟 sgRNA 配對的那條 DNA）的正電溝槽，把裡面 32 個帶正電的殘基逐一換成 alanine 掃過一遍，讓靜電握手變弱，再量專一性有沒有變好。最後定案的 eSpCas9(1.1) 是 K848A / K1003A / R1060A 三重突變（更早的 1.0 版是 K810A / K1003A / R1060A）。這不是隱喻，這就是字面上的設計邏輯。你的 workflow 如果要排序變體，而已發表的高保真變體有一整個家族是靠「刪電荷」設計出來的，那麼數鹽橋就是最直接對準這條設計邏輯的計算讀數，而且它幾乎不用錢，軌跡跑完就有。
用你熟悉的語言講兩件事。
第一，電荷的來源完全不同。在 DFT 裡電荷是自洽解出來的，兩個離子靠近時電子密度會重新分佈，極化與電荷轉移自動包含在裡面。古典力場正好相反：每個原子的部分電荷在你按下執行之前就已經被 RESP fit 決定死了，靠近之後一動也不動。力場裡的一座鹽橋 = 兩堆固定點電荷的庫倫項 + Lennard-Jones，缺掉的極化被折進一個有效介電常數裡。這件事有一個可預測的後果，等一下的 When NOT 會講。
第二，這個 metric 的統計性質，跟你在 VASP 結構裡用截斷距離數配位數是同一回事。它是一個閾值下的計數：某個 pair 在 3.9 Å 就算一座，在 4.1 Å 就不算，在閾值附近會來回跳。你在算配位數時對這種 step function 的脾氣應該很熟。差別只在於，MD 還多了時間軸，所以你得到的是一條會抖的計數時間序列，而不是一個靜態數字。
03What這是什麼
一座鹽橋（salt bridge，也叫 ion pair）是一對反號的可離子化基團靠到某個幾何判準之內。判準不只一種，這是所有麻煩的來源。
在 Cas9 這種蛋白與核酸的介面上，兩邊的成員幾乎是固定的：
- 陽離子側：Arg 的 NE / NH1 / NH2、Lys 的 NZ、質子化的 His 的 ND1 / NE2。
- 陰離子側：DNA 與 RNA 骨架磷酸的非橋接氧 OP1 / OP2（每個核苷酸一個磷酸，帶 −1 e），以及 Asp 的 OD1 / OD2、Glu 的 OE1 / OE2。
最常用的判準是重原子距離：
| 判準 | rcut | 出處與用法 | 脾氣 |
|---|---|---|---|
| 重原子 N–O 最短距離 | 4.0 Å | Barlow & Thornton 的 ion-pair 分析沿用下來的傳統值，最常見 | 寬鬆，會收進一部分幾乎沒有作用的 pair |
| 電荷基團質心距離 | 4.0 Å | Kumar / Nussinov 一系，要求帶電基團整體靠近 | 對側鏈方向較不敏感，數值較保守 |
| N–O 距離 + 氫鍵幾何 | 3.2 Å（VMD saltbr 預設） |
要求同時滿足氫鍵條件 | 最嚴格，數出來的量明顯較少 |
| 氫鍵 3.5 Å + 離子 4.0 Å 併用 | 兩段 | HBPLUS 一系 | 區分「有氫鍵的鹽橋」與「純靜電接近」 |
流通中的截斷值大致落在 3.2 到 4.5 Å 這個區間。同一條軌跡換一個判準，總數可以差好幾成。 所以「我們觀察到 N 座鹽橋」這句話，在沒有寫出判準之前是沒有意義的——它不是一個測量值，是一個測量值加上一個你沒公開的閾值。
還有第三種狀態，多數人的腳本會直接漏掉：solvent-separated ion pair（SSIP）。兩個電荷之間夾著一個水分子，重原子距離落在大約 5 到 6 Å。它在水裡是一個真實而且常常比直接接觸更有利的狀態，但任何 4 Å 的截斷都會把它記成零。
04Why為什麼重要
因為蛋白與核酸的介面在物理上就是靜電主導的。
做一下最簡單的算術：一段 20 bp 的 R-loop 區域，兩條股加起來大約 40 個磷酸，也就是大約 −40 e 的骨架電荷聚在一個很小的體積裡。Cas9 為了抓住它，演化出一條強正電的通道。這個介面不像蛋白與蛋白的介面那樣由疏水核心撐著，它主要是由一整排靜電握手撐著。鹽橋計數就是這件事最便宜的觀測量。
但這裡有一個對你的專案至關重要的不對稱性，值得單獨拉出來講：
DNA 骨架的磷酸是沿序列均勻排列、而且與鹼基身分無關的陰離子陣列。A、T、G、C 的磷酸長得一模一樣。所以蛋白與磷酸骨架之間的鹽橋，帶的序列資訊接近於零。
推論：鹽橋數量量得到「握得多緊」（affinity），量不到「握得多準」（specificity）。專一性的來源是 R-loop 內部的鹼基配對，以及 HNH 的構形檢查點，不是骨架的靜電。
這也正好解釋了為什麼 eSpCas9 的設計會有效，而且原作者的說法值得逐字記住：中和 non-target 股溝槽裡的正電荷，是為了削弱蛋白對 non-target 股的抓握，讓被頂開的那條 DNA 更容易跟 target 股重新配回去，於是要維持 R-loop 就得靠更嚴格的 Watson-Crick 配對。換句話說，它拆掉的不是「辨認鹼基」的作用，而是「把兩條 DNA 撐開」的那隻手；辨識工作被丟回給鹼基配對本身。鹽橋在這個故事裡的角色是設定一場競爭的門檻——R-loop 撐開 vs 雙股重新閉合——不是熱力學上的辨識器。
05How怎麼做
演算法只有三步，跟 contact map 幾乎一樣，差別在於原子選擇。
- 建兩個 AtomGroup。 陽離子側與陰離子側各一個，用明確的原子名而不是殘基名。MDAnalysis 的寫法大概長這樣：
cat = u.select_atoms("(resname ARG and name NE NH1 NH2) or "
"(resname LYS and name NZ) or "
"(resname HIP and name ND1 NE2)")
ani = u.select_atoms("(nucleic and name OP1 OP2) or "
"(resname ASP and name OD1 OD2) or "
"(resname GLU and name OE1 OE2)")
-
逐幀算距離矩陣、數落在截斷內的 pair。
MDAnalysis.analysis.distances.distance_array加上週期性邊界（box=u.dimensions）。注意去重：一個 Arg 有三個氮，同一座鹽橋不要數三次，正確做法是先對 residue pair 取最短距離，再判斷。這一步是最常見的 bug，忘記去重會讓 Arg 的貢獻被系統性放大。最容易漏掉的那一步如果實作是直接對整個距離矩陣取
(d < cutoff).sum()，而第二個軸是所有（陽離子 N × 磷酸 O）原子對，那數到的就不是「幾座鹽橋」，是「幾個落在截斷內的原子對」。三件事跟著改變：一，Arg（3 個氮）對 Lys（1 個氮）被系統性加權；二，一個磷酸有兩個非橋接氧，同一個接觸最多可以被數到 6 次；三，這個數字不能跟任何文獻的「salt bridge count」比。好消息是這不影響跨 replica 的 CV（同一個定義套在每一條軌跡上），也不影響 WT 對變體的相對比較。壞消息是 methods 必須誠實寫成「cationic-N／phosphate-O atom pairs within 0.4 nm」而不是「salt bridges」，否則就是一個會被抓到的用詞錯誤。要拿到真正的殘基級計數，就在算完距離之後對 (residue, residue) 取最小值再套閾值。
另外兩個一定要寫進 methods 的定義：陽離子側算不算 His——若只取 Arg 的 NE / NH1 / NH2 與 Lys 的 NZ、His 一刀切掉，那是個可辯護的選擇，但下面那段 His 質子化的討論就等於被跳過了，要說；陰離子側取不取 sgRNA 的磷酸氧——如果同一條流程裡 BSA 把 sgRNA 算進核酸側、鹽橋卻只取兩條 DNA 股，兩個定義就不一致，混著講會出事。
-
輸出兩種東西。 一是總數的時間序列（聚合量，穩）；二是逐 pair 的 occupancy，也就是每一對出現在多少比例的幀裡（逐點量，不穩，見 When NOT）。
現成工具：
- VMD 的
saltbr插件：預設 3.2 Å，會直接吐出每一對的時間序列檔案。方便，但預設值偏嚴，而且它的預設定義未必是你要寫進 methods 的那個。 - GROMACS 沒有原生的 salt bridge 工具。實務上是
gmx mindist逐對算，或用gmx hbond搭配自訂 index group，常見的判準是 0.32 nm。 - 自己寫 20 行 MDAnalysis 通常比找工具快，而且你會被迫把判準寫清楚，這反而是好事。
兩個一定要先處理好的前置條件：
- 質子化狀態。 Lys 與 Arg 在標準力場裡永遠是 +1，Asp 與 Glu 永遠是 −1，這沒得選。問題出在 His：pH 7 下大部分 His 是中性的，它算不算陽離子，完全取決於你當初
pdb2gmx/tleap給了哪個殘基名（HIS / HID / HIE / HIP）。這件事必須寫進 methods，而且如果某個 His 剛好在介面上，它會直接改變你的結論。 - 離子與 Mg²⁺。 Cas9 的催化需要 Mg²⁺，而古典力場的 Mg²⁺ 參數是出了名的難搞（水交換速率可以差好幾個數量級）。活性中心附近的鹽橋計數要格外保留。你的離子濃度與離子參數（例如 Joung-Cheatham 一系）也要寫出來，因為離子屏蔽會直接影響計數。
成本： 剝掉水之後，87 ns 的軌跡是分鐘等級。相對於產生軌跡的 GPU 時數，這是免費的。
06When什麼時候用
- 比較 WT 與電荷中和型變體（eSpCas9(1.1)、以及任何把 Arg / Lys 換掉的設計）在同一個底物上的握力差異。 這是鹽橋計數最站得住腳的用途，因為擾動的性質與 metric 的性質完全對齊。
- 定位蛋白抓在 R-loop 的哪一段。 把鹽橋按 DNA 位置作圖，可以看出 PAM 端與 PAM-distal 端的握力分布，這對理解錯配敏感區有直接幫助。
- 當作聚合打分軸的一個分量，和 BSA、總接觸數並列。三個量的失敗模式不同，一起看比單看一個穩。
- 抓建模產物。 一段建出來的 loop 如果掛著一座在晶體結構裡不存在的鹽橋，通常代表那段建錯了。這是很便宜的 sanity check。
- 檢查系統設定沒出包。 鹽橋總數在 equilibration 之後應該是平的。如果它單調上升或崩塌，先回去看 system setup，不要急著解釋生物學。
07When NOT什麼時候別用
-
不要拿鹽橋總數當專一性 ranker。 為什麼會壞：磷酸骨架的電荷與序列無關，一個帶三個錯配的 off-target 雙股，呈現給蛋白的磷酸陣列跟 on-target 幾乎一模一樣。這個 metric 在物理上就看不見鹼基身分。怎麼看出它壞了：把 on-target 與 mismatched 底物各跑一組，如果兩者的計數差落在你的跨 replica 誤差之內（你的 CV 是 9.5%），那它對這個效應量沒有解析力。先算出「最小可偵測效應」再決定要不要跑，不要跑完才發現。
-
不要從少數 replica 報逐 pair 的 occupancy。 為什麼會壞：逐 pair occupancy 是逐點量，而你自己的可靠性稽核已經把逐點量判死了（逐位點距離 CV 38 到 46%、groove 接觸計數 145.5%）。同一對殘基在一條 replica 佔 80%、在另一條佔 10%，是完全正常的。怎麼看出它壞了：把 occupancy 逐 replica 算出來看散布；只要有任何一對的跨 replica 全距超過它的平均值，整張 occupancy 表就不能當結論用，只能當敘事。
-
不要把鹽橋計數當結合自由能的代理。 為什麼會壞：它只數了分子項，完全沒算脫水的代價。要形成一個接觸離子對，兩個電荷都得先把自己的水合殼剝掉，這個代價在水裡非常大，以至於一個暴露在溶劑中的鹽橋，淨 ΔG 常常接近零甚至不利。接觸多不等於結合強。 怎麼看出它壞了：你的計數會跟 MM-PBSA 或實驗 Kd 幾乎不相關，而且會系統性地把高電荷介面排得太前面。
-
不要跨截斷值、跨力場、跨離子強度比較計數。 為什麼會壞：計數是一個任意閾值的階梯函數；而且固定電荷力場系統性地高估接觸離子對的穩定度（沒有極化、沒有電荷轉移，庫倫項在近距離太硬）。這正是 NBFIX / CUFIX 這類 pair-specific 修正，以及 ECC（electronic continuum correction）做法要處理的問題。ECC 這裡要特別講清楚，因為很容易看錯：它把離子電荷乘上 1/√εel，其中 εel 是水的高頻（電子）介電常數約 1.78，不是靜態的 78——所以縮放因子大約是 0.75，不是 0.11。這個修正的物理意義是把力場缺掉的電子極化用一個均勻介電背景補回來。換一套修正，計數會系統性地移動。怎麼看出它壞了：把同一條軌跡用 3.2 / 4.0 / 4.5 Å 各算一遍；如果你的結論在三個值之間翻面，那個結論是截斷假象，不是物理。
-
不要忽略水橋介入的 pair。 為什麼會壞：一對電荷從直接接觸（約 3.5 Å）滑到 solvent-separated（約 5.5 Å）時，作用其實大部分還在，但 4 Å 的截斷會把它記成「消失了」。一個變體如果只是把介面稍微推開、讓水擠進來，你的計數會誇張地報告「失去了五座鹽橋」，而真實的親和力可能幾乎沒變。怎麼看出它壞了：畫該 pair 的距離直方圖。如果出現雙峰、峰位大約在 3.5 Å 與 5.5 Å，那就是直接接觸與水橋兩態並存的特徵，此時單一截斷的計數必然低估。
-
不要在 on-target 底物上測高保真變體的鹽橋差異。 為什麼會壞：這類變體本來就設計成在 on-target 上跟 WT 一樣。你的專案量到的情況與這個預期一致：eSpCas9(1.1) 的三個 grip 位點，在 on-target 底物上只有 K848 一直貼著 non-target 股（4 條 replica 平均 0.7 nm，全距 0.3 到 1.0），K1003 與 R1060 平均在一到兩奈米外（1.6 與 1.5 nm），所以在 on-target 上看不到差異是預期中的事。但這一條要照本站的規矩讀：這些逐位點距離正是你自己判為 noise tier 的量（跨 replica CV 38 到 46%），那個平均值不是一個穩的數字——R1060 在其中一條 replica 就貼到 0.5 nm，在另一條卻是 2.0 nm。所以不可以從這裡推出「隔那麼遠所以不可能形成鹽橋」；能推的只有「這三個位點的握持狀態在 replica 之間根本不一致」。它足以支持「不要在 on-target 上做這個對照」，不足以支持「K1003 與 R1060 不重要」——按 Slaymaker 的設計邏輯，它們該起作用的時機本來就在 R-loop 還沒完全撐開、或底物有錯配的時候。怎麼看出它壞了：先把 WT 軌跡裡那些殘基到 DNA 的距離逐 replica畫出來；如果同一個位點在不同 replica 之間從 0.5 nm 跳到 2 nm，那就不要用單一平均值下結論，換 mismatched 底物再做。
你的鹽橋計數跨 4 條 replica 的 CV 是 9.5%，落在「reproducible（CV < 10%）」這一級，但它是四個綠燈裡最差的一個（BSA 2.1%、總 protein–DNA 接觸 4.1%、target-DNA 接觸 9.4%、鹽橋 9.5%）。
更該注意的是 CV 本身的不確定度。n = 4 時，樣本標準差的相對誤差大約是 1/√(2(n−1)) = 1/√6 ≈ 41%（常態近似下的經驗法則）。所以 9.5% 應該讀成「大約一成，可能低到 6%，也可能高到 13%」，而不是一個精確數字。它隨時可能因為多跑一條 replica 就掉到 marginal 那一級。
實務上的結論：可以當次要打分軸，不能當唯一的軸，而且不要宣稱能解析小於一到兩成的效應。
08Project Lens跟我的 project 多相關
為什麼是四星而不是五星： 因為它的可靠性剛好卡在門檻線上，而且它在物理上看不見專一性。
- 9.5% 是四個綠燈裡最差的一個，而且如上所述，這個 CV 在 n = 4 下本身就有大約四成的相對誤差。把它放在打分層的主軸上，等於把整條 funnel 的排序架在一個隨時可能掉成黃燈的量上。
- 它與序列無關。 這是結構性的限制，不是取樣不夠可以解決的問題。磷酸骨架長得都一樣，所以鹽橋計數對 on-target 與 off-target 的辨別力先天就很弱。你的整個 workflow 的立場是「專一性由動力學檢查點決定，不由平衡結合決定」，而鹽橋計數量的正是結合那一側。它可以支持那個立場（說明握力如何設定門檻），但它自己不是那條軸。
為什麼不是更低：
- 它是唯一與已發表變體設計邏輯直接對齊的 metric。 eSpCas9 家族就是靠中和正電殘基做出來的。當你要向 reviewer 解釋「為什麼我的 workflow 會挑出這一類變體」，鹽橋計數是一句話就能講完的機制敘事，不像 PCA 的第幾主成分那樣需要三段解釋。
- 成本近乎為零。 在「便宜」是硬約束的專案裡（單張消費級 GPU，加上共用叢集上一份同時只跑得動個位數 job、磁碟配額很緊的小額額度），任何軌跡跑完就順手拿到、不用額外 GPU 時數的量都值得保留。
- 它是聚合量。 你的專案最重要的方法學結論是「打分要用聚合量，不要用逐殘基／逐位點」，鹽橋總數站在對的那一邊。
具體怎麼用： 放在 Gate 2 的次要打分軸，與 BSA 與總接觸數並列，只報總數、不報逐 pair occupancy，而且一定要在 mismatched 底物上跑。動手前先確認一件事：計數的單位是原子對還是殘基對（見 How 第 2 步），在寫進論文之前要嘛改成殘基級去重、要嘛把名字誠實改掉——這件事早改只要幾行，等到圖都畫完再改就會全部重跑。方法段必須寫死判準（建議 4.0 Å 重原子最短距離，也就是 0.4 nm，並附上 3.2 Å 的敏感度檢查）、His 質子化狀態、離子參數與濃度。
09Reviewer Thinkingreviewer 會問什麼
你的鹽橋判準是什麼？換一個判準結論還成立嗎？點開看參考答案
這是第一個會被問的，因為流通中的截斷值從 3.2 到 4.5 Å 都有。標準答案：明寫主判準（例如「陽離子氮與陰離子氧的最短重原子距離 ≤ 4.0 Å，以 residue pair 為單位去重」），並且在補充材料放一張敏感度圖，把 3.2 / 4.0 / 4.5 Å 三個值的結果並列。如果你的變體排序在三個值之間不變，這一問就結束了；如果會變，你必須自己先講。
你報 9.5% 的 CV，但那是 4 條 replica 算出來的。這個 CV 本身的誤差是多少？點開看參考答案
會問這一題的 reviewer 是統計比較硬的那種，而且這一題非常致命，因為多數人答不出來。誠實的答案：n = 4 時樣本標準差的相對誤差約 1/√(2(n−1)) ≈ 41%，所以 9.5% 的合理區間大約是 6% 到 13%。正確的立場是把 9.5% 當成「大約一成」的量級陳述，並且明說任何小於這個量級的效應都不在解析範圍內。先自己說出來，比被問出來好。
你怎麼處理水橋介入的離子對？點開看參考答案
直接切入：4 Å 截斷會把 solvent-separated ion pair 記成零，所以純計數會低估。做法是對每一對關鍵殘基畫距離直方圖，報告雙峰的相對權重，而不是只報一個計數。如果你的比較涉及「介面被推開」這類擾動，這一步是必要的，否則你會把「水擠進來」誤讀成「作用消失」。
固定電荷力場高估接觸離子對，這對你的 WT 與變體比較有什麼影響？點開看參考答案
這是最有深度的一問。回答要分兩層：偏差的絕對值很大（沒有極化、沒有電荷轉移），所以絕對的鹽橋壽命與計數都不可信；但如果 WT 與變體用完全相同的力場、水模型、離子參數與濃度，偏差在很大程度上會抵消，相對排序仍然可用。要更保險，就補跑一組帶 NBFIX / CUFIX 修正或離子電荷縮放的對照，看排序會不會翻。要注意這個抵消是假設，不是定理，該講的時候要講成假設。
磷酸骨架的電荷與序列無關，這個 metric 怎麼可能分辨 on-target 與 off-target？點開看參考答案
它不能，而且你應該主動承認。正確的定位是：鹽橋計數量的是握力，握力設定的是「錯配底物能不能撐到 HNH docking」的門檻，所以它是機制鏈上的一環，不是終點的判別器。判別要靠 HNH 距離 這種對準構形檢查點的量。把這一段講清楚，反而會讓 reviewer 相信你知道自己在量什麼。
10Common Mistakes最常犯的錯
- 一個 Arg 的三個氮各數一次。 後果：Arg 的貢獻被系統性放大，Arg 密集的區域看起來假性重要；而且數字不能跟文獻的鹽橋數比。這是直接對距離矩陣取 sum 時最容易發生的預設行為（見 How 第 2 步）。正解：先對 residue pair 取最短距離再判斷，每對殘基最多一座；在改之前，methods 就寫「atom pairs within 0.4 nm」而不是「salt bridges」。
- 報了計數卻沒寫判準。 後果：數字不可重現，跨 paper 完全不能比。正解：判準、去重規則、週期性處理三件事都寫進 methods。
- His 的質子化狀態沒交代。 後果：介面上一個 His 算不算陽離子，會直接改變計數與結論，而且 reviewer 一眼就會看到。正解：列出所有非預設質子化的殘基與依據（pKa 預測工具或結構環境）。如果做法是整批不算 His，那就把這句話直接寫進 methods，並補一句有沒有 His 落在介面上。
- 報逐 pair occupancy 表當作主要結果。 後果：那是逐點量，跨 replica 會亂翻，等於把雜訊排版成表格。正解：只把 occupancy 當定性敘事，主張建立在總數上。
- 用 4 Å 截斷比較「介面被推開」的變體。 後果：把直接接觸滑成水橋誤讀成作用消失，誇大變體效應。正解：畫距離直方圖，報雙峰權重。
- 在 on-target 底物上測高保真變體。 後果：得到 null，然後浪費時間解釋 null。正解：先逐 replica 確認被突變的殘基在 WT 軌跡裡真的在接觸範圍內（K1003 與 R1060 平均在一到兩奈米外，而且逐 replica 從 0.5 跳到 2.3 nm，光看平均值判斷不了），再換 mismatched 底物。
- 把鹽橋數當作結合強度直接排序變體。 後果：忽略脫水代價，與實驗幾乎不相關。正解：要談自由能就去算自由能，鹽橋計數只當機制敘事與次要打分分量。
11Further Reading讀哪幾篇
12Summary三句話
鹽橋是帶正電殘基的氮與 DNA 磷酸氧靠到截斷距離之內的計數，它是蛋白與核酸介面上承重的靜電作用，也正是 eSpCas9 這一類高保真變體被設計出來的地方。在你的 4 條 replica 上它的 CV 是 9.5%，剛好卡在 reproducible 的門檻線上，而且 n = 4 讓這個 CV 本身有大約四成的相對誤差，所以它可以當次要打分軸，但解析不了小於一到兩成的效應。它最大的限制是結構性的：磷酸骨架的電荷與序列無關，所以鹽橋數量得到握得多緊，量不到握得多準。