ggClusterNet---一条代码完成全内容微生物网络

写在前面

今年的毕业生应该是6月底离校，但是这位朋友到了九月份终于完善了手上的工作，对这几年来了一个完善的结束，感谢送给我的一方印章，也祝铠鸣大展宏图。再有就是ggclusterNet功能也完善了90%，剩下的慢慢补充上来，尤其是代码的优化要抽出些时间多和大佬交流，优化。要知道这部分实在是也很头疼。原计划今年将微生物组分析的四个R包够写好，想来想去似乎目前只有easystat和ggClusterNet可以顺利完成，我也会抓紧时间在今年结束的时候至少完成一个包（phyloBiota）的功能。

导入所需要的R包

#--导入所需R包#-------
library(phyloseq)
library(igraph)
library(network)
library(sna)
library(ggplot2)
library(ggClusterNet)

数据导入

我在ggClusterNet包中加入了一个phyloseq封装的微生物组数据，这份数据来自刘永鑫老师的示例文件。这里直接使用data()提取即可

data(ps)

网络分析一网打尽

为什么说是一网打尽呢，我使用了network函数打包了ggclusterNet封装的各种功能，一般情况下网络所需要的内容都基本上囊括了，一条代码即可分析完成网络，并且可以按照分组同时分析完成多组网络。完成多组网络的共同比对和输出。

network 函数

  • 计算相关矩阵生成边和节点文件保存为Gephi格式；

  • 使用sne中十多种layout方法并用ggplotv出图，组图

  • 计算节点性质

  • 计算网络性质

  • 计算随机网络和样本网络比对，检测

  • 计算zipi

  • 分组网络计算自动保存先看看出来的结果

    函数测试

    这个例子表示的是提取微生物相对丰度在千分之一以上的OTU，，然后对其做 皮尔森相关（推荐使用斯皮尔曼相关），相关阈值选择大于0.6的，并且显著性阈值低于0.05的相关关系，并且使用ggrapg可视化网络，通过igraph计算网络节点和整体属性，后计算zipi模块化信息，将全部结果保存在path中。

    ```{R}path = “./result”dir.create(path)result = network(ps = ps,N = 0.001,r.threshold=0.6,p.threshold=0.05,label = FALSE,path = path ,zipi = TRUE)

除了保存的结果，这里还有全部分组的网络可视化输出，共比较，和全部网络全局参数一起。
```{R}
# 全部样本的网络比对
p = result[[1]]
p
# 全部样本网络参数比对
data = result[[2]]
plotname1 = paste(path,"/network_all.jpg",sep = "")
ggsave(plotname1, p,width = 22,height = 22)
plotname1 = paste(path,"/network_all.pdf",sep = "")
ggsave(plotname1, p,width = 22,height = 22)

tablename <- paste(path,"/co-occurrence_Grobel_net",".csv",sep = "")
write.csv(data,tablename)

结果展示-网络可视化

结果展示-Gephi输入的两个文件

上面的可视化如果你不喜欢，就可以在gephi中出图，使用下面的这两个文件。

  • KO_Gephi_edge.csv

  • KO_Gephi_node.csv

  • ···

    结果展示-网络性质同随机网络比对

  • KO_net_VS_erdos_properties.csv

  • WT_net_VS_erdos_properties.csv

  • ···

    结果展示-节点属性

  • KO_node_properties.csv

  • WT_node_properties.csv

  • ···

    结果展示-zipi计算结果和图形

  • KOZiPi.csv

  • KO_ZiPi.pdf

    结果展示-幂律分布

    这里看到不同网络的幂律分布是不同的，和随机网络之间的差距也不是很相同，这表明了不同分组的网络不一定都是无标度网络。

    二分网络一网打尽

    二分网络的介绍在之前我已经讲过了，参见推送：

    点击跳转。

    或者在专辑中查看。

    微生物结合环境因子等其他指标网络

    这里是我构造的一组模拟环境因子数据环境因子的格式和vegan中的一样，也就是样本作为列，指标作为行。

—导入多个环境变量及其分组

data1 <- read.delim(“..//ori_data/bios_network_normal//env.txt”) head(data1) row.names(data1) <- data1SampleID = NULL data1 <- as.data.frame(t(data1)) head(data1) data1$ID = row.names(data1) data1 <- data1 %>% select(ID,everything())
Gru <- read.delim(“..//ori_data/bios_network_normal//env_group.txt”)
细菌真菌phyloseq对象我们已经准备好了，这里可以参照学习， 如果不会构建phyloseq对象，可以参见之前的推文。```{R}# 导入细菌ps，通过相对丰度value1来过滤otu表格ps16 = readRDS("..//ori_data/bios-network_micro/16S/ps_16S.rds")#导入真菌ps，通过相对丰度value1来过滤otu表格psIT = readRDS("..//ori_data/bios-network_micro/ITS/ps_ITS.rds")ps.merge <- merge16S_ITS (ps16,psIT,N16s = 0.005,NITS = 0.005)ps.merge


同样是一条函数即可计算得到类似于上面的全部内容，由于zipi的计算不够快，所以我这里设置：zipi = FALSE。下面是结果文件的总览，少了zipi的部分：

library(ggpubr)
path = "./result_bio/"
dir.create(path)

result <- corBionetwork(ps = ps.merge,
N = 0.001,
r.threshold = 0.9, # 相关阈值
p.threshold = 0.05,
label = FALSE,
group = "Group",
env = data1, # 环境指标表格
envGroup = Gru,# 环境因子分组文件表格
# layout = "fruchtermanreingold",
path = path,# 结果文件存储路径
fill = "group", # 出图点填充颜色用什么值
size = "igraph.degree", # 出图点大小用什么数据
scale = TRUE, # 是否要进行相对丰度标准化
bio = TRUE, # 是否做二分网络
zipi = FALSE, # 是否计算ZIPI
step = 100, # 随机网络抽样的次数
width = 12,
height = 10
)p = result[[1]]
p
# 全部样本网络参数比对
data = result[[2]]
plotname1 = paste(path,"/network_all.jpg",sep = "")
ggsave(plotname1, p,width = 22,height = 22)
plotname1 = paste(path,"/network_all.pdf",sep = "")
ggsave(plotname1, p,width = 22,height = 22)
tablename <- paste(path,"/co-occurrence_Grobel_net",".csv",sep = "")
write.csv(data,tablename)

结果展示

这里我就简单展示一下出图，使用的ggclusterNet的函数出图，展示不同模块之间的关系：

赞 (0)

相关推荐