还是老习惯，给出官网教程，至于你是看还是不看，它就在那里，等着你的深入研究～
https://horvath.genetics.ucla.edu/html/CoexpressionNetwork/Rpackages/WGCNA/Tutorials/
如有理解错误，还请各位大神批评指正！
WGCNA分析图文详解专题中要解释的第二张图：
image.png
官方注释：
Figure 5：Clustering dendrogram of genes, with dissimilarity based on topological overlap, together with assigned merged module colors and the original module colors.
这张图有三个部分：
- 1，聚类树：那么你对聚类方法又了解多少？
- 2，Dynamic Tree Cut：根据某一height，聚类树分成了多少类（也即模块）
- 3，Merged dynamic：合并后的聚类模块，那么又是根据什么指标要将Dynamic Tree Cut中的类别进行合并得到新的聚类模块呢？
#前期代码
softPower = 6;
adjacency = adjacency(datExpr, power = softPower);
TOM = TOMsimilarity(adjacency);
dissTOM = 1-TOM
根据最佳软阈值和表达谱数据首先计算得到一个邻接矩阵adjacency，然后将邻接矩阵转换为TOM矩阵再计算一个距离矩阵dissTOM，这个距离矩阵就是后面用来画聚类树的图。
物以类聚，人以群分。物体与物体之间进行聚类总是根据某一相似度指标来进行刻画的，那么这里的距离矩阵dissTOM就是这个指标，只不过是距离（非相似度）。
而聚类的方法非常多，官方教程中使用的是层次聚类。
We now use hierarchical clustering to produce a hierarchical clustering tree (dendrogram) of genes. Note that we use the function hclust that provides a much faster hierarchical clustering routine than the standard hclust function.
使用R中的函数hclust()
geneTree = hclust(as.dist(dissTOM), method = "average");
采用的聚类方法为平均距离（method = "average"），
关于这个method，查看帮助文档，其实还有很多选项，大家可以根据结果进行调整。画出来的第一部分是这个样子的：
image.png
图片中，每一个树枝代表一个gene，纵坐标的Heigth为聚类距离。
设定一个聚类距离阈值，就可以将gene聚成module，这里又有很多方法设定这个阈值。而官网教程成给定的是dynamicTreeCut包中的函数cutreeDynamic，并且设定了最小模块中包含的基因个数为30个基因。
# We like large modules, so we set the minimum module size relatively high:
minModuleSize = 30;
# Module identification using dynamic tree cut:
dynamicMods = cutreeDynamic(dendro = geneTree, distM = dissTOM,
deepSplit = 2, pamRespectsDendro = FALSE,
minClusterSize = minModuleSize);
table(dynamicMods)
image.png
这里可以看到，总共得到了22个模块（1-22），以及每个模块中基因的数量。0为没有聚到任何模块中的基因。
由于Dynamic Tree Cut方法可能会产生两个模块，但是这两个模块的基因的表达谱是高度相似的，因此有必要对模块进行合并。为了定量每一个模块的共表达相似性，这里计算了每一个模块的eigengenes，即每个模块的特征值（这也是WGCNA分析中一个非常重要的概念，具体是什么，下次再进行说明～）。
# Calculate eigengenes
MEList = moduleEigengenes(datExpr, colors = dynamicColors)
MEs = MEList$eigengenes
# Calculate dissimilarity of module eigengenes
MEDiss = 1-cor(MEs);
# Cluster module eigengenes
METree = hclust(as.dist(MEDiss), method = "average");
# Plot the result
sizeGrWindow(7, 6)
plot(METree, main = "Clustering of module eigengenes",
xlab = "", sub = "")
这里也是直接对22个模块和一个灰色模块（就是没有聚到任何类别的基因的集合）进行了聚类树的绘制，然后设置聚类树的高度为0.25（下图中的红色阈值线），将22个模块合并成了18个模块。这里对模块进行合并采用的是模块特征值之间的皮尔森相关系数计算的距离矩阵，即两个模块之间的相关性大于0.75就合并为一个新的模块。
image.png
最后将聚类树，原始的模块和新的合并后的模块画在一张图里面就是最开始那张图了。
参考资料：
1,https://horvath.genetics.ucla.edu/html/CoexpressionNetwork/Rpackages/WGCNA/Tutorials/
2,A General Framework for Weighted Gene Co-Expression Network Analysis, Stat Appl Genet Mol Biol. 2005;4:Article17. Epub 2005 Aug 12
网友评论