library(readxl)
Random_Effects_Models <- read_excel("D:/Descargas/Random Effects Models.xlsx", 
    sheet = "Modelo de Efectos Fijos")
Rend= Random_Effects_Models$Rendimiento
Gen= Random_Effects_Models$Genotipos; as.factor(Gen)
Amb= Random_Effects_Models$Ambientes; as.factor(Amb)
Rep=Random_Effects_Models$Repeticiones; as.factor(Rep)

df<-Random_Effects_Models
library(collapsibleTree)
collapsibleTreeSummary(Random_Effects_Models, hierarchy=c ("Ambientes", "Genotipos", "Repeticiones", "Rendimiento"))
boxplot( Rend~Gen, df)

boxplot( Rend~Amb, df)

medias<-tapply(Rend, list(Gen, Amb), mean)
Desv<-tapply(Rend, list(Gen, Amb), sd)
cv=100*Desv/medias

Análisis de Varianza

Mod1<-aov(Rend~Gen*Amb,df)
summary(Mod1)
##              Df Sum Sq Mean Sq F value Pr(>F)
## Gen          29  33878  1168.2   1.010  0.455
## Amb           9   8699   966.6   0.836  0.583
## Gen:Amb     261 347738  1332.3   1.152  0.117
## Residuals   300 346840  1156.1

Revisión de Supuestos

#Comparación de medias

library(agricolae)
HSD_test = with(df, HSD.test (Rend, Gen, DFerror = 300,MSerror = 202.4))
HSD_test
## $statistics
##   MSerror  Df     Mean       CV      MSD
##     202.4 300 210.5131 6.758125 17.03265
## 
## $parameters
##    test name.t ntr StudentizedRange alpha
##   Tukey    Gen  30         5.354168  0.05
## 
## $means
##            Rend      std  r      Min      Max      Q25      Q50      Q75
## Col-11 215.4756 38.88270 20 155.7863 269.9634 186.9253 209.8975 246.8758
## Col-13 215.8136 36.69735 20 154.6437 267.2387 185.2718 215.5080 244.9687
## Col-14 205.2197 31.63515 20 150.8240 249.1620 182.6148 206.9311 231.1814
## Col-15 211.6978 37.74128 20 159.7964 268.0041 176.1071 216.6012 245.7341
## Col-16 207.1711 41.62277 20 155.9694 268.3959 172.3817 198.6215 252.7546
## Col-17 202.0186 38.46938 20 151.8860 263.6460 163.3067 203.5582 231.3956
## Col-19 209.1760 27.55284 20 160.1663 247.9827 192.5404 208.5846 230.4865
## Col-20 208.8964 37.54105 20 150.0183 256.7940 171.6611 222.2391 237.6305
## Col-21 213.8787 42.35520 20 150.2673 269.6558 171.9944 224.7166 253.7086
## Col-23 204.2602 38.02769 20 152.0545 268.6230 168.7551 201.5860 231.7490
## Col-27 209.8436 36.60541 20 157.3537 267.7551 177.4804 198.6654 236.0100
## Col-3  210.9220 30.80369 20 171.3764 268.1323 184.7233 204.8051 236.0988
## Col-30 203.9498 27.90003 20 151.4942 250.2387 186.7925 208.5974 223.1941
## Col-31 214.4432 34.59560 20 154.5924 264.1734 183.9103 211.0126 246.3988
## Col-32 215.7028 31.83354 20 168.0438 262.3423 184.6583 215.5702 243.2832
## Col-33 207.7607 33.71974 20 157.0425 262.3643 180.0650 203.1663 235.7418
## Col-34 226.2848 31.64556 20 155.9621 268.4216 204.4380 232.1400 253.5365
## Col-35 208.7629 36.99964 20 161.3932 257.5558 174.8582 201.6849 245.6069
## Col-36 203.0101 32.24556 20 156.4821 251.2018 175.7015 198.6544 233.5508
## Col-37 228.9787 31.98095 20 157.5369 269.6118 209.8572 235.4122 250.7157
## Col-38 216.9059 34.54156 20 159.4559 267.8356 182.6844 227.6794 242.7998
## Col-4  203.8810 26.70744 20 150.0146 250.9235 186.7064 204.0800 224.5344
## Col-40 192.9516 31.69027 20 151.7579 256.6219 167.2472 184.0513 219.3057
## Col-41 224.7566 35.83632 20 151.0657 269.1394 203.8365 230.6366 256.8810
## Col-42 214.7102 39.72763 20 154.1969 268.4069 176.3305 219.3789 249.6701
## Col-5  212.6157 38.67325 20 156.4711 268.8830 177.8677 206.3598 255.5828
## Col-6  200.8454 35.30256 20 151.3550 258.2699 169.9426 192.5843 240.3992
## Col-7  206.5744 35.94414 20 150.4688 263.4263 174.3959 213.4498 229.4473
## Col-8  210.2855 28.58007 20 160.5509 264.2796 193.0247 212.5123 224.9062
## Col-9  208.5994 40.47984 20 150.7471 267.0812 172.7442 223.6088 237.6168
## 
## $comparison
## NULL
## 
## $groups
##            Rend groups
## Col-37 228.9787      a
## Col-34 226.2848     ab
## Col-41 224.7566    abc
## Col-38 216.9059   abcd
## Col-13 215.8136   abcd
## Col-32 215.7028   abcd
## Col-11 215.4756   abcd
## Col-42 214.7102   abcd
## Col-31 214.4432   abcd
## Col-21 213.8787   abcd
## Col-5  212.6157   abcd
## Col-15 211.6978    bcd
## Col-3  210.9220    bcd
## Col-8  210.2855    bcd
## Col-27 209.8436   bcde
## Col-19 209.1760    cde
## Col-20 208.8964    cde
## Col-35 208.7629    cde
## Col-9  208.5994    cde
## Col-33 207.7607    cde
## Col-16 207.1711     de
## Col-7  206.5744     de
## Col-14 205.2197     de
## Col-23 204.2602     de
## Col-30 203.9498     de
## Col-4  203.8810     de
## Col-36 203.0101     de
## Col-17 202.0186     de
## Col-6  200.8454     de
## Col-40 192.9516      e
## 
## attr(,"class")
## [1] "group"

