> ## Content Index
> Fetch the complete content index at: https://taiwanding.com/llms.txt
> Use this file to discover other available public pages before exploring further.

# 從 DMS 到機器學習：預測藥物是否能拯救突變蛋白的完整 Pipeline
- URL: https://taiwanding.com/cong-dms-dao-ji-qi-xue-xi-yu-ce-yao-wu-shi-fou-neng-zheng-jiu-tu-bian-dan-bai-de-wan-zheng-pipeline/
- Published: 2026-04-06T10:17:19.000Z
- Updated: 2026-04-06T10:31:01.000Z
- Author: Kevin Chen
- Tags: computational biology

## 這篇文章在講什麼？

想像你有一個蛋白質，因為基因突變而「摺疊失敗」——它無法正確摺疊，所以被細胞當作垃圾丟掉，無法到達它應該工作的位置（例如細胞表面），這就是很多罕見疾病的成因。

現在有一種叫做 **藥理伴護蛋白（Pharmacological Chaperone, PC）** 的小分子藥物，它可以跟蛋白質結合、幫助它正確摺疊、讓它順利被運送到該去的地方，但問題是：**不是每個突變都能被救回來**。

那麼，能不能用機器學習來預測「哪些突變可以被藥物拯救」？

這就是這篇論文在做的事，在 V2R 內部交叉驗證中，最佳模型 AUROC 為 0.751；將模型凍結後直接遷移到另一個蛋白質 Rhodopsin 時，AUROC 達到 0.843。

## Part 1：背景知識-給剛入門的你

### 什麼是 GPCR？

GPCR（G protein-coupled receptor，G 蛋白偶聯受體）是人體中最大的受體蛋白家族，大約有 800 多個成員，它們橫跨在細胞膜上，負責接收外界訊號（荷爾蒙、神經傳導物質、光等等），然後把訊號傳進細胞內部。

你可能聽過的 GPCR 包括：多巴胺受體（精神科藥物的標靶）、腎上腺素受體（氣喘藥物的標靶）、以及本研究的主角——**V2R（血管加壓素 2 型受體）** 和 **Rhodopsin（視紫紅質）**。

V2R 負責腎臟的水分回收，V2R 突變會導致腎性尿崩症（nephrogenic diabetes insipidus）——患者的腎臟無法正常濃縮尿液。

Rhodopsin 是視網膜中負責感光的受體，Rhodopsin 突變會導致色素性視網膜炎（retinitis pigmentosa），最終可能失明。

### 什麼是 Deep Mutational Scanning（DMS）？

傳統上，研究一個突變的效果需要一個一個做實驗，但 DMS 技術可以**一次測量一個蛋白質上幾乎所有可能的單點突變**的效果。

以 V2R 為例：V2R 有 371 個胺基酸位置，每個位置可以突變成 19 種其他胺基酸，所以理論上有 371 × 19 = 7,049 個可能的錯義突變，DMS 實驗可以一次測量所有這些突變的表面表現量，不過經過原始資料的品質控制（排除 nonsense 突變、同義突變、以及缺少基線或藥物處理數值的變異），最終納入分析的是 6,649 個錯義突變。

2025 年的兩篇重要研究提供了本論文的資料基礎：

- **Mighell & Lehner (2025)**：測量了 V2R 所有突變在有/無 tolvaptan（一種 PC）處理下的表面表現量
- **Manian et al. (2025)**：測量了 Rhodopsin 所有突變在有/無 YC-001（另一種 PC）處理下的運輸效率

### 什麼是 Pharmacological Chaperone（藥理伴護蛋白）？

想像蛋白質摺疊就像摺紙：正常的蛋白質可以順利摺成正確的形狀，但突變可能讓某些步驟出錯，導致蛋白質被困在錯誤的中間狀態。

PC 就像一個「模具」——它跟蛋白質結合後，穩定了正確的構象，幫助蛋白質完成摺疊並通過細胞的品質管制系統。

