当韦恩图和upset都不能满足我的可视化要求

1

类韦恩图绘制

我们经常的需求是不同对象之间做韦恩图来查看不同组OTU共有和特有情况。这项功能我们使用VennDiagram包就可以解决。其次upset也为我们带来了全新的解决方式，可以再分组数量超过5组还可以很详细和清晰的带给我们庞大的信息。实际上以上两种图都充斥的大量的信息，反而展现共有和特有的otu，让我们感觉感觉更加利索。

ven图展示：

upset图形展示：

类韦恩图展示：

2

类韦恩图实战

如果你就是想做一张类韦恩图的话，下面这几条就足够了。（为避免占用资源，出图函数我已经放到推送的末尾，大家复制粘贴即可。） flower_plot函数中几个有意义的参数设置:

• sample：设置特有环的标签；
• value：设置特有环的数量信息；
• start：设置旋转角度，90即可；
• a：设置椭圆宽度；
• b：设置标签外围大小；
• cer_label：设置中心圆标签；
• col_c：设置中心园标签颜色；
• ellipse_col：设置椭圆颜色；
• circle_col：设置中心圆颜色；
• circle_text_cex：设置中心圆标签字体大小。

sinID = paste("A",1:10,sep = "_")
sin_num = 1:10

pdf("filenamecs.pdf", width = 6, height = 6)

flower_plot(sinID,
sin_num,
ellipse_col = mi[6],
circle_col = mi[2],
90, 0.4, 2,
cer_label = paste("all",":",all_num,sep = ""), col_c = mi[1]
)

dev.off()


3

微生物-生态学等组学中运用

简单的图形我们熟悉了用法和基本调参，这里我们将其运用到基因丰都，这些基因可以是我们扩增子得到细菌真菌群落矩阵，可以是不同功能基因的矩阵等，这里大家各展所长。

phyloseq已经免试好几年了，作为封装微生物组下有分析数据的优秀工具，我已经大部分整合在扩增子流程中了，如果不够了解，大家可以查看之前的推送。补充一下基础知识。因为下面的开始我直接使用的是封装好的ps文件来进行输入，这样大家可以避免输入多个文件造成不一致或者输入错误。整体归一也将更加简单方便。

完整的分析脚本和测试数据都已经存到github上了，大家可以自行下载。后台回复：类韦恩图


library("phyloseq")
library(tidyverse)
library(reshape2)

导入和提取数据

ps = readRDS("./ps_OTU_.ps")

mapping = as.data.frame(sample_data(ps))
ps1 = ps
ps1
vegan_otu <- function(physeq){
OTU <- otu_table(physeq)
if(taxa_are_rows(OTU)){
OTU <- t(OTU)
}
return(as(OTU,"matrix"))
}
aa = vegan_otu(ps1)
otu_table = as.data.frame(t(aa))
count = aa
countA = count
dim(count)
sub_design <- as.data.frame(sample_data(ps1))

# levels(sub_design$SampleType)[1]
# name1 = paste("name",levels(sub_design$SampleType)[1],sep = "")
sub_design $SampleType
levels(sub_design $SampleType)


数据转化为0和1，方便提取共有和特有成分

##########这里的操作为提取三个分组
pick_val_num <- num*2/3
count[count > 0] <- 1###这个函数只能用于0,1 的数据，所以我这么转换

count2 = as.data.frame(count)

library("tibble")
#数据分组
iris.split <- split(count2,as.factor(sub_design$SampleType))
#数据分组计算平均值
iris.apply <- lapply(iris.split,function(x)colSums(x[]))
# 组合结果
iris.combine <- do.call(rbind,iris.apply)
ven2 = t(iris.combine)
head(ven2)
ven2[ven2 < pick_val_num] = 0
ven2[ven2 >=pick_val_num] = 1
ven2 = as.data.frame(ven2)

#########更加高级的设置在这里可以查看#https://mp.weixin.qq.com/s/6l7gftKQfiyxNH66i19YtA
ven3 = as.list(ven2)
ven2 = as.data.frame(ven2)
head(ven2)

