No output estimating LD between SNPs in VCF using PLINK

Question

[GabrielMontenegro]

I want to estimate the LD (in r2) between two snps present in a VCF file. However, even if the variants are present in the VCF (checked this using grep) my output is empty (it only has a header).

Command:

plink --vcf myvcf.vcf --r2 --ld-snps 1:756848,1:757000 --out 2nps

Am I doing something wrong? Maybe having the SNP IDs with a ":" is causing this?

Thanks.

Answer 1

[Kevin Blighe]

The colons in your variant IDs should not matter. Here is a test on chr22 variants from 1000 Genomes:

Regular rs IDs:

plink1.90/plink --vcf test.vcf --r2 --ld-snps rs367963583,rs188945759 --out 2nps

cat 2nps.ld 
 CHR_A         BP_A         SNP_A  CHR_B         BP_B         SNP_B           R2 
    22     16050922   rs367963583     22     16050922   rs367963583            1 
    22     16050984   rs188945759     22     16050984   rs188945759            1

Modified IDs:

bcftools annotate -Ov -x ID -I +'%CHROM:%POS:%REF:%ALT' test.vcf > test2.vcf

plink1.90/plink --vcf test2.vcf --r2 --ld-snps 22:16051249:T:C-22:16052639:C:T --out 2nps

cat 2nps.ld 
 CHR_A         BP_A                 SNP_A  CHR_B         BP_B                 SNP_B           R2 
    22     16051249       22:16051249:T:C     22     16051249       22:16051249:T:C            1 
    22     16051249       22:16051249:T:C     22     16051453       22:16051453:A:C     0.780694 
    22     16051453       22:16051453:A:C     22     16051249       22:16051249:T:C     0.780694 
    22     16051453       22:16051453:A:C     22     16051453       22:16051453:A:C            1 
    22     16051453       22:16051453:A:C     22     16052962       22:16052962:C:T     0.625686 
    22     16051722      22:16051722:TA:T     22     16051722      22:16051722:TA:T            1 
    22     16052080       22:16052080:G:A     22     16052080       22:16052080:G:A            1 
    22     16052167   22:16052167:A:AAAAC     22     16052167   22:16052167:A:AAAAC            1 
    22     16052240       22:16052240:C:G     22     16052240       22:16052240:C:G            1 
    22     16052271       22:16052271:G:A     22     16052271       22:16052271:G:A            1 
    22     16052639       22:16052639:C:T     22     16052639       22:16052639:C:T            1

It will not work if your variant IDs are duplicated or if you have only a single sample in your VCF.

Kevin