在已上市的藥物中，**Migalastat**（治療法布瑞氏症，GLA 突變）是最典型的 pharmacological chaperone 範例，在更廣義的蛋白質穩定化與 corrector 脈絡中，也常提到 **Elexacaftor/Tezacaftor**（CFTR correctors，治療囊性纖維化）與 **Tafamidis**（transthyretin stabilizer，治療類澱粉心肌病），這些藥物的作用機制不完全相同，但都涉及透過小分子穩定蛋白質構象來恢復功能。

## Part 2：把生物問題轉成 ML 問題

### 定義標籤：什麼是「有缺陷」？什麼是「被拯救」？

這是整個 pipeline 最關鍵的設計決策之一，我們需要把連續的實驗數值轉換成二元標籤。

**Step 1：定義「有缺陷」（defective）的突變**

V2R 的 DMS 資料中，每個突變都有一個 baseline surface expression score，野生型（正常）的分數被標準化為 1.0，完全失去功能的 nonsense 突變約為 0。

我選擇 **0.7** 作為閾值：baseline < 0.7 的突變被定義為「有缺陷」。

這個閾值的來源：V2R 同義突變（synonymous variants，不改變蛋白質序列的突變）的中位數是 1.001，0.7 大約是中位數減去 2 個標準差，這代表表現量顯著低於正常的突變。

**Step 2：定義「被拯救」（rescued）**

在有缺陷的突變中，如果加了藥物後表現量回到 0.7 以上，就標記為「rescued」（label = 1），否則就是「not rescued」（label = 0）。

```
V2R：6,649 個錯義突變
  → 2,133 個有缺陷（baseline < 0.7）
    → 1,477 個被拯救（69.2%）
    → 656 個未被拯救（30.8%）

```

> 閾值選擇會影響結果嗎？我在論文中做了 threshold sensitivity analysis，把閾值從 0.5 測到 0.8，AUROC 範圍是 0.734–0.753，結論不變，這種 robustness check 在投稿時幾乎是必做的。

### 選擇 Features：我們知道哪些跟「能不能被救」有關的資訊？

Feature engineering 是這個專案的核心，我把 features 分成兩層：

**Layer 1（蛋白質無關的，可跨蛋白質遷移）：**

| Feature             | 它在衡量什麼         | 來源                |
| ------------------- | -------------- | ----------------- |
| AlphaMissense score | 這個突變有多「致病」     | DeepMind 的預訓練模型   |
| RaSP ΔΔG            | 預測的熱力學穩定性變化    | 結構預測工具            |
| ThermoMPNN ΔΔG      | 另一個穩定性預測       | 深度學習工具            |
| ESM1b score         | 蛋白質語言模型的突變效果評分 | Meta 的蛋白質語言模型     |
| delta\_hydro        | 疏水性變化          | Kyte-Doolittle 量表 |
| Baseline expression | 未加藥時的表現量       | DMS 實驗數據          |

**Layer 2（蛋白質特異的，只用於機制分析）：**

| Feature       | 它在衡量什麼       |
| ------------- | ------------ |
| 到 AVP 結合位點的距離 | 突變離藥物結合口袋有多遠 |
| 到 Gs 偶聯界面的距離  | 突變離訊號傳導界面有多遠 |

> **為什麼分兩層？** 因為 Layer 2 的 features 是 V2R 特有的——Rhodopsin 用的是完全不同的藥物和結合位點，如果我們想做跨蛋白質預測，就只能用 Layer 1 的 features，這個「可遷移性」的設計從一開始就是刻意的。

## Part 3：模型訓練與評估

### Progressive Baseline：一步一步加 feature

很多論文只報告最終模型的效能，但這樣你無法理解「每個 feature 到底貢獻了多少」，我採用 **progressive feature addition** 的方式：

```
B0: 隨機猜測                    → AUROC 0.500
B1: AlphaMissense               → AUROC 0.575  (+0.075)
B2: + RaSP, ThermoMPNN          → AUROC 0.648  (+0.073)
B3: + ESM1b                     → AUROC 0.654  (+0.006)
B4: + delta_hydro               → AUROC 0.676  (+0.022)
B5: + baseline expression       → AUROC 0.751  (+0.075)

```

**關鍵觀察：**

