论文
A telomere-to-telomere gap-free reference genome of watermelon and its mutation library provide important resources for gene discovery and breeding
今天的推文我们重复一下论文中的Figure1a
之前的推文介绍了如何绘制外圈的圈图,但是当时不知道如何实现图例,后来想了一下,图例其实可以用ggplot2来做,然后把ggplot2的内容叠加的circlize的图上。搜索了一下这种思路的实现代码,找到链接
https://cran.r-project.org/web/packages/ggplotify/vignettes/ggplotify.html
按照这个思路来构建图例
图例里就是各种点 线段 文字的组合,所以先构造这些元素的坐标,然后用ggplot2的代码依次叠加,这种方式我觉得稍微有些麻烦,但是目前还想不到比较好的办法
图例构造数据和画图的代码
library(ggplot2)
text.df01<-data.frame(
x=0.0125,
y=c(1,3,5,7,9),
label=c("InDels","SNPs",
"TE density","Gene density",
"GC content")
)
text.df02<-data.frame(
x=-0.0125,
y=c(3,7),
label=c("Inner","Outer")
)
segment.df01<-data.frame(
x=0.025,
xend=0.025,
y=c(0.5,2.5,8.5),
yend=c(1.5,3.5,9.5)
)
segment.df02<-data.frame(
x=0.025,xend=0.027,
y=c(0.5,1.5,2.5,3.5,8.5,9.5),
yend=c(0.5,1.5,2.5,3.5,8.5,9.5)
)
text.df03<-data.frame(
x=0.030,
y=c(0.5,1.5,2.5,3.5,8.5,9.5),
label=c(0,446,0,3270,0.25,0.35)
)
point.df01<-data.frame(
x=c(0.025,0.035),
y=c(5)
)
point.df02<-data.frame(
x=c(0.025,0.03,0.035),
y=c(7)
)
text.df04<-data.frame(
x=c(0.025,0.025,0.035,0.035),
y=c(5,7,5,7),
label=c("Low","Low","High","High")
)
ggplot()+
geom_segment(aes(x=0,y=10,xend=0,yend=0),
arrow = arrow(angle=20,
type = "closed",
length = unit(5,'mm')),
size=1)+
geom_segment(aes(x=0,y=10,xend=0,yend=1),
size=1)+
geom_text(data=text.df01,
aes(x=x,y=y,label=label),
size=2)+
geom_text(data=text.df02,aes(x=x,y=y,label=label),
size=2)+
geom_segment(data=segment.df01,
aes(x=x,y=y,xend=xend,yend=yend),
color=c("#d3026e","#1e4793",'#860a3b'),
size=0.5)+
geom_segment(data=segment.df02,
aes(x=x,y=y,xend=xend,yend=yend),
color=c("#d3026e","#d3026e",
"#1e4793","#1e4793",
'#860a3b','#860a3b'),
size=0.5)+
geom_text(data=text.df03,
aes(x=x,y=y,label=label),
size=2,hjust=1)+
geom_point(data=point.df01,
aes(x,y),
shape=15,
color=c("#0037ff","#fd0200"),
size=2)+
geom_point(data=point.df02,
aes(x,y),
shape=15,
color=c("#1efc05","#000000","#fc0003"),
size=2)+
geom_text(data=text.df04,
aes(x=x,y=y,label=label),
vjust=-2,
size=2)+
theme_void() +
xlim(-0.02,0.04) -> p1
p1
然后是圈图叠加图例的代码
最终结果
library(circlize)
library(grid)
library(ggplotify)
brk <- seq(0,40,by=2)*10^6
brk.label<-c()
for (i in brk){
ifelse(i%%10^7==0,brk.label<-append(brk.label,
paste0(i/10^7,"0M")),
brk.label<-append(brk.label,""))
}
brk.label[1]<-"0M"
brk.label
#circos.par(points.overflow.warning = FALSE)
pdf(file = "mp.pdf",
width = 10,height = 10)
circos.par(start.degree =86,clock.wise = T,
#track.height=0.00001,
track.margin = c(.001,0.001))
circos.initialize(factors = df$chr_id,
xlim = matrix(c(rep(0,11),df$chr_len),ncol=2))
circos.trackPlotRegion(df$chr_id,
ylim = c(0, 10),
track.height = 0.1,
bg.border = NA,
#ylim=CELL_META$ylim,
panel.fun = function(x, y) {
circos.text(mean(CELL_META$xlim), 12,
get.cell.meta.data("sector.index"))
})
circos.trackPlotRegion(df$chr_id,
ylim = c(0, 100),
track.height = 0.1,
bg.col = '#EEEEEE6E',
bg.border = NA)
for (chromosome in df$chr_id){
circos.axis(sector.index = chromosome,
h = 110,
major.at = brk,
minor.ticks = 0,
labels = brk.label,
labels.facing="clockwise",
labels.cex = 0.6)
}
circos.trackLines(gc$chr_id,gc$bin_start,gc$gc,
area = TRUE,
col = "red",
border="transparent")
circos.trackPlotRegion(df$chr_id,
ylim = c(0, 10),
track.height = 0.1,
bg.col = '#EEEEEE6E',
bg.border = NA)
for (chromosome in df$chr_id){
circos.segments(sector.index = chromosome,
x0=genedensity[genedensity$chr_id==chromosome,]$bin_start,
x1=genedensity[genedensity$chr_id==chromosome,]$bin_start,
y0=0,y1=10,col=genedensity[genedensity$chr_id==chromosome,]$group)
}
circos.trackPlotRegion(df$chr_id,
ylim = c(0, 10),
track.height = 0.1,
bg.col = '#EEEEEE6E',
bg.border = NA)
for (chromosome in df$chr_id){
circos.segments(sector.index = chromosome,
x0=genedensity[genedensity$chr_id==chromosome,]$bin_start,
x1=genedensity[genedensity$chr_id==chromosome,]$bin_start,
y0=0,y1=10,col=genedensity[genedensity$chr_id==chromosome,]$group)
}
circos.trackPlotRegion(df$chr_id,
ylim = c(0, 100),
track.height = 0.1,
bg.col = '#EEEEEE6E',
bg.border = NA)
for (chromosome in df$chr_id){
circos.barplot(sector.index = chromosome,
value = genedensity[genedensity$chr_id==chromosome,]$gene_count,
pos = genedensity[genedensity$chr_id==chromosome,]$bin_start,
col = "red",
bar_width = 500000,
border="transparent")
}
circos.trackPlotRegion(df$chr_id,
ylim = c(0, 100),
track.height = 0.1,
bg.col = '#EEEEEE6E',
bg.border = NA)
for (chromosome in df$chr_id){
circos.barplot(sector.index = chromosome,
value = genedensity[genedensity$chr_id==chromosome,]$gene_count,
pos = genedensity[genedensity$chr_id==chromosome,]$bin_start,
col = "blue",
bar_width = 500000,
border="transparent")
}
p2 <- as.grob(p1)
vp = viewport(x=.5, y=.5, width=.2, height=.2)
pushViewport(vp)
grid.draw(p2)
upViewport()
circos.clear()
dev.off()
有些文字的位置不是太合适,需要出图后调整
这里遇到一个问题是,我想把中间的空白区间刘大一点,我想到的是把每一圈的轨道宽度设置小一点,但是暂时没有找到参数进行设置,找到一个track.height好像不起作用。
构造图例的时候需要反复调整每个元素的位置坐标,还挺耗时间的
示例数据和代码可以给公众号推文点赞,点击在看,最后留言获取
欢迎大家关注我的公众号
小明的数据分析笔记本
小明的数据分析笔记本 公众号 主要分享:1、R语言和python做数据分析和数据可视化的简单小例子;2、园艺植物相关转录组学、基因组学、群体遗传学文献阅读笔记;3、生物信息学入门学习资料及自己的学习笔记!