下面提取共有和特有成分


all_num = dim(ven2[rowSums(ven2) == 4,])[1]

ven2[,1] == 1

A = rep("A",length(colnames(ven2)))
B = rep(1,length(colnames(ven2)))

i = 1
for (i in 1:length(colnames(ven2))) {
B[i] = length(ven2[rowSums(ven2) == 1,][,i][ven2[rowSums(ven2) == 1,][,i] == 1])

A[i] = colnames(ven2)[i]
}
A
B


出图


pdf("filename.pdf", width = 6, height = 6)

flower_plot(sinID,
sin_num,
ellipse_col = mi[4],
circle_col = mi[8],
90, 0.4, 2,
cer_label = paste("all",":",all_num,sep = ""), col_c = mi[1]
)

dev.off()


下面我构造了一个十个分组的文件来做十个花瓣的类韦恩图。

下面的函数其实之前就写过了，只是不够完善，我在函数中添加了中心标签参数和标签颜色修改参数。这是用基础包编写的，如果有时间我想在使用ggplot编写一份，这样边更加方便了。

sinID = paste("A",1:10,sep = "_")
sin_num = 1:10

pdf("filenamecs.pdf", width = 6, height = 6)

flower_plot(sinID,
sin_num,
ellipse_col = mi[6],
circle_col = mi[2],
90, 0.4, 2,
cer_label = paste("all",":",all_num,sep = ""), col_c = mi[1]
)

dev.off()

作图函数和配色工具

library("plotrix")

flower_plot <- function(sample, value, start, a, b, cer_label,col_c,
ellipse_col = rgb(135, 206, 235, 150, max = 255),
circle_col = rgb(0, 162, 214, max = 255),
circle_text_cex = 1.5
) {
par( bty = "n", ann = F, xaxt = "n", yaxt = "n", mar = c(1,1,1,1))
plot(c(0,10),c(0,10),type="n")
n <- length(sample)
deg <- 360 / n
res <- lapply(1:n, function(t){
draw.ellipse(x = 5 + cos((start + deg * (t - 1)) * pi / 180),
y = 5 + sin((start + deg * (t - 1)) * pi / 180),
col = ellipse_col,
border = ellipse_col,
a = a, b = b, angle = deg * (t - 1))
text(x = 5 + 2.5 * cos((start + deg * (t - 1)) * pi / 180),
y = 5 + 2.5 * sin((start + deg * (t - 1)) * pi / 180),
value[t]
)

if (deg * (t - 1) < 180 && deg * (t - 1) > 0 ) {
text(x = 5 + 3.3 * cos((start + deg * (t - 1)) * pi / 180),
y = 5 + 3.3 * sin((start + deg * (t - 1)) * pi / 180),
sample[t],
srt = deg * (t - 1) - start,
adj = 1,
cex = circle_text_cex
)

} else {
text(x = 5 + 3.3 * cos((start + deg * (t - 1)) * pi / 180),
y = 5 + 3.3 * sin((start + deg * (t - 1)) * pi / 180),
sample[t],
srt = deg * (t - 1) + start,
adj = 0,
cex = circle_text_cex
)
}
})
draw.circle(x = 5, y = 5, r = 1.3, col = circle_col, border = circle_col)

text(x = 5, y = 5,cer_label,cex = 2,col = col_c)
}

library(RColorBrewer)#调色板调用包

# #调用所有这个包中的调色板
# display.brewer.all()
# #提取特定个数的调色板颜色，会出图显示
# display.brewer.pal(9,"Oranges")
# display.brewer.pal(8,"Set1")
# ##仅仅只显示色号,我们要在颜色上并显示色号才好
mi = brewer.pal(9,"Blues")
mi = brewer.pal(12,"Set1")

library("scales")
show_col(mi)

Reference

https://www.researchgate.net/figure/The-pan-genome-of-Sinorhizobium-The-flower-plots-and-Venn-diagrams-illustrate-the-number_fig3_235681265