WGCNA-最短教程掌握高级分析
- 2019 年 12 月 30 日
- 筆記
WGCNA的目的
先给大家讲讲WGCNA的精髓，其实就一句话：关联表型和基因。WGCNA通过将基因进行分组（module），把基因模块和表型进行关联，实现了快速锁定核心基因的目的。
WGCNA的步骤
WGCNA需要读入的数据有两个：
- 基因表达谱数据
- 样本表型数据
我们采用WGCNA的示例数据来给大家进行解读。整体来讲，WGCNA的分析流程是很繁琐的，一个全套的分析可能会涉及8-9个步骤，但是，这里面有很多步骤其实无关紧要，跟分析的主线，也就是“筛选与表型相关的核心基因”是脱离的。所以，在精简之后分为以下几个步骤：
01
数据预处理
这部分内容包括以下部分：
- 读取和过滤基因表达数据
- 读取样本表型数据
- 可视化样本聚类和表型数据
官方的示例数据是一个小鼠的芯片表达谱数据，包含了135个雌性小鼠的数据，在提供的表达谱数据中，除了探针ID和样本表达量之外，还有额外的探针注释信息，在读取原始数据时，需要把多余注释信息去除。
对于基因的表达量数据，需要进行过滤，对于基因而言，可以过滤缺失值或者低表达的基因，对于样本而言，如果该样本中基因缺失值很多，也需要过滤，WGCNA内置了一个检验基因和样本的函数，通过该函数可以进行一个基本过滤。基础过滤之后，还可以看下是否存在离群值的样本，通过样本的聚类树进行判断，生成的图片如下。
从图上可以看出，F2_221 这个样本和其他样本差距很大，可以将该样本过滤掉。表型数据中也包含了不需要的列，而且其样本比表达谱的样本多，需要根据表达谱的样本提取对应的表型数据。表达谱数据和表型数据准备好之后，可以绘制样本聚类树和表型的热图，生成的图片如下。
上半部分为样本的聚类树，下半部分为样本对应的表型的热图，顺序和聚类树中的顺序一致，表达量从低到高，颜色从白色过渡到红色，灰色代表缺失值。
02
识别modules
在构建共表达网络时，将基因间的相关系数进行乘方运算来表征其相关性，首先需要确定乘方的值，也就是软阈值或者我们叫power值、beta值。这个power值至关重要！至关重要！至关重要！那怎么选呢？在WGCNA作者的官方介绍中，一般选6就可以。
不过因为大家用的表达值都不一样，每次计算出来的阈值也不一样，所以，为了挑选出最佳的阈值，在示例代码中，我们将从1到20所有的整数作为阈值进行计算，去测试最优的阈值。具体怎么测试呢？其中powerEstimate就是最佳的power值，fitIndices保存了每个power对应的网络的特征，也就是每个power构建的相关性网络中的连接度的统计值，k就是连接度值，可以看到，对于每个power值，提供了max, median, max3种连接度的统计量，这里对连接度的均值进行可视化，生成的图片如下。
从这个图片中怎么选择power值？最直接的方法：看y轴等于0.9的时候，与曲线的交点，在这个示例数据中，刚好等于6，那我们就选6作为power值。切记：在分析自己数据的时候，power值肯定在6的附近，4-10之内都是可以的，如果离的太远，那就有问题了。之后，就可以用6作为beta值去构建基因module了，生成的图片如下。
上方为基因的聚类树，聚类时的距离为1-TOM值，下方为基因对应的modules。在WGCNA中，module的概念我们已经讲过了，其实就可以理解成是一个gene group的意思。在WGCNA当中，module的命名是有特定的规则的：按照颜色进行命名。比如blue可以是一个module，brown是另外一个module。
03
筛选和表型相关的module
本质上是计算module的ME值与表型的相关系数，可以对module和表型间的系数的结果进行可视化，生成的图片如下。
指定一个我们感兴趣的表型，可以得到与其相关性最高的module。在这个图中，红色表示正相关，绿色表示负相关，颜色越深表示相关性越强。比如，以上结果中，和weight这个表型最相关的为module为red和brown，当然也可以自己指定一个阈值，筛选出多个候选的modules。
04
生成网络，筛选核心基因
找到和表型关联的module之后，就可以导出指定modules对应的基因共表达网络，方便可视化。所谓WGCNA（共表达网络）分析，说到底还是要构建一个网络的。那么，既然说到网络，就不得不说最厉害的网络图分析软件Cytoscape。它能够帮我们去实现网络图的可视化和深层次的挖掘。在WGCNA中，也可以直接对接Cytoscape，生成的文件可以导入cytoscape进行绘图。比如我们选定brown这个模块，就可以直接生成如下两个文件：
CytoscapeInput-edges-brown-red.txt
CytoscapeInput-nodes-brown-red.txt
在Cytoscape中直接导入这两个文件就可以做图了，然后通过Cytohubba、MCODE等一系列插件，可以筛选我们感兴趣的核心基因了。至于Cytoscape的操作，我们这里就不再展开，在之前的文章（从网络图探寻基因互作的蛛丝马迹）中我们详细的讲过Cytoscape的使用。
那么，到这里为止，一个WGCNA的分析就完成了。难吗？其实也就那么回事。导入数据->选择power值->划分module->挑选module->绘制网络，结束！当你把一个分析真正弄懂了之后，你就会知道里面的精髓和不可或缺的部分。就像一个经验丰富的外科主刀医生，下刀“快准狠”，绝不会有多余的动作，也绝不会少必须的动作。做数据分析也是一样，抓住关键问题才是真理。
