Question: Population-wise decay of linkage disequilibrium (LD)
0
3.1 years ago by
akang90
akang90 wrote:

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)
``````
plot ld plink • 1.6k views
modified 3.1 years ago • written 3.1 years ago by akang90

May be you can solve it on the level of R... How does the plink.ld file looks like?

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
``````

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.