scRNA-seq 高變異基因識別 (Feature Selection) 指南

1. 為什麼需要特徵選擇？

此區塊介紹特徵選擇在單細胞數據處理中的三大核心目的。透過剔除無效資訊，我們能聚焦於真正定義細胞多樣性的生物學信號。請將游標移至下方卡片以查看詳細說明。

⬋

克服維度災難

單細胞數據通常包含 20,000+ 個基因，但大部分基因在多數細胞中不表達或表達極低。篩選 2,000-3,000 個高變異基因能大幅降低維度，顯著提升計算效率並減少記憶體消耗。

⚙

過濾技術噪音

測序過程中存在大量隨機性。許多低表達基因的微小變異，主要來自技術噪音（如 Dropout 效應，即基因明明有表達卻沒被測到）。排除這些基因可避免演算法被雜訊誤導。

🔍

增強生物信號

管家基因（Housekeeping genes）在所有細胞中穩定高表達，它們雖然表達量大，但對方差貢獻低，無助於區分細胞亞型。特徵選擇能精準過濾這些「高表達但無差異」的干擾項。

2. 核心挑戰：均值與變異的關係

本區塊透過互動式散佈圖，展示測序數據中最棘手的現象：**表達量越高的基因，其觀測到的方差通常也越大**。若僅挑選「方差最大」的基因，只會選出一堆管家基因。點擊下方按鈕，觀察算法如何透過「標準化」找出真正的異常變異基因。

圖表說明：每個點代表一個基因。X軸為平均表達量，Y軸為方差（均為對數尺度模擬數據）。

3. 主流算法解析

了解數據特徵後，我們來看看目前主流的單細胞分析工具套件是如何在數學上處理這個挑戰的。此區塊整理了三種最廣泛被使用的特徵選擇算法機制。

Variance Stabilizing Transformation (VST)

這是 R 語言套件 Seurat 目前的預設方法，核心思想是基於回歸模型來平滑均值與方差的關係。

✓ 1. 對所有基因計算 $log(mean)$ 和 $log(variance)$。
✓ 2. 使用 Loess 回歸 擬合出一條趨勢線，這條線代表了在特定表達水平下的「預期技術變異」。
✓ 3. 計算每個基因觀測到的變異與趨勢線之間的差異，即殘差（Residual）。
✓ 4. 根據殘差數值進行排名，選取殘差最大的前 N 個基因（通常為 2000 個）。

Dispersion-based (分散度法)

這是 Python 套件 Scanpy 以及早期 Seurat 使用的經典方法，透過分箱（Binning）技術來克服均值偏差。

✓ 1. 計算每個基因的分散度（Dispersion = 方差 / 均值）。
✓ 2. 為了消除高表達基因分散度自然偏高的影響，將所有基因依據其均值劃分成不同的箱子（Bins）。
✓ 3. 在每個獨立的箱子內部，計算分散度的 Z-score（標準化分數）。
✓ 4. 這種分箱策略確保了無論在低、中、高表達區間，都能公平地選出相對變異程度較大的基因。

Pearson Residuals (皮爾森殘差)

通常與 SCTransform 流程結合使用，是目前處理測序深度偏差被認為最嚴謹的統計方法。

✓ 1. 建立一個正則化的負二項分布（Negative Binomial）模型。
✓ 2. 使用該模型預測每個基因在每一個單一細胞中的「預期計數」（Expected Counts）。
✓ 3. 計算實際觀測值與預期值之間的皮爾森殘差。
✓ 4. 殘差的方差即代表了該基因排除了測序深度影響後的純粹生物學變異。選取殘差方差最大的基因作為 HVGs。

4. 代碼實現與分析建議

此區塊提供實務操作的代碼範例以及分析過程中應注意的關鍵事項。您可以根據您使用的程式語言切換查看相應的代碼，並參照建議調整您的分析流程。

              # 讀取數據並標準化後執行
              pbmc <- FindVariableFeatures(pbmc, 
                                                                             selection.method = "vst", 
                                                                             nfeatures = 2000)
            
            # 可視化變異基因散佈圖
              plot1 <- VariableFeaturePlot(pbmc)
            
            # 標示出變異最大的前 10 個基因名稱
              LabelPoints(plot = plot1, 
                                                            points = head(VariableFeatures(pbmc), 10), 
                                                            repel = TRUE)
            

              import scanpy as sc
            
            # 假設 adata 是已經過濾與對數標準化的 AnnData 物件
              # 執行特徵選擇 (預設為 dispersion 法，這裡指定 seurat flavor)
              sc.pp.highly_variable_genes(adata, 
                                          n_top_genes=2000, 
                                          flavor='seurat')
            
            # 生成均值-方差可視化圖表
              sc.pl.highly_variable_genes(adata)
            

💡 實務分析建議

1. 決定基因數量 (n_features)

通常選擇 2,000 個基因是標準做法。但如果您的樣本細胞異質性極高（例如：全腦圖譜、發育中的多組織胚胎樣本），增加到 3,000-5,000 個往往會有更好的聚類解析度。
2. 處理批次效應 (Batch Effects)

若數據混雜了不同批次（如不同測序儀、不同天實驗），特徵選擇極易選出「批次特異性基因」。強烈建議在各個批次中分別尋找 HVGs，最後取交集（最嚴謹）或聯集（保留最多資訊）作為最終特徵集。
3. 剔除干擾基因清單

線粒體基因（mt-）或核糖體基因（RPS/RPL）常因細胞壓力或代謝狀態被識別為高變異。如果您的研究重點是細胞類型（如免疫分群）而非代謝狀態，務必在進行 PCA 降維前，將這些基因從 HVG 清單中剔除。