Population-wise decay of linkage disequilibrium (LD)
0
0
Entering edit mode
8.5 years ago
akang ▴ 110

I am trying to calculate the extent of LD in 10 populations. I have a plink file containing all the data in one ped file. Below is the command i used to calculate the ld and then plot the results in R. But the plot shows everything as one. I dont get the population wise plot.I am certainly doing alot wrong here. Do i need to specify family id here? How to go about it?Ill appreciate any help.

--noweb
--ped Input.txt
--map Input.map
--no-fid
--no-parents
--no-sex
--no-pheno
--compound-genotypes
--missing-genotype 9
--allow-no-sex
--maf 0.01
--hwe 0.0001
--geno 0.1
--mind 0.1
--cow
--r2
--ld-window 99999
--ld-window-r2 0
R
a<-(ld$BP_B-ld$BP_A)
plot(a,ld$R2)  Plink plot LD • 4.3k views ADD COMMENT 0 Entering edit mode May be you can solve it on the level of R... How does the plink.ld file looks like? ADD REPLY 0 Entering edit mode It has 7 columns like below and 175288 rows  CHR_A BP_A SNP_A CHR_B BP_B SNP_B R2 1 244698 SNP1 1 524506 SNP2 0.026080900 1 244698 SNP1 1 643558 SNP3 0.029830100  ADD REPLY 0 Entering edit mode It is not too helpful. Anyway, if you have some identification of populations in your file (lets say that CHR_A), than you can just use subset function ld <- read.table("plink.ld",header=T) ld$a <- (ld$BP_B-ld$BP_A) #I recomed to use meaningful names, not "a" for(family in levels(ld$CHR_A)){ pop_ld <- subset(ld, family == ld$CHR_A); pdf(paste('ld',family,'.pdf', sep = '')); plot(pop_ld$a,pop_ld$R2); dev.off() }

P.S. I am sorry, I do not know how to format code here... Just assume, that ; is a newline.

0
Entering edit mode

Thanks @kamiljaron but there are no population identifiers here. I think i need to subset the data right there in Plink while calculating the LD. So does 10 populations mean 10 different plink input .ped files?