Modelo de efectos aleatorios
Sebastian Castro
29/8/2021
#modelo de efectos fijos
library(readxl)
Modelo_de_efectos_ale <- read_excel("C:/Users/juanc/OneDrive/Escritorio/2021-1/Suelos/disenio/Modelo de efectos ale.xlsx",
sheet = "Fijos")
df=Modelo_de_efectos_ale
gp=Modelo_de_efectos_ale$Genotipo
rto=Modelo_de_efectos_ale$Rendimiento
am=Modelo_de_efectos_ale$Ambiente
rep= Modelo_de_efectos_ale$Repeticionlibrary(collapsibleTree)
collapsibleTreeSummary(Modelo_de_efectos_ale, hierarchy=c ("Ambiente", "Genotipo", "Repeticion", "Rendimiento"))library(lattice)
bwplot(rto~gp|am, df)
Mod1<-aov(rto~gp*am,df)
summary(Mod1)## Df Sum Sq Mean Sq F value Pr(>F)
## gp 29 1188 40.96 0.809 0.7487
## am 9 854 94.90 1.875 0.0553 .
## gp:am 261 11747 45.01 0.889 0.8356
## Residuals 300 15186 50.62
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1res=Mod1$residuals TukeyHSD(Mod1, "am" )## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = rto ~ gp * am, data = df)
##
## $am
## diff lwr upr p adj
## A10-A1 -0.78870859 -4.9294798 3.3520626 0.9998455
## A2-A1 -0.06710292 -4.2078741 4.0736683 1.0000000
## A3-A1 1.56194045 -2.5788307 5.7027116 0.9715723
## A4-A1 -0.12332170 -4.2640929 4.0174495 1.0000000
## A5-A1 -0.16710307 -4.3078743 3.9736681 1.0000000
## A6-A1 2.46128822 -1.6794830 6.6020594 0.6726347
## A7-A1 1.24199704 -2.8987742 5.3827682 0.9943232
## A8-A1 1.76818427 -2.3725869 5.9089555 0.9378510
## A9-A1 -1.51880753 -5.6595787 2.6219637 0.9764024
## A2-A10 0.72160568 -3.4191655 4.8623769 0.9999267
## A3-A10 2.35064905 -1.7901221 6.4914202 0.7287395
## A4-A10 0.66538689 -3.4753843 4.8061581 0.9999632
## A5-A10 0.62160552 -3.5191657 4.7623767 0.9999794
## A6-A10 3.24999682 -0.8907744 7.3907680 0.2716435
## A7-A10 2.03070564 -2.1100656 6.1714768 0.8645961
## A8-A10 2.55689287 -1.5838783 6.6976641 0.6219031
## A9-A10 -0.73009894 -4.8708701 3.4106723 0.9999191
## A3-A2 1.62904337 -2.5117278 5.7698146 0.9626363
## A4-A2 -0.05621879 -4.1969900 4.0845524 1.0000000
## A5-A2 -0.10000016 -4.2407714 4.0407710 1.0000000
## A6-A2 2.52839114 -1.6123801 6.6691623 0.6371875
## A7-A2 1.30909996 -2.8316712 5.4498712 0.9916493
## A8-A2 1.83528719 -2.3054840 5.9760584 0.9225546
## A9-A2 -1.45170461 -5.5924758 2.6890666 0.9826444
## A4-A3 -1.68526215 -5.8260333 2.4555090 0.9536933
## A5-A3 -1.72904352 -5.8698147 2.4117277 0.9457374
## A6-A3 0.89934777 -3.2414234 5.0401190 0.9995442
## A7-A3 -0.31994341 -4.4607146 3.8208278 0.9999999
## A8-A3 0.20624382 -3.9345274 4.3470150 1.0000000
## A9-A3 -3.08074798 -7.2215192 1.0600232 0.3467272
## A5-A4 -0.04378137 -4.1845526 4.0969898 1.0000000
## A6-A4 2.58460993 -1.5561613 6.7253811 0.6069427
## A7-A4 1.36531874 -2.7754525 5.5060899 0.9887049
## A8-A4 1.89150597 -2.2492652 6.0322772 0.9079513
## A9-A4 -1.39548583 -5.5362570 2.7452854 0.9868167
## A6-A5 2.62839130 -1.5123799 6.7691625 0.5831720
## A7-A5 1.40910011 -2.7316711 5.5498713 0.9858873
## A8-A5 1.93528734 -2.2054839 6.0760585 0.8954240
## A9-A5 -1.35170446 -5.4924757 2.7890667 0.9894837
## A7-A6 -1.21929118 -5.3600624 2.9214800 0.9950515
## A8-A6 -0.69310395 -4.8338751 3.4476672 0.9999479
## A9-A6 -3.98009576 -8.1208670 0.1606754 0.0711475
## A8-A7 0.52618723 -3.6145840 4.6669584 0.9999951
## A9-A7 -2.76080457 -6.9015758 1.3799666 0.5110663
## A9-A8 -3.28699180 -7.4277630 0.8537794 0.2566353trat=interaction(df$Ambiente, df$Genotipo)
bartlett.test(res,trat)##
## Bartlett test of homogeneity of variances
##
## data: res and trat
## Bartlett's K-squared = 282.83, df = 299, p-value = 0.7409shapiro.test(res)##
## Shapiro-Wilk normality test
##
## data: res
## W = 0.99926, p-value = 0.9985library(agricolae)## Warning: package 'agricolae' was built under R version 4.0.5HSD_test = with(df, HSD.test (rto,gp, DFerror = 300,MSerror =50.62 ))
HSD_test## $statistics
## MSerror Df Mean CV MSD
## 50.62 300 57.08848 12.46271 8.518008
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey gp 30 5.354168 0.05
##
## $means
## rto std r Min Max Q25 Q50 Q75
## G1 58.35899 7.207985 20 45.05842 70.94239 53.71168 58.00608 63.89135
## G10 55.70308 7.897205 20 36.70814 65.16678 50.90288 57.51116 61.98903
## G11 58.04696 6.673184 20 47.47905 77.04758 54.44748 58.35469 60.75634
## G12 57.25749 6.342240 20 48.46189 70.75955 51.36867 57.24249 60.47036
## G13 56.62283 7.035601 20 45.35577 69.84443 51.68183 55.83278 62.14081
## G14 56.44755 6.154186 20 42.70095 67.87662 53.48398 55.76298 59.49921
## G15 57.43690 8.117622 20 46.58297 77.65315 51.72482 55.16652 62.35256
## G16 57.78743 6.153433 20 48.30204 66.19857 52.45112 58.27253 63.62769
## G17 56.02427 7.257473 20 43.50510 70.75605 51.41274 55.09266 60.02712
## G18 58.31929 5.998491 20 47.77372 71.64907 54.08621 58.03184 60.95759
## G19 56.59429 6.867614 20 44.69799 68.60224 52.12413 58.08700 60.47657
## G2 57.53389 6.521398 20 45.55058 72.05589 54.81996 57.52891 58.48868
## G20 59.20955 8.408406 20 44.84559 76.73944 52.54745 57.89273 62.42598
## G21 56.63794 6.968888 20 47.70700 73.92101 51.01855 54.64944 62.62744
## G22 57.31181 7.865283 20 42.00656 68.36092 54.38791 58.59302 61.79965
## G23 57.68668 7.090801 20 43.16641 73.19179 53.36171 56.94719 62.77165
## G24 57.33594 5.480070 20 47.64349 66.08835 53.07716 58.47754 61.20624
## G25 57.76414 5.577331 20 42.01217 66.20714 54.44717 58.90349 61.20355
## G26 58.60604 8.486411 20 42.42955 70.94239 54.60162 59.27524 64.18655
## G27 55.87573 5.708256 20 44.68213 66.39034 52.81415 56.43196 59.45212
## G28 58.00699 6.958119 20 47.24766 69.34638 52.06523 57.95164 63.63498
## G29 56.76300 6.237932 20 49.21477 69.59528 51.18458 56.03837 60.02907
## G3 53.81733 7.590046 20 43.59831 68.82169 48.70935 51.41853 57.92458
## G30 56.57924 6.554179 20 43.58411 64.58179 52.74687 58.53536 61.32241
## G4 53.86848 7.566430 20 39.66982 64.47528 46.35471 55.01609 60.83103
## G5 54.98845 7.078077 20 42.23015 69.23220 50.09500 54.52636 61.35426
## G6 56.09069 7.651011 20 45.42026 76.56118 51.69196 55.29948 57.15114
## G7 57.66413 6.891836 20 40.93420 67.87150 53.37411 57.60743 62.53205
## G8 60.67447 7.239185 20 39.30133 73.51393 57.51044 61.70765 63.10049
## G9 57.64075 6.562881 20 46.40534 67.23692 52.93595 58.12550 62.73768
##
## $comparison
## NULL
##
## $groups
## rto groups
## G8 60.67447 a
## G20 59.20955 a
## G26 58.60604 a
## G1 58.35899 a
## G18 58.31929 a
## G11 58.04696 a
## G28 58.00699 a
## G16 57.78743 a
## G25 57.76414 a
## G23 57.68668 a
## G7 57.66413 a
## G9 57.64075 a
## G2 57.53389 a
## G15 57.43690 a
## G24 57.33594 a
## G22 57.31181 a
## G12 57.25749 a
## G29 56.76300 a
## G21 56.63794 a
## G13 56.62283 a
## G19 56.59429 a
## G30 56.57924 a
## G14 56.44755 a
## G6 56.09069 a
## G17 56.02427 a
## G27 55.87573 a
## G10 55.70308 a
## G5 54.98845 a
## G4 53.86848 a
## G3 53.81733 a
##
## attr(,"class")
## [1] "group"#Modelo de efectos mixtos
library(readxl)
Modelo_de_efectos_ale <- read_excel("C:/Users/juanc/OneDrive/Escritorio/2021-1/Suelos/disenio/Modelo de efectos ale.xlsx",
sheet = "Mixto")
df2= Modelo_de_efectos_ale
gp2=Modelo_de_efectos_ale$Genotipo
rto2=Modelo_de_efectos_ale$Rendimiento
am2=Modelo_de_efectos_ale$Ambiente
rep2= Modelo_de_efectos_ale$Repeticionlibrary(collapsibleTree)
collapsibleTreeSummary(Modelo_de_efectos_ale, hierarchy=c ("Ambiente", "Genotipo", "Repeticion", "Rendimiento"))library(lattice)
bwplot(rto2~gp2|am2, df2)
library(agricolae)
HSD_test = with(df2, HSD.test (rto2,am2, DFerror = 50,MSerror =39.38 ))
HSD_test## $statistics
## MSerror Df Mean CV MSD
## 39.38 50 57.03252 11.00311 9.290013
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey am2 10 4.681429 0.05
##
## $means
## rto2 std r Min Max Q25 Q50 Q75
## A1 54.91696 7.137799 10 42.01217 65.16678 52.07033 54.44584 60.80293
## A10 54.88788 6.585155 10 40.93420 63.20073 51.32055 55.55619 59.03405
## A2 59.03379 7.111374 10 46.64848 69.34638 53.96715 59.40725 64.76456
## A3 57.81941 5.029460 10 48.66816 66.16564 55.21142 57.97772 60.18572
## A4 55.01679 9.163595 10 40.09230 67.24593 51.23017 53.62719 62.81421
## A5 53.90801 9.041834 10 36.70814 65.33908 48.12927 54.29483 60.43193
## A6 59.55730 4.655292 10 52.10906 68.89239 56.76565 59.93862 62.04888
## A7 58.40209 6.113024 10 47.99314 67.74321 54.36918 59.92152 62.49634
## A8 59.94005 6.299767 10 49.58148 67.87150 54.98928 62.13899 63.33731
## A9 56.84292 6.128902 10 50.21001 70.75605 51.81814 57.40030 58.37172
##
## $comparison
## NULL
##
## $groups
## rto2 groups
## A8 59.94005 a
## A6 59.55730 a
## A2 59.03379 a
## A7 58.40209 a
## A3 57.81941 a
## A9 56.84292 a
## A4 55.01679 a
## A1 54.91696 a
## A10 54.88788 a
## A5 53.90801 a
##
## attr(,"class")
## [1] "group"library(lme4)## Warning: package 'lme4' was built under R version 4.0.5## Loading required package: Matrix## Warning: package 'Matrix' was built under R version 4.0.5Mod2 <- lmer(rto2 ~ am2 + (1 | gp2) + (1 | gp2:am2), data = df2)## boundary (singular) fit: see ?isSingularsummary(Mod2)## Linear mixed model fit by REML ['lmerMod']
## Formula: rto2 ~ am2 + (1 | gp2) + (1 | gp2:am2)
## Data: df2
##
## REML criterion at convergence: 623.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2297 -0.6079 0.0331 0.6404 1.8603
##
## Random effects:
## Groups Name Variance Std.Dev.
## gp2:am2 (Intercept) 8.814 2.969
## gp2 (Intercept) 0.000 0.000
## Residual 39.385 6.276
## Number of obs: 100, groups: gp2:am2, 50; gp2, 5
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 54.91696 2.38775 22.999
## am2A10 -0.02909 3.37678 -0.009
## am2A2 4.11683 3.37678 1.219
## am2A3 2.90245 3.37678 0.860
## am2A4 0.09982 3.37678 0.030
## am2A5 -1.00895 3.37678 -0.299
## am2A6 4.64033 3.37678 1.374
## am2A7 3.48512 3.37678 1.032
## am2A8 5.02308 3.37678 1.488
## am2A9 1.92596 3.37678 0.570
##
## Correlation of Fixed Effects:
## (Intr) am2A10 am2A2 am2A3 am2A4 am2A5 am2A6 am2A7 am2A8
## am2A10 -0.707
## am2A2 -0.707 0.500
## am2A3 -0.707 0.500 0.500
## am2A4 -0.707 0.500 0.500 0.500
## am2A5 -0.707 0.500 0.500 0.500 0.500
## am2A6 -0.707 0.500 0.500 0.500 0.500 0.500
## am2A7 -0.707 0.500 0.500 0.500 0.500 0.500 0.500
## am2A8 -0.707 0.500 0.500 0.500 0.500 0.500 0.500 0.500
## am2A9 -0.707 0.500 0.500 0.500 0.500 0.500 0.500 0.500 0.500
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see ?isSingularvar_am2_gp2 = 2.969**2
var_am2 = 0.000**2
var_gp2 = 1.242**2
var_res2 = 6.276**2
var_tot2 = var_am2_gp2 + var_am2 + var_gp2 +var_res2; var_tot2## [1] 49.7457100 * var_am2/var_tot2## [1] 0100 * var_gp2/var_tot2## [1] 3.100899#Modelo de efectos aleatorios
library(readxl)
Modelo_de_efectos_ale <- read_excel("C:/Users/juanc/OneDrive/Escritorio/2021-1/Suelos/disenio/Modelo de efectos ale.xlsx",
sheet = "Aleatorio")
df3= Modelo_de_efectos_ale
gp3=Modelo_de_efectos_ale$Genotipo
rto3=Modelo_de_efectos_ale$Rendimiento
am3=Modelo_de_efectos_ale$Ambiente
rep3= Modelo_de_efectos_ale$Repeticionlibrary(collapsibleTree)
collapsibleTreeSummary(Modelo_de_efectos_ale, hierarchy=c ("Ambiente", "Genotipo", "Repeticion", "Rendimiento"))library(lattice)
bwplot(rto2~gp2|am2, df2)
library(lme4)
Mod3 <- lmer(rto3 ~ (1 | gp3) + (1 | am3) + (1 | am3:gp3), data = df3)
summary(Mod3)## Linear mixed model fit by REML ['lmerMod']
## Formula: rto3 ~ (1 | gp3) + (1 | am3) + (1 | am3:gp3)
## Data: df3
##
## REML criterion at convergence: 267.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.97209 -0.57484 0.04893 0.62361 1.41452
##
## Random effects:
## Groups Name Variance Std.Dev.
## am3:gp3 (Intercept) 1.369e+01 3.700633
## gp3 (Intercept) 1.292e+00 1.136845
## am3 (Intercept) 2.332e-06 0.001527
## Residual 3.896e+01 6.242091
## Number of obs: 40, groups: am3:gp3, 20; gp3, 5; am3, 4
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 56.835 1.385 41.05library(agricolae)
HSD_test = with(df3, HSD.test (rto3,gp3, DFerror = 50,MSerror =39.38 ))
HSD_test## $statistics
## MSerror Df Mean CV MSD
## 39.38 50 56.83514 11.04132 8.879013
##
## $parameters
## test name.t ntr StudentizedRange alpha
## Tukey gp3 5 4.001952 0.05
##
## $means
## rto3 std r Min Max Q25 Q50 Q75
## G10 56.78451 8.771276 8 40.09230 64.79919 50.90288 60.00878 63.36261
## G17 52.57095 7.126517 8 43.50510 64.82336 47.65697 51.01275 57.04667
## G25 57.87483 5.099136 8 51.79051 66.20714 53.79919 56.84185 61.22682
## G28 61.09131 6.399892 8 50.30719 69.34638 58.04580 61.94887 65.03527
## G7 55.85409 7.670296 8 40.93420 67.74321 53.37411 55.98064 59.07018
##
## $comparison
## NULL
##
## $groups
## rto3 groups
## G28 61.09131 a
## G25 57.87483 a
## G10 56.78451 a
## G7 55.85409 a
## G17 52.57095 a
##
## attr(,"class")
## [1] "group"