高手题目一答案：nature高级组合图形绘制教程 / 开普饭

写在前面

开学了，大量工作来袭。我在成长，希望可以帮助大家一起成长。按照咱们的约定，在看数目超过50，我应该放出高手题目一的代码。

如果需要我重复高手题目二，请继续点赞，达到100个在看，献上高手题目二代码。点击这里查看

准备包和数据

使用内置数据，所以别再问我要了。各位大哥大姐。

library(ggplot2)# 用于绘图 library(psych)# 用于相关计算 dim(mtcars)

准备参数

由于矩形内是相关，所以需要指定两个阈值。由于是矩形展示，所以需要提供矩形长宽，这里我为了简化，直接作为正方形，所以只需要指定边长。

r.threshold=0.9 p.threshold=0.05


# 设置矩形半径

r = 5

#--坐标轴标度长度 leng = 0.2

矩形展示需要四个分组

模拟四个分组

#---模拟四个分组 netClu = data.frame(ID = colnames(mtcars2),group =rep(1:4,length(colnames(mtcars2)))[1:length(colnames(mtcars2))] ) netClu$group = as.factor(netClu$group) head(netClu) row.names(netClu) = netClu$ID> head(netClu) ID group 1 Mazda RX4 1 2 Mazda RX4 Wag 2 3 Datsun 710 3 4 Hornet 4 Drive 4 5 Hornet Sportabout 1 6 Valiant 2

添加展示丰度柱状图值

netClu = merge(netClu,as.data.frame(colSums(mtcars2)),by = "row.names",all= F)


colnames(netClu)[4] = "abu"

netClu$abu = netClu$abu/sum(netClu$abu)*30head(netClu) Row.names ID group abu 1 AMC Javelin AMC Javelin 3 1.0889636 2 Cadillac Fleetwood Cadillac Fleetwood 3 1.5676720 3 Camaro Z28 Camaro Z28 4 1.3906268 4 Chrysler Imperial Chrysler Imperial 1 1.5615073 5 Datsun 710 Datsun 710 3 0.5585488 6 Dodge Challenger Dodge Challenger 2 1.1181519

四条边分别准备数据

##提取一条边数据 gr = dplyr::filter(netClu, group == as.vector(as.factor(names(table(netClu$group)[1]))[1]) ) # 根据该组内元素多少准备坐标 x1 = seq(from=0, to=r, length.out=table(netClu$group)[1])

#--相关点坐标数据构建 axis1 = data.frame(row.names = gr$ID,x = x1,y = 0,element = gr$ID) #--丰度坐标构建 abu1 = data.frame(row.names = gr$ID,x = x1,y = -gr$abu,element = gr$ID )

其他三个边坐标构造

gr = dplyr::filter(netClu, group == as.vector(as.factor(names(table(netClu$group)[3]))[1]) ) x2 = seq(from=0, to=r, length.out=table(netClu$group)[3]) axis2 = data.frame(row.names = gr$ID,x = x2,y = r +2,element = gr$ID) abu2 = data.frame(row.names = gr$ID,x = x2,y = gr$abu + r +2 ,element = gr$ID )


gr = dplyr::filter(netClu, group == as.vector(as.factor(names(table(netClu$group)[2]))[1]) )

y1 = seq(from=0, to=r, length.out=table(netClu$group)[2])

axis3 = data.frame(row.names = gr$ID,x = -1,y = y1 + 1,element = gr$ID)

abu3 = data.frame(row.names = gr$ID,x = -1 -gr$abu,y = y1 + 1 ,element = gr$ID )

gr = dplyr::filter(netClu, group == as.vector(as.factor(names(table(netClu$group)[4]))[1]) ) y2 = seq(from=0, to=r, length.out=table(netClu$group)[4]) axis4 = data.frame(row.names = gr$ID,x = r +1,y = y2 +1,element = gr$ID) abu4 = data.frame(row.names = gr$ID,x = r +1 + gr$abu,y = y2 + 1 ,element = gr$ID )

图形坐标标度构建

a = data.frame(x = 6,y = seq(from=0, to=-max(abs(abu1$y)), length.out=3)) b = data.frame(x = 6 +leng,y = seq(from=0, to=-max(abs(abu1$y)), length.out=3))


a1 = data.frame(x = -1,y = seq(from=r +2, to=max(abs(abu2$y)), length.out=3))

a1

b1 = data.frame(x = -1- leng,y =seq(from=r +2, to=max(abs(abu2$y)), length.out=3))

b1
a3 = data.frame(x = seq(from=-1, to=-max(abs(abu3$x)), length.out=3),y = 0)

a3

b3 = data.frame(x = seq(from=-1, to=-max(abs(abu3$x)), length.out=3),y = 0-leng)

b3

a4 = data.frame(x = seq(from=6, to=max(abs(abu4$x)), length.out=3),y = r+2) a4 b4 = data.frame(x = seq(from=6, to=max(abs(abu4$x)), length.out=3),y = r+2+leng) b4

组合矩形坐标和柱状图坐标

point = rbind(axis1,axis2,axis3,axis4) head(point)


adun = rbind(abu1,abu2,abu3,abu4)

head(adun)

colnames(adun) = c("x2","y2","element2") abunpro = merge(point,adun,by = "row.names",all = F) head(abunpro)

可视化效果

ggplot() +geom_point(data = adun,aes(x =x2,y = y2 )) +geom_point(data = point,aes(x =x,y = y ))

ggplot() + geom_segment(aes(x = x, y = y, xend = x2, yend = y2), data = abunpro, size = 5)

矩形为网络结构构造网络边和节点

