Post date: Nov 23, 2015 10:17:27 PM
I extended the analysis described here to compare striped vs. unstriped.
Here are the results for LD for G.S. In contrast to the melanistic vs. green analyses, there is no evidence that LD between morphs relative to within is greater for the high HMM regions than expected by chance for the LG:
FHA: obs = -0.0093, null = -0.0012, 95th quantile = 0.00596, p = 0.993
N1: obs = -0.00046, null = -0.0015, 95th quantile = 0.00437, p = 0.374
In fact, for FHA, one could argue that LD is significantly reduced between vs. within morhps for the high HMM regions (not true in N1, and could mean nothing anyway).
Here is the final R script (testLd.R):
load("hmm2s_FHA_RxGISpar7.rwsp")
fstfha<-fit[[1]]
ldfha<-read.table("BurrowD_binnedSNPs_morphs_FHA.txt",header=TRUE)
L<-dim(ldfha)[1]
fstld<-rep(NA,L)
for(i in 1:L){
x<-which(nei[,1]==8 & nei[,2]==ldfha[i,2] & nei[,4]==ldfha[i,4])
if(length(x)==1){
fstld[i]<-fstfha[x]
}
}
x<-which(fstld==1)
## melanistic vs. green
obs<-mean(ldfha$D0_M.GS[x]-(ldfha$D0_M[x] + ldfha$D0_G.S[x])/2,na.rm=TRUE)
null<-rep(NA,1000)
for(i in 1:1000){
st<-sample(1:L,1)
y<-c(fstld[st:L],fstld[-c(st:L)])
x<-which(y==1)
null[i]<-mean(ldfha$D0_M.GS[x]-(ldfha$D0_M[x] + ldfha$D0_G.S[x])/2,na.rm=TRUE)
}
## striped vs. unstriped
obs<-mean(ldfha$D0_G.S[x]-(ldfha$D0_G[x] + ldfha$D0_S[x])/2,na.rm=TRUE)
null<-rep(NA,1000)
for(i in 1:1000){
st<-sample(1:L,1)
y<-c(fstld[st:L],fstld[-c(st:L)])
x<-which(y==1)
null[i]<-mean(ldfha$D0_G.S[x]-(ldfha$D0_G[x] + ldfha$D0_S[x])/2,na.rm=TRUE)
}
load("hmm2s_N1_RxGISpar12.rwsp")
fstN1<-fit[[1]]
ldN1<-read.table("BurrowD_binnedSNPs_morphs_N1.txt",header=TRUE)
L<-dim(ldN1)[1]
fstld<-rep(NA,L)
for(i in 1:L){
x<-which(nei[,1]==8 & nei[,2]==ldfha[i,2] & nei[,4]==ldN1[i,4])
if(length(x)==1){
fstld[i]<-fstN1[x]
}
}
x<-which(fstld==1)
## melanistic vs. green
obs<-mean(ldN1$D0_M.GS[x]-(ldN1$D0_M[x] + ldN1$D0_G.S[x])/2,na.rm=TRUE)
null<-rep(NA,1000)
for(i in 1:1000){
st<-sample(1:L,1)
y<-c(fstld[st:L],fstld[-c(st:L)])
x<-which(y==1)
null[i]<-mean(ldN1$D0_M.GS[x]-(ldN1$D0_M[x] + ldN1$D0_G.S[x])/2,na.rm=TRUE)
}
## striped vs. un-striped
obs<-mean(ldN1$D0_G.S[x]-(ldN1$D0_G[x] + ldN1$D0_S[x])/2,na.rm=TRUE)
null<-rep(NA,1000)
for(i in 1:1000){
st<-sample(1:L,1)
y<-c(fstld[st:L],fstld[-c(st:L)])
x<-which(y==1)
null[i]<-mean(ldN1$D0_G.S[x]-(ldN1$D0_G[x] + ldN1$D0_S[x])/2,na.rm=TRUE)
}