计算 RNA seq 和 ChIP 芯片数据集之间基因列表重叠的概率

机器算法验证 r 遗传学 生物信息学 微阵列 生物统计学
2022-02-01 04:10:33

希望这些论坛上的人可以帮助我解决基因表达研究中的这个基本问题。

我对实验组织和对照组织进行了深度测序。然后,我获得了实验样本中基因的倍数富集值,而不是对照。参考基因组有约 15,000 个基因。与对照相比,在我感兴趣的样本中,15,000 个基因中有 3,000 个富集高于某个截止值。

所以:A = 总基因种群 = 15,000 B = RNA-Seq 富集亚群 = 3,000。

在之前的 ChIP-chip 实验中,我发现了 ChIP-chip 富集的 400 个基因。在 400 个 ChIP 芯片基因中,100 个基因位于 3,000 个丰富的 RNA-Seq 转录本中。

所以:C = ChIP-chip 富集基因的总数 = 400。

我的 100 个 ChIP 芯片基因仅靠 RNA-Seq 偶然富集的概率是多少?换句话说,如果我观察到的 B 和 C 之间的重叠(100 个基因)比仅凭偶然获得的重叠更好,那么最谨慎的计算方法是什么?从我目前所读到的,最好的测试方法是使用超几何分布。

我使用在线计算器 (stattrek.com) 设置了具有以下参数的超几何分布测试: - 人口规模 = 15,000 - 人口成功数 = 3,000 - 样本规模 = 400, - 样本成功数 = 100。我得到以下超几何概率 P(x=100)= 0.00224050636447747

B 和 C 之间重叠的实际基因数 = 100。这是否比仅靠偶然更好?如果任何一个基因被富集的机会是 1:5(15,000 个中的 3,000 个),看起来就不是这样。这就是为什么我不明白我上面计算的 P(x=100) 是 0.0022 的原因。这相当于偶然发生重叠的几率为 0.2%。这不应该更高吗?

如果我从 15,000 个大列表中随机抽取 400 个基因,那么这些基因中的任何 80 个都将被认为是偶然富集的(1:5)。实际重叠的基因数量是 100,所以这比偶然发生的要好一些。

我还尝试使用 R 中的 dhyper 或 phyper 函数提出一个解决方案(使用我在另一篇文章中看到的):A=基因组中的所有基因(15,000)B=RNA-Seq 富集基因(3,000)C=ChIP -芯片富集基因(400)这是R输入/输出(改编自之前的stackexchange帖子):

> totalpop <- 15000    
> sample1 <- 3000    
> sample2 <- 400    
> dhyper(0:2, sample1, totalpop-sample1, sample2)    
[1] 4.431784e-40 4.584209e-38 2.364018e-36    
> phyper(-1:2, sample1, totalpop-sample1, sample2)    
[1] 0.000000e+00 4.431784e-40 4.628526e-38 2.410304e-36    

我不确定如何解释这些数字。我相信 2.36e-36 是偶然在 B 和 C 之间完全重叠的概率吗?但这没有任何意义,因为这个概率更接近 1:5。如果我从 15,000 个基因开始,将丰富 3,000 个。同样,如果我从 400 个 ChIP 芯片基因开始,其中 80 个应该仅在 RNA-Seq 中富集,因为在该数据集中富集的机会是 1:5。

根据超几何分布计算 B 和 C 重叠的 p 值的正确方法是什么?

2个回答

你很接近,你使用dhyperand phyper,但我不明白0:2-1:2来自哪里。

您想要的 p 值是从具有 3000 个白球和 12000 个黑球的瓮中获得 100 个或更多白球的概率,样本大小为 400。这里有四种计算方法。

sum(dhyper(100:400, 3000, 12000, 400))
1 - sum(dhyper(0:99, 3000, 12000, 400))
phyper(99, 3000, 12000, 400, lower.tail=FALSE)
1-phyper(99, 3000, 12000, 400)

这些给出 0.0078。

dhyper(x, m, n, k)给出准确绘制的概率x. 在第一行,我们总结了 100 – 400 的概率;在第二行中,我们用 1 减去 0 - 99 的概率之和。

phyper(x, m, n, k)给出得到x或更少的概率,所以phyper(x, m, n, k)与 相同sum(dhyper(0:x, m, n, k))

lower.tail=FALSE有点令人困惑。 phyper(x, m, n, k, lower.tail=FALSE)与 相同,或更多1-phyper(x, m, n, k)的概率也是如此。x+1[我从不记得这个,所以总是要仔细检查。]

在那个stattrek.com 网站上,您想查看最后一行,“累积概率:P(X100)”而不是第一行“超几何概率:P(X = 100)”。

您绘制的任何特定数字的概率都很小(实际上,max(dhyper(0:400, 3000, 12000, 400))给出0.050),得到 101 或 102 或任何更大的数字甚至比 100 更有趣,并且 p 值是如果原假设为真,得到结果与观察到的结果一样有趣或更有趣的概率。

这是这种情况下超几何分布的图片。您可以看到它以 80 为中心(400 的 20%),而 100 在右尾非常远。 在此处输入图像描述

这样看..如果您认为它是二项式,这可能不正确,但它应该是相当近似的..您的 sigma^2 是 .8*.2*400= 64,然后 sigma = 8。所以从 80 到 100,你有 2.5 个标准差。这非常重要。它应该有一个小的 p 值。