#--节点坐标匹配 cor = occor.r plotcord = point #----使用相关矩阵构建网络--network包构建网络#----- g <- network::network(cor, directed=FALSE) #---转化为0-1的相关矩阵 # origm <- network::as.matrix.network.adjacency(g) # get sociomatrix


edglist <- network::as.matrix.network.edgelist(g)

edglist = as.data.frame(edglist)

head(edglist)
#t添加权重#---

network::set.edge.value(g,"weigt",cor)

# g

# ?get.edge.attribute
# row.names(plotcord) = NULL
edglist$weight = as.numeric(network::get.edge.attribute(g,"weigt"))

head(edglist)

edges <- data.frame(plotcord[edglist[, 1], ], plotcord[edglist[, 2], ])

head(edges)

edges$weight = as.numeric(network::get.edge.attribute(g,"weigt"))
##这里将边权重根据正负分为两类

aaa = rep("a",length(edges$weight))

for (i in 1:length(edges$weight)) {

if (edges$weight[i]> 0) {

aaa[i] = "+"

}

if (edges$weight[i]< 0) {

aaa[i] = "-"

}

}

#添加到edges中

edges$wei_label = as.factor(aaa)

colnames(edges) <- c("X1", "Y1","OTU_1", "X2", "Y2","OTU_2","weight","wei_label")

edges$midX <- (edges$X1 + edges$X2)/2

edges$midY <- (edges$Y1 + edges$Y2)/2

head(edges)
#--构造边

head(netClu)
#--返回第一个在第二个中的位置
point[,"group"] = netClu[match(point$element,netClu$ID),][3]
nodes = point

head(nodes)
head(abunpro)
head(netClu)
library(dplyr)

abunpro = abunpro %>% inner_join(netClu, by = "Row.names")

角度调整

angle = data.frame(group = as.data.frame(table(netClu$group))$Var,c(90,0,90,0))

abunpro = abunpro %>% left_join(angle, by = "group") colnames(abunpro)[dim(abunpro)[2]] = "angle" abunpro$hjust = NULL hjust = data.frame(group = as.data.frame(table(netClu$group))$Var,c(1,1,0,0)) hjust abunpro = abunpro %>% left_join(hjust, by = "group") colnames(abunpro)[dim(abunpro)[2]] = "hjust"

出图

pnet <- ggplot() + geom_segment(aes(x = X1, y = Y1, xend = X2, yend = Y2,color = as.factor(wei_label)),color = "green", data = edges, size = 0.5) + geom_segment(aes(x = x, y = y, xend = x2, yend = y2,color = group), data = abunpro, size = 10) + # geom_point(aes(x = x, y = y,fill = group,group = group),pch = 21, data = nodes,size = 3) + geom_text(aes(x= x2, y = y2,label = element2,angle = angle,hjust =hjust), data = abunpro)+ 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 = pnet + geom_segment(aes(x = c(6), y = c(0), xend = c(6), yend = c(-max(abs(abu1$y)))), size = 1) + geom_segment(aes(x = a$x, y = a$y, xend = b$x, yend = b$y),size = 1) + geom_text(aes(x = b$x, y = b$y, label = abs(round(a$y,2))))+ geom_segment(aes(x = c(-1), y = c(r+2), xend = c(-1), yend = c(max(abs(abu2$y)))), size = 1) + geom_segment(aes(x = a1$x, y = a1$y, xend = b1$x, yend = b1$y),size = 1) + geom_text(aes(x = b1$x, y = b1$y, label = abs(round(a1$y,2))-r-2)) + geom_segment(aes(x = c(-1), y = c(0), xend = c(-max(abs(abu3$x))), yend = c(0)), size = 1) + geom_segment(aes(x = a3$x, y = a3$y, xend = b3$x, yend = b3$y),size = 1) + geom_text(aes(x = b3$x, y = b3$y, label = abs(round(a3$x,2))-1)) + geom_segment(aes(x = c(6), y = c(r+2), xend = c(max(abs(abu4$x))), yend = c(r+2)), size = 1) + geom_segment(aes(x = a4$x, y = a4$y, xend = b4$x, yend = b4$y),size = 1) + geom_text(aes(x = b4$x, y = b4$y, label = abs(round(a4$x,2))-r-1))

ggsave("./cs3.pdf",pnet,width =15,height = 10)

右侧相关

#-------构造右侧右侧数据和另外一组数据的相关。 data2 = data.frame(ID = paste("A",1:8,sep = ""), wei = c(1,-1,1,-1,1,-1,1,-1), x1= axis4$x + 5, y1 = axis4$y[8:1] )


head(data2 )

data3 = cbind(axis4,data2)

head(data3)

data3$wei = runif(length(data3$x), 0, 0.5)

pnet = pnet + geom_point(aes(x = x + 3,y = y),data = axis4,color = "blue",size = 4) +

geom_point(aes(x = x1 ,y = y1),data = data2,color = "blue",size = 4) +

geom_segment(aes(x = x +3, y = y, xend = x1, yend = y1,size = wei),data = data3,color = "yellow") +

geom_text(aes(x= x1,y = y1,label = ID),data = data2,hjust = -1)

ggsave("./cs4.pdf",pnet,width =20,height = 10)

pnet

欢迎加入微生信生物

快来微生信生物

微生信生物

轻松一刻 ◆ ◆

二师兄的日常

二师兄，何许人！小弟亲师兄也，硕士毕业于2018年，你就看着他，就有数不清的意思。在枯燥乏味的科研生活中有着独特的光芒，让我膜拜。如果你感到无力，请关注二师兄，看看他能带给你多少意思。

高手题目一答案：nature高级组合图形绘制教程