WGCNA中的eigengene有什么重要意义呢?

经过软阈值确定,TOM矩阵计算,层次聚类和动态切割后,
我们获得了一些基因模块,每个模块里面都可以看成是有共同目的的团伙。

为了方便后面纷至沓来的相关性分析,尤其是跟性状的联合分析,现在需要选取出每个模块中的代表人物。
这个过程就是寻找每个模块的eigengene(读音: 爱根金)

这样,模块信息就变成了,行是eigengene,列是样本的矩阵了,十分简化。

什么是eigengene

按照作者的说法,就是每个模块基因和样本矩阵进行PCA分析后的第一主成分。

WGCNA中提供了简单的方法来计算。
就是moduleEigengenes这个函数,名字也比较好记,就是模块的Eigengene。
输入表达量数据和模块的信息就可以。

moduleEigengenes(datExpr, moduleColors)

返回的是个列表,可以把我们需要的提取出来。

MEs0 = moduleEigengenes(datExpr, moduleColors)$eigengenes

行是样本,列是模块的eigengene

有了这个就可以跟性状信息求相关性,得到WGCNA最重要的热图了。

现在的问题是,这个值是怎么算出来的。
我们以第一列黑色模块来举例。
我们重复一下,是主成分分析的第1个主成分。

那就先把黑色模块提取出来,提取模块也是个重要操作。
因为后面确定了模块,肯定就需要提取模块

table(moduleColors)

黑色模块含有211个基因,提取表达矩阵。

datablack <- datExpr[,moduleColors == "black"]

行是样本,列是基因

有一些缺失值,先用knn技术把缺失值填充一下。

### 缺失值填充
### 要求就是行是基因,列是样本
datModule = t(as.matrix(datablack))
datModule = impute::impute.knn(datModule, k = min(10, nrow(datModule)-1))
## 提取数据,行是样本,列是基因
datModule_knn = t(datModule$data)

主成分分析很简单

pc <- prcomp(datModule_knn)

提取第一个主成分

pc$x[, 1]

跟官方计算的比较一下,把MEs0的第一列提取出来看一下

[1]  0.0138038628  0.0668404876  0.0667271646 -0.0643323004  0.0637495925 -0.0008963038

发现不怎么一样。
但是你如果计算一下相关性,会发现新的结果。

cor(MEs0[,1],pc$x[, 1])

相关性极其高

           [,1]
[1,] -0.9994166

说明两种有不可告人的关系。

奇异值分解

进过阅读文献,查看源码得知,作者在计算eigengene的时候用的是奇异值分解法。
这里是让我十分沮丧的,因为奇异值分解用到了高数矩阵的知识,我看了很多资料,大脑就是缺失一环。
所以,当前我并不能深入浅出的讲清楚这个事情,搞懂是迟早的,只是当前优先度不够。

但是,我可以做出来。
首先奇异值分解,数据是需要缩放的,就用scale函数
scale函数的理解,可以看看这一篇。
z-score的标准化究竟怎么弄?

datModule_scale= scale(datModule_knn)

那么复杂的技能奇异值分解,就是一个单词的事情

svd = svd(datModule_scale)

而结果里面的左边奇异值的第一列就是我们要的结果

svd$u[,1]

数据证明是数值是对的,只是符号不一样。

而主成分得到的结果和分解奇异值得到的结果高度相关,只是数值不同

cor(pc$x[, 1],svd$u[,1])
          [,1]
[1,] 0.9994166

但是现在的问题是,我们算出来的结果跟官方的有出入,正负相反。
他们处理的方法是,计算每个样本中,基因的平均值,用这个平均值跟主成分求相关性
如果相关性是负的,那么就把主成分的结果乘以负一。

### 计算平均值
averExpr = rowMeans(datModule_scale, na.rm = TRUE)   
### 求相关性
corAve = cor(averExpr,svd$u[,1], use = "p")
corAve 

确实是负数

          [,1]
[1,] -0.998023

明明可以用主成分算的事情,为什么要用奇异值分解?
又是算法层面的解释:

PCA当样本非常多的时候,计算协方差矩阵,还要进行特征值分解,计算量很大,而有些奇异值分解方法可以跳过这一步直接得到结果。

eigengene的意义在哪?

第一,性状关联。
模块里面的基因,经过压缩,变成了一个eigengene,这样跟性状信息可以求相关性。
有了相关性,我们就可以得到性状和模块的相关性热图,这也是WGCNA里面最重要的一张图。

通过这个图,可以选取感兴趣的模块进一步研究。
那模块中的基因怎么提取呢?我们已经展示了。
提取出的基因可以画热图,GO分析,KEGG分析。

第二,hubgene
eigengene的存在,让找每个模块中的是hub gene变得简单。
hubgene,通俗的解释就是

买个模块中跟eigengene 相关性最强的基因。

那些cytoscape画的网络图都是次要的,因为找hub gene不需要画图,直接计算就行了。

第三,迭代WGCNA

这里面来了一个新的概念就是迭代WGCNA,我们后面会介绍。
eigengene的出现,我们筛选模块里面高价值的基因就有了标准,我们可以把模块中的每个基因都跟eigengene计算相关性,把小于0.8的全部剔除(人为定义)
每个模块都这么做,现在就到了一个新的行是基因,列是样本的矩阵。
重新来做一次WGCNA,再来筛选一次,直到无法筛选基因为止。
此时,我们就得到了都比较纯净模块。而有些文章,模块基本上都是模糊的,胆子大,心也大。

 

生物信息学

ArchR如何克服scATAC-seq数据的稀疏性?

2020-6-26 20:56:30

生物信息学

生存分析为什么不需要矫正P值

2020-6-28 23:09:21

0 条回复 A文章作者 M管理员
    暂无讨论,说说你的看法吧
个人中心
购物车
优惠劵
今日签到
有新私信 私信列表
搜索