RNA-seq 詳細教程:假設檢驗和多重檢驗(8)
學習目標
- 瞭解模型擬合的過程
- 比較兩種假設檢驗方法(Wald test vs. LRT)
- 瞭解多重測試校正的重要性
- 瞭解用於多重測試校正的不同方法
1. 模型擬合和假設檢驗
DESeq2
工作流程的最後一步是對每個基因進行計數並將其擬合到模型中並測試差異表達。
2. 廣義線性模型
如前所述,RNA-seq
生成的計數資料表現出過度分散(方差 > 均值),用於對計數建模的統計分佈需要考慮到這一點。因此,DESeq2
使用負二項分佈通過以下公式對 RNA-seq
計數進行建模:
所需的兩個引數是size factor和dispersion estimate。接下來使用廣義線性模型 (GLM) 來擬合數據。建模是一種數學形式化的方法,用於在給定一組引數的情況下估算資料的方式。模型擬合後,將估計每個樣本組的係數及其標準誤差。係數是 log2
3. 假設檢驗
假設檢驗的第一步是為每個基因建立一個零假設。在我們的例子中,原假設是兩個樣本組之間沒有差異表達 (LFC == 0)。然後,我們使用統計檢驗來根據觀察到的資料確定原假設是否為真。
3.1. Wald test
在 DESeq2
中,Wald 檢驗是比較兩組時用於假設檢驗的預設值。 Wald 檢驗是通常對已通過最大似然估計的引數執行的檢驗。在我們的案例中,我們正在測試每個基因模型係數 (LFC),這些係數是使用分散等引數得出的,這些引數是使用最大似然估計的。
DESeq2
通過以下方式實施 Wald 測試:
- 取 LFC 並將其除以標準誤差,得到 z 統計量
- 將 z 統計量與標準正態分佈進行比較,並計算 p 值,報告隨機選擇至少與觀察值一樣極端的 z 統計量的概率
- 如果 p 值很小,我們拒絕零假設並宣告有證據反對零假設(即基因差異表達)
模型擬合和 Wald 檢驗先前已作為 DESeq()
函式的一部分執行:
# 以下僅作示例,上一個教程已經執行
dds <- DESeqDataSetFromTximport(txi, colData = meta, design = ~ sampletype)
dds <- DESeq(dds)
3.2. 似然比檢驗
當比較兩個以上的樣本類別時,DESeq2
還提供似然比檢驗 (LRT) 替代假設檢驗。LRT
- 這與 Wald 檢驗相比如何?
Wald 檢驗(預設)僅估計每個基因一個模型並評估 LFC == 0 的原假設。
對於似然比檢驗,還對已通過最大似然估計的引數執行。對於這個測試,每個基因估計兩個模型;將一個模型的擬合度與另一個模型的擬合度進行比較。
- m1 是簡化模型(即刪除主要因素項的設計公式)
- m2 是完整模型(即您在建立
dds
物件時提供的完整設計公式)
在這裡,我們正在評估完整模型與簡化模型一樣適合的原假設。如果我們拒絕零假設,這表明完整模型(以及我們感興趣的主要因素)解釋了大量變異,因此該基因在不同水平上差異表達。 DESeq2
通過使用偏差分析 (ANODEV) 來比較兩個模型擬合來實現 LRT。結果表明,LR 服從卡方分佈,這可用於計算和關聯的 p 值。
要使用 LRT
,我們使用 DESeq()
函式,但這次新增兩個引數:
- 指定我們要使用 LRT 測試
- “簡化”模型
# Likelihood ratio test
dds_lrt <- DESeq(dds, test="LRT", reduced = ~ 1)
由於我們的“完整”模型只有一個因素(樣本型別),“簡化”模型(刪除該因素)在我們的設計公式中沒有留下任何東西。 DESeq2
無法擬合設計公式中沒有任何內容的模型,因此在沒有其他協變數的情況下,截距使用語法 ~ 1 建模。
4. Multiple test correction
無論我們使用 Wald 檢驗還是 LRT,每個經過檢驗的基因都會與一個 p 值相關聯。我們正是用這個結果來確定哪些基因被認為是顯著差異表達的。但是,我們不能直接使用 p 值。
4.1. p-value
顯著性截斷值 p < 0.05 的基因意味著它有 5% 的機率是假陽性。例如,如果我們測試 20,000 個基因的差異表達,在 p < 0.05 時,我們預計會偶然發現 1,000 個基因。如果我們發現總共有 3000 個基因存在差異表達,那麼大約三分之一的基因是假陽性!
由於每個 p 值都是單個測試(單個基因)的結果。我們測試的基因越多,我們的假陽性率就越大。這就是多重檢驗問題。
4.2. 校正
多重檢驗有幾種常見的校正方法:
- Bonferroni:調整後的 p 值計算方式為:p 值 * m(m = 測試總數)。這是一種非常保守的方法,假陰性的可能性很高,因此通常不推薦使用。
- FDR/Benjamini-Hochberg: Benjamini 和 Hochberg (1995) 定義了錯誤發現率 (FDR) 的概念,並建立了一種演算法,以在給定獨立 p 值列表的情況下將預期 FDR 控制在指定水平以下。
- Q-value / Storey method: 稱該特徵為重要時可達到的最低 FDR。例如,如果基因 X 的 q 值為 0.013,則表示 p 值至少與基因 X 一樣小的基因中有 1.3% 是假陽性。
DESeq2
通過去除那些在測試前不太可能顯著 DE 的基因,例如那些具有低計數和異常樣本(基因級 QC)的基因,幫助減少測試的基因數量。但是,還實施了多重測試校正,以使用 Benjamini-Hochberg
程式的解釋來降低錯誤發現率。
4.3. FDR < 0.05
通過將 FDR
截止值設定為 < 0.05,我們是說我們預期差異表達基因中的假陽性比例為 5%。例如,如果您將 500 個基因稱為差異表達,FDR
截斷值為 0.05,您預計其中 25 個是假陽性。
歡迎Star -> 學習目錄
更多教程 -> 學習目錄
本文由mdnice多平臺釋出