什么是FDR?
错误发现率FDR(False discovery rate)是在所有结果显著的检验中,假阳(零假设H0为真时,拒绝H0的情况)所占的比率。如下表所示,N次假设检验中,FDR定义为V/R=V/(V+S)。
而经典的Benjamini-Hochberg (BH) 方法就是用于控制错误发现率FDR的一种方法,让FDR≤α。
Benjamini-Hochberg方法介绍
有N次假设检验,对每一次假设检验都计算其P值,然后将计算出的P值按照从小到大的方式排序,接着从最小的P值开始,按照P(k)≤α*k/N进行比较,然后可以找到最大的第K个满足上述不等式的P值,最终可以认为这K个P值是显著的,其余的P值不显著。
我们来看一个具体例子,假设需要检验的总体均值为6%,重复进行了30次抽样检验,最小的6个P值如下表所示,如果使用5%的显著性水平,仅考虑P值大小来评估的话,那么我们将会拒绝最小的5个P值所对应的检验(P值=0.0625>5%),但使用Benjamini-Hochberg方法修正后,只会拒绝一个P值(下表中第一个)。
Benjamini-Hochberg方法原理
我们将10000次假设检验分为2组:
1. 9000次检验的零假设H0:真;
2. 1000次检验的零假设H0:假。
然后,可以看到这两组检验的P值分布情况如下图所示:
H0为真时,P值均匀分布在0%-100%之间,为什么会是均匀分布呢?是因为在零假设H0条件下,P值有5%的可能性小于5%,有10%的可能性小于10%,有20%的可能性小于20%,以此类推,可以很直观理解P值的均匀分布。而上图之所以不是完全的均匀分布,是因为样本数量还不够大(当样本数量越大,P值也就越接近于均匀分布)。
H0为假时,P值就不再是均匀分布了,而是集中在0%附近,其他区间基本没有出现。这也比较好理解:H0是假,假设检验的功效越大,检验出H0为假的能力就越好,也就意味着P值越小,拒绝H0的证据越明显。
如果将所有P值合在一起统计的话,就如下图所示。0%附近的第一个直方块高度为1400次,而后面均匀分布的方块平均高度为453次(如下图红线所示),因此使用这个直方图,我们可以大致估计出零假设H0是假,应该被拒绝的数是:1400-453=947(和上图中的真实值1000非常接近)。
进一步思考,0%附近第一个直方块中包含的1400个P值,按照常规假设检验,都需要拒绝掉么?通过上述分析,一个合理的拒绝数量应该为947个。
实际上,按照Benjamini-Hochberg方法,从小到大的顺序对p值进行排序,按照P(k)≤α*k/N进行比较,拒绝最小的P值,其中有116个H0为真的情况,也就意味着错误发现率FDR=116/947=0.12。
如果要控制FDR≤α=0.1,则可重新使用Benjamini-Hochberg方法,这次我们从更加图形可视化的角度来理解这个过程。
如下图所示,横坐标是假设检验次数,纵坐标是P值。在坐标系中我们先以α/N为斜率画一条红线(P=α*k/N函数),然后将所有假设检验的P值分布在坐标系中,拒绝掉所有在红线下的P值(也就是≤α*k/N的P值)。
具体而言,α/N为斜率的红线和所有假设检验的P值相交的最大点对应频次为883,在883个最小P值中实际上有83个零假设H0为真的情况,也就意味着实际的FDR=83/883=0.094<0.1。
为什么Benjamini-Hochberg方法能保证FDR≤α呢?本质是什么?
我们将上面的具体问题提炼总结一下,如下图所示:横坐标是假设检验次数,纵坐标是P值(0%-100%之间),先画一条以α/N为斜率的红线,假设L是红线和所有假设检验的P值相交的最大点(L以下都是需要拒绝的P值)。而基于上文,我们知道零假设H0为真时的P值是均匀分布,也就是说P值落在任意[0%-100%]区间的期望值为:k%*H0为真的假设检验次数,而H0为真的假设检验次数N0一定是小于N的。
在L个被拒绝的P值中,其中H0为真的假设检验个数为:h* H0为真的假设检验次数=h*N0。
这就是为什么Benjamini-Hochberg方法能有效控制错误发现率FDR≤α的底层逻辑。