library(reshape2)
library(ggplot2)
library(vegan)
set.seed(1)

#### Females ####
load('clean_data/Figure3/Female.RData')

#renaming data..
names = Diet_Female_sample_name;
names$Sample_name = as.character(names$Sample_name)
u_dm = unweighted_unifrac_dm
u_pc1 = unweighted_unifrac_pc1[,-(1:2)]
w_dm = weighted_unifrac_dm
w_pc1 = weighted_unifrac_pc1[,-(1:2)]
rm(Diet_Female_sample_name,
   unweighted_unifrac_pc1,
   unweighted_unifrac_pc2,
   unweighted_unifrac_dm,
   weighted_unifrac_pc1,
   weighted_unifrac_pc2)


#sample name data
X = data.frame('P'=rep(0,107),
               'G'=rep(0,107),
               'W'=rep(0,107),
               'F'=rep(0,107),
               'Diet'=rep('',107))

X$P = substr(sapply(strsplit(names$Sample_name, 'P'), function(x) x[2]), 1, 1)
X$G = substr(sapply(strsplit(names$Sample_name, 'G'), function(x) x[2]), 1, 1)
X$W = substr(sapply(strsplit(names$Sample_name, 'W'), function(x) x[2]), 1, 1)
X$F = (sapply(strsplit(names$Sample_name, 'F'), function(x) x[2]))
X$Diet = names$Deit
X$names = names$Sample_name
X$matchmouse = paste0(X$G,X$F)

#remove mice that only have 1 observation.
ind = names(table(X$matchmouse))[which(table(X$matchmouse) != 2)]
X = X[!X$matchmouse %in% ind, ]

#match the order of mice id
u_dm = u_dm[match(X$names, rownames(u_dm)), match(X$names, colnames(u_dm))]
w_dm = w_dm[match(X$names, rownames(w_dm)), match(X$names, colnames(w_dm))]

#exploratory pca plot
u_svd = svd(u_dm)
plot(u_svd$v[,1] ~ u_svd$v[,2], col = X$Diet)
w_svd = svd(w_dm)
plot(w_svd$v[,1] ~ w_svd$v[,2], col = X$Diet)


##adonis## : total # of perm : 2^15
res_w = adonis(w_dm ~ Diet, X, strata = X$matchmouse, permutations = 1e4)
res_u = adonis(u_dm ~ Diet, X, strata = X$matchmouse, permutations = 1e4)