Self-Study

Benjamini-Hochberg 程序中錯誤發現率的證明/推導

  • November 5, 2020

Benjamini-Hochberg過程是一種校正多重比較的方法,錯誤發現率 (FDR) 等於 $ \alpha $ .

或者是家庭明智的錯誤率,FWER?我對此有點困惑。根據我下面的計算,它似乎是等於的 FWER $ \alpha $ 而不是羅斯福。

我們能證明這是真的嗎?

讓我們假設不同假設的多個 p 值是獨立的,並且 p 值的分佈(以零假設為真為條件)在 $ 0,1 $ .


我可以使用模擬來證明它已經接近了。使用以下數字 $ \alpha = 0.1 $ , 在這個模擬中我拒絕一個假設的次數是

$$ \begin{array}{rcl} \alpha& =& 0.1\ \text{observed FDR} &=& 0.100002 \pm 0.00030 \end{array} $$

有錯誤基於 $ \pm 2\sigma $ 在哪裡 $ \sigma = \sqrt{\frac{0.1 \cdot 0.9}{ n}} $

set.seed(1)
m <- 10^6
n <- 10
a <- 0.1
k <- 1:n

sample <- function( plotting = F) {
 p <- runif(n)
 p <- p[order(p)]
 counts <- max(0,which(p<k/n*a))
 if (plotting) {
   plot(k,p, ylim = c(0,1) )
   lines(k,k/n*a)
 }
 counts
}

x <- replicate(m, sample())
s <- sum(x>0)/m
err_s <- sqrt(s*(1-s)/m)
c(s-2*err_s,s,s+2*err_s) 

介紹

讓我們從一些符號開始:我們有 $ m $ 我們測試的簡單假設,每個空值都編號 $ H_{0,i} $ . 全局零假設可以寫成所有局部零的交集: $ H_0=\bigcap_{i=1}^{m}{H_{0,i}} $ . 接下來,我們假設每個假設 $ H_{0,i} $ 有檢驗統計量 $ t_i $ 我們可以計算 p 值 $ p_i $ . 更具體地說,我們假設對於每個 $ i $ 的分佈 $ p_i $ 在下面 $ H_{0,i} $ 是 $ U[0,1] $ .

例如(Efron的第 3 章),考慮比較兩組中的 6033 個基因: $ 1\le i\le 6033 $ , 我們有 $ H_{0,i}: $ “沒有區別 $ i^{th} $ 基因”和 $ H_{1,i}: $ “有區別 $ i^{th} $ 基因”。在這個問題中, $ m=6033 $ . 我們的 $ m $ 假設可分為 4 組:

我們不觀察 $ S,T,U,V $ 但我們確實觀察到 $ R $ .

更好的

一個經典的標準是家庭誤差,表示為 $ FWE=I{V\ge1} $ . 那麼家庭錯誤率是 $ FWER=E[I{V\ge1}]=P(V\ge1) $ . 比較方法控制水平- $ \alpha $ $ FWER $ 在強烈的意義上如果 $ FWER\le\alpha $ 對於任何組合 $ \tilde{H}\in\left{H_{0,i},H_{1,i}\right}_{i=1,…,m} $ ; 它控制水平- $ \alpha $ $ FWER $ 在弱的意義上如果 $ FWER\le\alpha $ 在全局 null 下。

示例 - Holm 方法

  1. 對 p 值進行排序 $ p_{(1)},…,p_{(m)} $ 然後分別是假設 $ H_{0,(1)},…,H_{0,(m)} $ .

  2. 我們一個接一個地檢查假設:

  3. 如果 $ p_{(1)}\le\frac{\alpha}{m} $ 我們拒絕 $ H_{0,(1)} $ 並繼續,否則我們停止。

  4. 如果 $ p_{(2)}\le\frac{\alpha}{m-1} $ 我們拒絕 $ H_{0,(2)} $ 並繼續,否則我們停止。

  5. 我們不斷拒絕 $ H_{0,(i)} $ 如果 $ p_{(i)}\le\frac{\alpha}{m-i+1} $

  6. 我們第一次得到時就停下來 $ p_{(i)}>\frac{\alpha}{m-i+1} $ 然後拒絕 $ H_{0,(1)},…,H_{0,(i-1)} $ .

Holm 方法示例:我們已經模擬 $ m=20 $ p 值,然後對它們進行排序。紅色圓圈表示被拒絕的假設,綠色圓圈沒有被拒絕,藍線是標準 $ \frac{\alpha}{m-i+1} $ : Holm 方法的可視化示例

被拒絕的假設是 $ H_{0,(1)},H_{0,(2)},H_{0,(3)} $ . 你可以看到 $ p_{(4)} $ 是大於標準的最小 p 值,因此我們拒絕所有 p 值較小的假設。

