R语言实现时序RNA-seq分析

提到RNA-Seq差异表达分析，大家首先想到的癌症与癌旁组织的表达差异分析。然而如果想探究不同时间下对目标产生的影响，此方法便失去作用，那么便出现了时序RNA-seq。今天我们为大家介绍一个可以做时序RNA-seq分析的R包maSigPro。

首先我们看下其安装还是需要借助bioconductor库进行安装，具体步骤参见以前的教程，我们呢不在此赘述了。

接下来我们看下此包的主要功能框架图：

然后我们首先看下如何实现RNA-seq的分析,我们先看下我们用到的数据：

data(edesign.abiotic, edesignCT)

data(data.abiotic)

其中Time指的取样的时间点；Replication指的数据集所属的组别；control指的控制组是哪些数据集；后面就是对应的各组织或者实验，以0，1代表所属的数据集。

1. make.design.matrix构建回归模型，生成数据矩阵。其主要的参数是degree代表我们选择模型为几次方程。此处一般数据默认第一列Time；第二列Replication。

a=make.design.matrix(edesign.abiotic)#二次 b=make.design.matrix(edesignCT, degree = 3) # 三次多项式

2. p.vector 计算差异的基因

示例如下：

fit <- p.vector(data.abiotic,a, Q = 0.05, min.obs = 20)

fit$i # 差异基因的数量

fit$p.vector # 获取各个基因的P值

fit$SELEC # 差异基因表达值矩阵

3. T.fit 利用递归法进一步确证每个基因模型的可靠性。

其和p.vector的用法一样，可以直接套用，同时也可以直接利用p.vector产生的数据结果。

示例如下：

tstep <- T.fit(fit, step.method ="backward", alfa = 0.05)

4. get.siggenes 获取差异基因

其中最主要的是rsq（R的平方），其值介于0到1，越接近1模型越好。

Vars参数可以是groups，all，each。

示例如下：

sigs <- get.siggenes(tstep, rsq = 0.6, vars ="groups")

sigs$summary

sigs <- get.siggenes(tstep, rsq = 0.6, vars = "all")

sigs$summary

sigs <- get.siggenes(tstep, rsq = 0.6, vars = "each")

sigs$summary

5. suma2Venn 数据的可视化，韦恩图。

suma2Venn(sigs$summary[, c(2:4)])

6. PlotGroups 对某个基因在所有组表达情况或者模型的可视化

示例如下：

STMDE66 <- data.abiotic[rownames(data.abiotic) =="STMDE66",]

PlotGroups(STMDE66, edesign = edesign.abiotic)

PlotGroups(STMDE66, edesign = edesign.abiotic, show.fit =T, dis = a$dis, groups.vector = a$groups.vector)

7. see.genes 基因的聚类分析

其中主要的参数cluster.method包括"hclust","kmeans"和"Mclust"。

示例如下：

see.genes(sigs$sig.genes$ColdvsControl, main ="ColdvsControl", show.fit = T, dis =a$dis,cluster.method="kmeans" ,cluster.data = 1, k = 9)#查看ColdvsControl差异基因的表达谱情况。