- 首选选择要计算的区域，这里我选择候选位点上下游1M的区域
grep -v "#" chr4.vcf |awk '{if($2>=68647209&&$2<=70647209){print$3}}'|grep -v "CM000096.5:69647209" >snp1.list
#提取目标区域的所有snp ID 并去掉我们用来与其他snp ID比较的位点
mkdir snp1
echo "4 68647209 70647209"|tr " " "\t">snp1.position
bgzip chr4.vcf
tabix -p vcf chr4.vcf.gz
bcftools filter -R snp1.position chr4.vcf.gz>snp1/snp1.vcf #这样会大大提高计算速度
cd snp1/
for i in `cat ../snp1.list`;do plink --vcf snp1.vcf --ld CM000096.5:69647209 $i;echo CM000096.5:69647209 $i `grep R-sq plink.log` >>snp1.ld;rm -f plink.log;done
