R Markdown

library(ggplot2)
library(lme4)
## Loading required package: Matrix
library(MuMIn)
library(vegan)
## Loading required package: permute
## Loading required package: lattice
## This is vegan 2.5-7
library(emmeans)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(lmerTest)
## 
## Attaching package: 'lmerTest'
## The following object is masked from 'package:lme4':
## 
##     lmer
## The following object is masked from 'package:stats':
## 
##     step
library(tibble)
library(plyr)
## ------------------------------------------------------------------------------
## You have loaded plyr after dplyr - this is likely to cause problems.
## If you need functions from both plyr and dplyr, please load plyr first, then dplyr:
## library(plyr); library(dplyr)
## ------------------------------------------------------------------------------
## 
## Attaching package: 'plyr'
## The following objects are masked from 'package:dplyr':
## 
##     arrange, count, desc, failwith, id, mutate, rename, summarise,
##     summarize
setwd("C:/Users/jc130714/OneDrive - James Cook University/Documents/Yogi Yasutake/Honours/Post-honours/Statistics and data/Final data")
getwd()
## [1] "C:/Users/jc130714/OneDrive - James Cook University/Documents/Yogi Yasutake/Honours/Post-honours/Statistics and data/Final data"

##PCA building

JuvPCA = read.csv('Juveniles at hatching_PCA.csv', strip.white=TRUE)
head(JuvPCA)
##   Juvenile.weight..mg. Juvenile.SL Juvenile.yolk.area..mm.2. Fulton.s.K
## 1                  3.3       4.476                     1.321      3.680
## 2                  3.5       5.010                     1.139      2.783
## 3                  3.5       5.122                     1.317      2.605
## 4                  3.6       5.137                     1.039      2.655
## 5                  3.7       5.137                     1.206      2.729
## 6                  3.7       5.148                     1.517      2.712
#standardising
Juv.std = stdize(JuvPCA, centre=TRUE, scale = TRUE)
head(Juv.std)
##   z.Juvenile.weight..mg. z.Juvenile.SL z.Juvenile.yolk.area..mm.2. z.Fulton.s.K
## 1             -1.5607845    -2.2627284                  -0.2470842    1.4593526
## 2             -1.1050896    -0.6188274                  -1.0218995   -0.1472788
## 3             -1.1050896    -0.2740392                  -0.2641131   -0.4660976
## 4             -0.8772421    -0.2278622                  -1.4476222   -0.3765417
## 5             -0.6493947    -0.2278622                  -0.7366653   -0.2439991
## 6             -0.6493947    -0.1939991                   0.5873323   -0.2744481
Juv.rda = rda(Juv.std)
summary(Juv.rda, display=NULL)
## 
## Call:
## rda(X = Juv.std) 
## 
## Partitioning of variance:
##               Inertia Proportion
## Total               4          1
## Unconstrained       4          1
## 
## Eigenvalues, and their contribution to the variance 
## 
## Importance of components:
##                          PC1    PC2    PC3      PC4
## Eigenvalue            2.1246 1.1869 0.6775 0.010922
## Proportion Explained  0.5312 0.2967 0.1694 0.002731
## Cumulative Proportion 0.5312 0.8279 0.9973 1.000000
## 
## Scaling 2 for species and site scores
## * Species are scaled proportional to eigenvalues
## * Sites are unscaled: weighted dispersion equal on all dimensions
## * General scaling constant of scores:
screeplot(Juv.rda)

Juv.sites.score = (scores(Juv.rda, display='sites'))
head(Juv.sites.score)
##              PC1        PC2
## sit1  0.53515229 -0.5144726
## sit2  0.16986147 -0.3325342
## sit3  0.00140987 -0.3457572
## sit4  0.12859378 -0.2394690
## sit5  0.07368103 -0.1900736
## sit6 -0.07089583 -0.2280509
Juv.species.score = as.data.frame(scores(Juv.rda, display='species'))
Juv.species.score$Species = rownames(Juv.species.score)
Juv.sites.score = data.frame(Juv.sites.score, JuvPCA)
#getting and adding groupings 
JuvGroup = read.csv('Juv NH groupings PCA.csv')
head(JuvGroup)
##   Tank.ID Treatment Dev.treatment Repro.treatment Male.family Male.SL..mm.
## 1     139      CC-C            CC               C           F        8.648
## 2     139      CC-C            CC               C           F        8.648
## 3     139      CC-C            CC               C           F        8.648
## 4     139      CC-C            CC               C           F        8.648
## 5     139      CC-C            CC               C           F        8.648
## 6     139      CC-C            CC               C           F        8.648
##   Male.weight..g. Male.condition Female.family Female.SL Female.weight
## 1           22.63       3.498948             A     9.144         31.27
## 2           22.63       3.498948             A     9.144         31.27
## 3           22.63       3.498948             A     9.144         31.27
## 4           22.63       3.498948             A     9.144         31.27
## 5           22.63       3.498948             A     9.144         31.27
## 6           22.63       3.498948             A     9.144         31.27
##   Condition...Female Clutch.number Avg.egg.area Fish.replicate
## 1           4.089962             1       4.9532              1
## 2           4.089962             1       4.9532              2
## 3           4.089962             1       4.9532              3
## 4           4.089962             1       4.9532              4
## 5           4.089962             1       4.9532              5
## 6           4.089962             1       4.9532              6
Juv.sites.score2 = cbind(Juv.sites.score, JuvGroup)
head(Juv.sites.score2)
##              PC1        PC2 Juvenile.weight..mg. Juvenile.SL
## sit1  0.53515229 -0.5144726                  3.3       4.476
## sit2  0.16986147 -0.3325342                  3.5       5.010
## sit3  0.00140987 -0.3457572                  3.5       5.122
## sit4  0.12859378 -0.2394690                  3.6       5.137
## sit5  0.07368103 -0.1900736                  3.7       5.137
## sit6 -0.07089583 -0.2280509                  3.7       5.148
##      Juvenile.yolk.area..mm.2. Fulton.s.K Tank.ID Treatment Dev.treatment
## sit1                     1.321      3.680     139      CC-C            CC
## sit2                     1.139      2.783     139      CC-C            CC
## sit3                     1.317      2.605     139      CC-C            CC
## sit4                     1.039      2.655     139      CC-C            CC
## sit5                     1.206      2.729     139      CC-C            CC
## sit6                     1.517      2.712     139      CC-C            CC
##      Repro.treatment Male.family Male.SL..mm. Male.weight..g. Male.condition
## sit1               C           F        8.648           22.63       3.498948
## sit2               C           F        8.648           22.63       3.498948
## sit3               C           F        8.648           22.63       3.498948
## sit4               C           F        8.648           22.63       3.498948
## sit5               C           F        8.648           22.63       3.498948
## sit6               C           F        8.648           22.63       3.498948
##      Female.family Female.SL Female.weight Condition...Female Clutch.number
## sit1             A     9.144         31.27           4.089962             1
## sit2             A     9.144         31.27           4.089962             1
## sit3             A     9.144         31.27           4.089962             1
## sit4             A     9.144         31.27           4.089962             1
## sit5             A     9.144         31.27           4.089962             1
## sit6             A     9.144         31.27           4.089962             1
##      Avg.egg.area Fish.replicate
## sit1       4.9532              1
## sit2       4.9532              2
## sit3       4.9532              3
## sit4       4.9532              4
## sit5       4.9532              5
## sit6       4.9532              6
#making a PCA figure 
#putting the sites on the graph
g = ggplot() + geom_point(data=Juv.sites.score2, aes(y=PC2, x=PC1, shape=Treatment, colour= Treatment),  show.legend = TRUE)   + theme_classic()

