生物富集在生物信息中有着重要的地位,做生物信息分析的时候总会遇到这样或者那样的富集分析,比如GO富集分析等。大多数情况下我们都是使用线上在线分析解决。是否会遇到这样一种情况,当我们不能使用在线分析的时候,如何对我们感兴趣的基因进行富集呢?
这似乎是一个很有趣的问题....
那我们今天就来追溯生物富集的前世今生....
比如我们现在遇到了这样一个问题:
在某次实验中,我们找到了一种全新的病毒X,通过病毒侵染细胞实验我们鉴定了6个病毒基因。通过查阅资料发现,这6个基因中有4个基因在一个特殊通路里,通过继续查阅资料我们发现,这个通路的全部基因有7个,而病毒的全部基因有15个,那么我们就产生了一个疑问,这个病毒是不是和这个通路有关系呢?
从数据来看,我们鉴定的6个基因中有有4个在这个通路里,显然这个病毒和这个通路有着特殊的联系,但似乎又觉得有一些不太对的地方。
显然我们需要一个指标来证明我们鉴定的6个基因确实通过特殊的富集才得到的...
我们将上面的问题简化为一个数学问题来看的话,这似乎就是一个摸球问题:
在一个瓶子里面有15个球,其中7个是黑球,我们从中取了6个球,其中4个是黑球,这是不是小概率事件?
那怎么计算这个概率呢?通过排列组合和概率论的知识可以得到
也就是说我们随机抽取一次,其中抽到4个黑球的概率是0.196
那我们将这个公式继续推广可以得到:
这似乎就是一个超几何分布模型
假设有限总体包含N个球,其中黑球为M个,则剩余的N-M个为白球,如果从该有限总体中抽取出n个球,其中有k个是黑球的概率。
回到我们的问题,从6个球中抽到了4个黑球这个事件到底是否具有显著性呢?
我们在前面计算出了抽到4个黑球的概率是0.196 这个值是不能拿来直接用的,必须要对其进行一个评估,因为我们必须要考虑到随机情况
在这里我们就用到了超几何分布检验( Hyper Geometric Test)
怎么理解超几何分布检验呢?
我们给一个假设,此次抽样与随机抽样没有差异(原假设)
P值就是当原假设为真时所得到的样本观察结果或更极端结果出现的概率。如果P值很小,说明这种情况的发生的概率很小,而如果出现了,根据小概率原理,我们就有理由拒绝原假设。
当前的样本观测结果是黑色的球为4个,更极端的结果就是k=5,k=6的情况
因此我们得到如下计算结果:
从而得到的p值结果为0.230769 因此我们按照95%的置信区间来看,0.230769>0.05
因此我们不认为这是一个小概率事件,因此我们认为此次抽球与随机抽取没有差别
回到我们最初的生物学问题,如果我们发现有5个基因在特殊通路中的话,那么
p=p(x=5)+p(x=6)=0.034965 <0.05
此时我们就拒绝原假设,认为这是一个小概率事件,也就是我们鉴定的基因和通路有比较强的联系。
如何对p值进行计算呢?
R语言
1-phyper(4-1, 7, 8, 6)
[1] 0.2307692
其中4 为抽取6个球中黑球的数目 ,7 为袋子黑球的数目,8为袋中白球的数目,6为所抽球的个数
python
from scipy.stats import hypergeom
hypergeom.sf(4-1,15,7,6)
out[1] 0.23076923076923062
其中4 为抽取6个球中黑球的数目,15为袋中黑白球总数,7 为袋子黑球的数目,6为所抽球的个数
解决了p值后还会遇到另一个问题就是,在我们对鉴定到的差异差异基因做通路富集后,通常会计算一个p值。当某个通路的p值小于0.05(5%)时,我们通常认为这个通路是通过富集得到的。但是仍旧有5%的概率就是这个通路是通过随机抽取得到的。那么我们就错误地否认了原假设,导致了假阳性的产生(犯错的概率为5%)。
如果检验一次,犯错的概率是5%;检测10000次,犯错的次数就是500次,即额外多出了500次差异的结论(即使实际没有差异)。
为了控制假阳性的次数,于是我们需要对p值进行多重检验校正,提高阈值。
具体细节可以移步浅谈多重检验校正FDR