先碼
WGCNA代碼1：WGCNA分析，簡單全面的最新教程
WGCNA代碼2：WGCNA實(shí)戰(zhàn)練習(xí)+離群點(diǎn)+隨機(jī)選取基因TOM作圖

含有離群點(diǎn)

clust = cutreeStatic(
  sampleTree,
  cutHeight = 15,
  minSize = 10)

keepSamples = (clust==1)
datExpr = datExpr0[keepSamples, ]
nGenes = ncol(datExpr)
nSamples = nrow(datExpr)
————————————————
版權(quán)聲明：本文為CSDN博主「weixin_43569478」的原創(chuàng)文章，遵循 CC 4.0 BY-SA 版權(quán)協(xié)議，轉(zhuǎn)載請附上原文出處鏈接及本聲明。
原文鏈接：https://blog.csdn.net/weixin_43569478/article/details/83747303

nSelect = 400
# For reproducibility, we set the random seed
set.seed(10);
select = sample(nGenes, size = nSelect);
selectTOM = dissTOM[select, select];
# There’s no simple way of restricting a clustering tree to a subset of genes, so we must re-cluster.
selectTree = hclust(as.dist(selectTOM), method = "average")
selectColors = moduleColors[select];
# Open a graphical window
sizeGrWindow(9,9)
# Taking the dissimilarity to a power, say 10, makes the plot more informative by effectively changing
# the color palette; setting the diagonal to NA also improves the clarity of the plot
plotDiss = selectTOM^7;
diag(plotDiss) = NA;
TOMplot(plotDiss, selectTree, selectColors, main = "Network heatmap plot, selected genes")

WGCNA代碼3：WGCNA相關(guān)性分析+模塊和性狀
WGCNA代碼4：WGCNA實(shí)例分析及解讀+檢驗(yàn)選定的β值下記憶網(wǎng)絡(luò)是否逼近 scale free+hub gene

k <- softConnectivity(dataExpr,power=sft$powerEstimate)
sizeGrWindow(10, 5)
par(mfrow=c(1,2))
hist(k)
scaleFreePlot(k,main="Check Scale free topology\n")
#可以看出k與p(k)成負(fù)相關(guān)(相關(guān)性系數(shù)0.87),說明選擇的β值能夠建立基因無尺度網(wǎng)絡(luò)

k

根據(jù)性狀與模塊特征向量基因的相關(guān)性及pvalue來挖掘與性狀相關(guān)的模塊

cor_ADR <- signif(WGCNA::cor(traitData,MEs,use="p",method="pearson"),5)
p.values <- corPvalueStudent(cor_ADR,nSamples=nrow(traitData))
Freq_MS_max_cor <- which.max(abs(cor_ADR["Insulin_ug_l",-which(colnames(cor_ADR) == "MEgrey")]))
Freq_MS_max_p <- which.min(p.values["Insulin_ug_l",-which(colnames(p.values) == "MEgrey")])

根據(jù)基因網(wǎng)絡(luò)顯著性，也就是性狀與每個基因表達(dá)量相關(guān)性在各個模塊的均值作為該性狀在該模塊的顯著性，顯著性最大的那個模塊與該性狀最相關(guān)：

GS1 <- as.numeric(WGCNA::cor(traitData[,'Insulin_ug_l'],dataExpr,use="p",method="pearson"))
GeneSignificance <- abs(GS1)
ModuleSignificance <- tapply(GeneSignificance,net$colors,mean,na.rm=T)

尋找與該性狀相關(guān)的樞紐基因（hub genes）,首先計(jì)算基因的內(nèi)部連接度和模塊身份，內(nèi)部連接度衡量的是基因在模塊內(nèi)部的地位，而模塊身份表明基因?qū)儆谀膫€模塊。
WGCNA代碼5：小張聊科研