主讲老师第六十二课：WGNCA
当年号称一个分析就可以发一篇 3-4 分文章的存在，WGCNA 分析（WGCNA.zip），又称为加权基因共表达网络分析，该分析方法旨在寻找协同表达的基因模块(module)，并探索基因网络与关注的表型之间的关联关系，以及网络中的核心基因。
WGCNA 分析适用于复杂的数据模式，根据官方推荐，一般需要 5 组或者 15 个样品以上的数据。同时，其适用的研究方向包括不同器官或组织类型发育调控、同一组织不同发育调控、病原菌侵染后不同时间点应答等等。下面，我们一起来看一下如何进行 WGCNA 分析
1.R 包安装与数据读取
1.1 R 包安装
关于 WGCNA 的分析方法，官方存在专门的 R 包，即 WGCNA 包，该包储存在 CRAN上，直接使用install.packages()命令即可完成安装。此外，由于在分析过程中涉及多线程分析，需要安装 R 的preprocessCore 包
1.2 数据读取
接下来，是输入数据的准备工作，对于 WGCNA 分析，测序数据和芯片数据均适用。如果是转录组数据，最好是 RPKM，TPM 或者其他归一化后的表达量，而如果是芯片数据的话，那么使用常规的归一化矩阵即可。在此，我们使用转录组测序的 TPM 数据来进行展示。
接下来，我们可以查看一下具体的数据结构，该数据集由 407 例样本，8840 个基因组成。
2.检查样本和基因
准备好输入数据以后，我们需要进一步对输入数据进行检查
2.1 数据格式整理
由于 WGCNA 分析是基于基因表达的聚类分析，所以在分析之前，我们需要对其进行转置。
2.2 检查缺失值
此外，我们需要判断一下基因和样品质量的好坏，使用 goodSamplesGenes() 函数，可以检查数据中缺失的条目，权重低于阈值的条目以及零方差基因，并返回具有最大缺失值或低权重值标准的样本和基因列表可以看到，经过判断，结果显示 TRUE，表明其中不包含有低质量的样品或基因名。当然，如果gsg$allOK 的结果为 FALSE，经过取反以后，!gsg$allOK 的结果为 TRUE，因此需要进一步分别筛选基因和样品质量。接下来，我们需要进一步判断和筛选质量不行的样品和基因，最后，得到 407 个样品，8840 个基因的表达矩阵。
3.样品聚类
3.1 聚类
随后，我们使用 dist()函数计算样品之间的欧式距离，并通过聚类查看是否存在异常样品。
3.2 画图
运算完毕后，将得到的聚类结果进行可视化。
可以看到，最左边明显有几个样品和别人不在一条枝叶上，因此，我们需要将那几个明显偏离的分支给剪去。根据聚类树的结果，我们需要将异常样品进行去除，在此，将剪切线的地方定于高度为 155 的位置。
以剪切线所在位置作为分界，异常样品为 1 例。提取适合样品的表达信息，最终得到了406 例样品纳入
随后的分析。
4.临床数据的准备
接下来，我们需要根据样品的分组信息制作临床信息表格，当然，除了自己制作外，也可以从工作目录中提取准备好的表格。
根据正常和肿瘤样品信息，得到了新的数据框储存相应的临床信息由于前面以及删除了 1 例异常样品，我们通过 intersect()函数获取共同的样品名，接着通过取子集的方式，最终得到了 406 例样品的表达和临床矩阵
5.重新聚类
随后，我们对新得到的样品重新聚类根据临床分组信息，使用 number2colors()函数为临床信息创建颜色集，其中红色表示阳性，白色表示阴性。接着，使用 plotDendroAndColors()函数在树状图中绘制了层次聚类树状图以及基于颜色注释的临床分组特征。
6.构建网络，识别模块
(1). 选择合适的软阈值（soft thresholding power），即 beta 值，选取 1 到 20 为候选的 power值，使用 pickSoftThreshold()函数对相应的 power 值进行计算。
根据结果，进一步绘制 power 值对应的散点图。
1 到 2，2 到 3，整个变化趋势都还很明显，而 3 往后，随着 power 值的增大，其变化趋势已经不明显了，因此，选取 3 或者 4 为最佳的 power 值。
(2). 选取最佳的 power 值建立邻近矩阵，根据连接度使我们的基因分布符合无尺度网络(3). 计算拓扑矩阵（TOM 矩阵）(4). 基因聚类使用基因表达得到的 TOM 矩阵，再次对基因继续聚类(5). 动态剪切模块识别设定最小模块基因数量，对基因聚类结果进行剪切，得到不同的基因模块。 使 用cutreeDynamic()函数，最终得到了 13 个不同的基因模块接着，我们通过 labels2colors()函数，对每一个基因模块进行颜色赋值最后，基于 plotDendroAndColors()函数，将剪切结果添加到聚类树上，并使用前面赋予的模块颜色来进行区分。
7.相似模块聚类
得到不同的基因模块后，我们需要对表达模式相似的模块进行聚类当然，如果对聚类结果还需要手工调整模块的话，我们可以根据 Height，设定MEDissTheres 的阈值，对模块进行人工的合并，根据这个结果，不同模块之间的距离相对可以，就不进行合并了。
最后，根据人工设定的剪切高度，对相似的基因模块进行合并至此，得到了放到文章里的第一张图
8.模块与性状数据热图
随后，我们需要计算不同的基因模块与临床性状之间的相关性系数（cor）和对应的 P值（Pvalue）。
然后，根据相关系数和 p 值，将其使用热图进行展示这是放到文章里的第二张图，结果显示了每一个基因模块与正常或肿瘤样品之间的相关性结果。我们可以根据相关系数和 p 值，挑选与临床性状之间最显著相关的模块进行后续分析。
9.不同模块与临床性状的具体分析
尽管计算得到了临床性状与基因模块之间的相关性，可以挑选最相关的那些模块来分析，但是模块本身仍然包含非常多的基因，还需进一步的寻找最重要的基因，所有的模块都可以跟基因算出相关系数，所有的连续型性状也可以跟基因的表达值算出相关系数(1). 计算模块与基因之间的相关性矩阵(2). 计算临床性状与基因之间的相关性矩阵(3). 批量输出性状与模块散点图
把两个相关性矩阵联合起来,指定感兴趣模块进行分析
10.输出每个模块的基因
最后，我们需要把每个模块里面的基因进行输出，以用于后续的分析。
逐个输出每个颜色模块的基因，根据与临床性状之间的相关性，选取最相关的基因模块，对这个基因模块进行后续的一系列分析。到此，WGCNA 分析的主线部分基本就完成了，当然，WGCNA 分析不只是如此，还可以指定基因模块，绘制表达与性状之间的热图等等，这个感兴趣的小伙伴在学完之后可以查找相关的资料进行拓展延伸。