五一是小麥人一年裡最忙碌的季節。在此小編代表小麥研究聯盟公眾號所有編輯,祝忙碌在科研一線的小麥工作者節日快樂,勞有所獲!
在節日裡,給大家帶了個壓箱底的東西,一幅圖的繪製方法。2018年中國春基因組在science發表時候,一幅圖非常驚豔。就是下面這個。
圖1. 中國春文章發表時候配圖--圈圈圖
這個圖是基因組測序文章的標配,很多人一看就知道是圈圈圖,還知道有一個有名的軟體Circos,是專門畫這種圖的。但是這個軟體是perl語言寫的,用過的人一看就知道,這圖的配色風格,絕對不是perl。所以這個圖應該不是perl,而可能是畫圖配色功能尤其強大的R。不錯,這個圖就是R的circlize包畫的。其實這個圖關鍵是配色方案,牛逼的像六月的小麥,金光閃閃。這個配色方案呢,其實是R另一個包viridis。可以網上查這兩個包,查看幫助文檔,學習更多。
下面就給你放上一段代碼,讓你能夠畫出這樣的圖,你需要注意的是,這裡面的數據都是模擬的數據,要畫自己的圖,你需要根據自己的生物學意義,按照格式進行數據的整理。
圖2. 本教程出圖
library(viridis)
library(RColorBrewer)
display.brewer.all()
n <- 30
track.h <- 0.08
cols = viridis(n)
library("circlize")
stand.inter <- function(x,n){
x.boder <- max(x) - min(x)
x.times <- n/x.boder
x.after <- round((x-min(x))*x.times)
}
setwd("C:\\Users\\.../")
chr <- read.csv(file="Chr.info.csv",head=T,as.is=T)
names(chr) <- c("chr","length","cent")
head(chr)
sample.data <- NULL
for(m in 1:7){
for (k in c("A","B","D")) {
start <- seq(1,500e6,5e6)
end <- start + 5e6
chr.at <- paste("Chr",paste(m,k,sep=""),sep="")
value <- -(start/250e6-1)^2
value[start > 240e6 & start < 260e6] <- 0.4
sub.data <- data.frame(chr=chr.at,start=start,end=end,value=value)
sample.data <- rbind(sample.data,sub.data)
}
}
sample.data$value.adjust <- stand.inter(sample.data$value,n)
pdf(file="小麥圈圈圖.pdf")
### 開始畫圖 ####
circos.clear()
circos.par(cell.padding=c(0,0,0,0),start.degree = 86, gap.degree =c(rep(2,20),2))
xlim <- matrix(NA,21,2);xlim[,1] <- 0;xlim[,2] <- chr$length/1000000
circos.initialize(factors = sample.data$chr, xlim=xlim,x=sample.data$end)
circos.track(factors = as.factor(chr$chr),x=chr$length/1000000,y=chr$cent/1000000,ylim=c(0,1),"track.height" = 0.04,
panel.fun = function(x, y) {
circos.rect(0, 0,y, 1,col = "grey80", border =NA)
circos.rect(y, 0,x, 1,col = "grey50", border =NA)
circos.text(CELL_META$xcenter, CELL_META$cell.ylim[2] + uy(6, "mm"),
CELL_META$sector.index)
circos.axis(labels.cex = 0.6)
})
circos.track(factors = sample.data$chr,ylim=c(0,1),"track.height" = track.h,y=sample.data$value.adjust,x=sample.data$end/1e6,
bg.border=NA,
panel.fun = function(x, y) {
circos.rect(x-5, rep(0,length(x)),x, rep(1,length(x)),
col = cols[y], border =NA)
})
dev.off()
染色體配置文件(實際畫圖時候,請按照實際內容進行修改)
再次祝大家節日快樂!