R语言中phyper做超几何检验

您所在的位置：网站首页 › r语言中求均值 › R语言中phyper做超几何检验

R语言中phyper做超几何检验

2023-04-13 12:22| 来源: 网络整理| 查看: 265

R语言中phyper做超几何检验 2023-04-11 Default Category 文章目录 phyper

基因集分析有两种，一种是GSEA（gene set enrichment analysis），需要根据所有基因logFC排序，根据rank来算enrichment score，还有一种是ORA（Over-representation analysis），看选出的显著的基因集是否和已知的基因集显著相关。对于ORA分析，常用超几何分布来检验。在R语言中，用的函数是phyper。

1 2 3 4 5 6 7 8 9 10 11 12 ############## 以基因集分析为例，超几何检验看通过差异筛选出的基因集genes是否和已知的某个基因集gene_set显著相关 # x: 研究者筛选出的基因集genes和待检验的已知基因集之间的重叠的基因数目gene_set----> length(intersect(genes, gene_set)) # m: 已知基因集内的基因数目---->length(gene_set) # n: 背景基因集减去m---->length(universe) - n # k: 研究者筛选出的差异基因的数目length(genes) ############## 以挑选红球和黑球为例，看挑选的球中是否富集红球 # 有一堆球，红球有n_red个，黑球有n_black个，总共为n个（n_red + n_black），拿出来n_out个球，其中有n_red_out个红球 # x: n_red_out拿出来的球中的红球数 # m: n_red红球数 # n: n - n_red = n_black黑球数 # k: n_out拿出来的球 phyper(x - 1, m, n, k, lower.tail = FALSE)

$$ phyper(…, lower.tail = FALSE)计算的是Pr(X>x)在计算p值时，x-1计算的才是Pr(X≥x) $$

case study

我有一个队列的数据，每个病人的分子分型和药物反应的信息如下表格。

Subtype1 Subtype2 Subtype3 Subtype4 not response 5 0 9 3 response 1 3 1 6

Q1：亚型3是否和not response相关？

1 2 3 4 # phyper(9 （not response中的subtype3的病人数）- 1, 10（subtype3的病人数）, 18（总病人数-subtype3的病人数）, 17（not response的病人数）, lower.tail = FALSE) > phyper(9 - 1, 10, 18, 17, lower.tail = FALSE) [1] 0.021859 # 显著相关

写到这里，很有意思的一件事，我想到了另外检验方法，用U（Z）检验来做test。首先把表格合并，得到

non-Subtype3 Subtype3 not response 8 9 response 10 1

零假设：not response在3型和非3型中的分布频率p1和p2是一致的

$$ p_1=8/18，p_2=9/10，加权平均为\overline{p}={(8+9)}/(8+10+9+1)=17/28 $$ $$ 频数的标准差为S_{p1,-p2}=\sqrt{p_1p_2({{1\over{n_1}} + {1\over{n_2}} })}=\sqrt{{8\over{18}}*{9\over{10}}({{1\over{18}} + {1\over{10}} })}=0.2494 $$ $$ 因为上面第一行的频数小于30，需要进行连续校正，得到u值为 $$ $$ u_c = {{{|p_1 - p_2|} - {0.5\over{n_1}} - {0.5\over{n_2}}}\over{S_{p1,-p2}}} = {|{{8\over{18}} - {9\over{10}} | - {0.5\over{18}} - {0.5\over{10}} }\over0.2494} = 1.514 $$ 经查表可知p值在0.0643和0.0655之间，我们取0.06445，这是双侧检验，如果是单侧检验的话，p值应该是0.06445/2=0.032，和我们用超几何检验的出来的类似。

因为样本数小于40，我们还可以用fisher精准检验来算P值，p值为0.0407

1 2 3 4 5 6 7 8 9 10 11 12 > fisher.test(matrix(c(8,10,9,1), nrow = 2)) Fisher's Exact Test for Count Data data: matrix(c(8, 10, 9, 1), nrow = 2) p-value = 0.04074 alternative hypothesis: true odds ratio is not equal to 1 95 percent confidence interval: 0.001853788 0.955674842 sample estimates: odds ratio 0.09663796

####################################################################

#Author: Jason

#####################################################################

文章作者 zzx

上次更新 2023-04-11

statistic Hypergeometric test

【本文地址】

R语言中phyper做超几何检验

R语言中phyper做超几何检验

今日新闻

推荐新闻