---
title: 'D. galeata & Fish Kairomones -part II: adaptive potential-'
author: "Verena Tams"
date: "November 2017"
output: html_document
---

In this second part of the analysis I test for the adaptive potential (response ~ treatment * pop + (1|pop/clone)). I use a linear mixed effect model and test wether a random intercept or random slope model is appropriate.

I define a significant interaction of presence of predator and genotype origin ("TreatmentxPopulation") as an indication for the possible presence of local adaptation for a predator environment.

For hypothesis testing (significant interaction between treatment and population) I use the function Anova(model, type=2). This test performs a Wald Chi-Square-Test with the underlying null hypothesis that there is no significant interaction between the two explanatory variables.

With a posthoc test I investigate which of the contrast is statistical significant correcting by bonferroni adjustments. I used the function lsmeans(model, pairwise ~ 'significant term', adjust="bon").

```{r SetWd, echo=FALSE, cache=FALSE, results='hide', message=FALSE}
#setwd()
```

```{r LoadData, echo=FALSE, cache=FALSE, results='hide', message=FALSE, warning=FALSE}
#load data
FKn15<- read.table("./input files/FKmaster.txt", header = T)
FKn15.no0<- FKn15[ which(FKn15$AFR > 5), ]
rates<- read.table("./input files/surv_repro_relfit.txt", header = T)
FKn15.no0$Indiv <- factor(1:nrow(FKn15.no0)) 
#load libraries (package)
library(lme4)
library(ggplot2)
library(RColorBrewer)
library(ordinal)
library(MASS)
library(car)
library(sjPlot)
library(lsmeans)
library(multcompView)
library(RVAideMemoire)
#library(lmerTest)
library(effects)

# Function to extract variances from mer objects
GetVariances <- function(mod) {
  summ <- summary(mod)
  if(exists("varcor", summ)) {
    res <- unlist(lapply(summ$varcor, function(lst) as.vector(lst)))
  } else {
    if(exists("ST", summ)) {
      res <- unlist(lapply(summ$ST, function(lst) as.vector(lst)))^2
    } else {
      stop("summary(mod) must have a varcor or ST item")
    }
  }
  if("glmerMod"%in%class(mod) | "lmerMod"%in%class(mod)) res <- c(res, residual=sigma(mod)^2)
  res
}
```

###Age at First Reproduction (AFR)
```{r ModelAFR, echo=FALSE, cache=FALSE, results='hide', fig.width=6, fig.height=4, message=FALSE, warning=FALSE}
#random intercept model, single fixed factor
AFR.model.gammap <- glmer(AFR ~ treatment * pop + (1|pop/clone), family=Gamma(link="log"), data=FKn15.no0)
AFR.null.gammap <- glmer(AFR ~ treatment + (1|pop/clone), family=Gamma(link="log"), data=FKn15.no0)
anova(AFR.null.gammap, AFR.model.gammap)
coef(AFR.model.gammap)#small intercept variation
#random slope & intercept model
AFR.model.gammap2 <- glmer(AFR ~ treatment *pop + (1+1|pop/clone), family=Gamma(link="log"), data=FKn15.no0)
anova(AFR.model.gammap, AFR.model.gammap2) #same, no difference
coef(AFR.model.gammap2)#small slope variation

#output ##################################################
summary(AFR.model.gammap) 
#sjt.glmer(AFR.model.gammap)
#1.Total variance
TotVar.AFR <- sd(FKn15.no0$AFR, na.rm = TRUE)^2 
#2.significant fixed effect?
Anova(AFR.model.gammap, type = 2)
#3.coefficients
AFR.model.gammap.coefs <- summary(AFR.model.gammap)$coefficients
fixef(AFR.model.gammap)
#4.effect size
AFR.model.gammap.Effs <- c(treatment=AFR.model.gammap.coefs["treatmentfish", "Estimate"], StdErr=AFR.model.gammap.coefs["treatmentfish", "Std. Error"], 
 GetVariances(AFR.model.gammap))
#5. pop effect % of total random variation
round(100*AFR.model.gammap.Effs["pop"]/TotVar.AFR,2) 
#6.clone effect % of total random variation
round(100*AFR.model.gammap.Effs["clone:pop"]/TotVar.AFR,2)  
#7.plot effects
plot(allEffects(AFR.model.gammap))
#8.posthoc lsmeans
#lsmeans(AFR.model.gammap, pairwise ~ treatment, adjust="bon")
lsmeans(AFR.model.gammap, pairwise ~ treatment|pop, adjust="bon")
lsmeans(AFR.model.gammap, pairwise ~ treatment*pop, adjust="bon")
```

