之前用TBtools做了一次圈图，颜值确实不低，circos早就折腾过，奈何没学过perl，用起来特别不顺手，顺便吐槽下circos的网站，跟着教程来，确实..一言难尽...
毕竟在接触过的这些语言中R还是最顺手的，早就听说circlize这个包挺不错，还没用过，今天花了一天时间，拿着cookbook现学现卖，基本100%复刻了之前那个TBtools的circos图。
不过折腾过circlize，circos最后得到一个结论，对于简单的圈图来说，TBtools真香！！
这里只放代码，有空了详细写下解读，不过我做的感觉已经够傻瓜的了..
还是国际惯例，先上图
###########################################
# Prj: Heterosis ==> circRNA
# Assignment: circos plot
# Author: Shawn Wang
# Date: Sep 11,2020
###########################################
library(circlize)
## 1st ==> need chr length to plot the chr track.
## import chrlen
chrlen <- read.table("~/02.Project/03.circRNA/CircRNAVersion2/02.data/04.other/Gh.Chr.Len.xls",
header = F,
sep = "\t",
stringsAsFactors = F)
head(chrlen)
## start from 0
chr = data.frame(Chr = chrlen$V1,
start = 0,
end = chrlen$V2)
## remove contig
chr = filter(chr,str_detect(chr$Chr,"^[A-D][0-9][0-9]"))
## 2nd ==> gene density plot. need chr bin with gene number
## import gene bin file
demo <- read.table("~/02.Project/03.circRNA/CircRNAVersion2/02.data/04.other/Gh.gene_density.txt",
header = F,
sep = "\t",
stringsAsFactors = F)
demo = data.frame(chr = demo$V1,
start = demo$V2 - demo[1,2],
end = demo$V2,
value = demo$V3)
## 3rd ==> circRNA and miRNA postion
## in this cycle, I'll plot miRNA as blue dot and circRNA as green dot in same layer.
# import miRNA postion file
mi <- read.delim("~/02.Project/03.circRNA/CircRNAVersion2/02.data/05.circos/tmp/04.mi-pos.xls",
header = F,
sep = "\t",
stringsAsFactors = F)
# add values
mi.pos = data.frame(Chr = mi$V2,
start = mi$V3,
end = mi$V4,
value = 10)
## adjust start and end position
for (i in 1:nrow(mi.pos)) {
start = mi.pos[i,2]
end = mi.pos[i,3]
if (start < end) {
mi.pos[i,2] = start
mi.pos[i,3] = end
} else{
mi.pos[i,2] = end
mi.pos[i,3] = start
}
}
## double check the start-end position
mi.pos[,2] < mi.pos[,3]
## import circ.pos
circ = read.delim("~/02.Project/03.circRNA/CircRNAVersion2/02.data/05.circos/02.workingDir/highlight/gh.circRNA.pos",
header = F,
sep = "\t",
stringsAsFactors = F)
circ.pos = data.frame(Chr = circ$V1,
start = circ$V2,
end = circ$V3,
value = 20)
head(circ.pos)
## merge two position file
bed1 = circ.pos
bed2 = mi.pos
bed_list = list(bed1,bed2)
## links
## At2At
At2At = read.table("~/02.Project/03.circRNA/CircRNAVersion2/02.data/05.circos/02.workingDir/links/AtvsAt.link",
header = F,sep = "\t",stringsAsFactors = F)
fromAt1= At2At[,c(1:3)]
toAt1 = At2At[,c(4:6)]
## At2Dt
At2Dt = read.table("~/02.Project/03.circRNA/CircRNAVersion2/02.data/05.circos/02.workingDir/links/AtvsDt.link",
header = F,sep = "\t",stringsAsFactors = F)
fromAt2= At2Dt[,c(1:3)]
toDt2 = At2Dt[,c(4:6)]
## Dt2Dt
Dt2Dt = read.table("~/02.Project/03.circRNA/CircRNAVersion2/02.data/05.circos/02.workingDir/links/DtvsDt.link",
header = F,sep = "\t",stringsAsFactors = F)
fromDt3= Dt2Dt[,c(1:3)]
toDt3 = Dt2Dt[,c(4:6)]
## Dt2At
Dt2At = read.table("~/02.Project/03.circRNA/CircRNAVersion2/02.data/05.circos/02.workingDir/links/DtvsAt.link",
header = F,sep = "\t",stringsAsFactors = F)
fromDt4= Dt2At[,c(1:3)]
toAt4 = Dt2At[,c(4:6)]
## plot
circos.clear() ## cleaning last work
## track
circos.genomicInitialize(chr,major.by = 30000000)
## colored At and Dt subgenome with different color
chrcolor = rep(c("turquoise","#FF6347"),each = 13)
circos.genomicTrackPlotRegion(chr,ylim=c(0,0.1),track.height=0.03,bg.col=chrcolor,cell.padding=c(0.01, 0.5, 0.01, 0.5),track.margin=c(0, 0))
## gene density heatmap
min1 = quantile(demo$value)[1]
quant25 = quantile(demo$value)[2]
mid = quantile(demo$value)[3]
quant75 = quantile(demo$value)[4]
max = quantile(demo$value)[5]
col_fun = colorRamp2(c(min1,quant75, max), c("green","black","red"))
circos.genomicTrack(demo, numeric.column = 4,track.height = 0.06,ylim = c(0,1.2),
panel.fun = function(region,value,...){
circos.genomicRect(region, value, col = col_fun(value[[1]]), border = NA)
})
circos.genomicTrack(bed_list, track.height = 0.06,
panel.fun = function(region, value, ...) {
i = getI(...)
circos.genomicPoints(region, value, pch = 16, cex = 0.3, col = i, ...)
})
## links
# At2Dt
circos.genomicLink(fromAt2, toDt2, col = "#1E90FF", border = "#1E90FF", lwd = 0.15)
# Dt2At
circos.genomicLink(fromDt4, toAt4, col = "#FFD700", border = "#FFD700", lwd = 0.15)
# At2At
circos.genomicLink(fromAt1, toAt1, col = "#7FFFAA", border = "#7FFFAA",lwd = 0.12)
# Dt2Dt
circos.genomicLink(fromDt3, toDt3, col = "#FF6347", border = "#FF6347", lwd = 0.12)