1. **AlphaMissense 單獨用效果很差（0.575）**， 這告訴我們一件重要的事：「致病性」≠「可拯救性」，一個突變可以是致病的，但我們不知道它能不能被藥物救回來。
2. **穩定性預測器（RaSP, ThermoMPNN）貢獻很大（+0.073）。** 蛋白質的熱力學穩定性確實跟能不能被 PC 拯救有關。
3. **ESM1b 幾乎沒有額外貢獻（+0.006）**， 這可能是因為 ESM1b 的訊號已經被 AlphaMissense 和穩定性預測器捕捉了。
4. **Baseline expression 是最強的單一 feature（+0.075）**， 越嚴重的缺陷越難被拯救，這很直覺。

> **AUROC 是什麼？**  
>  
> AUROC（Area Under the Receiver Operating Characteristic curve）衡量的是分類器區分兩個類別的能力，0.5 代表跟隨機猜一樣，1.0 代表完美區分。0.75 大約代表「還不錯，有實用價值」。  
>  
> 你可以這樣想像：如果隨機抽一個被拯救的突變和一個未被拯救的突變，模型正確判斷誰是誰的機率就是 AUROC，0.751 代表大約 75% 的機率能判斷正確。

### 為什麼用 Random Forest？

我同時訓練了 Random Forest（RF）和 Logistic Regression（LR）做比較，RF 在 4/5 個 feature set 上都優於 LR，特別是在純計算 features（B4）時差距最大（+0.100）。

這代表 features 之間存在**非線性交互作用**，例如：一個中度不穩定的突變如果在藥物結合位點附近，可能就救不回來；但同樣中度不穩定的突變如果在其他位置，就可能被救回來，RF 可以捕捉這種條件式的規則，LR 不行。

```python
# 核心訓練程式碼（簡化版）
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import StratifiedKFold, cross_val_predict

rf = RandomForestClassifier(
    n_estimators=200,
    class_weight='balanced',  # 處理類別不平衡
    random_state=42
)

cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=42)
y_pred = cross_val_predict(rf, X, y, cv=cv, method='predict_proba')[:, 1]  # 取正類（rescued）的機率

```

> **`class_weight='balanced'` 很重要。** 我們的資料是 69.2% rescued vs 30.8% not-rescued，如果不處理這個不平衡，模型會傾向於把所有東西都預測成 rescued，`balanced` 會自動根據類別比例調整權重。

### Positional Leakage Test：你的模型真的學到了什麼？

這是一個**很多人會忽略但非常重要的驗證步驟**。

問題是這樣的：蛋白質上不同位置有不同的結構角色，如果 cross-validation 的 fold 裡，訓練集和測試集包含了同一個位置的不同突變，模型可能只是學到了「哪些位置容易被救」，而不是真正學到了 variant-level 的特徵。

**Leave-positions-out cross-validation** 把整個 residue position 作為一個單位來分 fold——如果某個位置的突變在測試集中，那這個位置的所有突變都不會出現在訓練集中。

結果：

```
Standard 5-fold CV:          AUROC = 0.751
Leave-positions-out 10-fold: AUROC = 0.745
差距：0.006

```

差距幾乎為零，代表模型確實在學習 variant-level 的特徵，而不是在背位置資訊。

## Part 4：機制分析——為什麼有些突變救不回來？

這是我最喜歡的部分，因為這不只是數字，而是**可以用生物學解讀的分析**。

根據特徵分布的差異，656 個未被拯救的突變可以區分成兩組，它們的特徵模式支持兩類主要 failure mode 的解釋：

### 第一類：藥物結合位點被破壞（NR-binding, n=150）

這些突變位於 tolvaptan 的結合位點（由原始論文通過統計分析鑑定的 32 個位置）。

特徵：

- 中度的 baseline expression（0.184）
- 小的疏水性變化（-0.46）
- 低的不穩定性分數

**解釋：** 蛋白質的整體結構沒有被嚴重破壞，但藥物結合的口袋被改變了，藥物進不去，自然就救不了。

### 第二類：蛋白質太嚴重地被破壞（NR-non-binding, n=506）

這些突變不在結合位點，但就是救不回來。

特徵：

- 極低的 baseline expression（0.064）
- 大的疏水性變化（-3.09）
- 高的不穩定性分數