####Biological Meaning of AFR
The population effect is estimated to be zero. 
The clone effect is small to no effect (`r round(100*AFR.model.gammap.Effs["clone:pop"]/TotVar.AFR,2)` % of total random variation.

Looking at the treatment effect, we see it is estimated as `r round(AFR.model.gammap.Effs["treatment"],2)` with a standard error of `r round(AFR.model.gammap.Effs["StdErr"],2)`.
The treatment effect is statistically significant (Chisq =34.8918, Df=1, Pr(>Chisq)= 3.485e-09), as well as the interaction of TreatmentxPopulation (Chisq =19.9237, Df=3, Pr(>Chisq)= 0.000176), but not for population (Chisq =2.3979, Df=3, Pr(>Chisq)= 0.494031). The Age of First Reproduction (AFR) is significantly different between treatments for popJ (p=<.0001) and popG (p=0.0023) and between popG and popJ in control treatment (p=0.0347).

###Total Number of Broods (broods)

```{r ModelBroods, echo=FALSE, cache=FALSE, results='hide', message=FALSE, fig.cap="Predicted Number of Broods in Treatments and Controls: Daphnia in the treatment are more likely to have fewer broods.", fig.width=6, fig.height=4, warning=FALSE, fig.show='hide'}
# fit the model
FKn15.no0$broodsF <- factor(FKn15.no0$broods, levels = sort(unique(FKn15.no0$broods)), ordered=TRUE)

#random intercept model, single fixed factor
broods.model.clmmp <- clmm(broodsF ~ treatment * pop + (1|pop/clone), data=FKn15.no0, link="probit")
broods.null.clmmp <- clmm(broodsF ~ treatment + (1 |pop/clone), data=FKn15.no0, link="probit")
anova(broods.null.clmmp, broods.model.clmmp) 
coef(broods.model.clmmp)#small intercept variation

#random slope & intercept model
broods.model.clmmp2 <- clmm(broodsF ~ treatment*pop + (1 + 1|pop/clone), data=FKn15.no0, link="probit")

anova(broods.model.clmmp, broods.model.clmmp2) #same, no difference
coef(broods.model.clmmp2)#small slope variation

#output ##################################################
summary(broods.model.clmmp) 
#1.Total variance
TotVar.broods<- sd(FKn15.no0$broods, na.rm = TRUE)^2 
#2.significant fixed effect?
Anova(broods.model.clmmp, type = 2)
#3.coefficients
broods.model.clmmp.coefs <- summary(broods.model.clmmp)$coefficients
#4.effect size
broods.model.clmmp.Effs <- c(treatment=broods.model.clmmp.coefs["treatmentfish", "Estimate"], StdErr=broods.model.clmmp.coefs["treatmentfish", "Std. Error"], GetVariances(broods.model.clmmp), residual=1)
#5.pop effect % of total random variation
round(100*broods.model.clmmp.Effs["pop"]/TotVar.broods,2) 
#6.clone effect % of total random variation
round(100*broods.model.clmmp.Effs["clone:pop"]/TotVar.broods,2) 
#7. plot effects
#plot(allEffects(broods.model.clmm))#error
#8.posthoc lsmeans
#lsmeans(broods.model.clmmp, pairwise ~ treatment, adjust="bon")
lsmeans(broods.model.clmmp, pairwise ~ treatment|pop, adjust="bon")
lsmeans(broods.model.clmmp, pairwise ~ treatment*pop, adjust="bon")
```

####Biological Meaning of broods
The population effect  is estimated to be zero. 
The clone effect is small to no effect (`r round(100*broods.model.clmmp.Effs["clone:pop"]/TotVar.broods,2)` % of total random variation. Looking at the treatment effect, we see it is estimated as `r round(broods.model.clmmp.Effs ["treatment"],2)` with a standard error of `r round(broods.model.clmmp.Effs ["StdErr"],2)`. 
The treatment effect is statistically significant (Chisq =10.8998, Df=1, Pr(>Chisq)= 0.0009617), while the interaction of TreatmentxPopulation is not statistically significant (Chisq =3.9963, Df=3, Pr(>Chisq)= 0.2618684) as well as for population (Chisq =5.8112, Df=3, Pr(>Chisq)= 0.1211649).
Total number of broods differs significantly between popG and popJ in fish exposed water (p=0.0448) and between tretaments within popJ (p=0.0021) and popG (p=0.0492) resulting in a signifcant difference between popG-control and popJ-fish (p=0.0031).

###Total Number of Offspring (offspring)
```{r ModelOff, echo=FALSE, cache=FALSE, results='hide', fig.width=6, fig.height=4, message=FALSE, warning=FALSE}
#random intercept model, single fixed factor
offspring.model.glmmp <- glmer(offspring ~ treatment * pop + (1|pop/clone), data=FKn15.no0, family =poisson())
offspring.null.glmmp <- glmer(offspring ~ treatment + (1 |pop/clone), data=FKn15.no0, family = poisson())
anova(offspring.null.glmmp,offspring.model.glmmp) #
coef(offspring.null.glmmp)#small intercept variation

#random slope & intercept model
offspring.model.glmmp2 <- glmer(offspring ~ treatment * pop + (1 +1|pop/clone), data=FKn15.no0, family = poisson())
anova(offspring.model.glmmp, offspring.model.glmmp2) #same, no difference
coef(offspring.model.glmmp2)# small slope variation

#output ##################################################
summary(offspring.model.glmmp) 
#1.Total variance
TotVar.offspring <- sd(FKn15.no0$offspring, na.rm = TRUE)^2 
#2.significant fixed effect?
Anova(offspring.model.glmmp, type = 2)
#3.coefficients
offspring.model.glmmp.coefs <- summary(offspring.model.glmmp)$coefficients
fixef(offspring.model.glmmp)
#4.effect size
offspring.model.glmmp.Effs <- c(treatment=fixef(offspring.model.glmmp)["treatmentfish"],  StdErr=summary(offspring.model.glmmp)$coefficients["treatmentfish", "Std. Error"], GetVariances(offspring.model.glmmp), residual=1)
offspring.model.glmmp.Effs["residual"] <- offspring.model.glmmp.Effs["Indiv:(clone:pop)"]
offspring.model.glmmp.Effs <- offspring.model.glmmp.Effs[names(offspring.model.glmmp.Effs)!="Indiv:(clone:pop)"]
#5.pop effect % of total random variation
round(100*offspring.model.glmmp.Effs["pop"]/TotVar.offspring,2) 
#6.clone effect % of total random variation
round(100*offspring.model.glmmp.Effs["clone:pop"]/TotVar.offspring,2) 
#7. plot effects
plot(allEffects(offspring.model.glmmp))
#8.posthoc lsmeans
#lsmeans(offspring.model.glmmp, pairwise ~ pop, adjust="bon")
lsmeans(offspring.model.glmmp, pairwise ~ treatment|pop, adjust="bon")
lsmeans(offspring.model.glmmp, pairwise ~ treatment*pop, adjust="bon")
```

####Biological Meaning of offspring
The population and clone effect is estimated to be 0% of total random variation. 
Looking at the treatment effect, we see it is estimated as `r round(offspring.model.glmmp.Effs["treatment.treatmentfish"],2)` with a standard error of `r round(offspring.model.glmmp.Effs["StdErr"],2)`. 
The TreatmentxPopulation is statistical significant (Chisq =9.3807, Df=3, Pr(>Chisq)p=0.024653), as well as population (Chisq =14.0812, Df=3, Pr(>Chisq)= 0.002797) and treatment (Chisq =5.9065, Df=1, Pr(>Chisq)p=0.015085).
The total number of offspring differs significantly between popG and popJ (p=0.0012) in general and in detail in control(p=0.0198) and in fish (p=0.0023) resulting in a signifcant difference between popG-fish and popJ-control (p=0.0311). In addition there is a significant differences for popJ between treatments (p=0.0243).

###Total Number of Offspring of First Brood (brood1)
```{r ModelBrood1, echo=FALSE, cache=FALSE, results='hide', fig.width=6, fig.height=4, message=FALSE, warning=FALSE}
#random intercept model, single fixed factor
brood1.model.glmmp <- glmer(brood1 ~ treatment *pop + (1 |pop/clone), data=FKn15.no0, family = poisson())
brood1.null.glmmp <- glmer(brood1 ~ treatment + (1 |pop/clone), data=FKn15.no0, family = poisson())
anova(brood1.null.glmmp, brood1.model.glmmp) 
coef(brood1.null.glmmp)#small intercept variation

#random slope & intercept model
brood1.model.glmmp2 <- glmer(brood1 ~ treatment*pop + (1+1|pop/clone), family=poisson(), data=FKn15.no0)
anova(brood1.model.glmmp, brood1.model.glmmp2) #same, no difference
coef(brood1.model.glmmp2)#small slope variation

#output ##################################################
summary(brood1.model.glmmp) 
#1.Total variance
TotVar.brood1 <- sd(FKn15.no0$brood1, na.rm = TRUE)^2
#2.significant fixed effect?
Anova(brood1.model.glmmp, type = 2)
#3.coefficients
brood1.model.glmmp.coefs <- summary(brood1.model.glmmp)$coefficients
fixef(brood1.model.glmmp)
#4.effect size
brood1.model.glmmp.Effs <- c(treatment=fixef(brood1.model.glmmp)["treatmentfish"],  StdErr=summary(brood1.model.glmmp)$coefficients["treatmentfish", "Std. Error"], GetVariances(brood1.model.glmmp))
#5.pop effect % of total random variation
round(100*brood1.model.glmmp.Effs["pop"]/TotVar.brood1,2) 
#6.clone effect % of total random variation
round(100*brood1.model.glmmp.Effs["clone:pop"]/TotVar.brood1,2) 
#7. plot effects
plot(allEffects(brood1.model.glmmp))
#8.posthoc lsmeans
#lsmeans(brood1.model.glmmp, pairwise ~ pop, adjust="bon")
lsmeans(brood1.model.glmmp, pairwise ~ treatment|pop, adjust="bon")
lsmeans(brood1.model.glmmp, pairwise ~ treatment*pop, adjust="bon")
```

####Biological Meaning of brood1
The population effect and the clone effect is estimated to be 0% of total random variation. 
Looking at the treatment effect, we see it is estimated as `r round(brood1.model.glmmp.Effs ["treatment.treatmentfish"],2)` with a standard error of `r round(brood1.model.glmmp.Effs ["StdErr"],2)`. 
The interaction of TreatmentxPopulation is not statistical significant (Chisq =4.1949, Df=3, Pr(>Chisq)= 0.241176), as well as for treatment (Chisq =0.0287, Df=1, Pr(>Chisq)= 0.865369), while population is statistical significant (Chisq =11.6722, Df=3, Pr(>Chisq)= 0.008595).
The total number of offspring in first brood differs significantly between popG and popJ (p=0.0041) in general and in detail in control(p=0.0073) and in fish (p=0.0301)  resulting in a signifcant difference between popG-control and popJ-control (p

###Survival (survived_t14)
```{r ModelSurv, echo=FALSE, cache=FALSE, results='hide', message=FALSE, fig.width=6, fig.height=4, warning=FALSE}
#random intercept model, single fixed factor
surv.model.glmmp <- glmer(survived_t14 ~ treatment * pop + (1|pop/clone), data=FKn15, family = binomial(link = "logit"))
surv.null.glmmp <- glmer(survived_t14 ~ treatment + (1 |pop/clone), data=FKn15, family = binomial(link = "logit"))
anova(surv.null.glmmp,surv.model.glmmp) #no significance
coef(surv.model.glmmp)#small intercept variation

#random slope & intercept model
surv.model.glmmp2 <-glmer(survived_t14 ~ treatment *pop + (1 +1|pop/clone), data=FKn15, family = binomial(link = "logit"))
anova(surv.model.glmmp, surv.model.glmmp2) #same, no difference
coef(surv.model.glmmp2)#small slope variation

#output ##################################################
summary(surv.model.glmmp) 
#1.Total variance
TotVar.surv <- sd(FKn15.no0$survived_t14, na.rm = TRUE)^2
#2.significant fixed effect?
Anova(surv.model.glmmp, type = 2)
#3.coefficients
surv.model.glmmp.coefs <- summary(surv.model.glmmp)$coefficients
fixef(surv.model.glmmp)
#4.effect size
surv.model.glmmp.Effs <- c(treatment=fixef(surv.model.glmmp)["treatmentfish"], 
StdErr=summary(surv.model.glmmp)$coefficients["treatmentfish", "Std. Error"], GetVariances(surv.model.glmmp))
#5.pop effect % of total random variation
round(100*surv.model.glmmp.Effs["pop"]/TotVar.surv,2)
#6.clone effect % of total random variation
round(100*surv.model.glmmp.Effs["clone:pop"]/TotVar.surv,2) 
#7. plot effects
plot(allEffects(surv.model.glmmp))
```
####Biological Meaning of Survival
The population effect is estimated to be zero. The clone effect is estimated to be small (`r round(100*surv.model.glmmp.Effs["clone:pop"]/TotVar.surv,2)` % of total random variation).
Looking at the treatment effect, we see it is estimated as `r round(surv.model.glmmp.Effs ["treatment.treatmentfish"],2)` with a standard error of `r round(surv.model.glmmp.Effs ["StdErr"],2)`. 
There is no statistical significant main or interaction effect for treatment (Chisq=0.1582, Df=1, Pr(>Chisq)=0.6908), population (Chisq=1.0720, Df=3, Pr(>Chisq)=0.7838) or TreatmentxPopulation (Chisq=0.0584, Df=3, Pr(>Chisq)=0.69963).

###relative fitness within each population (relfit.nested)
```{r ModelRel.nest, echo=FALSE, cache=FALSE, results='hide', message=FALSE, warning=FALSE, fig.width=6, fig.height=4}
#random intercept model, single fixed factor
relnest.model.glmmp <- glmer(relfit.nested ~ treatment * pop + (1|pop/clone), data=rates, family = Gamma(link="log"))
relnest.null.glmmp <- glmer(relfit.nested ~ treatment+ (1 |pop/clone), data=rates, family = Gamma(link="log"))
anova(relnest.null.glmmp, relnest.model.glmmp)
coef(relnest.model.glmmp)#no intercept variation

#random slope & intercept model
relnest.model.glmmp2 <- glmer(relfit.nested ~ treatment *pop + (1+1 |pop/clone), data=rates, family = Gamma(link="log"))
anova(relnest.model.glmmp, relnest.model.glmmp2) #same, no difference
coef(relnest.model.glmmp2)#small slope variation

#output#############################################################################
summary(relnest.model.glmmp)
#1.Total variance
TotVar.relnest<- sd(rates$relfit.nested, na.rm = TRUE)^2
#2.significant fixed effect?
Anova(relnest.model.glmmp, type = 2)  
#3.coefficients
relnest.model.glmmp.coefs <- summary(relnest.model.glmmp)$coefficients
fixef(relnest.model.glmmp)
#4.effect size
relnest.model.glmmp.Effs <- c(treatment=fixef(relnest.model.glmmp)["treatmentfish"], StdErr=summary(relnest.model.glmmp)$coefficients["treatmentfish", "Std. Error"], GetVariances(relnest.model.glmmp))
#5.pop effect % of total random variation
round(100*relnest.model.glmmp.Effs["pop"]/TotVar.relnest,2) 
#6.clone effect % of total random variation
round(100*relnest.model.glmmp.Effs["clone:pop"]/TotVar.relnest,2) 
#7. plot effects
plot(allEffects(relnest.model.glmmp))
```
####Biological Meaning of relfit.nested
The clone effect within one population is 0% for each population.
Looking at the treatment effect, we see it is estimated as `r round(relnest.model.glmmp.Effs ["treatment.treatmentfish"],2)` with a standard error of `r round(relnest.model.glmmp.Effs ["StdErr"],2)`. 
There is no significant treatment effect (Chisq= 2.1050 , Df=1, Pr(<Chisq)=0.1468), while population (Chisq= 2.5030, Df=3, Pr(<Chisq)=0.4747) and TreatmentxPopulation (Chisq= 1.9149, Df=3, Pr(<Chisq)=0.5903) are not statistically significant.

###relative fitness among all populations (relfit.clone)
```{r ModelRel.clone, echo=FALSE, cache=FALSE, results='hide', message=FALSE, warning=FALSE, fig.width=6, fig.height=4}
#random intercept model, single fixed factor
relclone.model.glmmp <- glmer(relfit.clone ~ treatment * pop + (1|pop/clone), data=rates, family = Gamma(link="log"))
relclone.null.glmmp <- glmer(relfit.clone ~ treatment + (1 |pop/clone), data=rates, family = Gamma(link="log"))
anova(relclone.null.glmmp, relclone.model.glmmp)
coef(relclone.model.glmmp)#no intercept variation

#random slope & intercept model
relclone.model.glmmp2 <- glmer(relfit.clone ~ treatment *pop + (1+1 |pop/clone), data=rates, family = Gamma(link="log"))
anova(relclone.model.glmmp, relclone.model.glmmp2) #same, no difference
coef(relclone.model.glmmp2)#small slope variation

#output#############################################################################
summary(relclone.model.glmmp) 
#1.Total variance
TotVar.relclone<- sd(rates$relfit.clone, na.rm = TRUE)^2
#2.significant fixed effect?
Anova(relclone.model.glmmp, type = 2)
#3.coefficients
relclone.model.glmmp.coefs <- summary(relclone.model.glmmp)$coefficients
fixef(relclone.model.glmmp)
#4.effect size
relclone.model.glmmp.Effs <- c(treatment=fixef(relclone.model.glmmp)["treatmentfish"], 
 StdErr=summary(relclone.model.glmmp)$coefficients["treatmentfish", "Std. Error"], GetVariances(relclone.model.glmmp))
#5.pop effect % of total random variation
round(100*relclone.model.glmmp.Effs["pop"]/TotVar.relclone,2) 
#6.clone effect % of total random variation
round(100*relclone.model.glmmp.Effs["clone:pop"]/TotVar.relclone,2) 
#7. plot effects
plot(allEffects(relclone.model.glmmp))
```
####Biological Meaning of relfit.nested
The population effect and the clone effect is 100% since I compare all clones among all populations.
Looking at the treatment effect, we see it is estimated as `r round(relclone.model.glmmp.Effs ["treatment.treatmentfish"],2)` with a standard error of `r round(relclone.model.glmmp.Effs ["StdErr"],2)`. 
There is a significant treatment effect (Chisq= 24.0544, Df=1, Pr(<Chisq)=9.365e-07), while population (Chisq= 6.3916, Df=3, Pr(<Chisq)=0.09404) and TreatmentxPopulation (Chisq= 1.3568 Df=3, Pr(<Chisq)=0.71569) are not statistically significant.

###Somatic Growth rate (SGR)
```{r ModelSGR, echo=FALSE, cache=FALSE, results='hide', fig.width=6, fig.height=4, message=FALSE, warning=FALSE}
#problems wih NA, introduce new function to na.rm=T
summarySE <- function(data=NULL, measurevar, groupvars=NULL, na.rm=FALSE,
                      conf.interval=.95, .drop=TRUE) {
   library(plyr)

    # New version of length which can handle NA's: if na.rm==T, don't count them
    length2 <- function (x, na.rm=FALSE) {
        if (na.rm) sum(!is.na(x))
        else       length(x)
    }

    # This does the summary. For each group's data frame, return a vector with
    # N, mean, and sd
    datac <- ddply(data, groupvars, .drop=.drop,
      .fun = function(xx, col) {
        c(N    = length2(xx[[col]], na.rm=na.rm),
          mean = mean   (xx[[col]], na.rm=na.rm),
          sd   = sd     (xx[[col]], na.rm=na.rm)
        )
      },
      measurevar
    )

    # Rename the "mean" column    
    datac <- rename(datac, c("mean" = measurevar))

    datac$se <- datac$sd / sqrt(datac$N)  # Calculate standard error of the mean

    # Confidence interval multiplier for standard error
    # Calculate t-statistic for confidence interval: 
    # e.g., if conf.interval is .95, use .975 (above/below), and use df=N-1
    ciMult <- qt(conf.interval/2 + .5, datac$N-1)
    datac$ci <- datac$se * ciMult

    return(datac)
}

#shows me: n=, response, sd, se, ci
summarySE.SGR <- summarySE(FKn15, measurevar="SGR", groupvars=c("treatment", "clone", "pop"), na.rm=T)
##############################################################################
#gaussian
#random intercept model, single fixed factor
#SGR.model.glmm <- glmer(SGR ~ treatment + (1 |pop/clone), data=FKn15, family = gaussian(link="identity")) #same as 
SGR.model.glmmp <- lmer(SGR ~ treatment * pop + (1|pop/clone), data=FKn15)
SGR.null.glmmp <- lmer(SGR ~ treatment + (1 |pop/clone), data=FKn15)
anova(SGR.null.glmmp, SGR.model.glmmp)
coef(SGR.model.glmmp)#no intercept variation

#random slope & intercept model
SGR.model.glmmp2 <- lmer(SGR ~ treatment*pop + (1+1 |pop/clone), data=FKn15)
anova(SGR.model.glmmp, SGR.model.glmmp2) #same, no difference
coef(SGR.model.glmmp2)#small slope variation

#output#############################################################################
summary(SGR.model.glmmp)
#1.Total variance
TotVar.SGR<- sd(FKn15$SGR, na.rm = TRUE)^2
#2.significant fixed effect?
Anova(SGR.model.glmmp, type = 2)
#3.coefficients
SGR.model.glmmp.coefs <- summary(SGR.model.glmmp)$coefficients
fixef(SGR.model.glmmp)
#4.effect size
SGR.model.glmmp.Effs <- c(treatment=fixef(SGR.model.glmmp)["treatmentfish"], 
StdErr=summary(SGR.model.glmmp)$coefficients["treatmentfish", "Std. Error"], GetVariances(SGR.model.glmmp), na.rm=T)
#5.pop effect % of total random variation
round(100*SGR.model.glmmp.Effs["pop"]/TotVar.SGR,2) 
#6.clone effect % of total random variation
round(100*SGR.model.glmmp.Effs["clone:pop"]/TotVar.SGR,2)
#7. plot effects
plot(allEffects(SGR.model.glmmp))
#8. posthoc lsmeans
#lsmeans(SGR.model.glmmp, pairwise ~ pop, adjust="bon")
lsmeans(SGR.model.glmmp, pairwise ~ treatment|pop, adjust="bon")
lsmeans(SGR.model.glmmp, pairwise ~ treatment*pop, adjust="bon")
```

####Biological Meaning of SGR
The population effect is high (`r round(100*SGR.model.glmmp.Effs["pop"]/TotVar.SGR,2)`% of total random variation) and the clone effect is high (`r round(100*SGR.model.glmmp.Effs["clone:pop"]/TotVar.SGR,2)` % of total random variation).

Looking at the treatment effect, we see it is estimated as `r round(SGR.model.glmmp.Effs ["treatment.treatmentfish"],2)` with a standard error of `r round(SGR.model.glmmp.Effs ["StdErr"],2)`. 
There is a significant treatment effect (Chisq= 16.5043, Df=1, Pr(<Chisq)=4.854e-05) and interaction TreatmentxPopulation (Chisq= 21.8950, Df=3, Pr(<Chisq)=6.860e-05), while population has no sigificant effect on SGR (Chisq= 0.7885, Df=3, Pr(<Chisq)=0.8522).

The Somatic Growth Rate differs significantly for treatments within popG (p=<.0001) and popJ (p=0.0135) in general.

##Conclusion
There was a statistical significant interaction of TreatmentxPopulation for AFR (Chisq= 16.5043, Df=1, Pr(<Chisq)=4.854e-05), total number of offspring (Chisq =9.3807, Df=3, Pr(>Chisq)=0.024653) and SGR (Chisq= 21.8950, Df=3, Pr(<Chisq)=6.860e-05), meaning that for those traits a potential for local adaptation to fish exposure exists.