circlize这个包挺强大的,R语言里用来画圈图非常方便。
今天这篇文章记录用circlize这个包画圈图展示blast双序列比对结果的代码
植物线粒体基因组类的文章通常会分析细胞器基因组间基因转移情况,基本的分析方法就是blast比对。可视化展示可以选择用这个圈图来做
首先是使用blast建库比对
makeblastdb -in mt.fasta -dbtype nucl -out mt blastn -query cp.fasta -db mt -outfmt 6 > output.txt
然后是准备数据用来画最外圈用来展示两条序列的部分
df<-data.frame(chr=c(rep("chloroplast",2),rep("mitochondrial",2)), x=c(1,131478,1,444567), y=c(0,1,0,1)) > df chr x y 1 chloroplast 1 0 2 chloroplast 131478 1 3 mitochondrial 1 0 4 mitochondrial 444567 1
然后是读入blast的输出结果
df1<-read.csv("output6.txt",stringsAsFactors = F,header=F,sep="\t")
作图用到的代码
library(circlize) library(RColorBrewer) library(ComplexHeatmap) col<-RColorBrewer::brewer.pal(6,"Paired") circos.par("start.degree" = 130) circos.initialize(factors = df$chr,x=df$x) circos.trackPlotRegion(factors = df$chr,y=df$y, panel.fun = function(x,y){ circos.axis() },track.height = 0.1) highlight.sector(sector.index = "chloroplast",col=col[1]) highlight.sector(sector.index = "mitochondrial",col=col[2]) circos.text(x=70000,y=0.5, labels = "chloroplast", sector.index = "chloroplast") circos.text(x=220000,y=0.5, labels = "mitochondrial", sector.index = "mitochondrial", facing = "outside") col_fun = colorRamp2(c(70,90,100), c("green", "yellow", "red")) for (i in 1:13){ x<-sort(c(df1[i,8],df1[i,7])) y<-sort(c(df1[i,10],df1[i,9])) z<-df1[i,3] circos.link("chloroplast",x,"mitochondrial",y, col=add_transparency(col_fun(z))) } circos.clear() lgd_links = Legend(at = c(70, 80, 90, 100), col_fun = col_fun, title_position = "topleft", title = "identity(%)") lgd_list_vertical = packLegend(lgd_links) draw(lgd_list_vertical, x = unit(10, "mm"), y = unit(10, "mm"), just = c("left", "bottom"))
image.png
调整整体的角度
circos.par("start.degree" = 130)
调整用来表示染色体的外圈粗细
circos.trackPlotRegion(factors = df$chr,y=df$y, panel.fun = function(x,y){ circos.axis() },track.height = 0.1)
画图的时候可以加一个track.height参数
调整外圈的刻度,现在展示的有点多,我想增大间隔,减少展示的数字,暂时不知道如何实现。
https://jokergoo.github.io/circlize_book/book/legends.html
本文分享自微信公众号 - 小明的数据分析笔记本(gh_0c8895f349d3),作者:Punicagranatum
原文出处及转载信息见文内详细说明,如有侵权,请联系 yunjia_community@tencent.com 删除。
原始发表时间:2020-08-30
本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。
我来说两句