**解釋：** 蛋白質的摺疊被嚴重破壞了，破壞程度超過了藥物所能提供的穩定化能力。就像一棟已經塌了一半的房子，光靠加一根柱子是撐不起來的。

> **這個觀察為什麼值得注意？** 如果這個解釋成立，它預測了一件事：不同的藥物-蛋白質組合，這兩種失敗模式的比例會不同，如果藥物的結合位點很深、很保守，第一類失敗就會比較少；如果結合位點很淺、很靈活，第一類失敗就會比較多，這是可以在未來實驗中驗證的預測。

## Part 5：跨蛋白質遷移——外部驗證

### 為什麼這一步很重要？

如果模型只能預測 V2R 的突變，那它的價值有限——因為 V2R 的 DMS 資料已經有實驗結果了，模型的真正價值在於：**能不能預測我們還沒有實驗資料的蛋白質？**

我把 V2R 訓練好的模型**完全凍結（frozen）**，不做任何調整，直接拿去預測 Rhodopsin 的突變能不能被 YC-001 拯救。

注意這兩個系統有多不同：

- 不同的蛋白質（V2R vs. Rhodopsin）
- 不同的藥物（tolvaptan vs. YC-001）
- 不同的實驗方法（flow cytometry vs. transcriptional reporter）
- 不同的疾病（腎性尿崩症 vs. 色素性視網膜炎）

### 只用四個 feature 的原因

六個 Layer 1 features 中，有兩個（RaSP 和 ThermoMPNN）因為技術問題無法為 Rhodopsin 計算：

- RaSP 的 GitHub repository 已經不公開了
- ThermoMPNN 有無法解決的 dependency 衝突

所以跨蛋白質預測只用了四個 features：AlphaMissense、baseline expression、delta\_hydro、ESM1b。

> **提醒：** 工具不可用是真實研究中常見的問題，誠實報告這個限制比假裝它不存在更重要，這也是為什麼我在論文的 Limitations 部分詳細說明了這一點。

### 結果：

```
Model A (AlphaMissense only):                    AUROC = 0.677
Model B (AlphaMissense + baseline):              AUROC = 0.788
Model D (+ hydrophobicity + ESM1b):              AUROC = 0.843

```

0.843 的 AUROC 甚至高於 V2R 內部的 0.751，不過要注意：**兩個數字不能直接比較**，因為實驗系統不同、標籤的雜訊特性不同，Rhodopsin 使用的 composite score 是一個 meta-analytic 指標，可能比 V2R 的單一實驗量測更乾淨。

### Permutation Test：確認結果不是運氣

我做了 1,000 次 label permutation test：把 Rhodopsin 的標籤隨機打亂 1,000 次，每次都重新計算 AUROC。

null distribution 的平均值是 0.500（跟預期一致），而我們觀察到的 0.843 遠遠超過所有 permuted AUROC（p < 0.001）。

```python
# Permutation test（簡化版）
import numpy as np
from sklearn.metrics import roc_auc_score

# y_pred_rho 為 frozen V2R model 對 Rhodopsin 預測的正類（rescued）機率
observed_auroc = roc_auc_score(y_rho, y_pred_rho)
n_permutations = 1000
permuted_aurocs = []

for _ in range(n_permutations):
    y_perm = np.random.permutation(y_rho)
    perm_auroc = roc_auc_score(y_perm, y_pred_rho)
    permuted_aurocs.append(perm_auroc)

p_value = np.mean(np.array(permuted_aurocs) >= observed_auroc)

```

## 附錄 A：Feature Engineering 實戰細節（進階）

> 以下內容適合有 Python 和 ML 經驗、想深入了解或重現此研究的讀者，如果你只想了解研究的整體思路和結果，可以跳過。

### ESM1b Score 的計算

ESM1b 是 Meta 開發的蛋白質語言模型，類似 NLP 中的 BERT，但訓練在蛋白質序列上。它可以給每個位置的每個胺基酸一個「合理性」分數。