很容易證明這種方法控制了電平—— $ \alpha $ $ FWER $ 在強烈的意義上。一個反例是Simes 方法,它只控制水平- $ \alpha $ $ FWER $ 在弱的意義上。

羅斯福

錯誤發現比例 ( $ FDP $ ) 是比 $ FWE $ , 並定義為 $ FDP=\frac{V}{\max{1,R}}=\left{\begin{matrix} \frac{V}{R} & R\ge1\ 0 & o.w\end{matrix}\right. $ . 錯誤發現率為 $ FDR=E[FDP] $ . 控制水平- $ \alpha $ $ FDR $ 意味著如果我們多次重複實驗和拒絕標準是錯誤的, $ FDR=E[FDP]\le\alpha $ . 很容易證明 $ FWER\ge FDR $ : 我們首先聲明 $ I{V\ge1}\ge\left{\begin{matrix} \frac{V}{R} & R\ge1\ 0 & o.w\end{matrix}\right. $ .

如果 $ V=0 $ 然後 $ I{V\ge1}=0=\left{\begin{matrix} \frac{V}{R} & R\ge1\ 0 & o.w\end{matrix}\right. $ . 在上表中, $ R\ge V $ 因此對於 $ V>0 $ 我們得到 $ \frac{V}{\max{1,R}}\le1 $ 和 $ I{V\ge1}\ge\frac{V}{\max{1,R}} $ . 這也意味著 $ E\left[I{V\ge1}\right]\ge E\left[\frac{V}{\max{1,R}}\right] $ , 這正是 $ FWER\ge FDR $ .

如果控制水平,另一個非常簡單的主張 - $ \alpha $ $ FDR $ 意味著控制水平- $ \alpha $ $ FWER $ 弱意義上的,意思是在全局下 $ H_0 $ 我們得到 $ FWER=FDR $ : 在全球範圍內 $ H_0 $ 每個拒絕都是錯誤的拒絕,意思是 $ V=R $ , 所以$$ FDP=\left{\begin{matrix} \frac{V}{R} & R\ge1\ 0 & o.w\end{matrix}\right.=\left{\begin{matrix} \frac{V}{V} & V\ge1\ 0 & o.w\end{matrix}\right.=\left{\begin{matrix} 1 & V\ge1\ 0 & o.w\end{matrix}\right.=I{V\ge1} $$接著$$ FDR=E[FDP]=E[I{V\ge1}]=P(V\ge1)=FWER. $$

BH

Benjamini-Hochberg 方法如下:

  1. 對 p 值進行排序 $ p_{(1)},…,p_{(m)} $ 然後分別是假設 $ H_{0,(1)},…,H_{0,(m)} $
  2. 標記為 $ i_0 $ 最大的 $ i $ 為此 $ p_{(i)}\le \frac{i}{m}\alpha $
  3. 拒絕 $ H_{0,(1)},…,H_{0,(i_0)} $

與前面的說法相反,說明為什麼 BH 方法保持 $ FDR\le\alpha $ (更準確地說,它保持 $ FDR=\frac{m_0}{m}\alpha $ )。這也不是一個簡短的證明,他是一些高級統計課程材料(我在我的一門統計學碩士課程中看到過)。我確實認為,就這個問題而言,我們可以簡單地假設它控制著羅斯福。

BH 示例:我們再次模擬了 $ m=20 $ p 值,然後對它們進行排序。紅色圓圈表示被拒絕的假設,綠色圓圈沒有被拒絕,藍線是標準 $ \frac{i\cdot\alpha}{m} $ :

BH 方法的可視化示例

被拒絕的假設是 $ H_{0,(1)} $ 到 $ H_{0,(10)} $ . 你可以看到 $ p_{(11)} $ 是大於標準的最大 p 值,因此我們拒絕所有 p 值較小的假設——即使其中一些 ( $ p_{(6)},p_{(7)} $ ) 大於標準。比較這個(大於標準的最大p 值)和 Holm 方法(大於標準的**最小p 值)。

證明 $ FDR\le\alpha $ 對於 BH

為了 $ m_0=0 $ (這意味著每個 $ p_i $ 分佈在 $ H_{1,i} $ ) 我們得到 $ FDR\le\alpha $ 因為 $ V=0 $ , 所以假設 $ m_0\ge1 $ . 表示 $ V_i=I{H_{0,i}\text{ was rejected}} $ 和 $ \mathcal{H}0\subseteq{1,…,m} $ 一組假設 $ H{0,i} $ 是正確的,所以 $ FDP=\sum_{i\in\mathcal{H}0}{\frac{V_i}{\max{1,R}}}=\sum{i\in\mathcal{H}_0}{X_i} $ . 我們首先證明 $ i\in\mathcal{H}_0 $ 我們得到 $ E[X_i]=\frac{\alpha}{m} $ :

$$ X_i=\sum_{k=0}^{m}{\frac{V_i}{\max{1,R}}I{R=k}}=\sum_{k=1}^{m}{\frac{I{H_{0,i}\text{ was rejected}}}{k}I{R=k}}=\sum_{k=1}^{m}{\frac{I\left{ p_i\le\frac{k}{m}\alpha \right}}{k}I{R=k}} $$

讓 $ R(p_i) $ 我們得到的拒絕次數,如果 $ p_i=0 $ 其餘的 p 值不變。認為 $ R=k^* $ , 如果 $ p_i\le\frac{k^}{m}\alpha $ 然後 $ R=R(p_i)=k^ $ 所以在這種情況下, $ I\left{ p_i\le\frac{k^}{m}\alpha \right}\cdot I{R=k^}=I\left{ p_i\le\frac{k^}{m}\alpha \right}\cdot I{R(p_i)=k^} $ . 如果 $ p_i>\frac{k^}{m}\alpha $ 我們得到 $ I\left{ p_i\le\frac{k^}{m}\alpha \right}=0 $ 然後再次 $ I\left{ p_i\le\frac{k^}{m}\alpha \right}\cdot I{R=k^}=I\left{ p_i\le\frac{k^}{m}\alpha \right}\cdot I{R(p_i)=k^} $ , 所以總的來說我們可以推斷

$$ X_i=\sum_{k=1}^{m}{\frac{I\left{ p_i\le\frac{k}{m}\alpha \right}}{k}I{R=k}}=\sum_{k=1}^{m}{\frac{I\left{ p_i\le\frac{k}{m}\alpha \right}}{k}I{R(p_i)=k}}, $$

現在我們可以計算 $ E[X_i] $ 以所有 p 值為條件,除了 $ p_i $ . 在這種情況下, $ I{R(p_i)=k} $ 是確定性的,我們總體上得到:

$$ E\left[I\left{ p_i\le\frac{k}{m}\alpha \right}\cdot I{R(p_i)=k}\middle|p\backslash p_i\right]=E\left[I\left{ p_i\le\frac{k}{m}\alpha \right}\middle|p\backslash p_i\right]\cdot I{R(p_i)=k} $$ 因為 $ p_i $ 獨立於 $ p\backslash p_i $ 我們得到

$$ E\left[I\left{ p_i\le\frac{k}{m}\alpha \right}\middle|p\backslash p_i\right]\cdot I{R(p_i)=k}=E\left[I\left{ p_i\le\frac{k}{m}\alpha \right}\right]\cdot I{R(p_i)=k}\=P\left(p_i\le\frac{k}{m}\alpha\right)\cdot I{R(p_i)=k} $$

我們之前假設在 $ H_{0,i} $ , $ p_i\sim U[0,1] $ 所以最後一個表達式可以寫成 $ \frac{k}{m}\alpha\cdot I{R(p_i)=k} $ . 下一個,

$$ E[X_i|p\backslash p_i]=\sum_{k=1}^m\frac{E\left[I\left{ p_i\le\frac{k}{m}\alpha \right}\cdot I{R(p_i)=k}\middle|p\backslash p_i\right]}{k}=\sum_{k=1}^m\frac{\frac{k}{m}\alpha\cdot I{R(p_i)=k}}{k}\=\frac{\alpha}{m}\cdot\sum_{k=1}^m{I{R(p_i)=k}}=\frac{\alpha}{m}\cdot 1=\frac{\alpha}{m}. $$

使用總期望定律我們得到 $ E[X_i]=E[E[X_i|p\backslash p_i]]=E\left[\frac{\alpha}{m}\right]=\frac{\alpha}{m} $ . 我們之前已經獲得 $ FDP=\sum_{i\in\mathcal{H}_0}{X_i} $ 所以

$$ FDR=E[FDP]=E\left[\sum_{i\in\mathcal{H}0}{X_i}\right]=\sum{i\in\mathcal{H}0}{E[X_i]}=\sum{i\in\mathcal{H}_0}{\frac{\alpha}{m}}=\frac{m_0}{m}\alpha\le\alpha\qquad\blacksquare $$

概括

我們已經看到了強感和弱感的區別 $ FWER $ 以及 $ DFR $ . 我認為您現在可以發現自己的差異並理解原因 $ FDR\le\alpha $ 並不意味著 $ FWER\le\alpha $ 在強烈的意義上

引用自:https://stats.stackexchange.com/questions/495211

comments powered by Disqus