在ggClusterNet中仿造cytosccape添加多行的聚类布局 / 开普饭

写在前面

算法肯定不同，功能相近。

时间戳：现在2020年8月8日，我在上周便构造了这个函数，用于模仿cytoscape网络的矩阵布局。类似下图的样式。

但是我觉得这个算法我写的不够完善，后面应该还会改进，尤其是在多行布局的时候对于自定义不同模块之间距离的设置上。但是更新之后不会影响函数。

微生物网络

输入文件

#--导入所需R包#------- library(igraph) library(network) library(sna) library(ggplot2)# 用于出图 library(ggrepel)# 用于出图 library(ggClusterNet) library(tidyverse) library("ggalt")#用于将模块圈起来


vegan_tax <- function(physeq){

tax <- tax_table(physeq)
return(as(tax,"matrix"))

}

vegan_otu <- function(physeq){ OTU <- otu_table(physeq) if(taxa_are_rows(OTU)){ OTU <- t(OTU) } return(as(OTU,"matrix")) } #-----导入数据#------- ps = readRDS("../ori_data/ps_liu.rds") ps

corMicro函数用于计算相关

按照丰度过滤微生物表格，并却计算相关矩阵，按照指定的阈值挑选矩阵中展示的数值。

#-----微生物网络构建参数设置#----


#----------计算相关#----

result = corMicro (ps = ps,N = 0.02,r.threshold=0.6,p.threshold=0.05,method = "pearson")

#--提取相关矩阵 cor = result[[1]] # head(cor)

制作分组使用模块化程度作为分组

这是网络布局的基础，无论是什么聚类布局，都需要制作一个分组文件，这个文件有两列，一列是节点，一列是分组信息，这个分组信息名称为：group。

modulGroup使用的igraph包中的模块化分析，将节点按照模块化程度分了组，这里的分组其实并不是非要用这个函数，我们可以按照微生物分类水平分组，按照高丰度低丰度分组等，只有有一个两列的一个数据框文件，第一列是OTU，第二列是因子形式的分组列即可。

netClu = modulGroup( cor = cor,cut = NULL,method = "cluster_fast_greedy" )

head(netClu )

randomClusterG 根据分组，随机分布计算布局

计算布局的函数，本包中包含很多，可能都不太一样，核心是计算模块的中心位置，然后使用中心位置极端该模块内部节点的一个位置。

PolygonModsquareG：函数用按照行和列的整齐排布网络模块化模块，这个类型的函数需要以相关矩阵作为输入，所以大家不要怕，相关矩阵到处有的是。还有一个分组文件：对OTU或者网络节点的一个分组文件。其他参数就是不同布局各自有的一些参数。

这类布局函数额输出文件都是类似的：

第一个表格是节点和节点坐标文件，也是后续分析采用的重要文件。第二个表格往往使模块化的中心坐标位置。这用什么用？可以提取出来人工修改这个坐标然后使用ArtifClusterG函数重新根据指定的中心模块坐标来计算节点。

除了计算布局的函数外，我们要使用布局函数经计算的节点来构造边文件：

edgeBuild：需要一个相关矩阵和布局函数计算得到的节点坐标。当然我们也可以自己制定坐标，构造一个三列的数据框就好了，前两列是坐标第三列是接待名称。

#--------计算布局#--------- result2 <- PolygonModsquareG(cor = cor,nodeGroup =netClu,r1 = 1,N = 1.1,cut = 1)


node = result2[[1]]

head(node)
### nodeadd 节点注释的简单封装，便捷实用otu表格和分组文件进行注释
ps_net <- result[[3]]
tax_table = as.data.frame(vegan_tax(ps_net))

otu_table = as.data.frame(t(vegan_otu(ps_net)))
#---node节点注释#-----------

nodes = nodeadd(plotcord =node,otu_table = otu_table ,tax_table = tax_table)

head(nodes)
#-----计算边#--------

edge = edgeBuild(cor = cor,plotcord = node)

head(edge)

出图

pnet <- ggplot() + geom_segment(aes(x = X1, y = Y1, xend = X2, yend = Y2,color = as.factor(wei_label)), data = edge, size = 0.5,alpha = 0.05) + geom_point(aes(X1, X2,fill = Phylum,size = mean),pch = 21, data = nodes) + scale_colour_brewer(palette = "Set1") + scale_x_continuous(breaks = NULL) + scale_y_continuous(breaks = NULL) + # labs( title = paste(layout,"network",sep = "_"))+ # geom_text_repel(aes(X1, X2,label=Phylum),size=4, data = plotcord)+ # discard default grid + titles in ggplot2 theme(panel.background = element_blank()) + # theme(legend.position = "none") + 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())


pnet

# ggsave("1.png",pnet,width = 18,height = 2) # pnet <- pnet + geom_text_repel(aes(X1, X2,label=elements),size=4, data = nodes)

cut参数设置3，代表将全部分组分为三行展示，上面的就是一行展示，

#--------计算布局#--------- result2 <- PolygonModsquareG(cor = cor,nodeGroup =netClu,r1 = 1,N = 0.5,cut = 3) node = result2[[1]]


### nodeadd 节点注释的简单封装，便捷实用otu表格和分组文件进行注释

ps_net <- result[[3]]

tax_table = as.data.frame(vegan_tax(ps_net))

otu_table = as.data.frame(t(vegan_otu(ps_net)))

#---node节点注释#-----------

nodes = nodeadd(plotcord =node,otu_table = otu_table ,tax_table = tax_table)
#-----计算边#--------

edge = edgeBuild(cor = cor,plotcord = node)
### 出图

ggplot() + geom_segment(aes(x = X1, y = Y1, xend = X2, yend = Y2,color = as.factor(wei_label)),

data = edge, size = 0.5,alpha = 0.05) +

geom_point(aes(X1, X2,fill = Phylum,size = mean),pch = 21, data = nodes) +

scale_colour_brewer(palette = "Set1") +

scale_x_continuous(breaks = NULL) + scale_y_continuous(breaks = NULL) +

# labs( title = paste(layout,"network",sep = "_"))+

# geom_text_repel(aes(X1, X2,label=Phylum),size=4, data = plotcord)+

# discard default grid + titles in ggplot2

theme(panel.background = element_blank()) +

# theme(legend.position = "none") +

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())

# geom_text_repel(aes(X1, X2,label=elements),size=4, data = nodes)

我们可以将每个模块圈起来，用于明显的区分

上面我们的节点合并了物种注释信息，下面我们进一步合并网络模块化信息。

row.names(netClu) = netClu$ID nodeG = merge(nodes,netClu,by = "row.names",all =FALSE) dim(nodeG) head(nodeG)

这里使用到了geom_encircle函数。

ggplot() + geom_segment(aes(x = X1, y = Y1, xend = X2, yend = Y2,color = as.factor(wei_label),alpha = 0.05), data = edge, size = 0.5) + geom_point(aes(X1, X2,fill = group,size = mean),pch = 21, data = nodeG) + geom_encircle(aes(X1, X2,group = group,fill = group), linetype = 2,alpha = 0.1, data = nodeG) + scale_colour_brewer(palette = "Set1") + theme_void()

在ggClusterNet中仿造cytosccape添加多行的聚类布局