小云手把手教你用GS-MM散点图挖掘WGCNA的hub基因/SCI论文/科研/研究生/生信分析热点思

小云之间接触过WGCNA分析，但是浅尝辄止，只做到把基因分成模块，算出这些模块和基因之间的相关性及P值，就结束了，没有继续挖掘Hub基因，但小云也是个有求知欲的人，于是小云觉得还是要挖一下，于是就整理了下面这些代码。小云温馨提示，以下代码是在掌握了WGCNA分析模块基因的基础，即已经把基因分成模块并算出这些模块和基因之间的相关性及P值的前提下学习的。

library('WGCNA')

fpkm=read.table(“after_meger_mdd.txt”,header=T,row.names=1,sep=”\t”, comment.char="",check.names=F)

datExpr=as.data.frame(t(fpkm[,-1]))

names(datExpr)=fpkm$Tag #第一行第一列为Tag

rownames(datExpr)=names(fpkm[,-1]) #倒置表达矩阵，行为样本，列为基因

datExpr=read.table(“性状t”,header=T,row.names=1,sep=”\t”) #读取性状文件

nGenes = ncol(datExpr)

nSamples = nrow(datExpr)

# 指定datTrait中感兴趣的一个性状，这里选择SubB

SubB = as.data.frame(datTraits$SubB)

names(SubB) = "SubB"

#  各基因模块的名字（颜色）

modNames = substring(names(MEs), 3) #这里的MEs是在WGCNA分析区分模块之后形成的变量，里面记录了模块的信息。

# 计算MM的P值

geneModuleMembership = as.data.frame(cor(datExpr, MEs, use = "p"))

MMPvalue = as.data.frame(corPvalueStudent(as.matrix(geneModuleMembership

 ), nSamples))

names(geneModuleMembership) = paste("MM", modNames, sep="")

names(MMPvalue) = paste("p.MM", modNames, sep="")

# 计算性状和基因表达量之间的相关性（GS）

geneTraitSignificance = as.data.frame(cor(datExpr, SubB, use = "p"))

GSPvalue = as.data.frame(corPvalueStudent(as.matrix(geneTraitSignificance),

  nSamples))

names(geneTraitSignificance) = paste("GS.", names(SubB), sep="")

names(GSPvalue) = paste("p.GS.", names(SubB), sep="")

module = "green" #选择模块

column = match(module, modNames)

moduleGenes = moduleColors==module

green_module<-as.data.frame(dimnames(data.frame(datExpr))[[2]][moduleGenes])

names(green_module)="genename"

MM<-abs(geneModuleMembership[moduleGenes,column])

GS<-abs(geneTraitSignificance[moduleGenes, 1])

c<-as.data.frame(cbind(MM,GS)) #包含了MM和GS的数据，可以保留一下

rownames(c)=green_module$genename

green_hub <-abs(c$MM)>0.8&abs(c$GS)>0.2 #筛选hub基因

write.csv(green_hub, "hubgene_MMGS_green.csv")

verboseScatterplot(abs(geneModuleMembership[moduleGenes,

column]),abs(geneTraitSignificance[moduleGenes, 1]), xlab =

 paste("Module Membership in", module, "module"), ylab = "Gene

   significance for SubB", main = paste("Module membership

   vs. gene significance"), pch = 20,col="grey") #画散点图

abline(h=0.2,v=0.8,col="red",lwd=1.5) #添加参考线

    小伙伴们，看懂了没有，虽然看上去有点复杂，但一步一步的跟着小云走，其实也不难哦，小伙伴们有什么问题欢迎来和小云分享讨论哟。