查看: 46628|回复: 39

[软件使用] 如何绘制“上流”的共发生网络图?

  [复制链接]

管理员

Rank: 15Rank: 15Rank: 15Rank: 15Rank: 15

主题
672
注册时间
2020.6.16
在线时间
388 小时

发表于 2020.6.19 10:45:16 | 显示全部楼层 |阅读模式
在之前《如何绘制漂亮的网络图?》一文中我们已经为大家介绍过如何使用Gephi绘制下图这样的共发生网络图(Co-occurrence network),该文主要是介绍通过准备“边文件”和“点文件”来创建网络图。



但有时候绘图比较容易,难的是准备数据。前文也提过,最常见的共发生网络分析过程是先计算相关系数矩阵,接着使用igraph生成gml格式的网络图文件,最后导入Gephi对网络图进行个性化绘制。当然,如果你不愿使用R语言计算,也可以使用OmicShare Tools中的组内相关性分析工具计算得到相关系数矩阵,然后整理成“边文件”,使用Gephi或Cytoscape做个性化绘图。



工具网址:
https://www.omicshare.com/tools/

但如果数据量较大,手工整理显然非常费时费力,有OmicShare用户在OmicShare Forum中留言:“数据的准备太难了!”因此,下面我将从最初的OTU绝对丰度表开始,介绍如何使用R语言进行网络图数据的准备、网络图文件生成和网络图绘制。本文的重点内容是数据的过滤、相关系数的计算“graphml”格式网络图文件的生成

一、数据导入与过滤


为了便于大家练习使用,本文所用到的范例数据已上传到OmicShare论坛上,大家可前去下载。

下载链接:
https://www.omicshare.com/forum/thread-6116-1-1.html

[AppleScript] 纯文本查看 复制代码
#载入所需R包;
library(igraph)
library(Hmisc)
#工作目录设置;
setwd("C:\\Users\\MHY\\Desktop\\Co-occurrence_Network\\example")
#读入OTU绝对丰度表;
otu=read.table("otu_data.xls" ,header=T,row.names = 1,sep = "\t")
#转成矩阵,之后的相关性计算需要矩阵对象;
otu<-as.matrix(otu)
dim(otu)
head(otu)




[AppleScript] 纯文本查看 复制代码
#将丰度值大于1的值替换为1,便于计算不同otu的检测率;
dt<-otu
dt[dt>1]<-1
#将样本发现率低于20%的otu过滤掉;
no<-which(rowSums(dt)/ncol(dt)>0.2)
length(no)
otu<-otu[no,]


二、相关性系数计算

主要利用Hmisc包的rcorr()函数计算spearman相关系数,当然,也可计算pearson相关系数。rcorr()函数默认是对列进行两两计算,故须转置,计算otu间的相关性

[Python] 纯文本查看 复制代码
#计算相关性系数;
sp.cor<-rcorr(t(otu),type="spearman")

#提取r、p值矩阵;
r.cor<-sp.cor$r
p.cor<-sp.cor$P

#使用Benjamini-Hochberg("FDR-BH")法进行多重检验校正;
p.adj <- p.adjust(p.cor, method="BH")

#指定阈值;
r.cutoff=0.6
p.cutoff=0.001

#在只考虑“正相关”情况下;
r.matrix<-r.cor
p<-p.adj 
#将矩阵中不符合条件的r值替换为0;
r.matrix[which(r.cor <= r.cutoff)]=0
r.matrix[which(p.adj>p.cutoff)]=0

#删掉相关系数矩阵数据全都为0(对角线处的1不计)的行和列;
r.matrix<-r.matrix[which(rowSums(r.matrix)!=1),]
r.matrix<-r.matrix[,which(colSums(r.matrix)!=0)]

#查看过滤后的矩阵;
dim(r.matrix)
r.matrix[1:7,1:7]




三、生成网络图

从结构来看,对称矩阵是展示网络图非常合适的数据结构,比如相关系数矩阵是两个OTU的丰度相关系数,而网络图恰是Source和Target两个OTU结点的连线,相关系数可作为权重。如果是n维矩阵(即n个OTU),除去对角线,共有n(n-1)/2个关系对,当然用组合公式也可算得。