#adding the arrows
g = g + geom_segment(data=Juv.species.score, aes(y=0, x=0, yend=PC2, xend=PC1), arrow = arrow(length=unit(0.3, 'lines')), show.legend = FALSE)
hjust = ifelse(Juv.species.score$PC1>0,0,1)
vjust = ifelse(Juv.species.score$PC2>0,0,1)
g = g + geom_text(data=Juv.species.score, aes(y=PC2, x=PC1, label = Species))

#now adding cross hairs and axis lines 
g = g + geom_segment(data=NULL, aes(y=-Inf, x=-0, yend=Inf, xend=0))
g = g + geom_segment(data=NULL, aes(y=0, x=-Inf, yend=0, xend=Inf))
g

PCA_means = tbl_df(Juv.sites.score2)
## Warning: `tbl_df()` was deprecated in dplyr 1.0.0.
## Please use `tibble::as_tibble()` instead.
SLgroup_meanPC1 = ddply(PCA_means, "Female.SL", summarise, grp.mean=mean(PC1))
Wgroup_meanPC1 = ddply(PCA_means, "Female.weight", summarise, grp.mean=mean(PC1))
Condgroup_meanPC1 = ddply(PCA_means, "Condition...Female", summarise, grp.mean=mean(PC1))

SLgroup_meanPC2 = ddply(PCA_means, "Female.SL", summarise, grp.mean=mean(PC2))
Wgroup_meanPC2 = ddply(PCA_means, "Female.weight", summarise, grp.mean=mean(PC2))
Condgroup_meanPC2 = ddply(PCA_means, "Condition...Female", summarise, grp.mean=mean(PC2))
#Linear regression to investigate potential maternal effects
lm.1 = lm(grp.mean ~ Female.SL, data=SLgroup_meanPC1)
summary(lm.1)
## 
## Call:
## lm(formula = grp.mean ~ Female.SL, data = SLgroup_meanPC1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.48555 -0.18458  0.00557  0.15276  0.55865 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.049360   1.223401  -0.040    0.968
## Female.SL    0.004769   0.137501   0.035    0.973
## 
## Residual standard error: 0.2783 on 18 degrees of freedom
## Multiple R-squared:  6.684e-05,  Adjusted R-squared:  -0.05549 
## F-statistic: 0.001203 on 1 and 18 DF,  p-value: 0.9727
lm.2 = lm(grp.mean ~ Female.weight, data=Wgroup_meanPC1)
summary(lm.2)
## 
## Call:
## lm(formula = grp.mean ~ Female.weight, data = Wgroup_meanPC1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.42928 -0.18141 -0.00881  0.15065  0.53550 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)
## (Intercept)    0.257368   0.407890   0.631    0.536
## Female.weight -0.008932   0.013625  -0.656    0.520
## 
## Residual standard error: 0.2751 on 18 degrees of freedom
## Multiple R-squared:  0.02332,    Adjusted R-squared:  -0.03094 
## F-statistic: 0.4298 on 1 and 18 DF,  p-value: 0.5204
lm.3 = lm(grp.mean ~ Condition...Female, data=Condgroup_meanPC1)
summary(lm.3)
## 
## Call:
## lm(formula = grp.mean ~ Condition...Female, data = Condgroup_meanPC1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.45928 -0.19658 -0.02087  0.16694  0.50252 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)
## (Intercept)          0.5778     0.5005   1.155    0.263
## Condition...Female  -0.1388     0.1179  -1.177    0.255
## 
## Residual standard error: 0.2682 on 18 degrees of freedom
## Multiple R-squared:  0.07146,    Adjusted R-squared:  0.01987 
## F-statistic: 1.385 on 1 and 18 DF,  p-value: 0.2545
lm.1 = lm(grp.mean ~ Female.SL, data=SLgroup_meanPC2)
summary(lm.1)
## 
## Call:
## lm(formula = grp.mean ~ Female.SL, data = SLgroup_meanPC2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.34299 -0.17378 -0.00722  0.09401  0.43348 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept)   2.0056     1.0439   1.921   0.0707 .
## Female.SL    -0.2263     0.1173  -1.929   0.0697 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2375 on 18 degrees of freedom
## Multiple R-squared:  0.1713, Adjusted R-squared:  0.1252 
## F-statistic:  3.72 on 1 and 18 DF,  p-value: 0.06969
lm.2 = lm(grp.mean ~ Female.weight, data=Wgroup_meanPC2)
summary(lm.2)
## 
## Call:
## lm(formula = grp.mean ~ Female.weight, data = Wgroup_meanPC2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.43002 -0.15118 -0.06429  0.15202  0.47506 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)  
## (Intercept)    0.67021    0.35173   1.905   0.0728 .
## Female.weight -0.02282    0.01175  -1.942   0.0679 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2372 on 18 degrees of freedom
## Multiple R-squared:  0.1733, Adjusted R-squared:  0.1273 
## F-statistic: 3.773 on 1 and 18 DF,  p-value: 0.06791
lm.3 = lm(grp.mean ~ Condition...Female, data=Condgroup_meanPC2)
summary(lm.3)
## 
## Call:
## lm(formula = grp.mean ~ Condition...Female, data = Condgroup_meanPC2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.44961 -0.18598 -0.03389  0.15255  0.52515 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)
## (Intercept)         0.06049    0.48656   0.124    0.902
## Condition...Female -0.01558    0.11465  -0.136    0.893
## 
## Residual standard error: 0.2607 on 18 degrees of freedom
## Multiple R-squared:  0.001025,   Adjusted R-squared:  -0.05447 
## F-statistic: 0.01847 on 1 and 18 DF,  p-value: 0.8934

