首页
学习
活动
专区
工具
TVP
发布
精选内容/技术社群/优惠产品,尽在小程序
立即前往

差异分析③

  • 统计差异基因数目
代码语言:javascript
复制
tfit <- treat(vfit, lfc=1)
dt <- decideTests(tfit)
summary(dt)
        BasalvsLP BasalvsML LPvsML
Down        1417      1512    203
NotSig     11030     10895  13780
Up          1718      1758    182

一些研究需要不止一个调整后的p值cutoff值。 为了对重要性进行更严格的定义,可能需要log-fold-change(log-FC)超过最小值。 一般用来计算经验贝叶斯慢化t-统计的p值,并具有最小的log-FC要求。

  • 保存文件
代码语言:javascript
复制
de.common <- which(dt[,1]!=0 & dt[,2]!=0)
length(de.common)


vennDiagram(dt[,1:2], circle.col=c("turquoise", "salmon"))
write.fit(tfit, dt, file="results.txt")
#使用topTreat输出差异基因信息
#The top DE genes can be listed using topTreat for results using treat
# (or topTable for results using eBayes). 
#By default topTreat arranges genes from smallest to largest adjusted p-value with associated gene information, 
#log-FC, average log-CPM, moderated t-statistic, 
#raw and adjusted p-value for each gene. 
#The number of top genes displayed can be specified, where n=Inf includes all genes. 
basal.vs.lp <- topTreat(tfit, coef=1, n=Inf)
basal.vs.ml <- topTreat(tfit, coef=2, n=Inf)
head(basal.vs.lp)

维恩图显示仅比较基础与仅LP(左),基础与仅ML(右)之间比较基因DE的数量,以及两个比较(中心)中DE的基因数目。 任何比较中不是DE的基因的数目标记在右下角。

  • 差异基因可视化

为了总结目测所有基因的结果,可以使用plotMD函数生成显示来自线性模型的log-FC与平均对数-CPM值拟合的均值 - 差异图,其中突出显示差异表达的基因。

代码语言:javascript
复制
plotMD(tfit, column=1, status=dt[,1],
       main=colnames(tfit)[1], 
       xlim=c(-8,13))
  • 使用Glimma生成交互式均值差分图。

log-FC与log-CPM值显示在左侧面板中,与右侧面板中选定基因的每个样品的单个值相关。 结果表也显示在这些图下方,以及搜索栏以允许用户使用可用的注释信息来查找特定的基因。

代码语言:javascript
复制
glMDPlot(tfit, coef=1, status=dt,
         main=colnames(tfit)[1],
         side.main="ENTREZID",
         counts=x$counts,
         groups=group, launch=T)
  • 热图

使用来自gplots软件包的heatmap.2函数,从基础对比LP对比度的顶部100个DE基因(按调整的p值排列)创建热图。热图将样品按细胞类型正确聚类,并将基因重新排列成具有相似表达模式的区块。从热图中,我们观察到ML和LP样品的表达对于基础和LP之间的前100个DE基因非常相似。

代码语言:javascript
复制
library(gplots)
basal.vs.lp.topgenes <- basal.vs.lp$ENTREZID[1:100]
i <- which(v$genes$ENTREZID %in% basal.vs.lp.topgenes)
mycol <- colorpanel(1000,"blue","white","red")
heatmap.2(v$E[i,], scale="row",
          labRow=v$genes$SYMBOL[i], labCol=group,
          col=mycol, trace="none", density.info="none", 
          margin=c(8,6), lhei=c(2,10), dendrogram="column")
举报
领券