RNA-seq报告解读：差异表达基因分析

利用聚类方法计算样本和样本的距离，从而对样本之间的相似性进行考察。该分析结果能准确反映实验设计情况，属于同一设计的样本距离相近，优先聚在一起。根据基因表达情况获得测序样本聚类图如下图所示：

1、差异基因筛选

利用DESeq软件对各个样本基因的counts数目进行标准化处理（采用basemean值来估算表达量），计算差异倍数，并采用NB（负二项分布检验的方式）对reads数进行差异显著性检验，最终根据差异倍数及差异显著性检验结果来筛选差异蛋白编码基因。   在利用RNA-seq数据比较分析两个样品中同一个基因是否存在差异表达的时候，可以选取两个标准：一是FoldChange，就是两样品中同一个基因表达水平的变化倍数；二是pvalue或FDR（adjusted pvalue），FDR值的计算方法先要对每个基因进行p-value的计算，再用FDR错误控制法对p-value作多重假设检验校正。默认筛选差异的条件为p<0.05且差异倍数大于2。   

差异基因筛选结果文件各列说明见下表：

不同分组差异基因统计表见下表：

表格说明：
（1）Case：实验组名；
（2）Control：对照组名；
（3）Up_diff：显著性差异的上调基因数量；
（4）Down_diff：显著性差异的下调基因数量；
（5）Total_diff：显著性差异基因总数量。

差异表达基因统计柱状图

图片说明：横轴是各个比较组；纵轴是比较组的差异基因数量，其中Up是显著性差异的上调基因数量，Down是显著性差异的下调基因数量。

输出文件包括：（1）、各组差异筛选基因文件

（2）、差异未筛选基因文件

2、火山图

通过绘制火山图可以了解差异表达基因的整体分布情况。各组实验设计中筛选出来差异基因的火山图如下：

差异表达火山图

图片说明：将比较所产生的差异情况反映到火山图中，灰色为非显著性差异的基因，红色和绿色为显著性差异基因；X轴为log2 FoldChange的展示，Y轴方向为-log10 Pvalue的展示。

3、差异基因表达水平聚类分析

对差异表达基因进行非监督层次聚类。计算多个样品两两之间的距离，构成距离矩阵，合并距离最近的两类为一新类，计算新类与当前各类的距离，再合并、计算，直至只有一类为止，用挑选的差异基因的表达情况来计算样品直接的相关性，一般来说，同一类样品能通过聚类出现在同一个簇中，聚在同一个簇中的基因可能具有相似的生物学功能。

各组差异分组聚类分析结果如下图：

各差异分组聚类分析结果

图片说明：图中红色表示高表达基因，蓝色表示低表达蛋白编码基因。

4、差异基因GO富集分析

得到差异表达基因之后，我们对差异表达基因进行GO富集分析，对其功能进行描述（结合GO注释结果）。GO 功能富集分析的方法：将全部蛋白编码基因 / 转录本作为背景列表，差异蛋白编码基因 / 转录本列表作为从背景列表中筛选出来的候选列表，利用超几何分布检验计算代表 GO 功能集在差异蛋白编码基因 / 转录本列表中是否显著富集的 p 值，再对 p 值经Benjamini & Hochberg 多重检验纠正后得到 FDR。超几何分布检验计算 p 值的公式如下：

Enrichment score计算公式为：

其中，N为所有基因中具有GO注释的基因数目；n为N中差异表达基因中具有GO注释的基因数目；M为所有基因中注释为某特定GO term的基因数目；m为注释为某特定GO term的差异表达基因数目。可以根据GO分析的结果结合生物学意义从而挑选用于后续研究的基因。

GO富集性分析结果各列说明如下表所示：

GO富集分析top30 （筛选三种分类中对应差异基因数目大于2的GO条目，按照每个条目对应的-log10Pvalue由大到小排序的各10条），条形图展示如下：

GO富集分析结果展示

图片说明：图中X坐标为GO条目名称，Y坐标为-log10Pvalue。

