微精选
Hi大家好，我是老米！生信入门一个月，感谢生信技能树平台。这是我的第3篇Markdown。不知道多少人还记得我的前两个投稿：
原来一个星期真的可以零基础入门TCGA数据挖掘，甚至markdown写作公众号投稿
一个基因引发的血案
现在继续完成生信技能树的作业题：重复一篇WGCNA分析的文章。WGCNA学习教程：一文学会WGCNA分析。
复现的文献：A novel microglial subset plays a key role in myelinogenesis in developing brain。EMBO J，好杂志。
小胶质细胞是神经系统中的巨噬细胞，可介导脑稳态和神经炎症。作用至关重要，但在大脑发育中具体机制尚未完全阐明。这篇文章中，通过对比健康成人和炎症激活(构建的实验性自身免疫性髓鞘炎模型EAE)的细胞群，发现新生鼠的小胶质细胞对于髓鞘生成及神经形成起关键作用。其中CD11c +小胶质细胞亚群在发育中大脑的髓鞘区域高表达，介导神经和胶质细胞的存活，迁移和分化。这些细胞是IGF1的主要来源，而选择性耗竭会该亚群可导致性髓鞘形成受损。 引用Jimmy的话：
值得一提的是，在单细胞水平研究小神经胶质细胞（Microglia）动态发育和异质性已经有了不少研究。
波士顿儿童医院的研究者们分析了超过76,000个来自于发育、衰老和脑部感染后的小鼠脑部的小胶质细胞，结果表明至少有9种转录特异的小胶质细胞形态，它们可以表达特定的基因集，且位于特定的脑区。发表于免疫学杂志Immunity， doi:10.1016/j.immuni.2018.11.004 (2019).
斯坦福大学医学院的研究者采用高深度scRNA测序揭示了小胶质细胞和脑髓细胞的发育异质性，发表于Neuron，这些细胞取自于胚胎期、出生后早期和成年的小鼠不同脑区。我们发现大部分的成年小胶质细胞表达稳定的基因（homeostatic genes），且不同脑区间没有差异。相反，出生后早期的小胶质细胞异质性更高。doi:10.1016/j.neuron.2018.12.006 (2019).
德国弗莱堡大学医学院神经病理学研究所的研究者采用单细胞RNA测序揭示小鼠和人的小神经胶质细胞的空间和时间异质性，成果最近以Letter的形式发表于Nature杂志。doi:10.1038/s41586-019-0924-x (2019).
这些都是当前的脑发育研究的热点和前沿！（Jimmy老师，你知识怎么能这么渊博？我这个搞神经发育的都惭愧【捂脸】）
该文章对五个处理组，共17个老鼠：
orange represents neonatal CD11c+ microglia (n = 4),
green neonatal CD11c microglia (n = 4),
blue EAE CD11c+ microglia (n = 3),
purple EAE CD11c microglia (n = 3),
black adult microglia (n = 3).
其实就是 neonatal(新生的) 和 EAE的Microglia，还有CD11C阳性和阴性，然后和成年小鼠的Microglia进行比较。这些样本进行了RNASeq测序，数据在GEO可供下载：/geo/query/acc.cgi?acc=GSE78809。当然我偷了个懒，该文章Supplemental material提供了整理之后的csv矩阵，大概1万3千个基因。
作者进行了WGCNA分析，以此证明小胶质细胞的作用及寻找关键信号通路及基因（嗯，估计这个实验室后续还有大文章出来）。
组学及生信分析的篇幅占了整个文章的一半以上。在此，仅重复文章的Figure3。
其实我跌跌撞撞学了WGCNA半个月，大概了解了这个原理，现在还一知半解，对很多细节摸索了很久，头很大。下面进入正题，请大家一起跟着我头大。
结果如图：
这结果挺好的，符合实验设计，都聚类好了。但是！但是如果不对原始数值进行log2处理，会出现一个outlier老鼠，如图：
看到没，那个新生的NEO1号鼠有点调皮！他超脱于所有老鼠之上，统领其他所有老鼠！这要么是问题少年，要么是智商超群，比如科大少年班的那种！而我大概深究了一下，控制这个离群值的大概是500个基因，哈哈哈。。。
专业名词叫多维标度法（Multidimensional Scaling）MDS，是一种维数缩减方法，把高维的数据点映射到一个低维的流形上；同时也是一种可视化方法，实践中通常利用2D或3D的MDS 结果观察（投影后）点的分布和聚集来研究数据的性质。实际上你也可以把它看做是PCA分析，反正都是看样本间的相关性计算距离远近。
说人话叫做：通过基因表达log2组间差异，看看这17个老鼠是否能够分到一个类别，比如成年鼠归到一类，EAE归到一类，新生鼠归到一类。和上面的有点类似，又不完全一样。
参考学习了stat的视频，也参考了他提供的代码。如下：
结果如下：
嗯，，，结果表示完美！别问我为什么要先做这个图，因为它好看，我忍不住！
点评：老米在这里开始秀技术了，实际上大家随便跑一下PCA分析即可，暂时无需深入写这么长代码！
下面正式进入WGCNA分析。
这个没啥好讲的，代码如下：
4. 采用一步法构建加权共表达网络
(weight co-expression network)
结果如下：
这里需要停顿一下：SCI文章只取到8个模块，图注说用了12691个基因（被我发现全部加起来也才9999，没错！9999个基因）。事实上，能够分这么少模块，我不知道作者参数细节。有两个参数特别影响模块大小：
一个是minModuleSize，表示每个模块里面最少放多少个基因，这很好理解，设定越大，模块越少；
另一个是mergeCutHeight，这个参数表示你在哪里砍树。值设定越小，树枝越多，通常是0.25。
google一下，发现这么一个说法：
这里，我的mergeCutHeight = 0.28，最终取到13个模块！
5. 模块与基因与表型
这是关键，其实就是开始解读这些关联的基因群对老鼠表型的影响，比如敲除CDC11之后的老鼠，哪些信号通路激活或抑制，EAE的老鼠哪些基因激活？等等等。
通过模块与各种表型的相关系数，可发现自己感兴趣的模块。如图：
可以看到有些模块与成年鼠关联，有些与EAE相关性很高。和原SCI的图有点像。
就是上面的图转化成bar图。
部分结果如下：
截取部分如下：
5.4 关注感兴趣的模块，导出基因进行GO分析
GO结果分析如下：
通过与上面的几个图表配合，发现EAE鼠与炎症通路/免疫相关，而新生鼠的胶质细胞则是cell cycle/胶质细胞分化/神经细胞发育/轴树突形成通路相关。与原SCI分析相符，与预期相符。
如下：
5.6 提取指定模块的基因并做热图
都挺好的。
点评：这个热图有问题，具体代码是 n=scale(t(dat+1)) 里面少了一个转置！
如图：
6. 模块导出
感兴趣的模块导出到cytoscape，VisANT等软件进一步进行可视化分析。
至此，复现结束。如果你仔细的看到了这里，说明你是想学WGCNA的了。回到课题本身，是否对小胶质细胞的认识更近了一步？提个小小的科学问题：在PD中，是否有一个胶质细胞亚群可对MPP+诱导的黑质神经细胞死亡有抵抗和保护作用？哈哈，别说我是搞肿瘤的了。。。
这些内容，换个思路应该够一个硕士毕业了。。。
写完又是深夜【捂脸】。还在持续学习中fighting。。。
能跟下来这篇教程的前提是你学习过R语言，否则一切都是镜中花水中月！