#lmer testing PCA 1 and 2

lmer.JPC1a = lmer(PC1~ Dev.treatment * Repro.treatment  + (1|Male.family) + (1|Female.family) + (1|Fish.replicate) , data = Juv.sites.score2)
summary(lmer.JPC1a)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: PC1 ~ Dev.treatment * Repro.treatment + (1 | Male.family) + (1 |  
##     Female.family) + (1 | Fish.replicate)
##    Data: Juv.sites.score2
## 
## REML criterion at convergence: 79.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.9830 -0.7203 -0.0263  0.5483  3.2234 
## 
## Random effects:
##  Groups         Name        Variance Std.Dev.
##  Fish.replicate (Intercept) 0.01238  0.1113  
##  Female.family  (Intercept) 0.03566  0.1888  
##  Male.family    (Intercept) 0.01408  0.1187  
##  Residual                   0.05980  0.2445  
## Number of obs: 387, groups:  
## Fish.replicate, 21; Female.family, 5; Male.family, 5
## 
## Fixed effects:
##                                    Estimate Std. Error         df t value
## (Intercept)                        0.042390   0.107723   7.790611   0.394
## Dev.treatmentHH                   -0.007157   0.042207 348.794691  -0.170
## Repro.treatmentH                  -0.014795   0.046357 324.496807  -0.319
## Dev.treatmentHH:Repro.treatmentH   0.275558   0.067025 336.942982   4.111
##                                  Pr(>|t|)    
## (Intercept)                         0.704    
## Dev.treatmentHH                     0.865    
## Repro.treatmentH                    0.750    
## Dev.treatmentHH:Repro.treatmentH 4.95e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Dv.tHH Rpr.tH
## Dv.trtmntHH -0.148              
## Rpr.trtmntH -0.181  0.644       
## Dv.trHH:R.H  0.114 -0.720 -0.773
anova(lmer.JPC1a)
## Type III Analysis of Variance Table with Satterthwaite's method
##                               Sum Sq Mean Sq NumDF  DenDF F value    Pr(>F)    
## Dev.treatment                 1.1747  1.1747     1 361.79  19.643 1.239e-05 ***
## Repro.treatment               1.0382  1.0382     1 363.72  17.361 3.865e-05 ***
## Dev.treatment:Repro.treatment 1.0108  1.0108     1 336.94  16.903 4.948e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
hist(residuals(lmer.JPC1a), col="darkgray")

plot(fitted(lmer.JPC1a), residuals(lmer.JPC1a))

JPC1EMM = (emmeans(lmer.JPC1a, ~ Dev.treatment*Repro.treatment) %>% as.data.frame) 
JPC1EMM
##   Dev.treatment Repro.treatment     emmean        SE       df    lower.CL
## 1            CC               C 0.04239025 0.1080194 7.918628 -0.20714909
## 2            HH               C 0.03523346 0.1103176 8.432328 -0.21690691
## 3            CC               H 0.02759552 0.1099749 8.284918 -0.22449686
## 4            HH               H 0.29599701 0.1106439 8.683943  0.04430763
##    upper.CL
## 1 0.2919296
## 2 0.2873738
## 3 0.2796879
## 4 0.5476864
ggplot(JPC1EMM, aes(x=Repro.treatment, y=emmean, colour=Dev.treatment))+ geom_point() + geom_linerange(aes (ymin = emmean - SE, ymax = emmean + SE))  +  theme_classic() 

#including Female SL
lmer.JPC2d = lmer(PC2~ Dev.treatment * Repro.treatment  + (1|Male.family) + (1|Female.family) + (1|Fish.replicate ) + Female.SL , data = Juv.sites.score2)

summary(lmer.JPC2d)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: PC2 ~ Dev.treatment * Repro.treatment + (1 | Male.family) + (1 |  
##     Female.family) + (1 | Fish.replicate) + Female.SL
##    Data: Juv.sites.score2
## 
## REML criterion at convergence: -246.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.3996 -0.5718 -0.0644  0.6033  3.7799 
## 
## Random effects:
##  Groups         Name        Variance Std.Dev.
##  Fish.replicate (Intercept) 0.03088  0.1757  
##  Female.family  (Intercept) 0.04267  0.2066  
##  Male.family    (Intercept) 0.04762  0.2182  
##  Residual                   0.02235  0.1495  
## Number of obs: 387, groups:  
## Fish.replicate, 21; Female.family, 5; Male.family, 5
## 
## Fixed effects:
##                                    Estimate Std. Error         df t value
## (Intercept)                        4.532418   0.313146 148.303557  14.474
## Dev.treatmentHH                    0.026228   0.028501 360.970681   0.920
## Repro.treatmentH                   0.004432   0.029266 361.175532   0.151
## Female.SL                         -0.516717   0.030942 360.861122 -16.700
## Dev.treatmentHH:Repro.treatmentH  -0.040386   0.042188 361.117105  -0.957
##                                  Pr(>|t|)    
## (Intercept)                        <2e-16 ***
## Dev.treatmentHH                     0.358    
## Repro.treatmentH                    0.880    
## Female.SL                          <2e-16 ***
## Dev.treatmentHH:Repro.treatmentH    0.339    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Dv.tHH Rpr.tH Fml.SL
## Dv.trtmntHH -0.368                     
## Rpr.trtmntH  0.021  0.582              
## Female.SL   -0.892  0.380 -0.066       
## Dv.trHH:R.H -0.042 -0.642 -0.783  0.073
anova(lmer.JPC2d)
## Type III Analysis of Variance Table with Satterthwaite's method
##                               Sum Sq Mean Sq NumDF  DenDF  F value Pr(>F)    
## Dev.treatment                 0.0017  0.0017     1 358.24   0.0750 0.7843    
## Repro.treatment               0.0166  0.0166     1 359.71   0.7413 0.3898    
## Female.SL                     6.2328  6.2328     1 360.86 278.8768 <2e-16 ***
## Dev.treatment:Repro.treatment 0.0205  0.0205     1 361.12   0.9164 0.3391    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
hist(residuals(lmer.JPC2d), col="darkgray")

plot(fitted(lmer.JPC2d), residuals(lmer.JPC2d))

JPC2EMM = (emmeans(lmer.JPC2d, ~ Dev.treatment*Repro.treatment) %>% as.data.frame) 
JPC2EMM
##   Dev.treatment Repro.treatment      emmean        SE       df   lower.CL
## 1            CC               C -0.07019689 0.1414003 9.394089 -0.3880327
## 2            HH               C -0.04396911 0.1421929 9.593950 -0.3626218
## 3            CC               H -0.06576444 0.1420338 9.551068 -0.3842634
## 4            HH               H -0.07992224 0.1422787 9.630250 -0.3985966
##    upper.CL
## 1 0.2476389
## 2 0.2746836
## 3 0.2527345
## 4 0.2387521
ggplot(JPC2EMM, aes(x=Repro.treatment, y=emmean, colour=Dev.treatment))+ geom_point() + geom_linerange(aes (ymin = emmean - SE, ymax = emmean + SE))  +  theme_classic()

##Juvenile characteristics individually

Juv = read.csv('Juveniles at hatching_FINAL.csv', strip.white=TRUE)
head(Juv)
##   Tank.ID Treatment Dev.treatment Repro.treatment Pair.rep Female.ID Male.ID
## 1     139      CC-C            CC               C        1      ACC4     FCC
## 2     139      CC-C            CC               C        1      ACC4     FCC
## 3     139      CC-C            CC               C        1      ACC4     FCC
## 4     139      CC-C            CC               C        1      ACC4     FCC
## 5     139      CC-C            CC               C        1      ACC4     FCC
## 6     139      CC-C            CC               C        1      ACC4     FCC
##   Male.family Family.Number Male.SL..mm. Male.weight..g. Male.condition
## 1           F             6        8.648           22.63       3.498948
## 2           F             6        8.648           22.63       3.498948
## 3           F             6        8.648           22.63       3.498948
## 4           F             6        8.648           22.63       3.498948
## 5           F             6        8.648           22.63       3.498948
## 6           F             6        8.648           22.63       3.498948
##   Female.family Family.Number.1 Female.SL Female.weight Condition...Female
## 1             A               1     9.144         31.27           4.089962
## 2             A               1     9.144         31.27           4.089962
## 3             A               1     9.144         31.27           4.089962
## 4             A               1     9.144         31.27           4.089962
## 5             A               1     9.144         31.27           4.089962
## 6             A               1     9.144         31.27           4.089962
##   Clutch.number Avg.egg.area Fish.replicate Juvenile.weight..mg. Juvenile.SL
## 1             1       4.9532              1                  3.3       4.476
## 2             1       4.9532              2                  3.5       5.010
## 3             1       4.9532              3                  3.5       5.122
## 4             1       4.9532              4                  3.6       5.137
## 5             1       4.9532              5                  3.7       5.137
## 6             1       4.9532              6                  3.7       5.148
##   Juvenile.yolk.area..mm.2. Fulton.s.K
## 1                     1.321      3.680
## 2                     1.139      2.783
## 3                     1.317      2.605
## 4                     1.039      2.655
## 5                     1.206      2.729
## 6                     1.517      2.712
str(Juv)
## 'data.frame':    387 obs. of  24 variables:
##  $ Tank.ID                  : int  139 139 139 139 139 139 139 139 139 139 ...
##  $ Treatment                : chr  "CC-C" "CC-C" "CC-C" "CC-C" ...
##  $ Dev.treatment            : chr  "CC" "CC" "CC" "CC" ...
##  $ Repro.treatment          : chr  "C" "C" "C" "C" ...
##  $ Pair.rep                 : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ Female.ID                : chr  "ACC4" "ACC4" "ACC4" "ACC4" ...
##  $ Male.ID                  : chr  "FCC" "FCC" "FCC" "FCC" ...
##  $ Male.family              : chr  "F" "F" "F" "F" ...
##  $ Family.Number            : int  6 6 6 6 6 6 6 6 6 6 ...
##  $ Male.SL..mm.             : num  8.65 8.65 8.65 8.65 8.65 ...
##  $ Male.weight..g.          : num  22.6 22.6 22.6 22.6 22.6 ...
##  $ Male.condition           : num  3.5 3.5 3.5 3.5 3.5 ...
##  $ Female.family            : chr  "A" "A" "A" "A" ...
##  $ Family.Number.1          : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ Female.SL                : num  9.14 9.14 9.14 9.14 9.14 ...
##  $ Female.weight            : num  31.3 31.3 31.3 31.3 31.3 ...
##  $ Condition...Female       : num  4.09 4.09 4.09 4.09 4.09 ...
##  $ Clutch.number            : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ Avg.egg.area             : num  4.95 4.95 4.95 4.95 4.95 ...
##  $ Fish.replicate           : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ Juvenile.weight..mg.     : num  3.3 3.5 3.5 3.6 3.7 3.7 3.7 3.8 3.8 3.8 ...
##  $ Juvenile.SL              : num  4.48 5.01 5.12 5.14 5.14 ...
##  $ Juvenile.yolk.area..mm.2.: num  1.32 1.14 1.32 1.04 1.21 ...
##  $ Fulton.s.K               : num  3.68 2.78 2.6 2.65 2.73 ...
Juv$Pair.rep = factor (Juv$Pair.rep)
str(Juv)
## 'data.frame':    387 obs. of  24 variables:
##  $ Tank.ID                  : int  139 139 139 139 139 139 139 139 139 139 ...
##  $ Treatment                : chr  "CC-C" "CC-C" "CC-C" "CC-C" ...
##  $ Dev.treatment            : chr  "CC" "CC" "CC" "CC" ...
##  $ Repro.treatment          : chr  "C" "C" "C" "C" ...
##  $ Pair.rep                 : Factor w/ 6 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Female.ID                : chr  "ACC4" "ACC4" "ACC4" "ACC4" ...
##  $ Male.ID                  : chr  "FCC" "FCC" "FCC" "FCC" ...
##  $ Male.family              : chr  "F" "F" "F" "F" ...
##  $ Family.Number            : int  6 6 6 6 6 6 6 6 6 6 ...
##  $ Male.SL..mm.             : num  8.65 8.65 8.65 8.65 8.65 ...
##  $ Male.weight..g.          : num  22.6 22.6 22.6 22.6 22.6 ...
##  $ Male.condition           : num  3.5 3.5 3.5 3.5 3.5 ...
##  $ Female.family            : chr  "A" "A" "A" "A" ...
##  $ Family.Number.1          : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ Female.SL                : num  9.14 9.14 9.14 9.14 9.14 ...
##  $ Female.weight            : num  31.3 31.3 31.3 31.3 31.3 ...
##  $ Condition...Female       : num  4.09 4.09 4.09 4.09 4.09 ...
##  $ Clutch.number            : int  1 1 1 1 1 1 1 1 1 1 ...
##  $ Avg.egg.area             : num  4.95 4.95 4.95 4.95 4.95 ...
##  $ Fish.replicate           : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ Juvenile.weight..mg.     : num  3.3 3.5 3.5 3.6 3.7 3.7 3.7 3.8 3.8 3.8 ...
##  $ Juvenile.SL              : num  4.48 5.01 5.12 5.14 5.14 ...
##  $ Juvenile.yolk.area..mm.2.: num  1.32 1.14 1.32 1.04 1.21 ...
##  $ Fulton.s.K               : num  3.68 2.78 2.6 2.65 2.73 ...

##Standard length

SLgroup_meanJSL = ddply(Juv, "Female.SL", summarise, grp.mean=mean(Juvenile.SL))
Wgroup_meanJSL = ddply(Juv, "Female.weight", summarise, grp.mean=mean(Juvenile.SL))
Condgroup_meanJSL = ddply(Juv, "Condition...Female", summarise, grp.mean=mean(Juvenile.SL))
lm.1 = lm(grp.mean ~ Female.SL, data=SLgroup_meanJSL)
summary(lm.1)
## 
## Call:
## lm(formula = grp.mean ~ Female.SL, data = SLgroup_meanJSL)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.43293 -0.11326 -0.00078  0.06775  0.50561 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  5.34847    1.10091   4.858 0.000126 ***
## Female.SL   -0.01504    0.12373  -0.122 0.904606    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2505 on 18 degrees of freedom
## Multiple R-squared:  0.0008201,  Adjusted R-squared:  -0.05469 
## F-statistic: 0.01477 on 1 and 18 DF,  p-value: 0.9046
lm.2 = lm(grp.mean ~ Female.weight, data=Wgroup_meanJSL)
summary(lm.2)
## 
## Call:
## lm(formula = grp.mean ~ Female.weight, data = Wgroup_meanJSL)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.49689 -0.12054  0.00414  0.08445  0.47312 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   5.063710   0.369795  13.693 5.86e-11 ***
## Female.weight 0.005106   0.012353   0.413    0.684    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2494 on 18 degrees of freedom
## Multiple R-squared:  0.009405,   Adjusted R-squared:  -0.04563 
## F-statistic: 0.1709 on 1 and 18 DF,  p-value: 0.6842
lm.3 = lm(grp.mean ~ Condition...Female, data=Condgroup_meanJSL)
summary(lm.3)
## 
## Call:
## lm(formula = grp.mean ~ Condition...Female, data = Condgroup_meanJSL)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.43467 -0.08638 -0.00717  0.08737  0.43031 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)    
## (Intercept)         4.80481    0.45732  10.506 4.15e-09 ***
## Condition...Female  0.09732    0.10776   0.903    0.378    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2451 on 18 degrees of freedom
## Multiple R-squared:  0.04335,    Adjusted R-squared:  -0.009802 
## F-statistic: 0.8156 on 1 and 18 DF,  p-value: 0.3784
lmer.JSL = lmer(Juvenile.SL~ Dev.treatment * Repro.treatment  + (1|Male.family) + (1|Female.family) + (1|Fish.replicate) , data = Juv)

summary(lmer.JSL)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Juvenile.SL ~ Dev.treatment * Repro.treatment + (1 | Male.family) +  
##     (1 | Female.family) + (1 | Fish.replicate)
##    Data: Juv
## 
## REML criterion at convergence: -14.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.2098 -0.5852  0.0365  0.7091  2.5497 
## 
## Random effects:
##  Groups         Name        Variance Std.Dev.
##  Fish.replicate (Intercept) 0.03097  0.1760  
##  Female.family  (Intercept) 0.04198  0.2049  
##  Male.family    (Intercept) 0.02004  0.1416  
##  Residual                   0.04375  0.2092  
## Number of obs: 387, groups:  
## Fish.replicate, 21; Female.family, 5; Male.family, 5
## 
## Fixed effects:
##                                   Estimate Std. Error        df t value
## (Intercept)                        5.12130    0.12127   8.87615  42.230
## Dev.treatmentHH                    0.05653    0.03651 362.15069   1.548
## Repro.treatmentH                   0.03467    0.04027 355.53623   0.861
## Dev.treatmentHH:Repro.treatmentH  -0.28144    0.05812 359.87026  -4.842
##                                  Pr(>|t|)    
## (Intercept)                      1.53e-11 ***
## Dev.treatmentHH                     0.122    
## Repro.treatmentH                    0.390    
## Dev.treatmentHH:Repro.treatmentH 1.91e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Dv.tHH Rpr.tH
## Dv.trtmntHH -0.112              
## Rpr.trtmntH -0.138  0.651       
## Dv.trHH:R.H  0.085 -0.723 -0.777
lmer.JSLf = lmer(Juvenile.SL~ Dev.treatment * Repro.treatment + Condition...Female + (1|Male.family) + (1|Female.family) + (1|Fish.replicate), data = Juv)

summary(lmer.JSLf)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Juvenile.SL ~ Dev.treatment * Repro.treatment + Condition...Female +  
##     (1 | Male.family) + (1 | Female.family) + (1 | Fish.replicate)
##    Data: Juv
## 
## REML criterion at convergence: -22.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.2189 -0.5538  0.0841  0.6711  2.6321 
## 
## Random effects:
##  Groups         Name        Variance Std.Dev.
##  Fish.replicate (Intercept) 0.03185  0.1785  
##  Female.family  (Intercept) 0.04417  0.2102  
##  Male.family    (Intercept) 0.01973  0.1405  
##  Residual                   0.04226  0.2056  
## Number of obs: 387, groups:  
## Fish.replicate, 21; Female.family, 5; Male.family, 5
## 
## Fixed effects:
##                                   Estimate Std. Error        df t value
## (Intercept)                        4.62814    0.18390  40.00319  25.167
## Dev.treatmentHH                    0.09771    0.03767 360.34135   2.594
## Repro.treatmentH                   0.07935    0.04148 353.55095   1.913
## Condition...Female                 0.11288    0.03135 360.40027   3.601
## Dev.treatmentHH:Repro.treatmentH  -0.29524    0.05729 359.58351  -5.153
##                                  Pr(>|t|)    
## (Intercept)                       < 2e-16 ***
## Dev.treatmentHH                  0.009870 ** 
## Repro.treatmentH                 0.056531 .  
## Condition...Female               0.000362 ***
## Dev.treatmentHH:Repro.treatmentH 4.23e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Dv.tHH Rpr.tH Cn...F
## Dv.trtmntHH -0.294                     
## Rpr.trtmntH -0.305  0.683              
## Cndtn...Fml -0.744  0.302  0.296       
## Dv.trHH:R.H  0.102 -0.707 -0.760 -0.064
anova(lmer.JSL)
## Type III Analysis of Variance Table with Satterthwaite's method
##                                Sum Sq Mean Sq NumDF  DenDF F value    Pr(>F)
## Dev.treatment                 0.48176 0.48176     1 363.27  11.012 0.0009967
## Repro.treatment               0.76060 0.76060     1 364.94  17.386 3.815e-05
## Dev.treatment:Repro.treatment 1.02570 1.02570     1 359.87  23.446 1.911e-06
##                                  
## Dev.treatment                 ***
## Repro.treatment               ***
## Dev.treatment:Repro.treatment ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
hist(residuals(lmer.JSL), col="darkgray")

plot(fitted(lmer.JSL), residuals(lmer.JSL))

JSLEMM = (emmeans(lmer.JSL, ~ Dev.treatment*Repro.treatment) %>% as.data.frame) 
JSLEMM
##   Dev.treatment Repro.treatment   emmean        SE       df lower.CL upper.CL
## 1            CC               C 5.121298 0.1214039 9.092307 4.847088 5.395508
## 2            HH               C 5.177825 0.1229252 9.474479 4.901857 5.453793
## 3            CC               H 5.155965 0.1227008 9.385730 4.880126 5.431805
## 4            HH               H 4.931053 0.1230448 9.582998 4.655266 5.206839
#Plot
ggplot(JSLEMM, aes(x=Repro.treatment, y=emmean, colour=Dev.treatment))+ geom_point() + geom_linerange(aes (ymin = emmean-SE, ymax = emmean+SE))  +  theme_classic()

##Weight

SLgroup_meanJW = ddply(Juv, "Female.SL", summarise, grp.mean=mean(Juvenile.weight..mg.))
Wgroup_meanJW = ddply(Juv, "Female.weight", summarise, grp.mean=mean(Juvenile.weight..mg.))
Condgroup_meanJW = ddply(Juv, "Condition...Female", summarise, grp.mean=mean(Juvenile.weight..mg.))
lm.1 = lm(grp.mean ~ Female.SL, data=SLgroup_meanJW)
summary(lm.1)
## 
## Call:
## lm(formula = grp.mean ~ Female.SL, data = SLgroup_meanJW)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.47331 -0.25214 -0.01798  0.12886  0.58514 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   7.0053     1.4241   4.919 0.000111 ***
## Female.SL    -0.3406     0.1601  -2.128 0.047385 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.324 on 18 degrees of freedom
## Multiple R-squared:  0.201,  Adjusted R-squared:  0.1567 
## F-statistic:  4.53 on 1 and 18 DF,  p-value: 0.04739
lm.2 = lm(grp.mean ~ Female.weight, data=Wgroup_meanJW)
summary(lm.2)
## 
## Call:
## lm(formula = grp.mean ~ Female.weight, data = Wgroup_meanJW)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.6065 -0.1954 -0.1037  0.2387  0.6496 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    4.96712    0.48302  10.283  5.8e-09 ***
## Female.weight -0.03341    0.01613  -2.071    0.053 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3257 on 18 degrees of freedom
## Multiple R-squared:  0.1924, Adjusted R-squared:  0.1475 
## F-statistic: 4.288 on 1 and 18 DF,  p-value: 0.05305
lm.3 = lm(grp.mean ~ Condition...Female, data=Condgroup_meanJW)
summary(lm.3)
## 
## Call:
## lm(formula = grp.mean ~ Condition...Female, data = Condgroup_meanJW)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.62814 -0.26582 -0.02711  0.24200  0.71916 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)    
## (Intercept)         4.01708    0.67632   5.940 1.28e-05 ***
## Condition...Female -0.00918    0.15936  -0.058    0.955    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3624 on 18 degrees of freedom
## Multiple R-squared:  0.0001843,  Adjusted R-squared:  -0.05536 
## F-statistic: 0.003318 on 1 and 18 DF,  p-value: 0.9547
lmer.JWf = lmer(Juvenile.weight..mg.~ Dev.treatment * Repro.treatment  + Female.SL + (1|Male.family) + (1|Female.family) + (1|Fish.replicate ) , data = Juv)

summary(lmer.JWf)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Juvenile.weight..mg. ~ Dev.treatment * Repro.treatment + Female.SL +  
##     (1 | Male.family) + (1 | Female.family) + (1 | Fish.replicate)
##    Data: Juv
## 
## REML criterion at convergence: -11
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.6193 -0.5700 -0.0990  0.5809  3.5744 
## 
## Random effects:
##  Groups         Name        Variance Std.Dev.
##  Fish.replicate (Intercept) 0.05742  0.2396  
##  Female.family  (Intercept) 0.08289  0.2879  
##  Male.family    (Intercept) 0.08641  0.2939  
##  Residual                   0.04142  0.2035  
## Number of obs: 387, groups:  
## Fish.replicate, 21; Female.family, 5; Male.family, 5
## 
## Fixed effects:
##                                    Estimate Std. Error         df t value
## (Intercept)                       1.049e+01  4.269e-01  1.470e+02  24.569
## Dev.treatmentHH                   3.895e-02  3.881e-02  3.608e+02   1.004
## Repro.treatmentH                 -8.355e-04  3.985e-02  3.610e+02  -0.021
## Female.SL                        -7.401e-01  4.213e-02  3.607e+02 -17.565
## Dev.treatmentHH:Repro.treatmentH -5.536e-02  5.744e-02  3.609e+02  -0.964
##                                  Pr(>|t|)    
## (Intercept)                        <2e-16 ***
## Dev.treatmentHH                     0.316    
## Repro.treatmentH                    0.983    
## Female.SL                          <2e-16 ***
## Dev.treatmentHH:Repro.treatmentH    0.336    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Dv.tHH Rpr.tH Fml.SL
## Dv.trtmntHH -0.367                     
## Rpr.trtmntH  0.021  0.582              
## Female.SL   -0.891  0.380 -0.066       
## Dv.trHH:R.H -0.041 -0.642 -0.783  0.072
anova(lmer.JWf)
## Type III Analysis of Variance Table with Satterthwaite's method
##                                Sum Sq Mean Sq NumDF  DenDF  F value Pr(>F)    
## Dev.treatment                  0.0058  0.0058     1 358.05   0.1412 0.7073    
## Repro.treatment                0.0542  0.0542     1 359.47   1.3092 0.2533    
## Female.SL                     12.7805 12.7805     1 360.69 308.5450 <2e-16 ***
## Dev.treatment:Repro.treatment  0.0385  0.0385     1 360.90   0.9289 0.3358    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
hist(residuals(lmer.JWf), col="darkgray")

plot(fitted(lmer.JWf), residuals(lmer.JWf))

JWEMM = (emmeans(lmer.JWf, ~ Dev.treatment*Repro.treatment) %>% as.data.frame) 
JWEMM
##   Dev.treatment Repro.treatment   emmean        SE       df lower.CL upper.CL
## 1            CC               C 3.895528 0.1935436 9.402870 3.460545 4.330511
## 2            HH               C 3.934479 0.1946170 9.600867 3.498389 4.370569
## 3            CC               H 3.894692 0.1944014 9.558428 3.458811 4.330573
## 4            HH               H 3.878283 0.1947330 9.636570 3.442162 4.314403
#Plot
ggplot(JWEMM, aes(x=Repro.treatment, y=emmean, colour=Dev.treatment))+ geom_point() + geom_linerange(aes (ymin = emmean-SE, ymax = emmean+SE))  +  theme_classic()

##Yolk area

SLgroup_meanJYa = ddply(Juv, "Female.SL", summarise, grp.mean=mean(Juvenile.yolk.area..mm.2.))
Wgroup_meanJYa = ddply(Juv, "Female.weight", summarise, grp.mean=mean(Juvenile.yolk.area..mm.2.))
Condgroup_meanJYa = ddply(Juv, "Condition...Female", summarise, grp.mean=mean(Juvenile.yolk.area..mm.2.))
lm.1 = lm(grp.mean ~ Female.SL, data=SLgroup_meanJYa)
summary(lm.1)
## 
## Call:
## lm(formula = grp.mean ~ Female.SL, data = SLgroup_meanJYa)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.33393 -0.14351  0.01688  0.10027  0.34179 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept)  2.57726    0.79116   3.258  0.00437 **
## Female.SL   -0.13420    0.08892  -1.509  0.14859   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.18 on 18 degrees of freedom
## Multiple R-squared:  0.1123, Adjusted R-squared:  0.06301 
## F-statistic: 2.278 on 1 and 18 DF,  p-value: 0.1486
lm.2 = lm(grp.mean ~ Female.weight, data=Wgroup_meanJYa)
summary(lm.2)
## 
## Call:
## lm(formula = grp.mean ~ Female.weight, data = Wgroup_meanJYa)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.34399 -0.16740  0.03817  0.14471  0.33017 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    1.496589   0.282027   5.307 4.81e-05 ***
## Female.weight -0.003779   0.009421  -0.401    0.693    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1902 on 18 degrees of freedom
## Multiple R-squared:  0.008859,   Adjusted R-squared:  -0.0462 
## F-statistic: 0.1609 on 1 and 18 DF,  p-value: 0.6931
lm.3 = lm(grp.mean ~ Condition...Female, data=Condgroup_meanJYa)
summary(lm.3)
## 
## Call:
## lm(formula = grp.mean ~ Condition...Female, data = Condgroup_meanJYa)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.36281 -0.06290 -0.01345  0.13687  0.28644 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)  
## (Intercept)         0.83418    0.33167   2.515   0.0216 *
## Condition...Female  0.13068    0.07815   1.672   0.1118  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1777 on 18 degrees of freedom
## Multiple R-squared:  0.1344, Adjusted R-squared:  0.08636 
## F-statistic: 2.796 on 1 and 18 DF,  p-value: 0.1118
lmer.JYa = lmer(Juvenile.yolk.area..mm.2.~ Dev.treatment * Repro.treatment  + (1|Male.family) + (1|Female.family) + (1|Fish.replicate), data = Juv) 
## boundary (singular) fit: see ?isSingular
summary(lmer.JYa)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Juvenile.yolk.area..mm.2. ~ Dev.treatment * Repro.treatment +  
##     (1 | Male.family) + (1 | Female.family) + (1 | Fish.replicate)
##    Data: Juv
## 
## REML criterion at convergence: -97.8
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -2.62263 -0.73338 -0.03937  0.70882  2.57065 
## 
## Random effects:
##  Groups         Name        Variance Std.Dev.
##  Fish.replicate (Intercept) 0.000000 0.00000 
##  Female.family  (Intercept) 0.016885 0.12994 
##  Male.family    (Intercept) 0.007539 0.08683 
##  Residual                   0.040989 0.20246 
## Number of obs: 387, groups:  
## Fish.replicate, 21; Female.family, 5; Male.family, 5
## 
## Fixed effects:
##                                   Estimate Std. Error        df t value
## (Intercept)                        1.37227    0.07471   7.30682  18.369
## Dev.treatmentHH                    0.02268    0.03465 349.21430   0.654
## Repro.treatmentH                  -0.06178    0.03794 304.78198  -1.629
## Dev.treatmentHH:Repro.treatmentH  -0.06982    0.05492 328.33539  -1.271
##                                  Pr(>|t|)    
## (Intercept)                      2.21e-07 ***
## Dev.treatmentHH                     0.513    
## Repro.treatmentH                    0.104    
## Dev.treatmentHH:Repro.treatmentH    0.205    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Dv.tHH Rpr.tH
## Dv.trtmntHH -0.178              
## Rpr.trtmntH -0.216  0.640       
## Dv.trHH:R.H  0.139 -0.717 -0.770
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular
anova(lmer.JYa)
## Type III Analysis of Variance Table with Satterthwaite's method
##                                Sum Sq Mean Sq NumDF  DenDF F value    Pr(>F)
## Dev.treatment                 0.01040 0.01040     1 376.02  0.2536    0.6148
## Repro.treatment               0.64959 0.64959     1 376.89 15.8480 8.236e-05
## Dev.treatment:Repro.treatment 0.06624 0.06624     1 328.34  1.6160    0.2046
##                                  
## Dev.treatment                    
## Repro.treatment               ***
## Dev.treatment:Repro.treatment    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
hist(residuals(lmer.JYa), col="darkgray")

plot(fitted(lmer.JYa), residuals(lmer.JYa))

JYaEMM = (emmeans(lmer.JYa, ~ Dev.treatment*Repro.treatment) %>% as.data.frame) 
JYaEMM
##   Dev.treatment Repro.treatment   emmean         SE       df lower.CL upper.CL
## 1            CC               C 1.372274 0.07506646 7.431329 1.196837 1.547710
## 2            HH               C 1.394950 0.07726475 8.042569 1.216941 1.572959
## 3            CC               H 1.310489 0.07694533 7.846340 1.132446 1.488532
## 4            HH               H 1.263347 0.07772558 8.480500 1.085864 1.440831
#Plot
ggplot(JYaEMM, aes(x=Repro.treatment, y=emmean, colour=Dev.treatment))+ geom_point() + geom_linerange(aes (ymin = emmean-SE, ymax = emmean+SE))  +  theme_classic()

#Fulton’s K condition

SLgroup_meanJFK = ddply(Juv, "Female.SL", summarise, grp.mean=mean(Fulton.s.K))
Wgroup_meanJFK = ddply(Juv, "Female.weight", summarise, grp.mean=mean(Fulton.s.K))
Condgroup_meanJFK = ddply(Juv, "Condition...Female", summarise, grp.mean=mean(Fulton.s.K))
lm.1 = lm(grp.mean ~ Female.SL, data=SLgroup_meanJFK)
summary(lm.1)
## 
## Call:
## lm(formula = grp.mean ~ Female.SL, data = SLgroup_meanJFK)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.6924 -0.3236 -0.1239  0.3432  1.0251 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept)   5.0510     2.1818   2.315   0.0326 *
## Female.SL    -0.2473     0.2452  -1.009   0.3266  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4963 on 18 degrees of freedom
## Multiple R-squared:  0.05349,    Adjusted R-squared:  0.0009037 
## F-statistic: 1.017 on 1 and 18 DF,  p-value: 0.3266
lm.2 = lm(grp.mean ~ Female.weight, data=Wgroup_meanJFK)
summary(lm.2)
## 
## Call:
## lm(formula = grp.mean ~ Female.weight, data = Wgroup_meanJFK)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.6890 -0.3180 -0.1461  0.3467  0.9032 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    3.87737    0.71606   5.415 3.82e-05 ***
## Female.weight -0.03460    0.02392  -1.446    0.165    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4829 on 18 degrees of freedom
## Multiple R-squared:  0.1041, Adjusted R-squared:  0.05437 
## F-statistic: 2.092 on 1 and 18 DF,  p-value: 0.1652
lm.3 = lm(grp.mean ~ Condition...Female, data=Condgroup_meanJFK)
summary(lm.3)
## 
## Call:
## lm(formula = grp.mean ~ Condition...Female, data = Condgroup_meanJFK)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.7117 -0.3761 -0.0637  0.3482  1.1111 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)   
## (Intercept)          3.5553     0.9374   3.793  0.00133 **
## Condition...Female  -0.1666     0.2209  -0.754  0.46051   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5023 on 18 degrees of freedom
## Multiple R-squared:  0.03063,    Adjusted R-squared:  -0.02323 
## F-statistic: 0.5687 on 1 and 18 DF,  p-value: 0.4605
lmer.JFK = lmer(Fulton.s.K~ Dev.treatment * Repro.treatment  + (1|Male.family) + (1|Female.family) + (1|Fish.replicate ) , data = Juv) 

summary(lmer.JFK)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Fulton.s.K ~ Dev.treatment * Repro.treatment + (1 | Male.family) +  
##     (1 | Female.family) + (1 | Fish.replicate)
##    Data: Juv
## 
## REML criterion at convergence: 557.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.5218 -0.7871 -0.1323  0.4151  3.8853 
## 
## Random effects:
##  Groups         Name        Variance Std.Dev.
##  Fish.replicate (Intercept) 0.01051  0.1025  
##  Female.family  (Intercept) 0.04446  0.2109  
##  Male.family    (Intercept) 0.03248  0.1802  
##  Residual                   0.22147  0.4706  
## Number of obs: 387, groups:  
## Fish.replicate, 21; Female.family, 5; Male.family, 5
## 
## Fixed effects:
##                                   Estimate Std. Error        df t value
## (Intercept)                        2.87212    0.13971   8.84836  20.558
## Dev.treatmentHH                    0.10101    0.07939 306.65530   1.272
## Repro.treatmentH                  -0.11480    0.08612 231.57682  -1.333
## Dev.treatmentHH:Repro.treatmentH   0.51316    0.12547 278.67265   4.090
##                                  Pr(>|t|)    
## (Intercept)                      8.94e-09 ***
## Dev.treatmentHH                     0.204    
## Repro.treatmentH                    0.184    
## Dev.treatmentHH:Repro.treatmentH 5.65e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Dv.tHH Rpr.tH
## Dv.trtmntHH -0.223              
## Rpr.trtmntH -0.267  0.630       
## Dv.trHH:R.H  0.171 -0.712 -0.761
anova(lmer.JFK)
## Type III Analysis of Variance Table with Satterthwaite's method
##                               Sum Sq Mean Sq NumDF  DenDF F value    Pr(>F)    
## Dev.treatment                 9.0092  9.0092     1 352.79 40.6800 5.627e-10 ***
## Repro.treatment               1.4243  1.4243     1 344.28  6.4313   0.01165 *  
## Dev.treatment:Repro.treatment 3.7046  3.7046     1 278.67 16.7278 5.652e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
hist(residuals(lmer.JFK), col="darkgray")

plot(fitted(lmer.JFK), residuals(lmer.JFK))

JFKEMM = (emmeans(lmer.JFK, ~ Dev.treatment*Repro.treatment) %>% as.data.frame) 
JFKEMM
##   Dev.treatment Repro.treatment   emmean        SE        df lower.CL upper.CL
## 1            CC               C 2.872121 0.1411734  8.735380 2.551284 3.192958
## 2            HH               C 2.973130 0.1472613  9.550212 2.642904 3.303355
## 3            CC               H 2.757324 0.1462223  9.168547 2.427470 3.087177
## 4            HH               H 3.371497 0.1486572 10.560094 3.042635 3.700360
#Plot
ggplot(JFKEMM, aes(x=Repro.treatment, y=emmean, colour=Dev.treatment))+ geom_point() + geom_linerange(aes (ymin = emmean - SE, ymax = emmean + SE))  +  theme_classic()