net|R语言中的网络可视化( 二 )


如果跳到生物信息学领域 , 有一个 WGCNA 包用的特别多 , 但很多写教程的人都没搞清楚原理与模型假设 。WGCNA 包是基于巴拉巴西的无尺度网络构建的 , 基本原理就是先对所有基因构建两两相关性矩阵 , 然后从相关性矩阵中探索出共表达的基因模块 , 相当于把几千维的基因降维到不到十维且最好能联系上生物学意义例如某个通路啥的 。 具体计算则是每个模块进行主成分分析(其实是 SVD 分解) , 然后用第一个主成分作为这个模块的代表对你的研究分组进行差异分析 , 找出哪个模块有影响然后解释 。
这里面核心步骤里有两个坑 , 第一个在相关性矩阵到基因模块这里 , 第二个在主成分分析那边 。 第一个坑是因为其探索模块用了无尺度网络的假设 , 首先得去选一个幂级数来计算邻接矩阵 , 这个幂级数是拟合无尺度网络的度分布搞出来的 , 很多数据本身不符合无尺度网络的度分布假设 , 所以硬套这个假设是不合适的 。 第二个坑跟第一个有关系 , 只有模块内部第一个主成分可解释方差很高才能这么用 , 但由于第一个坑很多人用了默认值 , 第二个坑也就只用了第一个主成分 , 很多时候方差解释连三分之一都不到 , 虽然能讲故事 , 但明显是有偏的 。 当然这个包里也是涵盖了很多对于用户而言天书级的概念 , 很多人不求甚解套默认值也把文章给发了 , 完全就当神奇降维盒子在用了 。
其实说白了网络分析是另一层意义上的因子分析 , 起一个降维作用 , 只是降维方式不是简单的线性组合而是引入了图论的一些统计量罢了 , 但我看到很多人用起来套代码 , 解释上完全就是意会 , 特别是代谢组学分析里会出现套基因组学的分析方法而不验证假设盲目追求自动化的默认值 。 不过我也看到了很多基于图论的生物信息学文章 , 很多想法非常超前但引用非常少且真正生产实验数据的人基本看不懂或不会去读 , 这就还不如心理测量学那边研究人员的学科内科普做得好 。
这里说个题外话 , 非统计与计算机科学背景开发者写的包其实非常多 , 多数都不在 CRAN 上 。 他们一般只会用基础R包函数来实现自己想要的功能 , 如果没有就会去依赖其他包来调用函数 , 很少用 Rcpp, 不开并行计算 , 基本不关心速度 , 用 S3 对象而不是 S4, 这倒是科研编程的日常状态 。 有时候读他们的源码有种见字如面的感觉:有的人明显是其他语言转过来的 , 有次读一个包的代码怎么看怎么别扭 , 后来发现这个开发者的母语是 java, 很多定义方式都是那边继承过来的 。 有的人注释掉的代码中其实有另一重意思 , 跟最终版的实现方式完全不同 , 源码里保留了很多进化遗迹 , 类似化石 。 有的人严谨 , 每个函数都写测试 , 文档明显打磨过语言 。 有的人飘逸 , 通篇找不到注释 , 很多编码风格都不一致 , 感觉是爆栈网复制过来的 。 有的人很明显是 tidyverse 风格出现后才开始学的 R, 对基础函数用法非常不熟 。 非统计与计算机科学开发者的代码通常存在很多不严谨的地方 , 没有经过软件工程的训练 , 更多是为了解决特定目的而快速实现的 , 不过很多代码展示出了丰富多彩的想象力 。
下面生成几个网络并做下基础可视化 , 这里不会用最常见的那种点对点数据结构 , 因为这个东西是需要从原始数据生成的 , 很少有原始数据本身就是这种关系结构 , 而且此处也不涉及 ggplot2 风格的绘图包 , 用基础绘图系统来做 , 函数统一为 plot, 当然不同的对象类型会有不同的绘图参数 , 此处也不会区分有向图与无向图 , 按默认值来 。 理解清楚了这些最基础的代码实现过程 , 自定义出图上就会很轻松了 。
network 版

set.seed(110) library(network) # 生成一个3节点网络 net <- network.initialize(3) # 画出来 plot(net,vertex.cex=10)


推荐阅读