R: p value矫正

Python011

R: p value矫正,第1张

p.adjust 提供多种矫正方法:

Adjust P-values for Multiple Comparisons

"holm", "hochberg", "hommel", "bonferroni", "BH", "BY","fdr"

1 模拟p value

The Normal Distribution

rnorm(n, mean = 0, sd = 1) random 函数,生成随机数,底物

dnorm(x, mean = 0, sd = 1) density 函数

pnorm(q, mean = 0, sd = 1) probability 函数

qnorm(p, mean = 0, sd = 1) quantile 函数

2 Bonferroni矫正

3 FDR矫正

更多

R语言中dnorm, pnorm, qnorm与rnorm以及随机数

setwd("E:/GSE25066")#环境设置

library(limma)#加载差异分析包limma

#将分组文件加载到环境中,分组信息第一列为样本名,第二列为分组信息如“high”“low”

targets<-read.csv("group.csv")

#将表达矩阵加载到环境中,行为基因,列为样本,这里应该注意去除重复项。

eset<-read.csv("expreset-basal1.csv",row.names = "symbol")

targets$Target=gsub("_",".",targets$Target)##若数据中存在特殊符号,将"_"替换成“.”,也可以不替换

##该数据集中实际存在不符合R的命名原则,所以在没个分类前加一个“F”,具体自己定

targets$Target=c(paste0("F",c(targets$Target),collapse = NULL,sep=""))

colnames(targets)=c("FileName","Target")#更改列名,为了和limma包中的一致

lev<-unique(targets$Target)##使用unique()函数进行去重

f <- factor(targets$Target, levels=lev)

design <- model.matrix(~0+f)

colnames(design) <- lev

cont.wt <- makeContrasts("high-low",

                      +  levels=design)

fit <- lmFit(eset, design)#前面矩阵的row.name=“symbol”

fit2 <- contrasts.fit(fit, cont.wt)

fit2 <- eBayes(fit2)

tT=topTable(fit2, adjust="BH",sort.by="logFC",n=Inf)

tT = subset(tT, select=c("adj.P.Val","P.Value","logFC"))

colnames(tT)=c("FDR","P.Value","logFC")

write.csv(tT,"DEGbasal.csv")

#最后的tT就是得到的差异基因,其中包含基因,P.Value和logFC

最新版的r语言没有方差分析表格的解决方法如下:

ANOVA对各疗法的F检验表明,4种药品用于缓解术后疼痛的疗效不同,但是并不能得出哪种药品疗法与其他不同。

多重比较可以解决这个问题.e.g. TukeyHSD()函数提供了对各组均值差异的成对检验;multcomp包中的glht()函数提供了多重均值比较更为全面的方法,既适用于线性模型,也适用于广义线性模型;多重t检验方法针对每组数据进行t检验。代码如下: TukeyHSD(medicine.aov) #par()函数旋转轴标签,增大左边界面积,使标签摆放更美观。 par(las = 2) par(mar = c(5, 8, 4, 2)) plot(TukeyHSD(medicine.aov))

图形中置信区间包含0的药品对比,说明差异不显著。 library(multcomp) #为适合字母阵列摆放,par语句用来增大顶部边界面积 par(mar = c(5, 4, 6, 2)) tuk <- glht(medicine.aov, linfct = mcp(Treatment = "Tukey")) #cld()函数中level选项为设置的显著性水平(这里的0.05对应95%置信区间) plot(cld(tuk, level = 0.05), col = "lightgrey")

有相同字母的组(用箱线图表示)说明均值差异不显著。

多次重复使用t检验会增大犯第一类错误的概率,为了克服这一缺点,需要调整p-值。R软件调整p-值用的是p.adjust()函数,函数使用的不同参数代表不同的调整方法。 attach(medicine) #求数据在各水平下的均值 mu<-c(mean(Response[Treatment==1]), mean(Response[Treatment==2]), mean(Response[Treatment==3]),mean(Response[Treatment==4]))mu #作多重t检验。这里用到的pairwise.t.test()函数用来得到多重比较的p值 pairwise.t.test(Response, Treatment, p.adjust.method = "none")

#观察两个作调整后的p值的情况。p.adjust.method()函数的参数也可换为"hochberg","hommel","bonferroni","BH","BY","fdr"等。 pairwise.t.test(Response, Treatment, p.adjust.method = "holm") #绘制箱线图 plot(medicine$Response~medicine$Treatment)