0

我正在使用WGCNA包进行网络分析，步骤如下：

  1. 数据输入
  2. 生成模块
  3. 获取基因id
  4. 表型 x 模块相关性

我想使用该程序包将表型数据与基因表达矩阵一起包含在内，以查找哪些基因与表型分组。然后，我想获取感兴趣的模块并制作网络图并检查哪些基因与表型相关。

我生成了如下模块：

library(WGCNA)
options(stringsAsFactors = FALSE)
enableWGCNAThreads()
lnames = load(file = "dataInput.RData");

# Choose a set of soft-thresholding powers
powers = c(c(1:10), seq(from = 12, to=40, by=2))
# Call the network topology analysis function
sft = pickSoftThreshold(datExpr, powerVector = powers, verbose = 5)
# Plot the results:
sizeGrWindow(9, 5)
par(mfrow = c(1,2));
cex1 = 0.9;
# Scale-free topology fit index as a function of the soft-thresholding power
plot(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");
# this line corresponds to using an R^2 cut-off of h
abline(h=0.90,col="red")
# Mean connectivity as a function of the soft-thresholding power
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")


##Constructing the gene network and identifying modules is now a simple function call:
net_unsigned = blockwiseModules(datExpr, power = 6,
                                  TOMType = "unsigned", minModuleSize = 30, maxBlockSize = 300,
                                  reassignThreshold = 0, mergeCutHeight = 0.25,
                                  numericLabels = TRUE, pamRespectsDendro = FALSE,
                                  saveTOMs = TRUE,
                                  saveTOMFileBase = "PopulusTOM_signed",
                                  verbose = 5)
##maxBlockSize = The total number of genes you have in your gene expression matrix that passed the filter from Data_Input scrip
##Plotting graph
pdf("Dendogram_Modules_signed.pdf", width = 30, height = 30);
##Convert labels to colors for plotting
mergedColors = labels2colors(net_unsigned$colors)
##Plot the dendrogram and the module colors underneath
plotDendroAndColors(net_unsigned$dendrograms[[1]], mergedColors[net_unsigned$blockGenes[[1]]],
                    "Module colors",
                    dendroLabels = FALSE, hang = 0.03,
                    addGuide = TRUE, guideHang = 0.05)
dev.off()
##Save
moduleLabels = net_unsigned$colors
moduleColors = labels2colors(net_unsigned$colors)
MEs = net_unsigned$MEs;
geneTree = net_unsigned$dendrograms[[1]];
save(MEs, moduleLabels, moduleColors, geneTree,
     file = "unsigned-networkConstruction-auto.RData")

这会生成模块，然后我将模块与一种表型相关联。如何将表型数据包含在基因表达中？谢谢！

4

0 回答 0