您好,登录后才能下订单哦!
这篇文章主要介绍“WGCNA怎么挖掘潜在的共表达基因”,在日常操作中,相信很多人在WGCNA怎么挖掘潜在的共表达基因问题上存在疑惑,小编查阅了各式资料,整理出简单好用的操作方法,希望对大家解答”WGCNA怎么挖掘潜在的共表达基因”的疑惑有所帮助!接下来,请跟着小编一起来学习吧!
共表达基因指的是表达量具有协同变化趋势的基因集合,通常认为这些基因参与相同的生物学过程,比如参与同一个代谢通路,正是由于功能上的协同作用,导致表达量呈现出高度相关性。
在WGCNA中,对传统的相关系数进行乘方运算,用最终得到的值来表征基因间的相关性。在计算出这样的相关性统计量值之后,如何确定哪些基因是共表达的呢?
WGCNA的做法是聚类分析,聚类分析属于一种非监督的机器学习算法,通过聚类树,可以观察到哪些基因在聚类树中属于同一分支,属于同一分支的基因可以归为一类。实际操作中,考虑到基因数目较多等情况,肯定需要算法来自动化的进行分类,WGCNA采用的是dynamicTreeCut
这个R包。
对于聚类算法而言,需要输入基因间的距离矩阵,首先就需要将基因间的邻接矩阵转换为距离矩阵,对相关系数进行乘方运算,可以计算出邻接矩阵,但是这个值本质上反映的是基因间的相似度,并不是距离。在计算距离矩阵时,WGCNA采用了TOM
这种统计量,该统计量可以表征网络中节点的相似性,计算公式如下
对于两个基因i
和j
而言,a
表示两个基因邻接矩阵中对应的值,就是相关系数的乘方,K
代表的每个基因的连接度, 公式如下
对于加权网络而言,就是该节点的边对应数值的总和,比如在网络中基因A与3个基因相连,基因A的连接度就是对应3条边的数值之和。两个基因间的l
值代表的是两个基因所有边的权重乘积的总和,公式如下
公式只是帮助我们理解计算的过程,其实只需要理解TOM
是表征节点的相似度就行,我们要的是距离,所以直接用1减去相似度即可,公式如下
借助TOM
值,将基因间的相关系数转换为了距离,然后就可以用该距离矩阵进行聚类。上述的计算方法在WGCNA中都有对应的公式,代码如下
# 确定乘方运算中power的最佳取值 powers <- c(c(1:10), seq(from = 12, to=20, by=2)) sft <- pickSoftThreshold( datExpr, powerVector = powers, verbose = 5) softPower <- sft$powerEstimate # 计算邻接矩阵 adjacency <- adjacency(datExpr, power = softPower) # 计算TOM相似度矩阵 TOM <- TOMsimilarity(adjacency) # 计算距离矩阵 dissTOM <- 1-TOM # 聚类 geneTree <- hclust(as.dist(dissTOM), method = "average")
根据聚类结果和距离矩阵,就可以调用dynamicTreeCut
的算法来识别modules, 代码如下
# 指定每个module中基因数目的最小值 minModuleSize <- 30 # 识别modules dynamicMods <- cutreeDynamic( dendro = geneTree, distM = dissTOM, deepSplit = 2, pamRespectsDendro = FALSE, minClusterSize = minModuleSize)
通过table
函数可以查看modules的结果,用法如下
> table(dynamicMods) dynamicMods 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 88 614 316 311 257 235 225 212 158 153 121 106 102 100 94 91 78 76 65 58 20 21 22 58 48 34
可以看到,识别出22个modules, 0
代表那些没有归入任何modules的基因。通过plotDendroAndColors
函数可视化聚类树对应的modules, 代码如下
dynamicColors = labels2colors(dynamicMods) plotDendroAndColors( geneTree, dynamicColors, "Dynamic Tree Cut", dendroLabels = FALSE, hang = 0.03, addGuide = TRUE, guideHang = 0.05, main = "Gene dendrogram and module colors")
生成的图片如下
整个图片分为两个部分,上方为基因的聚类树,下方为识别到的modules, 不同的modules对应不同的颜色,其中灰色对应那些没有归入任何modules的基因。
通过dynamicTreeCut
识别到modules之后,还会结合每个modules的基因表达量数据,来识别相关性很高的modules, 从而进行合并,其原理是对modules进行聚类,每个module下的基因表达量是一个二维矩阵,做相关性分析我们只需要一个一维向量就可以了,可以利用PCA分析提取第一主成分来表征原始的矩阵,在WGCNA中,把每个module的表达谱数据对应的一维向量称之为Module eigengene E
。获取一维向量之后,就可以计算相关性,直接用1减去相关性作为距离,进行聚类,代码如下
MEList <- moduleEigengenes( datExpr, colors = dynamicColors) MEs <- MEList$eigengenes MEDiss <- 1-cor(MEs) METree <- hclust(as.dist(MEDiss), method = "average") plot(METree, main = "Clustering of module eigengenes", xlab = "", sub = "")
modules的聚类树示意如下
每个modules的名字用对应的颜色表示,在该聚类数中,分支长度为1减去两个module间的相关系数,在合并modules时,将高相关性的合并为一类,可以指定一个阈值,比如将相关系数大于0.8的合并为一类,在该聚类树中,对应的就是height小于0.2的modules, 对应下图红色的线
可以看到有8个modules都满足条件,在合并时,会将原本属于同一分支的modules直接合并为一个,从图上可以看出,合并后会减少4个modules。合并的代码如下
MEDissThres <- 0.2 merge <- mergeCloseModules( datExpr, dynamicColors, cutHeight = MEDissThres, verbose = 3) mergedColors <- merge$colors mergedMEs <- merge$newMEs plotDendroAndColors( geneTree, cbind(dynamicColors, mergedColors), c("Dynamic Tree Cut", "Merged dynamic"), dendroLabels = FALSE, hang = 0.03, addGuide = TRUE, guideHang = 0.05)
合并之后的modules 对应的图片如下
到此,关于“WGCNA怎么挖掘潜在的共表达基因”的学习就结束了,希望能够解决大家的疑惑。理论与实践的搭配能更好的帮助大家学习,快去试试吧!若想继续学习更多相关知识,请继续关注亿速云网站,小编会继续努力为大家带来更多实用的文章!
免责声明:本站发布的内容(图片、视频和文字)以原创、转载和分享为主,文章观点不代表本网站立场,如果涉及侵权请联系站长邮箱:is@yisu.com进行举报,并提供相关证据,一经查实,将立刻删除涉嫌侵权内容。