GEO芯⽚分析中的⼤坑，差异基因完全相反！
接着讲了，差异分析的另外⼀种⽅式，以及多组如何同时做差异分析。
GEO芯⽚如果超过了两组，也可以⼀次搞定差异分析
今天我们讲⼀下，配对样如何做差异分析。该部分内容取材于这个诚意GEO帖⼦。
最有诚意的GEO数据库教程
如果不配对，分析应该是这样的
library(limma)
## 如果不是配对的
group 'before',18),rep('after',18))
group group,levels = c("before","after"),ordered = F)
design group)
fit fit2
最终可以通过添加p.value直接提取差异基因
allDiff1=topTable(fit2,adjust='fdr',coef=2,number=Inf ,p.value=0.05)
有6796个差异基因，部分数据如下
如果要是配对样本，这个例⼦中是1:18和19:36依次配对。差异分析这样做
第⼀种形式
pairinfo = factor(rep(1:18,2))
group 'before',18),rep('after',18))
group group,levels = c("before","after"),ordered = F)
design group+pairinfo)
fit fit2
通过添加配对信息pairinfo到model.matrix函数中实现。
最终获取差异基因有两种⽅式：
allDiff2=topTable(fit2,adjust='fdr',coef=2,number=Inf ,p.value=0.05)
allDiff21=topTable(fit2,adjust='fdr',coef="groupafter",number=Inf ,p.value=0.05)怎么做数据分析
差异基因是7635，多出来的基因就是因为加了配对信息才显⽰出差异的，我们后⾯作图来展⽰。第⼆种形式
如果采取昨天提到的第⼆种⽅法，可以这样做
pairinfo = factor(rep(1:18,2))
group 'before',18),rep('after',18))
group group)
design 0 + group +pairinfo)
colnames(design)[1:length(levels(group))] group)
contrast.matrix ## 线性拟合
fit fit fit2 ## 提取差异结果，注意这⾥的coef是1
allDiff3=topTable(fit2,adjust='fdr',coef=1,number=Inf ,p.value=0.05)
最终获得的结果是⼀样的。
那么现在的问题是，这些配对的数据如何展⽰，最好的办法是每个配对的样本间给出连线。
我们可以通过ggplot2来实现。
⾸先，创建```ggplot2```喜欢的数据格式：
pairinfo = factor(rep(1:18,2))
group 'before',18),rep('after',18))
group group,levels = c("before","after"),ordered = F)
data group=group,t(exprSet))
使⽤ggplot2来作图
⾼表达的
library(ggplot2)
## ⾼表达
ggplot(data, aes(group,CYP1B1,fill=group)) +
geom_boxplot() +
geom_point(size=2, alpha=0.5) +
geom_line(aes(group=pairinfo), colour="black", linetype="11") +
xlab("") +
ylab(paste("Expression of ","CYP1B1"))+
theme_classic()+
theme(legend.position = "none")
低表达的
library(ggplot2)
ggplot(data, aes(group,CH25H,fill=group)) +
geom_boxplot() +
geom_point(size=2, alpha=0.5) +
geom_line(aes(group=pairinfo), colour="black", linetype="11") +
xlab("") +
ylab(paste("Expression of ","CH25H"))+
theme_classic()+
theme(legend.position = "none")
这些基因配不配对，都有差异，所以我们出那个相差的1000多个墙头草基因来看⼀下。genelist
选前4个基因来分别普通作图和配对作图⽐较⼀下
## 第⼀个⽅式
plotlist <- lapply(genelist[1:4], function(x){data =data.frame(data[,c("pairinfo","group")],gene=data[,x])ggplot(data, aes(group,gene,col=group)) +geom_jitter(size=2 })library(cowplot)p1 <- plot_grid(plotlist=plotlist, ncol=4,labels = LETTERS[1:4])
## 第⼆种⽅式plotlist <- lapply(genelist[1:4], function(x){data =data.frame(data[,c("pairinfo","group")],gene=data[,x])ggplot(data, aes(group,gene,fill=group)) +geom })library(cowplot)p2 <- plot_grid(plotlist=plotlist, ncol=4,labels = LETTERS[1:4])plot_grid(p1,p2,ncol=1)
很明显，差异分析的配对作图效果要好很多啊。
好了，明天我们要扩展⼀下这个贴⼦中，对因⼦的level进⾏排序的那个技能，到底还有什么实⽤的功能。
GEO芯⽚分析中的⼤坑，差异基因完全相反！
这是果⼦开启的长达6个⽉连续写作计划的第4天，
我们明天见。
发布评论