RNA-Seq基因富集的简单解决方法

RNA-Seq已经成为一种对于转录组研究很重要的技术，因为他能够同时的定量基因表达数据同时还能对每个基因的可变剪切进行识别。所以就有越来越多的方法及工具针对于RNA-Seq数据的差异表达分析以及功能分析，旨在发现一些新的生物学及生物医疗方面的进展。但是这些方法都主要针对于RNA-Seq数据的差异表达分析，而忽略了可变剪切这块的信息。

SeqGSEA就是一个综合了基因差异表达及可变剪切的分析方法的R包，同时还提供了这些基因的GSEA功能分析。

由于RNA-Seqcount数据本身的特性，R包使用的是负二项分布作为背景分布（如DESeqR包）。然后在将两个测度结合起来共同评价基因的差异情况，并提供GSEA分析。

R包的使用方法如下：

首先建立newReadCountSet，利用到的数据就是count数据，exonID数据，geneID数据，其中count、exonID、geneID数据形式如下：

例子数据，其中表示G1基因有4个外显子分别是E1-E4，G2基因有6个外显子E1-E6。如果是自己的数据，就可以利用这样的方式建立newReadCountSet。下面使用R包自带数据进行演示。

以RCS_example为例，其中一共包含20个样本，5000个特征，也就是5000个exon的read count信息。

首先需要将不表达或者表达量太少的EXON和GENE去除掉，使用下列语句：

将count数据进行过滤后，就可以直接使用下列语句进行直接的基因表达差异分析及可变剪切差异分析：

则可以获得如下分析结果：

这是经过负二项分布标准化后的每个特征的表达，而我们还需要确定结果是否具有统计学意义，所以需要利用随机扰动来确定P值。使用如下命令：

将结果进行随机1000次。

通过随机的分布估计出标准化Factor，在通过normFactor将count值进行标准化。

通过刚刚构建的随机矩阵计算DS的显著性，如上图所示，第三列为显著性P值，第四列为FDR矫正后P值。

而对于DE的计算则几乎沿用了DESeq这个R包的计算方法：

与DS计算几乎相似的步骤就能计算出DE差异显著性。获得了DS和DE差异，下一步则是将两个得分进行综合，如下图所示：

利用线性组合的方式，DE权重为0.3，我们可以对比出随机的结果与真实的结果的差异如上图所示。根据不同的结果组合方式，将会有不一样的结果，如下为rank组合方式：

有了这个综合分数，我们就可以进行GSEA的分析，如下图所示：

这样我们就能获得GSEA的分析结果，如上图所示，可以知道与哪些基因集相关。