主编寄语
写在前面
Wen, T., et al., The best practice for microbiome analysis using R. Protein & Cell, 2023. 14(10): p. 713-725.
数据代码获取地址
实例介绍
网络分析:主要基于ggClusterNet R包具体可以查看这篇文章:https://onlinelibrary.wiley.com/doi/full/10.1002/imt2.32#/
实践操作
网络稳定性:指微生物组中物种相互作用的网络能够抵抗外界干扰并维持其结构和功能的能力。
抵抗力:网络在面对干扰(如环境变化、物种入侵或移除)时保持不变的能力。
恢复力:网络在遭受干扰后恢复到原有状态的能力。
网络稳定性分析可以做什么:
生态系统评估:通过分析微生物组网络的稳定性,可以评估生态系统的健康状况和生物多样性的丰富程度。
环境变化影响:预测环境压力(如气候变化、污染)对微生物群落结构和功能的潜在影响。
关键物种识别:识别对网络稳定性贡献大的关键物种,为保护和管理生态系统提供依据。
微生物间研究:增进对微生物组内相互作用、物种共存机制及群落构建过程的理解。
复杂性:考虑微生物组内复杂的物种相互作用,包括正相互作用(如共生、互利共生)和负相互作用(如竞争、捕食)。
动态性:分析微生物组网络随时间和环境条件变化的动态特性。
模块化:识别网络中的模块,理解不同模块间的相互作用及其对整体网络稳定性的影响。
多尺度性:可以在多个尺度上进行分析。
应用场景:
环境监测:评估环境压力对微生物群落的影响。
农业土壤管理:理解土壤微生物组对作物健康的影响。
生物修复:设计微生物组以提高污染物的生物降解效率。
network stability analysis(网络稳定性分析)
在网络稳定性计算中,负相关(NegativeCorrelation)通常指的是网络中两个或多个节点之间的一种关系,其中当一个节点的状态或行为向某个方向变化时,另一个节点的状态或行为则向相反的方向变化。
library(ggClusterNet)
library(phyloseq)
library(tidyverse)
library(igraph)
library(tidyfst)
netpath = paste(otupath,"/network_stab/",sep = "")
dir.create(netpath)
ps = ps0
res = module.compare.m(
ps = ps,
Top = 200,
degree = TRUE,
zipi = FALSE,
r.threshold= 0.8,
p.threshold=0.05,
method = "spearman",
padj = F,
n = 3)
p = res[[1]]
p
dat = res[[2]]
head(dat)
dat2 = res[[3]]
head(dat2)
dat2$m1 = dat2$module1 %>% strsplit("model") %>%
sapply(`[`, 1)
dat2$m2 = dat2$module2 %>% strsplit("model") %>%
sapply(`[`, 1)
dat2$cross = paste(dat2$m1,dat2$m2,sep = "_Vs_")
head(dat2)
p2 = ggplot(dat2) + geom_bar(aes(x = cross,fill = cross)) +
labs(x = "",
y = "numbers.of.similar.modules"
)+ theme_classic()
p2
getwd()
library(ggplot2)
#--发现分组1和分组3网络更相似一些
FileName <- paste(netpath,"module.compare.groups.pdf", sep = "")
ggsave(FileName, p, width = 10, height = 10)
FileName <- paste(netpath,"numbers.of.similar.modules.pdf", sep = "")
ggsave(FileName, p2, width = 8, height = 8)
FileName <- paste(netpath,"module.otu.csv", sep = "")
write.csv(dat,FileName, quote = F)
FileName <- paste(netpath,"module.compare.groups.csv", sep = "")
write.csv(dat2,FileName, quote = F)
```
### 鲁棒性
```{r}
#--随即取出任意比例节点-网络鲁棒性#---------
res = Robustness.Random.removal(ps = ps,
Top = 500,
r.threshold= 0.8,
p.threshold=0.05,
method = "spearman"
)
p = res[[1]]
p
dat = res[[2]]
head(dat)
# dir.create("./Robustness_Random_removal/")
path = paste(netpath,"/Robustness_Random_removal/",sep = "")
fs::dir_create(path)
write.csv(dat,
paste(path,"random_removal_network.csv",sep = ""))
ggsave(paste(path,"random_removal_network.pdf",sep = ""), p,width = 8,height = 4)
res= Robustness.Targeted.removal(ps = ps,
Top = 500,
degree = TRUE,
zipi = FALSE,
r.threshold= 0.8,
p.threshold=0.05,
method = "spearman")
p = res[[1]]
p
dat = res[[2]]
path = paste(netpath,"/Robustness_Targeted_removal/",sep = "")
fs::dir_create(path)
write.csv(dat,
paste(path,"Robustness_Targeted_removal_network.csv",sep = ""))
ggsave(paste(path,"Robustness_Targeted_removal_network.pdf",sep = ""), p,width = 8,height = 4)
网络易损性
res = Vulnerability.micro(ps = ps,
Top = 500,
degree = TRUE,
zipi = FALSE,
r.threshold= 0.8,
p.threshold=0.05,
method = "spearman")
p = res[[1]] + theme_bw()
p
dat = res[[2]]
path = paste(netpath,"/Vulnerability/",sep = "")
fs::dir_create(path)
write.csv(dat,
paste(path,"Vulnerability_network.csv",sep = ""))
ggsave(paste(path,"Vulnerability_network.pdf",sep = ""), p,width = 4,height = 4)
计算负相关的比例
id <- sample_data(ps)$Group %>% unique()
res = negative.correlation.ratio(ps = ps,
Top = 500,
degree = TRUE,
zipi = FALSE,
r.threshold= 0.8,
p.threshold=0.05,
method = "spearman",
corg = NULL)
p = res[[1]]
p
dat = res[[2]]
head(dat)
path = paste(netpath,"/Vulnerability/",sep = "")
fs::dir_create(path)
write.csv(dat,
paste(path,"Vnegative.correlation.ratio_network.csv",sep = ""))
ggsave(paste(path,"negative.correlation.ratio_network.pdf",sep = ""), p,width = 4,height = 4)
网络抗毁性
# library(tidyfst)
library("pulsar")
# library(ggClusterNet)
# library(phyloseq)
# library(tidyverse)
id = sample_data(ps)$Group %>% unique()
res = natural.con.microp (
ps = ps,
Top = 500,
r.threshold= 0.6,
p.threshold=0.05,
method = "spearman",
norm = F,
end = 150,# 小于网络包含的节点数量
start = 0,
con.method = "pulsar"
)
p = res[[1]]
p
dat = res[[2]]
path = paste(netpath,"/Natural_connectivity/",sep = "")
fs::dir_create(path)
write.csv(dat,
paste(path,"/Natural_connectivity.csv",sep = ""))
ggsave(paste(path,"/Natural_connectivity.pdf",sep = ""), p,width = 5,height = 4)
network modularity analysis- Module eigenvectors(网络模块化分析-模块特征向量)
id = sample_data(ps)$Group %>% unique()
id
i = 1
netpath = paste(otupath,"/network3_MEs/",sep = "")
dir.create(netpath)
netpath;library(ggClusterNet)
library(igraph)
for (i in 1:length(id)) {
ps.1 = phyloseq::subset_samples(
ps,Group %in% c(id[i])
)
result = network.2(ps = ps.1, N = 500,
big = TRUE,
select_layout = TRUE,
layout_net = "model_maptree",
r.threshold=0.8,
p.threshold=0.05,
label = FALSE,
path = netpath,
zipi = F,
ncol = 1,
nrow = 1,
# method = "sparcc",
fill = "Phylum"
)
# 节点的模块化信息的输出 group列为模块化信息
tem <- ggClusterNet::model_maptree(cor =result[[4]],
method = "cluster_fast_greedy",
seed = 12
)
node_model = tem[[2]]
head(node_model)
tablename <- paste(netpath,"/node_model",".csv",sep = "")
write.csv(node_model,tablename)
head(node_model)
otu = ps.1 %>%
filter_OTU_ps(500) %>%
vegan_otu() %>%
as.data.frame()
node_model = node_model[match(colnames(otu),node_model$ID),]
MEList = WGCNA::moduleEigengenes(otu, colors = node_model$group)
MEs = MEList$eigengenes %>% as.data.frame()
tablename <- paste(netpath,"/",id[i],"node_characteristic_variables",".csv",sep = "")
write.csv(MEs,tablename)
}
detach("package:igraph")
Comparison and visualization of distances between groups(分组之间距离比较和可视化)
机器学习:在聚类分析中,比较不同群组的中心点距离。
复杂网络:在网络科学中,比较不同网络或网络中不同社区的连接距离
group = "Group"
map = as.data.frame(sample_data(ps))
alppath = paste(otupath,"/distance/",sep = "")
dir.create(alppath)
gro = map[,group] %>% unique()
colnames(gro) = "group"
conbgroup = combn(gro$group,2)
# 计算包括终点均值的所有样品bray距离
bray_curtis = vegan::vegdist(vegan_otu(ps), method = "bray")
bray_curtis = as.matrix(bray_curtis)
for (i in 1:dim(conbgroup)[2]) {
a = conbgroup[,i]
map = as.data.frame(sample_data(ps))
# head(map)
chose1 = map[as.matrix(map[,group]) %>% as.vector() == a[1],] %>% row.names()
chose2 = map[as.matrix(map[,group]) %>% as.vector() == a[2],] %>% row.names()
dat = data.frame(group = paste(a[1],a[2],sep = "_VS_"), Distance =bray_curtis[chose1,chose2] %>% as.dist() %>% as.vector() )
head(dat)
if (i == 1) {
table = dat
}
if (i != 1) {
table = rbind(table,dat)
}
}
head(table)
table$id = 1:dim(table)[1]
data <- table %>% dplyr::select(id,everything())
result = EasyStat::MuiKwWlx(data = data,num = c(3))
FileName <- paste(alppath,"/distance_label.csv", sep = "")
write.csv(result,FileName,sep = "")
FileName <- paste(alppath,"/distance_index.csv", sep = "")
write.csv(data,FileName,sep = "")
result1 = EasyStat::FacetMuiPlotresultBox(data = data,num = c(3),result = result,sig_show ="abc",ncol = 1 )
p1_1 = result1[[1]] +
# scale_x_discrete(limits = axis_order) +
mytheme2 +
guides(fill = guide_legend(title = NULL)) +
scale_fill_manual(values = colset4)
p1_1
res = EasyStat::FacetMuiPlotresultBar(data = data,num = c(3),result = result,sig_show ="abc",ncol = 1)
p1_2 = res[[1]]+
# scale_x_discrete(limits = axis_order) +
guides(color = FALSE) +
mytheme2+
guides(fill = guide_legend(title = NULL))+
scale_fill_manual(values = colset4)
p1_2
res = EasyStat::FacetMuiPlotReBoxBar(data = data,num = c(3),result = result,sig_show ="abc",ncol = 1)
p1_3 = res[[1]]+
# scale_x_discrete(limits = axis_order) +
mytheme2 +
guides(fill = guide_legend(title = NULL))+
scale_fill_manual(values = colset4)
p1_3
FileName <- paste(alppath,"distance_box", ".pdf", sep = "")
ggsave(FileName, p1_1, width = ((4+ gnum) ), height =8,limitsize = FALSE)
FileName <- paste(alppath,"distance_bar", ".pdf", sep = "")
ggsave(FileName, p1_2, width = ((4+gnum) ), height = 8,limitsize = FALSE)
FileName <- paste(alppath,"distance_boxbar", ".pdf", sep = "")
ggsave(FileName, p1_3, width = ((4+gnum) ), height = 8,limitsize = FALSE)
FileName <- paste(alppath,"distance_box", ".jpg", sep = "")
ggsave(FileName, p1_1, width = (( 4+gnum) ), height =8,limitsize = FALSE)
FileName <- paste(alppath,"distance_bar", ".jpg", sep = "")
ggsave(FileName, p1_2, width = ((4+ gnum) ), height = 8,limitsize = FALSE)
FileName <- paste(alppath,"distance_boxbar", ".jpg", sep = "")
ggsave(FileName, p1_3, width = ((4+ gnum) ), height = 8,limitsize = FALSE)
运行时的错误
Error in full_join(., igraph.degree, na_matches = "never") :
could not find function "full_join"
library(dplyr)
就可以解决
Error in ggplot() : could not find function "ggplot"
library(ggplot2)
就可以解决
Error in as.tibble(.) : could not find function "as.tibble"
library(tibble)可以解决
Error in (function (cond) :
在为'taxa_are_rows'函数选择方法时评估'physeq'参数出了错: object of type 'closure' is not subsettable
根际互作生物学研究室 简介
根际互作生物学研究室是沈其荣院士土壤微生物与有机肥团队下的一个关注于根际互作的研究小组。本小组由袁军教授带领,主要关注:1.植物和微生物互作在抗病过程中的作用;2 环境微生物大数据整合研究;3 环境代谢组及其与微生物过程研究体系开发和应用。团队在过去三年中在 Nature Communications,ISME J,Microbiome,SCLS,New Phytologist,iMeta,Fundamental Research, PCE,SBB,JAFC(封面),Horticulture Research,SEL(封面),BMC plant biology等期刊上发表了多篇文章。欢迎关注 微生信生物 公众号对本研究小组进行了解。
撰写:牛国庆 文浩仰
修改:文涛
排版:杨雯儀
审核:袁军 团队工作及其成果 (点击查看)
了解 交流 合作
小组负责人邮箱 袁军:junyuan@njau.edu.cn;
小组成员文涛:taowen@njau.edu.cn等
团队公众号:微生信生物 添加主编微信,或者后台留言。
1.仅限相关专业或研究方向人员添加,必须实名,不实名则默认忽略。
2.非相关专业的其他人员及推广宣传人员禁止添加。
3.添加主编微信需和简单聊一聊专业相关问题,等待主编判断后,可拉群。
4 微生信生物VIP微信群不受限制,给微生信生物供稿一次即可加入(群里发送推文代码+高效协助解决推文运行等问题+日常问题咨询回复)。
加主编微信 加入群聊
目前营销人员过多,为了维护微生信生物3年来维护的超5500人群聊,目前更新进群要求: