
R包安裝# 1-包的安裝------- remotes::install_github("taowenmicro/ggClusterNet",force = TRUE) library(ggClusterNet)
R包導(dǎo)入library(ggClusterNet) library(phyloseq) library(tidyverse) library(ggtree) library(igraph) library(sna) library(network)
查看幫助文檔# 查看幫助文檔 ?network.2()
多組網(wǎng)絡(luò)流程微生物網(wǎng)絡(luò)-基于多組pipeline流程#-2-微生物網(wǎng)絡(luò)分析流程#-------
path = "./result_big_igraph/" dir.create(path)
result = network.2(ps = ps16s, N = 800, big = TRUE, maxnode = 4, select_layout = TRUE, layout_net = "model_igraph", r.threshold=0.6, p.threshold=0.01, label = FALSE, path = path, zipi = FALSE)
# 多組網(wǎng)絡(luò)繪制到一個面板 p = result[[1]] # 全部樣本網(wǎng)絡(luò)參數(shù)比對 data = result[[2]] num= 3 # plotname1 = paste(path,"/network_all.jpg",sep = "") # ggsave(plotname1, p,width = 16*num,height = 16,dpi = 72)
plotname1 = paste(path,"/network_all.pdf",sep = "") ggsave(plotname1, p,width = 16*num,height = 16,limitsize = FALSE)
tablename <- paste(path,"/co-occurrence_Grobel_net",".csv",sep = "") write.csv(data,tablename)


跨域網(wǎng)絡(luò)(細(xì)菌-真菌)-pipeline#2-跨域網(wǎng)絡(luò)-16s_ITS#--------
#--創(chuàng)建文件夾用于存儲文件 Envnetplot<- paste("./16S_ITS_network",sep = "") dir.create(Envnetplot)
library(ggClusterNet) library(tidyverse) library(sna) library(igraph) library(phyloseq) library(WGCNA)
envRDA = env row.names(envRDA) = env$ID envRDA$ID = NULL head(envRDA)
#--細(xì)菌和真菌ps對象中的分組信息要一致 ps.merge <- ggClusterNet::merge16S_ITS(ps16s = ps16s, psITS = psITS, N16s = 500, NITS = 500 )
ps.merge
map = phyloseq::sample_data(ps.merge) head(map)
#依據(jù)實(shí)際情況決定是否修改分組 map$Group = "one" phyloseq::sample_data(ps.merge) <- map
# 可以選定節(jié)點(diǎn)進(jìn)行標(biāo)識 data = data.frame(ID = c("fun_ASV_205","fun_ASV_316","fun_ASV_118"), c("Coremode","Coremode","Coremode"))
result <- corBionetwork(ps = ps.merge, N = 0, lab = data,# 提供data用于特定節(jié)點(diǎn)的注釋 r.threshold = 0.6, # 相關(guān)閾值 p.threshold = 0.05, group = "Group", # env = data1, # 環(huán)境指標(biāo)表格 # envGroup = Gru,# 環(huán)境因子分組文件表格 # layout = "fruchtermanreingold", path = Envnetplot,# 結(jié)果文件存儲路徑 fill = "Phylum", # 出圖點(diǎn)填充顏色用什么值 size = "igraph.degree", # 出圖點(diǎn)大小用什么數(shù)據(jù) scale = TRUE, # 是否要進(jìn)行相對豐度標(biāo)準(zhǔn)化 bio = TRUE, # 是否做二分網(wǎng)絡(luò) zipi = F, # 是否計(jì)算ZIPI step = 100, # 隨機(jī)網(wǎng)絡(luò)抽樣的次數(shù) width = 12, label = TRUE, height = 10, big = TRUE, select_layout = TRUE, layout_net = "model_maptree", clu_method = "cluster_fast_greedy"
)
p = result[[1]] p
# 全部樣本網(wǎng)絡(luò)參數(shù)比對 data = result[[2]]
plotname1 = paste(Envnetplot,"/network_all.pdf",sep = "") ggsave(plotname1, p,width = 10,height = 8) tablename <- paste(Envnetplot,"/co-occurrence_Grobel_net",".csv",sep = "") write.csv(data,tablename)
tablename <- paste(Envnetplot,"/nodeG_plot",".csv",sep = "") write.csv(result[[4]],tablename) tablename <- paste(Envnetplot,"/edge_plot",".csv",sep = "") write.csv(result[[3]],tablename) tablename <- paste(Envnetplot,"/cor_matrix",".csv",sep = "") write.csv(result[[5]],tablename)

跨域網(wǎng)絡(luò)(細(xì)菌-真菌-環(huán)境因子)-pipeline
#2-跨域網(wǎng)絡(luò)16s-ITS-環(huán)境因子的網(wǎng)絡(luò)#--------
Envnetplot<- paste("./16s_ITS_Env_network",sep = "") dir.create(Envnetplot)
ps16s =ps16s%>% ggClusterNet::scale_micro() psITS = psITS%>% ggClusterNet::scale_micro()
library(phyloseq) #--細(xì)菌和真菌ps對象中的map文件要一樣 ps.merge <- ggClusterNet::merge16S_ITS(ps16s = ps16s, psITS= psITS, NITS = 200, N16s = 200)
map = phyloseq::sample_data(ps.merge) head(map) map$Group = "one" phyloseq::sample_data(ps.merge) <- map
#--環(huán)境因子導(dǎo)入 data1 = env envRDA = data.frame(row.names = env$ID,env[,-1])
envRDA.s = vegan::decostand(envRDA,"hellinger") data1[,-1] = envRDA.s
Gru = data.frame(ID = colnames(env)[-1],group = "env" ) head(Gru)
library(sna) library(ggClusterNet) library(igraph)
result <- ggClusterNet::corBionetwork(ps = ps.merge, N = 0, r.threshold = 0.4, # 相關(guān)閾值 p.threshold = 0.05, big = T, group = "Group", env = data1, # 環(huán)境指標(biāo)表格 envGroup = Gru,# 環(huán)境因子分組文件表格 # layout = "fruchtermanreingold", path = Envnetplot,# 結(jié)果文件存儲路徑 fill = "Phylum", # 出圖點(diǎn)填充顏色用什么值 size = "igraph.degree", # 出圖點(diǎn)大小用什么數(shù)據(jù) scale = TRUE, # 是否要進(jìn)行相對豐度標(biāo)準(zhǔn)化 bio = TRUE, # 是否做二分網(wǎng)絡(luò) zipi = F, # 是否計(jì)算ZIPI step = 100, # 隨機(jī)網(wǎng)絡(luò)抽樣的次數(shù) width = 18, label = TRUE, height = 10 )
p = result[[1]] p # 全部樣本網(wǎng)絡(luò)參數(shù)比對 data = result[[2]] plotname1 = paste(Envnetplot,"/network_all.jpg",sep = "") ggsave(plotname1, p,width = 15,height = 12,dpi = 72) # plotname1 = paste(Envnetplot,"/network_all.png",sep = "") # ggsave(plotname1, p,width = 10,height = 8,dpi = 72) plotname1 = paste(Envnetplot,"/network_all.pdf",sep = "") ggsave(plotname1, p,width = 15,height = 12) tablename <- paste(Envnetplot,"/co-occurrence_Grobel_net",".csv",sep = "") write.csv(data,tablename)

ggClusterNet中的網(wǎng)絡(luò)布局算法使用
布局算法使用1-PolygonClusterG#3-布局算法使用1-PolygonClusterG#-------
#-計(jì)算相關(guān)性 result = corMicro (ps = ps, N = 200, r.threshold=0.8, p.threshold=0.05 ) cor = result[[1]] dim(cor) ps_net = result[[3]] otu_table = as.data.frame(t(vegan_otu(ps_net))) tax_table = as.data.frame(vegan_tax(ps_net))
#-模擬一個分組 netClu = data.frame(ID = row.names(tax_table),group =rep(1:5,length(row.names(tax_table)))[1:length(row.names(tax_table))] ) netClu$group = as.factor(netClu$group) set.seed(12) result2 = PolygonClusterG(cor = cor,nodeGroup =netClu,zoom = 0.8,zoom2 = 0.8 ) node = result2[[1]] head(node) # ---node節(jié)點(diǎn)注釋 nodes = nodeadd(plotcord =node,otu_table = otu_table,tax_table = tax_table) #-----計(jì)算邊 edge = edgeBuild(cor = cor,node = node) head(edge) # colnames(edge)[8] = "cor" p1 <- ggplot() + geom_segment(aes(x = X1, y = Y1, xend = X2, yend = Y2,color = cor), data = edge, size = 0.1,alpha = 0.6) + geom_point(aes(X1, X2,fill = Phylum,size = mean),pch = 21, data = nodes) + # geom_text_repel(aes(X1, X2,label = elements),pch = 21, data = nodes) + # geom_text(aes(X1, X2,label = elements),pch = 21, data = nodes) + scale_colour_manual(values = c("#377EB8","#E41A1C")) + scale_size(range = c(4, 14)) + scale_x_continuous(breaks = NULL) + scale_y_continuous(breaks = NULL) + theme(panel.background = element_blank()) + theme(axis.title.x = element_blank(), axis.title.y = element_blank()) + theme(legend.background = element_rect(colour = NA)) + theme(panel.background = element_rect(fill = "white", colour = NA)) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank()) p1

布局算法使用2-model_igraph
#3-布局算法使用2-model_igraph#----- p.threshold = 0.05 r.threshold = 0.6 x = ps %>% # scale_micro(method = "TMM") %>% vegan_otu() %>% t() %>% as.data.frame() occor<-WGCNA::corAndPvalue(t(x)/colSums(x),method = 'pearson') mtadj<-multtest::mt.rawp2adjp(unlist(occor$p),proc='BH') adpcor<-mtadj$adjp[order(mtadj$index),2] occor.p<-matrix(adpcor,dim(t(x)/colSums(x))[2]) ## R value occor.r<-occor$cor diag(occor.r) <- 0
# occor.r[occor.p > 0.05|abs(occor.r)<0.4] = 0
occor.r[occor.p > p.threshold | abs(occor.r)< r.threshold] = 0 occor.r[is.na(occor.r)]=0 cor = occor.r
result2 <- model_igraph(cor = cor, method = "cluster_fast_greedy", seed = 12 ) node = result2[[1]] head(node)
dat = result2[[2]] head(dat) tem = data.frame(mod = dat$model,col = dat$color) %>% dplyr::distinct( mod, .keep_all = TRUE) col = tem$col names(col) = tem$mod
#---node節(jié)點(diǎn)注釋 otu_table = as.data.frame(t(vegan_otu(ps))) tax_table = as.data.frame(vegan_tax(ps)) nodes = nodeadd(plotcord =node,otu_table = otu_table,tax_table = tax_table) head(nodes) #-----計(jì)算邊 edge = edgeBuild(cor = cor,node = node) colnames(edge)[8] = "cor" head(edge)
tem2 = dat %>% dplyr::select(OTU,model,color) %>% dplyr::right_join(edge,by =c("OTU" = "OTU_1" ) ) %>% dplyr::rename(OTU_1 = OTU,model1 = model,color1 = color) head(tem2)
tem3 = dat %>% dplyr::select(OTU,model,color) %>% dplyr::right_join(edge,by =c("OTU" = "OTU_2" ) ) %>% dplyr::rename(OTU_2 = OTU,model2 = model,color2 = color) head(tem3)
tem4 = tem2 %>%inner_join(tem3) head(tem4)
edge2 = tem4 %>% mutate(color = ifelse(model1 == model2,as.character(model1),"across"), manual = ifelse(model1 == model2,as.character(color1),"#C1C1C1") ) head(edge2) col_edge = edge2 %>% dplyr::distinct(color, .keep_all = TRUE) %>% select(color,manual) col0 = col_edge$manual names(col0) = col_edge$color
library(ggnewscale) # p1 <- ggplot() + geom_segment(aes(x = X1, y = Y1, xend = X2, yend = Y2,color = color), # data = edge2, size = 1) + # scale_colour_manual(values = col0)
p2 = p1 + new_scale_color() + geom_point(aes(X1, X2,color =model), data = dat,size = 4) + scale_colour_manual(values = col) + scale_x_continuous(breaks = NULL) + scale_y_continuous(breaks = NULL) + theme(panel.background = element_blank()) + theme(axis.title.x = element_blank(), axis.title.y = element_blank()) + theme(legend.background = element_rect(colour = NA)) + theme(panel.background = element_rect(fill = "white", colour = NA)) + theme(panel.grid.minor = element_blank(), panel.grid.major = element_blank())
# ggsave("./cs2.pdf",p2,width = 16,height = 14)


根際互作生物學(xué)研究室 簡介根際互作生物學(xué)研究室是沈其榮院士土壤微生物與有機(jī)肥團(tuán)隊(duì)下的一個關(guān)注于根際互作的研究小組。本小組由袁軍副教授帶領(lǐng),主要關(guān)注:1.植物和微生物互作在抗病過程中的作用;2 環(huán)境微生物大數(shù)據(jù)整合研究;3 環(huán)境代謝組及其與微生物過程研究體系開發(fā)和應(yīng)用。團(tuán)隊(duì)在過去三年中在 isme J, Microbiome, PCE,SBB,Horticulture Research等期刊上發(fā)表了多篇文章。歡迎關(guān)注 微生信生物 公眾號對本研究小組進(jìn)行了解。
|