安装WGCNA
BioManager::install("impute")
BioManager::install("preprocessCore")
install.packages("WGCNA")
安装过程可能会不太顺利，大家可以自行查找解决方案。
导入表达数据
library(WGCNA)
exp0 <- read.delim(test1, row.names=1)
有很多基因是不表达或者表达了很低，所以我们过滤掉在所有样品中表达之和小于5的基因
exp <- exp0[rowSums(exp0) > 5, ]
我们的表达矩阵通常每一行是一个基因，每一列是一个样本；但WGCNA要求每一行是一个样本，每一列是一个基因，所以我们需要进行转置
Texp0 = t(exp0)
寻找最佳 β 值
# 开启多线程模式
enableWGCNAThreads(nThreads = 20)
# 通过对 power 的多次迭代，确定最佳 power
powers <- c(1:20)
x = pickSoftThreshold(Texp0, powerVector = powers, verbose = 5,
networkType = "signed")
# 计算出来的最佳 β 存放在
x$powerEstimate
# 画图
## 设置画布
sizeGrWindow(9, 5)
par(mfrow = c(1,2))
cex1 = 0.9
# Scale independence
plot(x$fitIndices[,1], -sign(x$fitIndices[,3])*x$fitIndices[,2],
xlab="Soft Threshold (power)",ylab="Scale Free Topology Model Fit,signed R^2",type="n",
main = paste("Scale independence"));
# 增加一条阈值线
abline(h=0.8,col="red")
# Mean connectivity
plot(x$fitIndices[,1], x$fitIndices[,5],
xlab="Soft Threshold (power)",ylab="Mean Connectivity", type="n",
main = paste("Mean connectivity"))