```python
# ESM1b log-likelihood ratio 計算
import torch
import esm

# 載入模型
model, alphabet = esm.pretrained.esm1b_t33_650M_UR50S()
batch_converter = alphabet.get_batch_converter()
model.eval()

# 準備序列
data = [("rhodopsin", wild_type_sequence)]
batch_labels, batch_strs, batch_tokens = batch_converter(data)

# 取得每個位置每個胺基酸的 log probability
with torch.no_grad():
    results = model(batch_tokens, repr_layers=[33])
    logits = results["logits"]  # shape: [1, seq_len, vocab_size]

# 計算 log-likelihood ratio
# LLR = log P(mutant | context) - log P(wildtype | context)
log_probs = torch.log_softmax(logits, dim=-1)

```

> **解釋：** 你可以把 ESM1b 想成一個「讀過」幾百萬條蛋白質序列的模型，當你問它「在這個位置放一個 Glycine 合理嗎？」，它會根據周圍的胺基酸序列給你一個分數，如果這個位置通常是疏水性胺基酸，你放一個親水性的胺基酸，分數就會很低。

### Hydrophobicity Change 的計算

這是最簡單的 feature，但貢獻不小。

```python
# Kyte-Doolittle hydrophobicity scale
kd_scale = {
    'A': 1.8, 'R': -4.5, 'N': -3.5, 'D': -3.5, 'C': 2.5,
    'Q': -3.5, 'E': -3.5, 'G': -0.4, 'H': -3.2, 'I': 4.5,
    'L': 3.8, 'K': -3.9, 'M': 1.9, 'F': 2.8, 'P': -1.6,
    'S': -0.8, 'T': -0.7, 'W': -0.9, 'Y': -1.3, 'V': 4.2
}

delta_hydro = kd_scale[mutant_aa] - kd_scale[wildtype_aa]

```

一個大的負 delta\_hydro 通常代表疏水性胺基酸被替換成親水性胺基酸——如果這發生在蛋白質的疏水核心，就可能嚴重破壞摺疊。

## 附錄 B：重現這個研究（進階）

### 環境設置：

```bash
# 建議使用 conda 環境
conda create -n rescue python=3.10
conda activate rescue

# 安裝核心套件
pip install pandas numpy scikit-learn openpyxl matplotlib seaborn

# 如果要計算 ESM1b features
pip install fair-esm torch

```

### 資料取得：

1. **V2R 資料：** Mighell & Lehner (2025) 的 Supplementary Table S1（DOI: 10.1038/s41594-025-01659-6）
2. **Rhodopsin 資料：** Manian et al. (2025) 的 GitHub（https://github.com/octantbio/rho-dms）
3. **分析程式碼：** https://github.com/pcc402-art/project-rescue

### 重現步驟：

```bash
git clone https://github.com/pcc402-art/project-rescue.git
cd project-rescue

# Phase 1: 基礎模型訓練
python phase1_baseline.py

# Phase 3: 完整重跑（含所有 robustness checks）
python phase3_v2_full_rerun.py

# Phase 5: 最終投稿版分析
python phase5_presubmission.py

# 生成圖表
python figures_v2.py

```

## 延伸閱讀

如果你對這個領域有興趣，以下資源可能有幫助：

- **AlphaMissense：** Cheng et al. (2023) Science. 蛋白質致病性預測的代表性工作。
- **ESM 蛋白質語言模型：** Brandes et al. (2023) Nature Genetics. 用語言模型預測突變效果。
- **DMS 入門：** Fowler & Fields (2014) Nature Methods. Deep Mutational Scanning 的經典綜述。
- **PC 藥物開發：** Beerepoot et al. (2017) Pharmacological Research. 藥理伴護蛋白的綜述。

---

## 關於作者

我是Po-Chun Chen，獨立研究者，研究方向為計算生物學與資訊安全，這是我以獨立研究者身份發表的第一篇論文。

如果你有問題、想討論，歡迎透過 pcc402@nyu.edu 聯繫我。

> **論文連結：** Chen, P.-C. (2026). Computational prediction of pharmacological chaperone rescuability in GPCR missense variants with cross-protein transfer from V2R to Rhodopsin. Zenodo. [https://doi.org/10.5281/zenodo.19434570](https://doi.org/10.5281/zenodo.19434570?ref=taiwanding.com)  
>  
> **程式碼：** [https://github.com/pcc402-art/project-rescue](https://github.com/pcc402-art/project-rescue?ref=taiwanding.com)