在GWAS文章中,我们经常会看到SNP连锁不平衡图,该图可以直观地将SNP间连锁不平衡程度展示出来。今天来教大家使用R包“LDheatmap”快速绘制SNP连锁不平衡图。
https://doi.org/10.1111/pbi.13126
在绘图前,我们需要准备两个输入文件:
第一个文件为SNP标记的基因型信息。第一行为SNP编号,每列为SNP标记在各个样本上的基因型。
第二个文件为SNP标记的位置信息。每一行的数字代表与之对应的SNP物理位置。
准备好这两个文件后我们就可以开始画图啦!
## 导入R包
library("LDheatmap")
library("genetics")
## 导入SNP标记基因型信息
SNPdata <- read.table("snp.geno.txt",header=T,sep="\t")
## 导入SNP标记位置信息
SNPpos <- read.table("snp.pos.txt",header=F,sep="\t")
## 将SNP基因型信息转换为genotype格式
num <- ncol(SNPdata)
for(i in 1:num){
SNPdata[,i]<-as.genotype(SNPdata[,i])
}
## 将SNP位置信息转换为vector
pos <- as.vector(unlist(SNPpos))
## 设置热图颜色
color.rgb <- colorRampPalette(rev(c("white","red")),space="rgb")
## 绘制连锁不平衡图
LDheatmap(SNPdata,pos,color=color.rgb(20),flip=TRUE)
我们的连锁不平衡图画好啦!
但是图中方块间分割不是非常明显,我们可以输入下面的代码给方块之间加分割线。
library("grid")##添加分割线grid.edit(gPath("ldheatmap","heatMap","heatmap"),gp=gpar(col="white",lwd=1.5))
这次我们的图就顺眼多了!
如果你想在图中标记某个SNP,可以输入下面的代码:
##添加SNP标记名称LDheatmap(SNPdata,pos,color=color.rgb(20),flip=TRUE,SNP.name=c("SNP1","SNP2"))##调整SNP标记名称的字体大小和颜色
grid.edit(gPath("ldheatmap", "geneMap","SNPnames"), gp = gpar(col="black",cex=1.5))
最终效果看起来还不错!
为了节省大家整理两个输入文件的时间,我写了一个python脚本,直接输入vcf文件和位置信息即可获得连锁不平衡图,用法如下:
##该脚本在Linux下使用,使用前需安装python、R及R包"LDheatmap"、"genetics"和"grid"
python ./LDplot.py -vcf ./snp.vcf -pos snp.pos.txt -chr chr_name -out ./out_prefix
-vcf 输入包含SNP基因型的vcf文件
-pos 输入需要作图的连锁标记的位置(与上文所讲位置文件一致)
-chr 输入需要作图的连锁标记的染色体名称
-out 输出文件名称的前缀
##输出文件包含连锁标记的基因型、位置及连锁不平衡图。
脚本存放在(https://github.com/biozhp/LDplot),欢迎大家使用!
参考资料:
https://sfustatgen.github.io/LDheatmap/
https://www.jianshu.com/p/0b92a3134a5c
https://www.biostars.org/p/105577/