天天看點

根據vcf檔案計算SNP密度并用circlize可視化結果

參考

收集vcftools所有用法

指令

vcftools --vcf snp.bialles.vcf --SNPdensity 100000 --out StatResults/SNPdensity
           

複制

100000 是指定視窗長度

--out

是輸出檔案的字首

使用R語言中的circlize包畫圖

參考

用circlize包繪制circos plot

代碼

df<-read.table("SNPdensity.snpden",sep="\t",header=T)
head(df)
df<-df[,c(1,2,4)]
colnames(df)<-c("Chr","X","Y")
head(df)
df$X<-df$X/1000000
options(scipen=999)
library(circlize)
library(RColorBrewer)
col<-RColorBrewer::brewer.pal(8,"Paired")
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.05)

for(i in 1:8){
  highlight.sector(sector.index = paste0("LG",i),col=col[i])
  circos.text(CELL_META$xcenter, CELL_META$ycenter,
              labels = paste0("LG",i),
              sector.index = paste0("LG",i),cex=0.5)
}
circos.trackPlotRegion(factors=df$Chr,y=df$Y)
circos.trackLines(df$Chr,df$X,df$Y,col=col)
circos.trackPlotRegion(factors=df$Chr,y=df$Y)
circos.trackPoints(df$Chr,df$X,df$Y,col=col,cex=0.5)
circos.trackPlotRegion(factors=df$Chr,y=df$Y)
circos.trackHist(df$Chr,df$X,col=col)
circos.clear()
           

複制

結果

根據vcf檔案計算SNP密度并用circlize可視化結果

image.png

還想實作的其他功能:

  • 想開一個口子用倆給每一圈添加文字标簽,如何用代碼實作還不太清楚,想到一個解決辦法是多加一條染色體,然後出圖後手工編輯将這條染色體删除掉,然後添加文字标簽
  • 如何填加一圈柱形圖呢?
  • vcftools計算SNP密度的時候是否可以設定滑動視窗?