#Normalidad de residuales

res<-Mod1$residuals
shapiro.test(res)
## 
##  Shapiro-Wilk normality test
## 
## data:  res
## W = 0.99289, p-value = 0.006006
plot(Mod1,2)

Modelo efectos Aleatorios

Random_Effects_Models2 <- read_excel("D:/Descargas/Random Effects Models.xlsx", 
    sheet = "Modelo de Efectos Aleatorios")
Rend2= Random_Effects_Models2$Rendimiento
Gen2= Random_Effects_Models2$Genotipos; as.factor(Gen2)
Amb2= Random_Effects_Models2$Ambientes; as.factor(Amb2)
Rep2=Random_Effects_Models2$Repeticiones; as.factor(Rep2)

df2<-Random_Effects_Models2; df2
collapsibleTreeSummary(Random_Effects_Models2, hierarchy=c ("Ambientes", "Genotipos", "Repeticiones", "Rendimiento"))
library(lattice)
boxplot( Rend2~Gen2, df2)

boxplot( Rend2~Amb2, df2)

bwplot(Rend2~Gen2|Amb2, df2)

library(Rcpp)
library(Matrix)
library(lme4)
Mod2 <- lmer(Rend2 ~ (1 | Gen2) + (1 | Amb2) + (1 | Amb2:Gen2), data = df2)
## boundary (singular) fit: see ?isSingular
summary(Mod2)
## Linear mixed model fit by REML ['lmerMod']
## Formula: Rend2 ~ (1 | Gen2) + (1 | Amb2) + (1 | Amb2:Gen2)
##    Data: df2
## 
## REML criterion at convergence: 391.6
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.54924 -0.83340  0.02417  0.84793  1.53691 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Amb2:Gen2 (Intercept)    0.000  0.000  
##  Gen2      (Intercept)    5.962  2.442  
##  Amb2      (Intercept)    0.000  0.000  
##  Residual              1217.327 34.890  
## Number of obs: 40, groups:  Amb2:Gen2, 20; Gen2, 5; Amb2, 4
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)  204.726      5.624    36.4
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingular
Mod2
## Linear mixed model fit by REML ['lmerMod']
## Formula: Rend2 ~ (1 | Gen2) + (1 | Amb2) + (1 | Amb2:Gen2)
##    Data: df2
## REML criterion at convergence: 391.5919
## Random effects:
##  Groups    Name        Std.Dev.
##  Amb2:Gen2 (Intercept)  0.000  
##  Gen2      (Intercept)  2.442  
##  Amb2      (Intercept)  0.000  
##  Residual              34.890  
## Number of obs: 40, groups:  Amb2:Gen2, 20; Gen2, 5; Amb2, 4
## Fixed Effects:
## (Intercept)  
##       204.7  
## optimizer (nloptwrap) convergence code: 0 (OK) ; 0 optimizer warnings; 1 lme4 warnings
var_amb_gen = 0.000**2
var_amb = 0.000**2
var_gen = 2.442**2
var_res = 34.890**2

