GEARS 完整教學

Repository: https://github.com/snap-stanford/GEARS Stars: 381 | Forks: 85 | Language: Python | License: MIT 論文: Roohani, Huang, Leskovec. Predicting transcriptional outcomes of novel multigene perturbations with GEARS. Nature Biotechnology, 2023.

1. 專案概覽 (Project Overview)

1.1 專案背景、研究團隊與動機

GEARS（Genetic Effects Analysis using Relational Structure，官方全稱為「Predicting transcriptional outcomes of novel multi-gene perturbations」）出自史丹佛大學 Jure Leskovec 實驗室（SNAP — Stanford Network Analysis Platform），作者為 Yusuf Roohani、Kexin Huang 與 Jure Leskovec。這篇論文於 2023 年發表在 Nature Biotechnology，是 geometric deep learning (GDL; 幾何深度學習) 應用於單細胞擾動預測 (single-cell perturbation prediction; 單細胞擾動預測) 領域的代表作。

背景脈絡是：CRISPR 技術讓實驗室可以做 Perturb-seq（結合 CRISPR 基因擾動與單細胞 RNA 定序 single-cell RNA sequencing (scRNA-seq; 單細胞 RNA 定序) 的高通量篩選技術），一次實驗可以敲除/激活數十到數百個基因並讀出每個細胞的全轉錄體 (transcriptome; 轉錄體) 反應。但問題在於：基因組合的排列組合是天文數字級——假設有 20,000 個基因，兩兩組合的擾動空間就高達 2 億種，實驗室不可能窮舉所有組合做濕實驗 (wet-lab experiment; 濕實驗)。GEARS 想解決的核心問題就是：能否從已知的單基因擾動實驗結果，計算式地預測從未做過的多基因（甚至全新）組合擾動會產生什麼轉錄反應？

這件事重要的原因在於藥物標靶發現 (drug target discovery; 藥物標靶發現)、合成致死 (synthetic lethality; 合成致死) 篩選、細胞重編程 (cell reprogramming; 細胞重編程) 等下游應用，都極度依賴「預先篩掉不值得做濕實驗的組合」，把濕實驗預算集中在計算模型認為最有希望的候選上。

1.2 解決什麼問題、為什麼重要

具體來說 GEARS 解決三類預測任務：

1. 未見過的單基因擾動 (unseen single-gene perturbation)：訓練資料中沒出現過的基因被單獨敲除時，轉錄體會如何變化？
2. 未見過的組合擾動 (unseen combinatorial perturbation)：兩個基因一起被敲除（其中一個或兩個都可能是訓練時沒單獨看過的基因），下游表現量如何變化？
3. 遺傳交互作用 (genetic interaction; GI) 預測：兩基因合併擾動的效果，是加成的 (additive)、協同的 (synergistic; 增強)，還是拮抗的 (epistatic/suppressive; 相互抵銷)？GEARS 提供 GI_predict API 直接輸出這類交互作用分數（例如 magnitude、model fit、equality of contribution 等指標）。

這件事的重要性在於，它把「基因擾動效果預測」從純粹統計外推 (extrapolation) 提升為結合基因調控關係先驗知識 (prior biological knowledge; 生物學先驗知識) 的歸納式學習——這也是它與純黑箱深度學習模型的本質差異。

1.3 在 snap-stanford 生態系中的定位

GEARS 是 SNAP 實驗室「Graph ML for Biology / Med」系列工具鏈中的重要一員。SNAP 的生態系包含：

• PyG (PyTorch Geometric)：GEARS 的底層 GNN 框架依賴，也是同源專案（Leskovec 團隊主導）。
• scGPT、GenePT 等後續 foundation model 研究常把 GEARS 當作 baseline 對照。
• gears_misc（yhr91/gears_misc）：論文圖表重現的配套 repo。
• 下游延伸：許多 perturbation prediction 社群工具（如 scPerturb 資料整理計畫）採用 GEARS 的資料格式與評估協議作為標準基準 (benchmark; 基準)。

GEARS 因此常被視為「Perturb-seq 預測」這個子領域的 參考實作與事實標準基準（de facto standard)，後續論文幾乎都會拿它做對照表格中的一行。

1.4 相關論文引用

1@article{roohani2023predicting,
2  title={Predicting transcriptional outcomes of novel multigene perturbations with gears},
3  author={Roohani, Yusuf and Huang, Kexin and Leskovec, Jure},
4  journal={Nature Biotechnology},
5  year={2023},
6  publisher={Nature Publishing Group US New York}
7}

論文連結：https://www.nature.com/articles/s41587-023-01905-6

驗證用資料集主要來自三篇 Perturb-seq 原始論文：

• Norman et al. 2019：131 個雙基因組合 + 105 個單基因擾動，K562 白血病細胞株。
• Adamson et al. 2016：Unfolded protein response (UPR; 未折疊蛋白反應) 相關的 87 個單基因擾動。
• Dixit et al. 2016：24 個轉錄因子/訊號分子擾動，資料量較小，常用於快速驗證。

2. 核心架構 (Core Architecture)

2.1 系統架構圖

GEARS 的核心洞察是：把「擾動」本身也表示成圖上的節點/嵌入 (embedding; 嵌入)，並用兩張圖分別編碼「基因彼此的共表現關係」與「擾動彼此的功能相似性」，再讓兩者在共享的隱空間 (latent space; 隱空間) 中融合，最後用 gene-specific decoder 逐基因解碼出表現量變化。