使用fisher算法分别对样品间差异基因进行CC，BP，MF富集分析， 并使用topGO[7]对富集到的Term 绘制有向无环图。topGO有向无环图能直观展示差异表达基因富集的GO节点（Term）及其层级关系，是差异表达基因 GO 富集分析的结果图形化展示，分支代表的包含关系，从上至下所定义的功能描述范围越来越具体。

差异基因topGO有向无环图展示

图片说明：对每个GO term进行富集，最显著的10个节点用矩形表示。矩形的颜色代表富集显著性，从黄色到红色显著性越来越高。每个节点的基本信息显示在相应的图形中，为GO ID和GO term。

差异基因和所有基因在GO Level2水平分布比较图如下：

差异表达基因及所有基因在GO Level2水平分布比较图

图片说明：蓝色表示所有基因富集的GO Level2条目，红色表示差异基因富集的GO Level2条目，横轴为条目名称，纵轴表示对应条目的基因数量和其百分比。

上调差异基因和下调差异基因在GO Level2水平分布比较图如下：

上调差异基因和下调差异基因在GO Level2水平分布比较图

图片说明：红色表示上调差异表达基因富集的GO Level2条目，绿色表示下调差异表达基因富集的GO Level2条目，横轴为条目名称，纵轴表示对应条目的基因数量和其百分比。

5、差异基因KEGG富集分析

KEGG是有关Pathway的主要公共数据库，利用KEGG数据库对差异蛋白编码基因进行Pathway分析（结合KEGG注释结果），并用超几何分布检验的方法计算每个Pathway条目中差异基因富集的显著性。计算的结果会返回一个富集显著性的p值，小的p值表示差异基因在该Pathway中出现了富集。相应的计算公式参见GO富集分析。Pathway分析对实验结果有提示的作用，通过差异基因的Pathway分析，可以找到富集差异基因的Pathway条目，寻找不同样品的差异蛋白编码基因可能和哪些细胞通路的改变有关。

KEGG富集分析top20（筛选对应差异基因数目大于2的pathway条目，按照每个条目对应的-log10Pvalue由大到小排序）气泡图如下：

KEGG富集top20气泡图

图片说明：图中X轴Enrichment Score为富集分值，气泡越大的条目包含的差异蛋白编码基因数目越多，气泡颜色由紫-蓝-绿-红变化，其富集pvalue值越小，显著程度越大。

下载KEGG通路图片，把差异基因标注在通路图上，使用红色表示上调基因，绿色表示下调基因，黄色表示既有上调又有下调的基因。

KEGG网络通路图示例图如下：

KEGG通路网络图示例图

图片说明：KEGG通路图网上红色表示上调基因，绿色表示下调基因，黄色表示对应的基因既有上调又有下调的，把鼠标放在相应的蛋白编码基因上会自动显示注释上该蛋白编码基因的基因及对应的FC值和上下调信息，点击该蛋白编码基因可以进入KEGG数据库，查看详细信息。

输出文件：（1）差异基因KEGG富集结果

（2）差异基因KEGG注释结果

差异表达基因及所有基因在KEGG Level2水平分布比较图如下：

差异表达基因及所有基因在KEGG Level2水平分布比较图

图片说明：横轴是注释到各Level2代谢通路的基因（差异表达基因）和所有注释到KEGG通路的基因（差异表达基因）总数的比值（%），纵轴表示Level2 pathway的名称，柱子右边数字代表注释到该Level2 pathway下的差异表达基因数量。

上调差异表达基因及下调差异表达基因在KEGG Level2水平分布图如下：

上调差异表达基因及下调差异表达基因在KEGG Level2水平分布图

图片说明：横轴是注释到各Level2代谢通路的上调（下调）差异表达基因和所有注释到KEGG通路的上调（下调）差异表达基因总数的比值（%），纵轴表示Level2 pathway的名称，柱子右边数字代表注释到该Level2 pathway下的上调（下调）差异表达基因数量。