在转录组数据分析过程中，我们最常做的是不同处理方式的样本之间的比较（Treated vs Control），这时候我们采用“DEG分析+pathway分析”的方式就可基本完成对数据的分析。但偶尔我们也会碰到一类特殊的数据，即同一种处理在不同的时间点可能呈现不同的生物学特征，因此我们需要分析不同基因随着时间的动态变化。这一类分析即所谓的“时序分析”。
那么如何对转录组数据进行时序分析呢？这里小编给大家介绍一款简单好用的R包——maSigPro包。
该包的分析流程如下图所示：
该包最初由Conesa等人于2006年开发，于2014年经Nueda等人进一步优化。maSigPro采用二步回归的策略对数据进行分析，并且可以找出不同处理组之间的差异基因。并且在第二次回归模型中可以将具有相似表达模式的基因进行聚类并可视化。
需要注意的是，该包不会对数据进行标准化处理，因此需要自己提前进行预处理哦。
其分析代码如下：
1. 加载maSigPro包及演示数据
library(maSigPro)# load maSigPro library
data(data.abiotic)
data(edesign.abiotic)
数据格式如下：
> edesign.abiotic
> data.abiotic
2. 定义回归模型
design <- make.design.matrix(edesign.abiotic, degree = 2)
其中degree参数与时间点的数量有关，时间点较多时degree数值应相应增加（该演示数据中为3个时间点）。design为列表形式，其中group.vector可以看回归参数是如何分配的。
> design$groups.vector
[1] "ColdvsControl" "HeatvsControl" "SaltvsControl" "Control" "ColdvsControl" "HeatvsControl" "SaltvsControl" "Control"
[9] "ColdvsControl" "HeatvsControl" "SaltvsControl"
3. 差异基因分析
通过p.vector()函数对每一个基因进行回归分析，同时给出每一个基因的p-value及p.adjust（MT.adjust、Q等参数）。
fit <- p.vector(data.abiotic, design, Q = 0.05, MT.adjust = "BH", min.obs = 20)
由此返回以下参数：
> fit$i # returns the number of significant genes
> fit$alfa # gives p-value at the Q false discovery control level
> fit$SELEC # is a matrix with the significant genes and their expression values
4. 差异分析
通过T.fit()函数可以进一步分析不同组别的差异基因之间的不同特征，该步采用逐步回归的方式（step.method参数）。
tstep <- T.fit(fit, step.method = "backward", alfa = 0.05)
5. 提取差异基因列表
sigs <- get.siggenes(tstep, rsq = 0.6, vars = "groups")
其中rsq为回归模型R-squared的cutt-off值，vars用于表示分组变量。
6. Venn图查看差异基因情况
suma2Venn(sigs$summary[, c(2:4)])
suma2Venn(sigs$summary[, c(1:4)])
7. see.genes()
采用see.genes()查看差异基因结果
> sigs$sig.genes$SaltvsControl$g
[1] 191
该函数可以将表达模式相似的基因进行聚类并可视化。其中k参数为cluster的数量，cluster.method为聚类的方式（hclust，kmeans，Mclust）
see.genes(sigs$sig.genes$ColdvsControl, show.fit = T, dis =design$dis, cluster.method="hclust" ,cluster.data = 1, k = 9)
8. PlotGroups()
使用PlotGroups()函数可以查看某一特定基因的表达情况。
STMDE66 <- data.abiotic[rownames(data.abiotic)=="STMDE66", ]
PlotGroups (STMDE66, edesign = edesign.abiotic)
并且我们可以通过以下方式添加回归曲线：
PlotGroups (STMDE66, edesign = edesign.abiotic, show.fit = T, dis = design$dis, groups.vector = design$groups.vector)
怎么样，是不是很简单呢？快来试一试吧。
网友评论