jamesbang

V1

2023/02/11阅读：21主题：雁栖湖

🤩 WGCNA | 值得你深入学习的生信分析方法！~（网状分析-第二步补充-大数据的网络构建与模块识别）

1写在前面

`WGCNA`的包内其实也提供了解决方案，基本思想是`分级聚类`。🧐

1️⃣ 首先，我们使用快速但相对粗糙的聚类方法，用于将基因预聚类成大小接近的`模块`，且不超过你所设定的基因`最大值`。😂

2️⃣ 然后我们分别在每个`模块`中执行完整的网络分析。🤠

3️⃣ 最后，合并特征基因高度相关的模块。😏

2用到的包

``rm(list = ls())library(WGCNA)library(tidyverse)``

3示例数据

``load("FemaleLiver-01-dataInput.RData")``

4软阈值

4.1 topology analysis

``powers <-  c(c(1:10), seq(from = 12, to=20, by=2))sft <-  pickSoftThreshold(datExpr, powerVector = powers, verbose = 5)``

4.2 可视化

``sizeGrWindow(9, 5)par(mfrow = c(1,2))cex1 = 0.9plot(sft\$fitIndices[,1], -sign(sft\$fitIndices[,3])*sft\$fitIndices[,2],     xlab = "Soft Threshold (power)",ylab="Scale Free Topology Model Fit,signed R^2",     type="n", main = paste("Scale independence"))text(sft\$fitIndices[,1], -sign(sft\$fitIndices[,3])*sft\$fitIndices[,2],labels=powers,cex=cex1,col="red")abline(h=0.90,col="red")plot(sft\$fitIndices[,1], sft\$fitIndices[,5],     xlab="Soft Threshold (power)",ylab="Mean Connectivity",      type="n", main = paste("Mean connectivity"))text(sft\$fitIndices[,1], sft\$fitIndices[,5], labels=powers, cex=cex1,col="red")``

5构建网络与模块识别

5.1 网络构建

``bwnet <-  blockwiseModules(  datExpr, maxBlockSize = 2000,  power = 6, TOMType = "unsigned", minModuleSize = 30,  reassignThreshold = 0, mergeCutHeight = 0.25,  numericLabels = T,  saveTOMs = T,  saveTOMFileBase = "femaleMouseTOM-blockwise",  verbose = 3)``

5.2 查看模块数

``table(bwnet\$colors)``

6对比结果

6.1 载入之前的结果

``load(file = "FemaleLiver-02-networkConstruction-auto.RData")``

6.2 匹配颜色

`colors``label`匹配起来，嘿嘿。🤨

``bwLabels <-  matchLabels(bwnet\$colors, moduleLabels);bwModuleColors <-  labels2colors(bwLabels)``

6.3 分次结果可视化

``sizeGrWindow(6,6)# Block 1plotDendroAndColors(bwnet\$dendrograms[[1]], bwModuleColors[bwnet\$blockGenes[[1]]],"Module colors", main = "Gene dendrogram and module colors in block 1",dendroLabels = F, hang = 0.03,addGuide = T, guideHang = 0.05)# Block 2plotDendroAndColors(bwnet\$dendrograms[[2]], bwModuleColors[bwnet\$blockGenes[[2]]],"Module colors", main = "Gene dendrogram and module colors in block 2",dendroLabels = F, hang = 0.03,addGuide = T, guideHang = 0.05)``

7对比两种方法的结果差异

7.1 对比一下

``sizeGrWindow(12,9)plotDendroAndColors(geneTree,                    cbind(moduleColors, bwModuleColors),                    c("Single block", "2 blocks"),                    main = "Single block gene dendrogram and module colors",                    dendroLabels = F, hang = 0.03,                    addGuide = T, guideHang = 0.05)``

7.2 对比eigengenes

``singleBlockMEs <-  moduleEigengenes(datExpr, moduleColors)\$eigengenesblockwiseMEs <-  moduleEigengenes(datExpr, bwModuleColors)\$eigengenes``

`match`之后看一下结果，嘿嘿。🫶

``single2blockwise <-  match(names(singleBlockMEs), names(blockwiseMEs))signif(diag(cor(blockwiseMEs[, single2blockwise], singleBlockMEs)), 3)``

8如何引用

📍
```Langfelder, P., Horvath, S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics 9, 559 (2008). https://doi.org/10.1186/1471-2105-9-559 ```

📍 往期精彩

jamesbang
V1

wx🔍: Grassssss 卷起来了