flowchart TB
    subgraph input["輸入層"]
        A["Perturb-seq 資料
AnnData (scanpy)"]
        B["對照組 (control) 平均表現量"]
    end

    subgraph graphs["雙圖結構 (Dual Graph Structure)"]
        C["基因共表現圖
G_coexpress
(Pearson correlation KNN)"]
        D["基因本體相似圖
G_go (Gene Ontology)
擾動彼此的功能相似性"]
    end

    subgraph embed["嵌入層 (Embedding Layer)"]
        E["gene_emb: 基因嵌入表"]
        F["pert_emb: 擾動嵌入表"]
        G["emb_pos: 基因位置嵌入
(經 SGConv 編碼共表現圖)"]
    end

    subgraph gnn["雙路 GNN 編碼器"]
        H["Gene Co-expression GNN
(SGConv × num_gene_gnn_layers)"]
        I["GO Similarity GNN
(SGConv × num_go_gnn_layers)"]
    end

    subgraph fuse["融合與解碼"]
        J["pert_fuse MLP
(擾動嵌入 + GO-GNN 輸出)"]
        K["cross_gene_state MLP
(跨基因狀態整合)"]
        L["gene-specific decoder
(indv_w1/b1, indv_w2/b2)"]
    end

    subgraph output["輸出層"]
        M["預測表現量變化
Δexpression per gene"]
        N["不確定性估計
(uncertainty_w, 選用)"]
    end

    A --> C
    A --> D
    B --> E
    C --> G
    D --> H
    D --> I
    E --> J
    F --> J
    G --> J
    I --> J
    J --> K
    K --> L
    L --> M
    L -.-> N

2.2 關鍵模組、類別與資料流說明

GEARS 套件（gears/）由六個核心檔案組成：

資料流的邏輯順序如下：

flowchart LR
    A["原始 AnnData
(.h5ad)"] --> B["PertData.load()
下載/讀取"]
    B --> C["PertData.prepare_split()
simulation/combo_seen 切分"]
    C --> D["PertData.get_dataloader()
建 PyG Data 物件"]
    D --> E["GEARS(pert_data)
初始化模型容器"]
    E --> F["model_initialize()
建 GEARS_Model + 兩張圖"]
    F --> G["train()
MSE + direction-aware loss"]
    G --> H["predict() / GI_predict()
推論新擾動組合"]
    H --> I["evaluate() / deeper_analysis()
Pearson delta, DE genes 分析"]

PertData（pertdata.py） 是整個框架的資料樞紐：它把 scanpy 的 AnnData 物件轉換成模型可用的圖資料結構。它會呼叫 get_similarity_network() 建構兩張圖：

1. G_coexpress：以對照組表現量計算基因兩兩的 Pearson correlation，取 KNN (K-nearest neighbor; K 近鄰) 建圖，代表「共表現關係」。
2. G_go：以 Gene Ontology 註解計算擾動基因彼此的功能相似度，代表「擾動的先驗生物學相似性」——這是 GEARS 相對於純資料驅動模型的關鍵差異化設計。

GEARS_Model（model.py） 是模型的神經網路本體，使用 torch_geometric.nn.SGConv（Simplified Graph Convolution，簡化圖卷積層）作為兩條 GNN 路徑的骨幹。它有兩層 embedding lookup：gene_emb（基因本身的表示）與 pert_emb（擾動事件的表示），並透過 pert_fuse MLP 把兩者融合，再用 cross_gene_state 讓所有基因間的資訊互相流通（近似一個簡化版的 self-attention），最後用 per-gene 的參數 indv_w1/b1、indv_w2/b2 做 gene-specific decoding——這代表每個基因都有專屬的最後一層解碼參數，而不是共享同一組權重，這對應到生物學上「不同基因對擾動的響應曲線本質不同」的假設。

GEARS（gears.py） 是使用者實際互動的頂層類別，封裝了 model_initialize()、train()、predict()、GI_predict()、save_model()/load_pretrained() 等方法，是典型的 Facade 模式 ——把底下複雜的圖建構、資料切分、模型初始化都藏在幾個方法呼叫背後。

2.3 設計哲學與技術選擇分析

GEARS 的設計哲學可以總結成三個關鍵抉擇：

1. 把「擾動」當成一等公民 (first-class citizen) 來嵌入，而非只是資料標籤。 多數傳統方法把「敲除基因 A」視為分類標籤或 one-hot 向量，GEARS 則賦予每個擾動一個可學習的嵌入向量，並讓這個嵌入透過 GO 相似圖與其他擾動共享資訊——這是模型能夠外推到全新基因組合的關鍵，因為即使某個組合從未在訓練資料出現過，只要組成的基因與訓練過的擾動在 GO 空間相近，模型就能合理內插。

2. 雙圖分工而非單一大圖。 共表現圖負責建模「基因-基因」關係（下游會表現出什麼），GO 相似圖負責建模「擾動-擾動」關係（做這件事跟做那件事有多像）。這種關注點分離 (separation of concerns; 關注點分離) 讓模型可以分別調整兩條路徑的層數（num_gene_gnn_layers vs num_go_gnn_layers），也讓生物學可解釋性更好。

3. SGConv 而非更複雜的 GAT/Transformer。 SGConv 是簡化圖卷積，本質上等同於先對鄰接矩陣做 K 次冪次傳播再接一個線性層，計算成本遠低於多頭注意力機制。這反映了團隊在「表達力」與「在數萬基因規模圖上仍可訓練」之間的務實取捨——這正好呼應本文件開頭 CLAUDE.md 提到的「簡化優先」原則：能用簡單機制解決就不用複雜機制。

3. 安裝與設定 (Installation & Setup)

3.1 環境需求

• Python：3.8 – 3.10（官方測試版本；3.11+ 可能有 torch_geometric 相容性問題，需自行確認）
• PyTorch：需先安裝好對應 CUDA 版本的 PyTorch
• PyTorch Geometric (PyG)：必須先手動安裝，pip install cell-gears 不會自動幫你裝好 PyG 的 CUDA 相依套件
• GPU：非必須，但訓練建議至少 8GB VRAM 的 GPU（README 範例用 cuda:8，代表官方是在多 GPU 叢集上跑的）；CPU 也可訓練小型資料集（如 Dixit），只是速度慢很多
• 其他相依套件（來自 requirements.txt）：

 1numpy==1.26.4
 2pandas==2.2.2
 3tqdm==4.66.5
 4scikit-learn==1.5.1
 5torch
 6torch_geometric
 7scanpy==1.10.2
 8networkx==3.3
 9dcor
10scipy==1.14.1

依照全域 CLAUDE.md 的 Python 工具鏈規範，本教學的安裝流程一律以 uv 建立隔離環境，禁止把 pip 當主要安裝方式（僅在 uv pip install 語境下使用 pip 相容介面）。

3.2 完整安裝步驟

 1# 1) 用 uv 建立乾淨的虛擬環境（Python 3.10）
 2uv venv gears-env --python 3.10
 3source gears-env/bin/activate
 4
 5# 2) 先裝 PyTorch（依你的 CUDA 版本調整，這裡以 CUDA 12.1 為例）
 6uv pip install torch --index-url https://download.pytorch.org/whl/cu121
 7
 8# 3) 安裝 PyTorch Geometric 及其相依 wheel（版本需對齊上一步的 torch/CUDA 版本）
 9uv pip install torch_geometric
10uv pip install pyg_lib torch_scatter torch_sparse torch_cluster torch_spline_conv \
11    -f https://data.pyg.org/whl/torch-2.3.0+cu121.html
12
13# 4) 安裝 GEARS 本體（PyPI 套件名為 cell-gears）
14uv pip install cell-gears
15
16# 5) 驗證安裝
17python -c "import gears; from gears import PertData, GEARS; print('GEARS OK')"

重要：第 3 步的 -f URL 中的 torch-X.Y.Z+cuXXX 字串必須跟你實際安裝的 torch 版本完全對應，否則會出現 ImportError: undefined symbol 這類 C extension 版本不匹配錯誤——這是 PyG 生態系最常見的安裝地雷。

3.3 資料集下載與準備

GEARS 內建三個論文資料集的自動下載器，透過 Harvard Dataverse 拉取：

1from gears import PertData
2
3# 建立資料容器，指定本地快取路徑
4pert_data = PertData('./data')
5
6# 下載並載入 Norman et al. 2019 資料集（自動快取，之後重複執行不會重下）
7pert_data.load(data_name='norman')   # 也可以是 'adamson' 或 'dixit'

第一次執行會自動下載壓縮檔並解壓到 ./data/norman/，包含：

• perturb_processed.h5ad：處理過的 AnnData 表現矩陣
• data_pyg/：預先計算好的 PyG 圖快取（v0.1.1 起針對 Norman、Replogle RPE1/K562 essential 提供預先切分好的 dataloader，省去重新建圖的時間）

若要使用自己的資料，需先準備一個 scanpy 的 AnnData 物件，滿足：

• adata.var 必須有 gene_name 欄位
• adata.obs 必須有 condition（記錄擾動基因名，多基因用 + 分隔，如 "CBL+CNN1"；對照組記為 "ctrl"）與 cell_type 欄位

1pert_data.new_data_process(dataset_name='XXX', adata=adata)
2pert_data.load(data_path='./data/XXX')

3.4 常見問題排解

4. 核心概念詳解 (Key Concepts)

4.1 Perturb-seq 資料的樣貌

想像一張巨大的試算表：每一列 (row) 是一個細胞，每一欄 (column) 是一個基因，格子裡的數字是這個基因在這個細胞裡的表現量（mRNA 讀數）。Perturb-seq 實驗會把細胞分成很多組，每組被 CRISPR 敲掉一個或兩個特定基因，然後量測全部細胞的表現矩陣。GEARS 要學的，就是「敲掉基因 X（和 Y）之後，表現矩陣的哪些欄位會系統性地上升或下降」。

4.2 為什麼要用「圖」而非直接用表格資料做迴歸

如果只把「敲掉哪個基因」當成一個 one-hot 分類標籤丟進神經網路，模型完全沒有辦法對「從未見過的基因」做出合理預測——因為 one-hot 向量之間彼此正交，模型學不到「敲掉基因 A 的效果」和「敲掉基因 B 的效果」之間有什麼關聯。

GEARS 的解法是引入兩張圖當作先驗知識的橋樑：

flowchart TB
    subgraph analogy["類比：圖書館推薦系統"]
        A["讀者 A 借過《深度學習》"]
        B["讀者 B 從沒借過任何書"]
        C["書籍相似度圖
(類比 GO 相似圖)"]
        A --> C
        C -->|"《深度學習》與《機器學習》相似度高"| D["預測讀者 B 若借書
大機率也會喜歡《機器學習》"]
        B -.->|"透過相似度圖
間接獲得推薦"| D
    end

同樣道理，即使某個基因從未在訓練資料中被單獨擾動過，只要它跟訓練過的擾動在「基因本體 (Gene Ontology; GO)」語意空間中相近（例如同屬於某個訊號傳導路徑），模型依然可以透過圖上的訊息傳遞 (message passing; 訊息傳遞) 得到一個合理的嵌入表示，進而做出預測。這就是 inductive learning（歸納式學習，能推廣到訓練時未見過的節點） 相對於 transductive learning（直推式學習，只能對訓練時就存在的節點做預測） 的差異，也是 GEARS 論文標題強調「novel」（新穎、未見過的）的技術核心。

4.3 GNN 訊息傳遞的直觀理解

SGConv（Simplified Graph Convolution）的運作可以想成「多輪傳話遊戲」：每一輪，每個節點把自己目前的嵌入向量傳給所有鄰居，並把收到的鄰居訊息加權平均後更新自己的嵌入。重複 K 輪之後，每個節點的嵌入就融合了「K 步以內鄰居」的資訊。GEARS 分別對兩張圖各做這件事：

flowchart LR
    subgraph coexpr["基因共表現圖上的傳播"]
        g1["基因 A"] <--> g2["基因 B
(表現量常一起變動)"]
        g2 <--> g3["基因 C"]
    end
    subgraph go["GO 相似圖上的傳播"]
        p1["擾動 X"] <--> p2["擾動 Y
(功能相近)"]
        p2 <--> p3["擾動 Z"]
    end
    coexpr -.->|"K 層 SGConv 後"| fused["融合嵌入"]
    go -.->|"K 層 SGConv 後"| fused
    fused --> decode["gene-specific decoder
逐基因輸出預測值"]

4.4 Gene-specific decoder：為什麼每個基因要有專屬解碼參數

GEARS 最後一層並非所有基因共用同一組線性層權重，而是 indv_w1（shape 為 [num_genes, hidden_size, 1]）與 indv_w2（shape 為 [1, num_genes, hidden_size+1]）——也就是每個基因都有自己專屬的一組解碼參數。這對應到一個生物學直覺：同樣受到「上游調控變化」的影響，不同基因的響應曲線形狀本質不同（有些基因對任何擾動都很敏感、有些很穩定），共用權重會抹平這種基因異質性 (gene heterogeneity; 基因異質性)。

4.5 資料切分策略：Simulation Split 的意義

PertData.prepare_split(split='simulation', seed=1) 是理解 GEARS 評估協議最容易被忽略、卻最重要的一環。一般機器學習的隨機切分（訓練/測試集隨機分割）在 perturbation prediction 場景下並不合適，原因是：如果訓練集和測試集中都包含「同一個基因」的擾動樣本（只是不同組合），模型很容易靠著記住這個基因本身的效應就得高分，卻沒有真正展現「外推到新組合」的能力。

simulation 切分策略的做法是：刻意把一部分基因/組合完全排除在訓練集之外，讓測試集包含模型從未見過的基因或基因組合，藉此模擬真實部署場景——研究者想預測的往往正是「還沒做過濕實驗的全新組合」。這個切分邏輯由 data_utils.py 中的 DataSplitter 類別實作，會依照使用者指定的 train_gene_set_size（訓練時使用的基因比例）決定要保留多少基因/組合作為「未見過」的測試對象。

flowchart LR
    subgraph all["完整基因組合空間"]
        direction TB
        seen["訓練集
(部分單基因 + 部分組合擾動)"]
        unseen0["測試集 A：
未見過的單基因擾動"]
        unseen1["測試集 B：
兩基因皆訓練過，
但這個組合沒訓練過"]
        unseen2["測試集 C：
其中一個基因
訓練時完全沒出現過"]
    end
    seen -.->|"simulation split
刻意排除"| unseen0
    seen -.->|"simulation split
刻意排除"| unseen1
    seen -.->|"simulation split
刻意排除"| unseen2

論文中把這三種測試情境分開報告指標，是因為它們的難度依序遞增：「兩基因都訓練過，只是沒看過這個組合」相對容易（模型只需要學會怎麼組合已知的兩個嵌入），「其中一個基因完全沒出現過」則需要模型真正利用 GO 相似圖做嵌入外推，難度高出許多。理解這個分層有助於正確解讀模型報告的效能數字，避免把「容易情境的高分」誤讀成「模型在所有情境下都表現優異」。

4.6 不確定性估計 (Uncertainty Estimation)

GEARS 提供選用的 uncertainty=True 模式，額外訓練一個 uncertainty_w MLP head 輸出每個預測值的 log-variance（對數變異數）。這讓使用者不只拿到「預測表現量會變成多少」，還能知道「這個預測有多不確定」——在藥物開發場景中，這對應到「哪些預測結果值得優先送去濕實驗驗證」的決策依據：不確定性低的預測可以直接信任，不確定性高的則需要人工複核或安排濕實驗確認。

5. 使用方式與程式碼範例 (Usage & Code Examples)

5.1 基礎用法：完整訓練與推論流程

 1from gears import PertData, GEARS
 2
 3# Step 1: 載入資料
 4pert_data = PertData('./data')
 5pert_data.load(data_name='norman')
 6
 7# Step 2: 資料切分（simulation split 模擬「訓練時沒看過某些組合」的情境）
 8pert_data.prepare_split(split='simulation', seed=1)
 9
10# Step 3: 建立 dataloader
11pert_data.get_dataloader(batch_size=32, test_batch_size=128)
12
13# Step 4: 初始化模型（GPU 若無則改用 'cpu'）
14gears_model = GEARS(pert_data, device='cuda:0')
15gears_model.model_initialize(hidden_size=64)
16
17# Step 5: 訓練
18gears_model.train(epochs=20)
19
20# Step 6: 儲存/載入模型（方便重複使用，不必每次重新訓練）
21gears_model.save_model('gears_norman_ckpt')
22gears_model.load_pretrained('gears_norman_ckpt')
23
24# Step 7: 預測全新的多基因擾動組合
25predictions = gears_model.predict([['CBL', 'CNN1'], ['FEV']])
26print(predictions)
27# 輸出格式：{'CBL_CNN1': array([...每個基因的預測表現量...]),
28#            'FEV':       array([...])}

5.2 進階用法：遺傳交互作用 (Genetic Interaction) 預測

 1from gears import PertData, GEARS
 2
 3pert_data = PertData('./data')
 4pert_data.load(data_name='norman')
 5pert_data.prepare_split(split='simulation', seed=1)
 6pert_data.get_dataloader(batch_size=32, test_batch_size=128)
 7
 8gears_model = GEARS(pert_data, device='cuda:0')
 9gears_model.load_pretrained('gears_norman_ckpt')
10
11# 預測 CBL + CNN1 雙基因擾動的交互作用類型
12gi_result = gears_model.GI_predict(['CBL', 'CNN1'], GI_genes_file=None)
13
14# gi_result 通常包含以下指標（依論文定義）：
15# - magnitude:  交互作用的整體幅度（越大代表雙基因擾動效果偏離"單純加總"越遠）
16# - model_fit:  加成模型 (additive model) 對觀測值的擬合優度
17# - dcor:       distance correlation，衡量兩基因擾動反應模式的非線性相關性
18# - equality_of_contribution: 兩基因各自貢獻是否對等
19
20print(gi_result)

這個結果可以幫助研究者快速篩選出「協同效應強烈」（可能代表合成致死或功能互補）或「拮抗抵銷」（可能代表兩基因在同一路徑上）的候選組合，優先安排濕實驗驗證，而非盲目地做全組合篩選。

5.3 進階用法：使用不確定性模式訓練

 1from gears import PertData, GEARS
 2
 3pert_data = PertData('./data')
 4pert_data.load(data_name='adamson')
 5pert_data.prepare_split(split='simulation', seed=1)
 6pert_data.get_dataloader(batch_size=32, test_batch_size=128)
 7
 8gears_model = GEARS(pert_data, device='cuda:0')
 9
10# 開啟不確定性估計 head
11gears_model.model_initialize(hidden_size=64, uncertainty=True)
12gears_model.train(epochs=20)
13
14predictions = gears_model.predict([['PERK']])
15# 在 uncertainty=True 模式下，predict() 額外提供每個基因預測的
16# log-variance，可用來排序「模型信心較低」的預測結果，
17# 交由人工複核或優先安排濕實驗確認。

5.4 使用自己的資料集（實際應用場景）

 1import scanpy as sc
 2from gears import PertData, GEARS
 3
 4# 假設你已經有一個 Perturb-seq 實驗的 AnnData
 5adata = sc.read_h5ad('my_perturbseq_experiment.h5ad')
 6
 7# 確保欄位符合 GEARS 規範
 8assert 'gene_name' in adata.var.columns
 9assert 'condition' in adata.obs.columns   # 例：'ctrl', 'TP53', 'TP53+MDM2'
10assert 'cell_type' in adata.obs.columns
11
12pert_data = PertData('./data')
13pert_data.new_data_process(dataset_name='my_experiment', adata=adata)
14pert_data.load(data_path='./data/my_experiment')
15pert_data.prepare_split(split='simulation', seed=1)
16pert_data.get_dataloader(batch_size=32, test_batch_size=128)
17
18gears_model = GEARS(pert_data, device='cuda:0')
19gears_model.model_initialize(hidden_size=64)
20gears_model.train(epochs=20)

5.5 評估模型效果：evaluate() 與 compute_metrics()

inference.py 提供了一套標準的評估流程，理解它有助於判斷模型訓練得好不好，也是論文報告指標的實際計算方式：

 1from gears.inference import evaluate, compute_metrics
 2
 3# 對測試集 dataloader 跑一輪推論，蒐集預測值與真實值
 4test_res = evaluate(
 5    loader=pert_data.dataloader['test_loader'],
 6    model=gears_model.model,
 7    uncertainty=False,
 8    device='cuda:0'
 9)
10
11# 計算整體與逐擾動 (per-perturbation) 的指標
12metrics, metrics_pert = compute_metrics(test_res)
13
14print("整體平均指標：")
15print(f"  MSE (全部基因):        {metrics['mse']:.4f}")
16print(f"  MSE (差異表現基因 DE): {metrics['mse_de']:.4f}")
17print(f"  Pearson (全部基因):     {metrics['pearson']:.4f}")
18print(f"  Pearson (DE 基因):      {metrics['pearson_de']:.4f}")
19
20# 查看特定擾動組合的表現
21print(metrics_pert['CBL+CNN1'])

這裡有兩個關鍵設計值得特別說明：

1. _de 後綴代表「只看差異表現基因 (differentially expressed genes; DE genes; 差異表現基因)」。GEARS 評估時不只看全部基因的整體吻合度，還特別挑出因為這個擾動而顯著變化的基因子集（通常是每個擾動的 top 20 DE genes）單獨計算指標——因為多數基因在任何擾動下都幾乎不變，如果只看全部基因的平均指標，容易被「一大堆沒變化、預測也沒變化」的基因稀釋掉真正重要的訊號，DE 基因指標才能反映模型是否抓到了「真正被這個擾動影響」的基因。這對應到論文demo tutorial_plot_top20_DE.ipynb 的可視化重點。

2. Pearson correlation 是用「跨細胞平均後」的表現量向量計算，而非逐細胞計算（見 results['pred'][p_idx].mean(0)）。原因是單一細胞層級的雜訊 (noise; 雜訊) 很大，取同一擾動下所有細胞的平均表現量再比較，才能凸顯系統性的生物學訊號而非隨機雜訊。

5.6 資料前處理細節：如何辨識單基因與組合擾動

utils.py 中的 parse_single_pert / parse_combo_pert / parse_any_pert 是理解 GEARS 資料格式慣例的關鍵函式：

1# GEARS 用 "+" 分隔多基因擾動，用 "ctrl" 代表對照組
2# 例如："CBL+ctrl" 代表只敲除 CBL（單基因擾動）
3#      "CBL+CNN1"  代表同時敲除 CBL 和 CNN1（組合擾動）
4#      "ctrl"      代表完全沒有擾動的對照組
5
6from gears.utils import parse_any_pert
7
8print(parse_any_pert('CBL+ctrl'))   # ['CBL']         單基因擾動
9print(parse_any_pert('CBL+CNN1'))   # ['CBL', 'CNN1'] 組合擾動

理解這個命名慣例非常重要——如果要把自己的資料餵給 GEARS，adata.obs['condition'] 欄位就必須遵守這套 "geneA+geneB" / "geneA+ctrl" / "ctrl" 的字串格式，否則 PertData.new_data_process() 在解析擾動類型時會出錯或誤判。

6. 進階應用場景 (Advanced Use Cases)

6.1 真實世界的應用案例

1. 藥物標靶優先排序：在藥廠的標靶發現流程中，先用 GEARS 對候選基因組合做計算篩選，把預測「表現變化幅度大」且「不確定性低」的組合排在濕實驗佇列前段，可大幅縮短濕實驗週期。
2. 合成致死篩選輔助：癌症研究中常見的合成致死篩選（找出「兩基因同時失活才致死，單獨失活不致死」的配對）可以用 GI_predict 的交互作用分數做初篩，再挑選 top candidate 做濕實驗確認。
3. 細胞重編程路徑設計：在幹細胞 (stem cell; SC; 幹細胞) 重編程研究中，需要找出能把細胞從一種狀態轉換到另一種狀態的轉錄因子組合，GEARS 可以模擬多個轉錄因子共同擾動後的表現體變化，輔助設計候選配方。

6.2 與其他 snap-stanford 工具的整合

• PyTorch Geometric (PyG)：GEARS 直接構築於 PyG 之上，SGConv、Data、DataLoader 皆來自 PyG；理解 PyG 的圖資料格式是修改 GEARS 內部圖結構（例如替換成自訂的基因調控網路）的前提。
• 與 gears_misc 配合：yhr91/gears_misc 提供論文圖表重現腳本，適合用來驗證自己的環境設定是否正確重現論文報告的指標數字。
• 與後續 foundation model（如 scGPT）比較：多篇後續論文把 GEARS 當作有先驗知識圖結構的代表性 baseline，與純 transformer-based foundation model 對照，是評估「圖結構先驗 vs. 大規模預訓練」孰優孰劣的重要參照點。

6.3 Distance Correlation 在遺傳交互作用判斷中的角色

GI_predict 底層使用 dcor 套件的 distance_correlation（距離相關性 distance correlation; dCor）而非傳統的 Pearson correlation 來衡量雙基因擾動反應模式的相似性。這個選擇的原因是 Pearson correlation 只能捕捉線性關係，而基因調控網路中的交互作用經常是非線性的（例如某基因的效果只有在另一個基因也被擾動時才會顯現，屬於「條件式」交互作用）。dcor 的優勢是它對任意類型的統計相依關係（不論線性與否）都能給出非零值，只有在兩變數真正獨立時才會趨近於零，因此比 Pearson 更適合捕捉複雜的遺傳交互作用模式。

此外，utils.py 中也使用了 TheilSenRegressor（Theil-Sen 穩健回歸 robust regression; 穩健回歸）來擬合「加成模型」的基準線：如果兩基因各自擾動的效果是簡單相加，那麼雙基因擾動的觀測值應該落在這條回歸線上；偏離越大，代表交互作用（協同或拮抗）越強。使用穩健回歸而非普通最小平方法 (OLS)，是為了避免少數離群值 (outlier; 離群值) 基因過度影響對「加成關係」基準線的估計。

6.4 效能調校與最佳實踐

• hidden_size：預設 64，資料集越大（基因/擾動數量越多）可以嘗試調高到 128 或 256，但要留意 GPU 記憶體會隨基因數量近似線性成長（因為 indv_w1/indv_w2 的參數量與 num_genes 成正比）。
• num_go_gnn_layers / num_gene_gnn_layers：層數過深容易發生 over-smoothing（訊息傳遞太多輪後所有節點嵌入趨於一致），論文預設值通常落在 1-2 層之間即可取得不錯效果，非必要不建議盲目加深。
• 資料切分策略：prepare_split 提供 simulation（隨機留出部分基因組合模擬"從未見過"情境）與自訂切分（v0.1.1 新增 custom split），建議依照下游應用場景選擇合適的切分方式，模擬實際部署時會遇到的分佈外 (out-of-distribution; OOD) 情境。
• CPU/GPU 混合開發流程：開發階段先用 Dixit（資料量最小）在 CPU 上跑通整條 pipeline，確認資料格式與程式邏輯無誤後，再切到 GPU 跑 Norman/Adamson 等大型資料集，避免在大資料集上反覆除錯浪費算力。

7. 與 AIKT 的關聯分析 (AIKT Integration Analysis)

7.1 此專案可作為 AIKT 的資料來源或分析工具嗎？

GEARS 本身是一個領域專用的科研模型（single-cell perturbation prediction），而非通用知識管理工具，因此它不會直接成為 AIKT 的核心 pipeline 元件，但有兩個明確的整合切入點：

1. 作為知識擷取對象：GEARS 的 GitHub repo、論文、demo notebooks 本身就是值得被 AIKT L1-L15 各層擷取、索引、教學化的「知識資產」——這正是本篇教學文件產出的情境（透過 gh-tutorial-qd 生成完整教學）。
2. 作為 tu-plan-generator (L19) 的技術背景素材：若使用者的藥物開發計畫涉及基因標靶篩選、multi-gene perturbation 相關的臨床前 (preclinical; 臨床前) 策略，GEARS 論文與方法論可以作為 L19 生成計畫時引用的技術背景文獻。

7.2 AIKT 的哪些 Layer 可以包裝或編排此專案？

flowchart TB
    subgraph L2["L2 ai-gh-save"]
        a1["gh: snap-stanford/GEARS
→ inbox/ 收錄 repo 摘要"]
    end
    subgraph L9["L9 paper-search"]
        a2["搜尋 Roohani 2023 Nat Biotech
+ 後續引用該論文的研究"]
    end
    subgraph L10["L10 paper-qa-lite"]
        a3["針對已下載的 GEARS 論文 PDF
做本地 RAG 問答"]
    end
    subgraph L12["L12 gh-tutorial-qd"]
        a4["本篇教學文件的產生流程
(GitHub repo → md → HTML)"]
    end
    subgraph L18["L18 research-pipeline-v2"]
        a5["多輪迭代研究：
GEARS 相關文獻 + 程式碼
整合成研究筆記"]
    end
    subgraph L19["L19 tu-plan-generator"]
        a6["藥物開發計畫中
引用 GEARS 作為
標靶篩選技術背景"]
    end

    a1 --> a4
    a2 --> a3
    a3 --> a5
    a4 --> a5
    a5 --> a6

具體來說：

• L2 ai-gh-save：gh: snap-stanford/GEARS 前綴可以快速把 repo metadata（README、star 數、更新時間）收錄進 inbox/，作為後續分析的起點——這正是本次任務的第一步。
• L9 paper-search + L10 paper-qa-lite：先用 paper-search 找到 GEARS 論文與其引用文獻，下載後用 paper-qa-lite 做本地問答（例如「GEARS 的 simulation split 具體怎麼定義？」），比人工翻閱論文快得多。
• L12 gh-tutorial-qd：本教學文件本身的產出方式——分析 repo 結構、README、關鍵程式碼，生成雙語教學 md 再編譯成 HTML，整條流程正是 L12 的標準作業。
• L18 research-pipeline-v2：若研究目標是「評估多個 perturbation prediction 模型的優劣」，可以用多輪迭代 pipeline 把 GEARS、scGPT、其他模型的文獻與程式碼一起彙整比較。
• L19 tu-plan-generator：在藥物開發計畫（尤其是標靶確認 target validation 階段）文件中，引用 GEARS 作為「計算式篩選降低濕實驗成本」的方法論佐證。

7.3 潛在的整合場景與價值

• 場景一：新藥標靶篩選文獻回顧自動化。分析師輸入「我想研究 XX 癌症的合成致死標靶」，AIKT 可以透過 L9 搜尋相關文獻、L2 收錄 GEARS 等相關工具 repo、L18 彙整成一份包含「計算工具清單 + 代表論文 + 方法論比較」的研究筆記，大幅縮短前期文獻回顧時間。
• 場景二：內部教學資料庫建置。生資團隊若需要培訓新進人員理解 perturbation prediction 領域，可以用 L12 gh-tutorial-qd 針對 GEARS、後續的 GenePT、scGPT 等工具各自產出教學文件，累積成內部知識庫，比逐一閱讀原始論文/程式碼有效率。
• 場景三：計畫書技術背景撰寫。L19 tu-plan-generator 在生成臨床前開發計畫時，若涉及基因標靶多重擾動策略，可以直接引用本篇教學文件整理過的技術要點（如雙圖架構、不確定性估計的臨床決策意義）作為方法論段落素材，避免每次都要重新從論文原文摘錄。

7.4 整合時需注意的邊界

GEARS 本身運算需要 GPU 與 PyG 環境，並非輕量 CLI 工具，因此不適合被包裝成即時互動的 AIKT script（例如不會出現 gears: predict TP53+MDM2 這種前綴指令）；比較務實的定位是「作為被研究、被教學化、被引用的對象」，而非「被編排執行的 pipeline 元件」。這個判斷本身也呼應全域 CLAUDE.md 的核心原則：先問「這是真問題還是想像的整合需求」，避免為了整合而整合，增加不必要的維護負擔。

8. 優缺點與生態系定位 (Strengths, Limitations & Ecosystem Position)

8.1 核心優勢與創新點

1. 首創把擾動本身表示成可學習嵌入，並用 GO 相似圖賦予歸納能力，讓模型能對訓練時從未見過的基因組合做出有意義的預測，這是相對於傳統迴歸/分類方法最本質的突破。
2. Gene-specific decoder 設計尊重了不同基因響應曲線的異質性，避免「一套權重打天下」造成的預測失真。
3. 提供不確定性估計，讓下游使用者（尤其是需要決定濕實驗優先順序的實驗室）能量化模型信心，而非只拿到一個數字就盲目相信。
4. API 設計簡潔（PertData + GEARS 兩個類別即可完成整條流程），降低了非深度學習背景的濕實驗研究者的使用門檻。
5. 在 Nature Biotechnology 發表並開源，具備高學術可信度，也讓其成為該子領域事實上的標準基準。

8.2 目前限制與改進空間

官方 README 明確列出了幾項使用限制，值得特別注意：

• 不支援跨細胞類型訓練或遷移：GEARS 目前設計上不處理「在細胞類型 A 訓練、預測細胞類型 B」的情境，這對於想要跨癌症類型泛化的應用是明顯限制。
• 未在 bulk RNA-seq 資料上驗證：只在單細胞資料上測試過，若想用於 bulk 定序資料需自行驗證有效性。
• 純粹從單基因資料無法可靠預測組合擾動：模型必須在訓練資料中看過至少部分組合擾動樣本，才能學會「組合效應」的模式，單基因資料本身不足以外推到多基因場景。
• 對稀疏資料集（每個擾動細胞數太少、擾動種類太少）效果不佳：小樣本場景下 GNN 訊息傳遞容易學不到穩定的圖結構訊號。
• 計算成本：對於基因數量龐大的資料集（如全轉錄體、非篩選過的基因集），indv_w1/indv_w2 這類 per-gene 參數會讓模型參數量與記憶體需求隨基因數線性成長，超大規模資料集的訓練成本不容小覷。

8.3 與同領域工具的比較

GEARS 的差異化定位在於：它是**少數明確結合圖結構生物學先驗（而非僅依賴大規模預訓練或純資料驅動）**的代表性方法，因此常被拿來驗證「加入結構化先驗是否真的比純黑箱模型更有效」這個研究問題。

8.4 常見誤區與實務提醒

在實務導入 GEARS 時，以下幾個誤區特別容易踩坑：

1. 誤把 simulation split 的高分當作「模型能預測任意新組合」。如第 4.5 節所述，測試集內不同難度層級的樣本表現差異很大，務必分開檢視「兩基因皆訓練過但組合未訓練過」與「基因本身未訓練過」這兩種情境的指標，不能只看整體平均值。
2. 忽略 DE genes 指標，只看全基因組平均指標。多數基因在任何擾動下變化都很小，全基因組平均的 Pearson/MSE 容易呈現「看起來很準」的假象，但模型可能根本沒抓到真正被擾動影響的少數關鍵基因。務必同時檢視 mse_de / pearson_de。
3. 在訓練集擾動種類太少（例如少於 20-30 種）的資料集上硬套用 GEARS。GNN 需要足夠的圖結構訊號才能學到有意義的相似度傳播模式，過少的擾動種類會讓 GO 相似圖過於稀疏，模型退化成近似記憶訓練樣本而非真正泛化。
4. 把 GI_predict 的交互作用分數直接當成濕實驗驗證的替代品。這些分數是計算篩選的排序依據，本質上仍是統計模型的推論結果，不能取代濕實驗驗證，尤其是用於後續臨床前開發決策時，必須明確標註「計算預測，待濕實驗確認」。
5. 忽略 PyG 版本與 CUDA 版本的嚴格對應關係。如 3.4 節所述，這是安裝階段最常見、也最耗時的除錯來源，建議在乾淨的 uv venv 中一步步照 3.2 節順序安裝，避免混用不同管道（conda + pip 混裝）造成的相依衝突。

8.5 適用場景建議

• 適合：已有一定規模 Perturb-seq 資料（至少涵蓋數十種擾動、每種擾動有足夠細胞數）、需要預測全新基因組合效果、且該生物系統有較完整 GO 註解的場景。
• 需謹慎評估：資料極度稀疏、需要跨細胞類型泛化、或處理的是 bulk 定序資料的場景，應先參考官方 README 列出的限制，必要時考慮混合其他方法或先做小規模驗證再投入生產環境。
• 建議搭配：訓練前務必先用小型資料集（如 Dixit）跑通整條 pipeline，確認資料格式與環境安裝正確後，再投入大型資料集訓練，避免在除錯階段浪費 GPU 資源。

9. 完整實戰演練：從資料檢視到結果視覺化

本節整合前面各節的知識，示範一個完整、可直接執行（假設已下載好 Norman 資料集）的端到端流程，對應 demo/ 目錄下四份官方 notebook 的核心邏輯，並加上中文說明。

9.1 檢視原始資料結構

在真正訓練模型前，先花時間檢視資料本身的樣貌，是任何資料科學工作流程都該有的第一步——這也呼應全域 CLAUDE.md「先看資料結構再看程式碼」的原則。

 1from gears import PertData
 2
 3pert_data = PertData('./data')
 4pert_data.load(data_name='norman')
 5
 6adata = pert_data.adata
 7print(f"細胞數量: {adata.n_obs}")
 8print(f"基因數量: {adata.n_vars}")
 9print(f"擾動條件種類數: {adata.obs['condition'].nunique()}")
10print(f"前 5 種擾動條件: {adata.obs['condition'].unique()[:5]}")
11
12# 區分單基因擾動與組合擾動的數量
13conditions = adata.obs['condition'].unique()
14single_pert = [c for c in conditions if c != 'ctrl' and c.endswith('+ctrl')]
15combo_pert = [c for c in conditions if c != 'ctrl' and not c.endswith('+ctrl')]
16print(f"單基因擾動種類: {len(single_pert)}")
17print(f"組合擾動種類: {len(combo_pert)}")

9.2 視覺化 Top 20 差異表現基因

demo/tutorial_plot_top20_DE.ipynb 展示了如何檢視模型對「真正重要的少數基因」的預測準確度，這是判斷模型是否真正學到生物學訊號（而非只是預測「幾乎不變」這個安全答案）的關鍵可視化步驟：

 1import matplotlib.pyplot as plt
 2import numpy as np
 3
 4from gears.inference import evaluate
 5
 6test_res = evaluate(
 7    loader=pert_data.dataloader['test_loader'],
 8    model=gears_model.model,
 9    uncertainty=False,
10    device='cuda:0'
11)
12
13# 挑選一個特定擾動組合，比較預測值與真實值在 top DE genes 上的分佈
14pert_name = 'CBL+CNN1'
15idx = np.where(test_res['pert_cat'] == pert_name)[0]
16
17pred_de = test_res['pred_de'][idx].mean(0)
18truth_de = test_res['truth_de'][idx].mean(0)
19
20fig, ax = plt.subplots(figsize=(6, 6))
21ax.scatter(truth_de, pred_de, alpha=0.6)
22lims = [min(truth_de.min(), pred_de.min()), max(truth_de.max(), pred_de.max())]
23ax.plot(lims, lims, 'r--', label='完美預測基準線 (y=x)')
24ax.set_xlabel('真實表現量變化 (Top 20 DE genes)')
25ax.set_ylabel('模型預測表現量變化')
26ax.set_title(f'{pert_name} 擾動：預測 vs. 真實')
27ax.legend()
28plt.tight_layout()
29plt.savefig('gears_de_scatter.png', dpi=150)

散點越貼近紅色對角線（y=x），代表模型在這個擾動組合的差異表現基因上預測越準確。這種視覺化方式比單一數字指標更容易讓濕實驗研究者直覺理解「這個預測值得信任的程度」。

9.3 深度分析：deeper_analysis 與 non_dropout_analysis

inference.py 除了 compute_metrics 之外，還提供 deeper_analysis 與 non_dropout_analysis 兩個函式，用於更細緻地拆解模型表現：

 1from gears.inference import deeper_analysis, non_dropout_analysis
 2
 3# deeper_analysis 額外計算了方向性一致度 (direction accuracy)：
 4# 模型預測的表現量變化方向（上升/下降）是否與真實方向一致，
 5# 這比單純的數值誤差更貼近下游應用「這個基因是被活化還是被抑制」的實際決策需求
 6deep_metrics = deeper_analysis(adata, test_res)
 7
 8# non_dropout_analysis 則排除了單細胞資料中常見的 dropout
 9# （技術性零值：基因其實有表現，但因定序深度不足而讀不到）造成的雜訊，
10# 只在「非零表現」的基因子集上計算指標，避免 dropout 雜訊拉低模型評分的公平性
11non_dropout_metrics = non_dropout_analysis(adata, test_res)

這兩個輔助分析函式反映了 GEARS 團隊對單細胞資料特有雜訊來源（dropout）的細緻處理，也是評估單細胞相關深度學習模型時必須考慮的領域知識——通用的迴歸評估指標若不排除 dropout 干擾，容易低估模型的真實預測能力。

9.4 端到端流程總覽圖

flowchart TB
    A["1. pert_data.load()
下載/讀取 Norman/Adamson/Dixit"] --> B["2. 檢視資料結構
(9.1 節)"]
    B --> C["3. prepare_split(simulation)
刻意排除部分基因/組合"]
    C --> D["4. get_dataloader()
建 PyG 圖資料"]
    D --> E["5. GEARS.model_initialize()
建雙圖 GNN 模型"]
    E --> F["6. train()
訓練 N epochs"]
    F --> G["7. evaluate() + compute_metrics()
(5.5 節)"]
    G --> H["8. deeper_analysis /
non_dropout_analysis
(9.3 節)"]
    H --> I["9. 視覺化 Top DE genes
(9.2 節)"]
    I --> J["10. predict() / GI_predict()
對全新組合做推論
(5.1-5.2 節)"]
    J --> K["11. 依信心排序
候選組合，交付濕實驗"]

這個流程完整呼應了本文件第 1 節提到的核心價值主張：用計算方法把上百萬種可能的基因組合，收斂成一小份「有科學根據、值得優先驗證」的候選清單，大幅降低濕實驗的時間與金錢成本。

附錄：關鍵檔案結構速查

 1GEARS/
 2├── gears/
 3│   ├── __init__.py
 4│   ├── pertdata.py      # PertData：資料載入、切分、dataloader 建構
 5│   ├── gears.py         # GEARS：頂層訓練/推論 API（Facade）
 6│   ├── model.py         # GEARS_Model, MLP：GNN 模型定義
 7│   ├── inference.py     # evaluate, compute_metrics, deeper_analysis
 8│   ├── data_utils.py    # DataSplitter, get_DE_genes
 9│   ├── utils.py         # get_similarity_network, loss_fct, GeneSimNetwork
10│   └── version.py
11├── demo/
12│   ├── data_tutorial.ipynb          # 資料載入與自訂資料教學
13│   ├── model_tutorial.ipynb         # 訓練流程教學
14│   ├── tutorial_plot_top20_DE.ipynb # 繪製 top 20 差異表現基因
15│   ├── tutorial_uncertainty.ipynb   # 不確定性模式教學
16│   └── tutorial_inference_Norman.ipynb
17├── requirements.txt
18├── setup.py / setup.cfg
19└── README.md

本教學文件由 AIKT (AI Knowledge Template) L12 gh-tutorial-qd layer 分析 snap-stanford/GEARS repository 生成，內容基於官方 README、原始碼結構（gears/model.py、gears/gears.py、gears/pertdata.py）與 Nature Biotechnology 2023 論文摘要整理而成。