### this a new script for the amovas and such 
### code adjusted for dealing with the new matrices 
rm(list=ls())

### packages - don't need all of them 
require(lattice)
require(ggplot2)
require(moments)
require(vegan)
require(reshape2)
require(phangorn)
library(plyr)
library(dplyr)
library(tidyr)
library(magrittr)

library(broom)
require(cowplot)
library(seqinr)
library(ape)

library(devtools)
#install_github('leffj/mctoolsr')
library(mctoolsr)

#may need in addition from bioconductor 
#source("http://bioconductor.org/biocLite.R")
#biocLite("Biostrings")
library(Biostrings)
library (poppr)
library (ade4)

setwd("~/Dropbox/Thally long fluc paper")
absnt<-read.csv("all_variants_fixed-only_using_5x_depth_threshold.csv", sep="\t", header=TRUE)
head(absnt)
str(absnt)
colnames(absnt)
#keep only the ones that we want to keep (chromosome info and treatment?, also just the estimates for now)
absnt_est<-absnt[c(6,7,11,15,19,23,27,31,35,39,43,47,51,55,59,63,67,71,75,79,83,87,91,95,99,103,107 ,111,115,119,123,127)]
head(absnt_est)
colnames(absnt_est)
absnt_est<-absnt_est[-1]
#absnt_est<-absnt_est[]*100 # in this step make absolute values 
str(absnt_est)
structures3<-data.frame(c("Anc","22","22","22","22","22","22","26","26","26","26","26","26","32","32","32","32","32","32","FL","FL","FL","FL","FL","FL","FS","FS","FS","FS","FS","FS"))
colnames(structures3)<-"Selection.Regime"
samples3<-absnt_est
head(samples3)
str(samples3)

distances3 <- dist(samples3) #these are the distances - we calculate them from the orginal data frame? 
amovatest3 <- amova(samples3,sqrt(distances3),structures3) #why does this take so long to run? 
amovatest3
#amovasig3<-randtest(amovatest3, nrepet=10)
#amovasig3
#plot(amovasig3)


####now just the FL ones  - but split
colnames(absnt)
absnt_estFL<-absnt[c(6,83,87,91,95,99,103)]
head(absnt_estFL)
colnames(absnt_estFL)
absnt_estFL<-absnt_estFL[-1]
#absnt_est<-absnt_est[]*100 # in this step make absolute values 
str(absnt_estFL)
structures4<-data.frame(c("FLb1","FLb1","FLb2","FLb2","FLb1","FLb2"))
colnames(structures4)<-"Selection.Regime"
samples4<-absnt_estFL
head(samples4)
str(samples4)

distances4 <- dist(samples4) 
amovatest4 <- amova(samples4,sqrt(distances4),structures4) 
amovatest4

### now all again, but this time with the FLs split up 

absnt<-read.csv("all_variants_fixed-only_using_5x_depth_threshold.csv", sep="\t", header=TRUE)
head(absnt)
str(absnt)
colnames(absnt)
#keep only the ones that we want to keep (chromosome info and treatment?, also just the estimates for now)
absnt_est<-absnt[c(6,7,11,15,19,23,27,31,35,39,43,47,51,55,59,63,67,71,75,79,83,87,91,95,99,103,107 ,111,115,119,123,127)]
head(absnt_est)
colnames(absnt_est)
absnt_est<-absnt_est[-1]
#absnt_est<-absnt_est[]*100 # in this step make absolute values 
str(absnt_est)
structures3<-data.frame(c("Anc","22","22","22","22","22","22","26","26","26","26","26","26","32","32","32","32","32","32","FL1","FL1","FL2","FL2","FL1","FL2","FS","FS","FS","FS","FS","FS"))
colnames(structures3)<-"Selection.Regime"
samples3<-absnt_est
head(samples3)
str(samples3)

distances3 <- dist(samples3) #these are the distances - we calculate them from the orginal data frame? 
amovatest3 <- amova(samples3,sqrt(distances3),structures3) #why does this take so long to run? 
amovatest3

#Aeutsignif   <- randtest(amovatest3, nrepet = 999)