var_tot = var_amb_gen + var_amb + var_gen + var_res; var_tot
## [1] 1223.275
100 * var_amb/var_tot
## [1] 0
100 * var_gen/var_tot
## [1] 0.4874915
## Plot data ####
with(df2, interaction.plot(x.factor = Gen2, trace.factor = Amb2, response = Rend2))

Modelo efectos Mixtos

Random_Effects_Models3 <- read_excel("Random Effects Models.xlsx", 
    sheet = "Modelo de Efectos Mixtos")

Rend3= Random_Effects_Models3$Rendimiento
Gen3= Random_Effects_Models3$Genotipos; as.factor(Gen3)
Amb3= Random_Effects_Models3$Ambientes; as.factor(Amb3)
Rep3=Random_Effects_Models3$Repeticiones; as.factor(Rep3)

df3<-Random_Effects_Models2; df3
collapsibleTreeSummary(Random_Effects_Models3, hierarchy=c ("Ambientes", "Genotipos", "Repeticiones", "Rendimiento"))
boxplot( Rend3~Gen3, df3)

boxplot( Rend3~Amb3, df3)

bwplot(Rend3~Gen3|Amb3, df3)

Mod3 <- lmer(Rend3 ~ Amb3 + (1 | Gen3) + (1 | Gen3:Amb3), data = df3)
summary(Mod3)
## Linear mixed model fit by REML ['lmerMod']
## Formula: Rend3 ~ Amb3 + (1 | Gen3) + (1 | Gen3:Amb3)
##    Data: df3
## 
## REML criterion at convergence: 917.5
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.80571 -0.66681  0.04136  0.59822  1.74151 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev.
##  Gen3:Amb3 (Intercept) 328.48   18.124  
##  Gen3      (Intercept)  15.24    3.903  
##  Residual              953.57   30.880  
## Number of obs: 100, groups:  Gen3:Amb3, 50; Gen3, 5
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)  205.998     12.810  16.081
## Amb3A10       27.942     17.947   1.557
## Amb3A2        15.274     17.947   0.851
## Amb3A3        -8.850     17.947  -0.493
## Amb3A4        11.630     17.947   0.648
## Amb3A5        -4.897     17.947  -0.273
## Amb3A6        -7.903     17.947  -0.440
## Amb3A7         8.708     17.947   0.485
## Amb3A8         6.769     17.947   0.377
## Amb3A9       -10.266     17.947  -0.572
## 
## Correlation of Fixed Effects:
##         (Intr) Am3A10 Amb3A2 Amb3A3 Amb3A4 Amb3A5 Amb3A6 Amb3A7 Amb3A8
## Amb3A10 -0.701                                                        
## Amb3A2  -0.701  0.500                                                 
## Amb3A3  -0.701  0.500  0.500                                          
## Amb3A4  -0.701  0.500  0.500  0.500                                   
## Amb3A5  -0.701  0.500  0.500  0.500  0.500                            
## Amb3A6  -0.701  0.500  0.500  0.500  0.500  0.500                     
## Amb3A7  -0.701  0.500  0.500  0.500  0.500  0.500  0.500              
## Amb3A8  -0.701  0.500  0.500  0.500  0.500  0.500  0.500  0.500       
## Amb3A9  -0.701  0.500  0.500  0.500  0.500  0.500  0.500  0.500  0.500
var_amb_gen2 = 0.000**2
var_amb2 = 0.000**2
var_gen2 = 2.442**2
var_res2 = 34.890**2

var_tot = var_amb_gen + var_amb + var_gen + var_res; var_tot
## [1] 1223.275
anova(Mod3)
## Analysis of Variance Table
##      npar Sum Sq Mean Sq F value
## Amb3    9 8251.3  916.81  0.9615