[Python] 纯文本查看 复制代码
#使用邻接矩阵(即相关系数矩阵)创建网络;
g1<-graph.adjacency(r.matrix,weight=T,mode="undirected")
#去掉冗余的边(multiple edges、loop edges);
g1<-simplify(g1)

#生成网络图的结点标签(OTU id)和degree属性;
V(g1)$label <- V(g1)$name
V(g1)$degree <- degree(g1)

#查看网络图的对象结构;
print(g1)




第一行:‘U’ 表示 undirected ,相反, ‘D’ 表示 directed 网络图;‘N’ 表示网络图具有name属性;‘W’表示weighted graphs;接下来的‘-’表示type属性不存在;接下来是网络共有240个结点和1068条边。

四、导出graph对象

[Python] 纯文本查看 复制代码
#将网络图导出为"graphml"、"gml"格式,方便导入Gephi中使用;
write_graph(g1, "g1.graphml", format = "graphml")
write_graph(g1, "g1.gml", format = "gml")
#支持的格式有"edgelist", "pajek", "ncol", "lgl","graphml", "dimacs", "gml", "dot", "leda"), ...


导出"graphml"格式的网络对象,之后可用Gephi和Cytoscape等网络图可视化软件直接打开,做进一步调整、美化。

五、使用igraph绘制网格图

虽然直接使用igraph绘制个性化网络图没有Gephi简单,igraph主要用于计算网络图的拓扑性质(topological properties),我这里仍然做了一些尝试。

[Python] 纯文本查看 复制代码
#计算群体结构(short random walks);
c <- cluster_walktrap(g1)
#使用默认颜色列表;
V(g1)$color <- c$membership+1
#绘制网络图;
layout <- layout.fruchterman.reingold
#也可以尝试其他的Layout;
layout <- layout_as_tree
layout <- layout_nicely
layout <- layout_on_sphere
#……
#绘制网络图;
plot.igraph(g1,vertex.size=4, 
            vertex.label=NA,
            edge.curved=F,
            edge.size=1.5,
            layout = layout.fruchterman.reingold)
#默认的结点大小为15,这里改为4,默认边的粗细为1,这里改为1.5;edge.curved若改为TURE,可实现类似文章开头那样的效果。


绘制的效果如下:



当然如果想做个性化调整,使用Gephi直接打开igraph导出的“g1.graphml”文件对网络图进行可视化,如下,具体的使用方法可参考《如何绘制漂亮的网络图?》一文。




使用Gephi导出网络图的效果如下,使用同样的数据和布局方式,我们可以看到使用igraph和Gephi得到网络图虽有差异,但在模块上还是比较一致的



这次就分享到这里啦!

参考资料:
https://igraph.org/r/

本文作者:基迪奥-莫北


               



本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
新的一天加油!
回复

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
4
注册时间
2016.4.19
在线时间
37 小时

发表于 2020.6.19 12:28:55 | 显示全部楼层
1.请问“将样本发现率低于20%的otu过滤掉”指的是什么意思?
回复 支持 反对

使用道具 举报

中华鲟

Rank: 5Rank: 5

主题
3
注册时间
2019.3.2
在线时间
226 小时

发表于 2020.6.19 13:06:52 | 显示全部楼层
学习学习
回复

使用道具 举报

帝王蝶

Rank: 4

主题
0
注册时间
2016.8.30
在线时间
176 小时

发表于 2020.6.19 15:33:50 | 显示全部楼层
这个好,学习了
新的一天加油!
回复 支持 反对

使用道具 举报

管理员

Rank: 15Rank: 15Rank: 15Rank: 15Rank: 15

主题
218
注册时间
2017.7.3
在线时间
715 小时

活跃会员荣誉管理


发表于 2020.6.19 17:52:56 | 显示全部楼层
chenye2015 发表于 2020.6.19 12:28
1.请问“将样本发现率低于20%的otu过滤掉”指的是什么意思?

例如,假设共有100个样本,如果otu1、otu2、otu3只在19个样本中能检测到,那么这三个otu被过滤掉,不参与后面网络图绘制
新的一天加油!
回复 支持 反对

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
5
注册时间
2020.4.11
在线时间
8 小时

发表于 2020.6.22 09:28:43 | 显示全部楼层
学习!!
加油
回复

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
4
注册时间
2016.4.19
在线时间
37 小时

