Post date: Nov 21, 2017 12:9:35 AM
Patrik and Clarissa suggested lateral color might be better and this lets us add T. podura. Thus, I re-ran the analyses from here, with lateral color. Below are the results and additional code.
Overall these actually make for an even nicer story, with one truly annoying exception. T. cristinae comes out as having more QTL/QTN, which appears to be mostly driven by RG. There wasn't really any trend in this direction for the dorsal color. Why are these different?
RG and GB:
mean median 2.5th_q 97.5th_q
Tcalifornicum 1.9696 2 0 5
Tchumash 8.1623 8 5 12
Tcristinae 4.0474 4 1 8
Tcuri 1.6515 1 0 4
Tknulli 1.4539 1 0 4
Tlandelsensis 1.7491 2 0 5
Tpetita 1.0797 1 0 3
Tpodura 0.9414 1 0 3
Tpoppensis 1.3528 1 0 4
GB only:
mean median 2.5th_q 97.5th_q
Tcalifornicum 1.0879 1 0 3
Tchumash 4.6752 5 2 7
Tcristinae 1.6236 1 0 4
Tcuri 0.9661 1 0 3
Tknulli 0.6602 0 0 3
Tlandelsensis 1.2354 1 0 4
Tpetita 0.4959 0 0 2
Tpodura 0.2762 0 0 2
Tpoppensis 0.7361 1 0 3
RG only:
mean median 2.5th_q 97.5th_q
Tcalifornicum 0.8824 1 0 3
Tchumash 3.7907 4 1 7
Tcristinae 2.4336 2 0 6
Tcuri 0.6865 0 0 3
Tknulli 0.7941 1 0 3
Tlandelsensis 0.5149 0 0 2
Tpetita 0.5845 0 0 2
Tpodura 0.6659 1 0 2
Tpoppensis 0.6175 0 0 3
######## Lateral color #####################
## read in list of files for RG and GB
files_GB<-read.table("latGB_files.txt",header=FALSE)
files_RG<-read.table("latRG_files.txt",header=FALSE)
### based on indel only
mat_ng<-matrix(NA,nrow=9,ncol=4)
mat_ngx<-matrix(NA,nrow=9,ncol=4)
mat_ngy<-matrix(NA,nrow=9,ncol=4)
## loop over 9 species
for(i in 1:9){
mapGB_all<-read.table(as.character(files_GB[i,]),header=TRUE)
mapRG_all<-read.table(as.character(files_RG[i,]),header=TRUE)
## subset indel SNPs
indel<-which(mapGB_all$CHR==8 & mapGB_all$PS >= 5000000 & mapGB_all$PS <= 6000000 & mapGB_all$scaffold == 128)
mapRG<-mapRG_all[indel,]
mapGB<-mapGB_all[indel,]
L<-length(indel)
## summarize posterior
N<-10000
ng<-rep(NA,N)
ngx<-rep(NA,N)
ngy<-rep(NA,N)
for(j in 1:N){
x<-rbinom(n=L,size=1,prob=mapGB$gamma)
y<-rbinom(n=L,size=1,prob=mapRG$gamma)
ng[j]<-sum(apply(cbind(x,y),1,max))
ngx[j]<-sum(x)
ngy[j]<-sum(y)
}
mat_ng[i,]<-c(mean(ng),quantile(ng,probs=c(0.5,0.025,0.975)))
mat_ngx[i,]<-c(mean(ngx),quantile(ngx,probs=c(0.5,0.025,0.975)))
mat_ngy[i,]<-c(mean(ngy),quantile(ngy,probs=c(0.5,0.025,0.975)))
}
ids<-read.table("ids.txt",header=FALSE)
ids<-c(as.character(ids[1:7,1]),"Tpodura",as.character(ids[8,1]))
rownames(mat_ng)<-ids
colnames(mat_ng)<-c("mean","median","2.5th_q","97.5th_q")
### based on mel--stripy
mat_ng_ms<-matrix(NA,nrow=9,ncol=4)
mat_ngx_ms<-matrix(NA,nrow=9,ncol=4)
mat_ngy_ms<-matrix(NA,nrow=9,ncol=4)
## loop over 9 species
for(i in 1:9){
mapGB_all<-read.table(as.character(files_GB[i,]),header=TRUE)
mapRG_all<-read.table(as.character(files_RG[i,]),header=TRUE)
## subset indel SNPs
mlb<-4139489
mub<-6414835
mel<-c(which(mapGB_all$scaffold==702.1 & mapGB_all$PS < mlb),
which(mapGB_all$scaffold==128 & mapGB_all$PS < mub))
mapRG<-mapRG_all[mel,]
mapGB<-mapGB_all[mel,]
L<-length(mel)
## summarize posterior
N<-10000
ng<-rep(NA,N)
ngx<-rep(NA,N)
ngy<-rep(NA,N)
for(j in 1:N){
x<-rbinom(n=L,size=1,prob=mapGB$gamma)
y<-rbinom(n=L,size=1,prob=mapRG$gamma)
ng[j]<-sum(apply(cbind(x,y),1,max))
ngx[j]<-sum(x)
ngy[j]<-sum(y)
}
mat_ng_ms[i,]<-c(mean(ng),quantile(ng,probs=c(0.5,0.025,0.975)))
mat_ngx_ms[i,]<-c(mean(ngx),quantile(ngx,probs=c(0.5,0.025,0.975)))
mat_ngy_ms[i,]<-c(mean(ngy),quantile(ngy,probs=c(0.5,0.025,0.975)))
}
rownames(mat_ng_ms)<-ids
colnames(mat_ng_ms)<-c("mean","median","2.5th_q","97.5th_q")
rownames(mat_ngx_ms)<-ids
colnames(mat_ngx_ms)<-c("mean","median","2.5th_q","97.5th_q")
rownames(mat_ngy_ms)<-ids
colnames(mat_ngy_ms)<-c("mean","median","2.5th_q","97.5th_q")