安装必要软件包 LDBlockShow
conda create -n ldblockshow
conda activate ldblockshow
conda install -c bioconda perl-svg
conda install -c bioconda ldblockshow
export PERL5LIB=/root/anaconda3/envs/ldblockshow/lib/site_perl/5.26.2:$PERL5LIB
提取单独的vcf(基于每个髙显著SNP前后0.05Mb的区间)
bcftools annotate -x INFO combine_all.pass.recode.vcf -O z -o combine_all.pass.clear.recode.vcf.gz
1. Region 1: chr09 (25.0 - 25.12 Mb)
包含 Rhdel09G0125900, Rhdel09G0126000
vcftools --gzvcf combine_all.pass.clear.recode.vcf.gz --chr chr09 --from-bp 25000000 --to-bp 25120000 --recode --recode-INFO-all --stdout | gzip -c > Reg1_chr09_25.0Mb.vcf.gz
2. Region 2: chr09 (30.8 - 30.93 Mb)
包含 Rhdel09G0159500, Rhdel09G0159600, Rhdel09G0159700
vcftools --gzvcf combine_all.pass.clear.recode.vcf.gz --chr chr09 --from-bp 30800000 --to-bp 30930000 --recode --recode-INFO-all --stdout | gzip -c > Reg2_chr09_30.8Mb.vcf.gz
3. Region 3: chr09 (15.1 - 15.9 Mb)
包含 Rhdel09G0097100, Rhdel09G0097200, Rhdel09G0100900
vcftools --gzvcf combine_all.pass.clear.recode.vcf.gz --chr chr09 --from-bp 15100000 --to-bp 15900000 --recode --recode-INFO-all --stdout | gzip -c > Reg3_chr09_15.5Mb.vcf.gz
4. Region 4: chr01 (58.01 - 58.04 Mb)
包含 Rhdel01G0327700, Rhdel01G0327800, Rhdel01G0327900
vcftools --gzvcf combine_all.pass.clear.recode.vcf.gz --chr chr01 --from-bp 58010000 --to-bp 58040000 --recode --recode-INFO-all --stdout | gzip -c > Reg4_chr01_58.0Mb.vcf.gz
5. Region 5: chr03 (8.94 - 8.99 Mb)
包含 Rhdel03G0056200, Rhdel03G0056300, Rhdel03G0056400
vcftools --gzvcf combine_all.pass.clear.recode.vcf.gz --chr chr03 --from-bp 8940000 --to-bp 8990000 --recode --recode-INFO-all --stdout | gzip -c > Reg5_chr03_8.9Mb.vcf.gz
6. Region 6: chr09 (21.1 - 21.2 Mb)
对应三个极强的连续信号 chr09: 21,163,565 / 21,163,594 / 21,163,579 以及 chr09: 20,323,347 / 20,323,348 chr09: 14,066,909
vcftools --gzvcf combine_all.pass.clear.recode.vcf.gz --chr chr09 --from-bp 21100000 --to-bp 21200000 --recode --recode-INFO-all --stdout | gzip -c > Reg6_chr09_21.1Mb.vcf.gz
7. Region 7: chr08 (12.1 - 12.2 Mb)
chr08: 12,145,370
vcftools --gzvcf combine_all.pass.clear.recode.vcf.gz --chr chr08 --from-bp 12100000 --to-bp 12200000 --recode --recode-INFO-all --stdout | gzip -c > Reg7_chr08_12.1Mb.vcf.gz
8. Region 8: chr13 (1.75 - 1.85 Mb)
chr13: 1,793,076
vcftools --gzvcf combine_all.pass.clear.recode.vcf.gz --chr chr13 --from-bp 1750000 --to-bp 1850000 --recode --recode-INFO-all --stdout | gzip -c > Reg8_chr13_1.8Mb.vcf.gz
需要准备的文件
VCF文件(基于上一步)
Reg1_chr09_25.0Mb.vcf.gz
Reg2_chr09_30.8Mb.vcf.gz
Reg3_chr09_15.5Mb.vcf.gz
Reg4_chr01_58.0Mb.vcf.gz
Reg5_chr03_8.9Mb.vcf.gz
Reg6_chr09_21.1Mb.vcf.gz
Reg7_chr08_12.1Mb.vcf.gz
Reg8_chr13_1.8Mb.vcf.gz
GWAS 结果文件(.assoc)
格式必须包含 Chr, Pos, P-value 三列。可以只保留该区域内的 SNP 以加快速度。
gwas_signals.assoc
GFF3 基因注释文件
Rhododendron_delavayi.gene.gff3
LDBlockShow绘制
Region 1: chr09 (25.06 Mb 附近)
LDBlockShow -InVCF Reg1_chr09_25.0Mb.vcf.gz -OutPut Reg1_Result -Region chr09:25050000:25080000 -InGWAS gwas_signals.assoc -InGFF Rhododendron_delavayi.gene.gff3 -SeleVar 3 -TopSite -Cutline 5.5
Region 2: chr09 (30.86 Mb 附近)
LDBlockShow -InVCF Reg2_chr09_30.8Mb.vcf.gz -OutPut Reg2_Result -Region chr09:30840000:30880000 -InGWAS gwas_signals.assoc -InGFF Rhododendron_delavayi.gene.gff3 -SeleVar 3 -TopSite -Cutline 5.5
Region 3: chr09 (15.2-15.8 Mb) - 由于此区域跨度大,建议先看 15.2Mb 附近
LDBlockShow -InVCF Reg3_chr09_15.5Mb.vcf.gz -OutPut Reg3_Result -Region chr09:15210000:15250000 -InGWAS gwas_signals.assoc -InGFF Rhododendron_delavayi.gene.gff3 -SeleVar 3 -TopSite -Cutline 5.5
Region 4: chr01 (58.02 Mb 附近)
LDBlockShow -InVCF Reg4_chr01_58.0Mb.vcf.gz -OutPut Reg4_Result -Region chr01:58010000:58040000 -InGWAS gwas_signals.assoc -InGFF Rhododendron_delavayi.gene.gff3 -SeleVar 3 -TopSite -Cutline 5.5
Region 5: chr03 (8.96 Mb 附近)
LDBlockShow -InVCF Reg5_chr03_8.9Mb.vcf.gz -OutPut Reg5_Result -Region chr03:8910000:8990000 -InGWAS gwas_signals.assoc -InGFF Rhododendron_delavayi.gene.gff3 -SeleVar 3 -TopSite -Cutline 5.5
Region 6: chr09 (21.16 Mb 附近) - 这一簇信号非常集中
LDBlockShow -InVCF Reg6_chr09_21.1Mb.vcf.gz -OutPut Reg6_Result -Region chr09:21150000:21180000 -InGWAS gwas_signals.assoc -InGFF Rhododendron_delavayi.gene.gff3 -SeleVar 3 -TopSite -Cutline 5.5
Region 7: chr08 (12.14 Mb 附近)
LDBlockShow -InVCF Reg7_chr08_12.1Mb.vcf.gz -OutPut Reg7_Result -Region chr08:12130000:12160000 -InGWAS gwas_signals.assoc -InGFF Rhododendron_delavayi.gene.gff3 -SeleVar 3 -TopSite -Cutline 5.5
Region 8: chr13 (1.79 Mb 附近)
LDBlockShow -InVCF Reg8_chr13_1.8Mb.vcf.gz -OutPut Reg8_Result -Region chr13:1780000:1810000 -InGWAS gwas_signals.assoc -InGFF Rhododendron_delavayi.gene.gff3 -SeleVar 3 -TopSite -Cutline 5.5