发表于 2020.6.22 20:33:16 | 显示全部楼层
谢谢,请问#在只考虑“正相关”情况下;是指如果OTU之间是负相关的话这个代码没法表现出来么?
回复 支持 1 反对 0

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
1
注册时间
2016.4.8
在线时间
1168 小时

发表于 2020.6.25 15:15:53 | 显示全部楼层
一直提示没有jiava
回复 支持 反对

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
1
注册时间
2016.4.8
在线时间
1168 小时

发表于 2020.6.25 15:19:28 | 显示全部楼层
bucuoshibucuo ,jiushibuhuiyong
新的一天加油!
回复 支持 反对

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
1
注册时间
2016.4.8
在线时间
1168 小时

发表于 2020.6.25 15:21:59 | 显示全部楼层
buhuiyong
回复

使用道具 举报

管理员

Rank: 15Rank: 15Rank: 15Rank: 15Rank: 15

主题
218
注册时间
2017.7.3
在线时间
715 小时

活跃会员荣誉管理


发表于 2020.6.28 09:37:47 | 显示全部楼层
chenye2015 发表于 2020.6.22 20:33
谢谢,请问#在只考虑“正相关”情况下;是指如果OTU之间是负相关的话这个代码没法表现出来么? ...

需要考虑一下r的阈值正负或者取下绝对值
新的一天加油!
回复 支持 反对

使用道具 举报

管理员

Rank: 15Rank: 15Rank: 15Rank: 15Rank: 15

主题
218
注册时间
2017.7.3
在线时间
715 小时

活跃会员荣誉管理


发表于 2020.6.28 09:40:10 | 显示全部楼层
qinyuan 发表于 2020.6.25 15:15
一直提示没有jiava

gephi和Cytoscape都需要安装JRE
新的一天加油!
回复 支持 反对

使用道具 举报

草履虫

Rank: 2

主题
0
注册时间
2018.11.20
在线时间
1 小时

发表于 2020.6.29 05:32:32 | 显示全部楼层
请问 为什么我的OTU表转成矩阵后样本名称都无法显示 > dim(otu) [1] 4373    0 > head(otu)                                                                                                                                                                                                                        OTU483,0,0,0,0,0,0,0,0,0,0,0,0,0,0,3,0,0,11,0,0,11,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0                                       OTU523,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,24,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,0,6,0,0,0,0,0,0,0,0,0,0,0,0,7,0,6,0,0,0,0,15,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0                                       OTU842,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,4,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0                                         OTU1498,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,13,28,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0                                      OTU1149,35,284,2,30,70,31,41,43,404,29,114,72,0,0,0,0,1,0,0,0,4,0,0,0,383,440,245,122,313,151,0,0,133,0,0,0,21,29,37,6,7,4,90,109,12,79,79,108,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 OTU836,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,50,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0                                 
回复 支持 反对

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
1
注册时间
2016.4.8
在线时间
1168 小时

发表于 2020.6.30 09:30:34 | 显示全部楼层
xiazaibuxialai
回复 支持 反对

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
1
注册时间
2018.7.8
在线时间
14 小时

发表于 2020.6.30 16:47:01 | 显示全部楼层
真高级
回复

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
0
注册时间
2018.9.19
在线时间
6 小时

发表于 2020.7.1 08:28:53 | 显示全部楼层
好棒好棒好棒
回复 支持 反对

使用道具 举报

帝王蝶

Rank: 4

主题
0
注册时间
2016.12.6
在线时间
220 小时

发表于 2020.7.7 11:22:17 | 显示全部楼层
代码真令人头大
回复 支持 反对

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
1
注册时间
2016.4.8
在线时间
1168 小时

发表于 2020.7.9 09:42:04 | 显示全部楼层
真希望我也能做出这样的图
回复 支持 反对

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
1
注册时间
2016.4.8
在线时间
1168 小时

发表于 2020.7.9 19:57:16 | 显示全部楼层
酷酷酷
回复

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
1
注册时间
2016.4.8
在线时间
1168 小时

发表于 2020.7.10 07:54:49 | 显示全部楼层
好棒好棒
回复

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

快速回复 返回顶部 返回列表