使用R语言包circlize可视化展示blast双序列比对结果

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