set.seed(123)
# RESPUESTA
diam_geom = c(rnorm(4, 1.8, 0.1), rnorm(4, 2.0, 0.12), rnorm(4, 1.9, 0.09))
# Factor
gen = gl(3,4, 12, paste0('g_', 1:3))
#Bloqueo
procedencia = gl(4, 1, 12, paste0('l_',1:4))
data = data.frame(gen, procedencia, diam_geom)
head(data)
## gen procedencia diam_geom
## 1 g_1 l_1 1.743952
## 2 g_1 l_2 1.776982
## 3 g_1 l_3 1.955871
## 4 g_1 l_4 1.807051
## 5 g_2 l_1 2.015515
## 6 g_2 l_2 2.205808
library(collapsibleTree)
## Warning: package 'collapsibleTree' was built under R version 4.1.3
collapsibleTreeSummary(data, c('procedencia', 'gen', 'diam_geom'), collapsed = FALSE)
library(ggplot2)
ggplot(data)+aes(gen, diam_geom)+geom_point(size=7, color='blue')+facet_wrap(~procedencia)+theme_dark()
\[H_0: \mu_{g_1}=\mu_{g_2}=\mu_{g_3}\]
mod1 = aov(diam_geom ~ procedencia + gen, data)
summary(mod1)
## Df Sum Sq Mean Sq F value Pr(>F)
## procedencia 3 0.04373 0.01458 1.334 0.3483
## gen 2 0.08908 0.04454 4.078 0.0762 .
## Residuals 6 0.06554 0.01092
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#sin bloques#
mod1 = aov(diam_geom ~ gen, data)
summary(mod1)
## Df Sum Sq Mean Sq F value Pr(>F)
## gen 2 0.08908 0.04454 3.669 0.0684 .
## Residuals 9 0.10927 0.01214
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Advertencia no interpretar el p-value de los bloques, solo de los factores
En este caso p-value de los genotipos > 5% (7.62%) se rechaza la hipotesis nula, por lo tanto podemos asumir estadisticamente que los genotipos son iguales
Eficiencia de bloqueo (H) corresponde a la pregunta ¿Valio la pena el bloqueo? ( Tener en cuenta la procedencia). En este caso H=1.334 (F-value de los bloques), cuando H>1 sugiere que si valio la pena bloquear.
Revision de supuestos
#1. extraer residuales del modelo (con bloques)
res_mod1 = mod1$residuals
#2. Prueba normalidad de residuales
shapiro.test(res_mod1)
##
## Shapiro-Wilk normality test
##
## data: res_mod1
## W = 0.96413, p-value = 0.8407
#3. Probando homocedasticidad (varianzas iguales)
bartlett.test(res_mod1, data$gen)
##
## Bartlett test of homogeneity of variances
##
## data: res_mod1 and data$gen
## Bartlett's K-squared = 1.1617, df = 2, p-value = 0.5594
En ambos casos los p-value de las pruebas fueron >5% por lo cual se cumplen los siguientes supuestos.
plot(data$diam_geom, res_mod1, pch=16)
Si los residuales muestran un patron se dice que estan autocorrelacionados y es un problema para el analisis de varianza .
Conclusion * Valio la pena bloquear * Estadisticamente no difieren los genotipos * Se cumplen los supuestos del anova
set.seed(123)
# RESPUESTA
diam_geom = c(rnorm(20, 1.8, 0.1), rnorm(20, 2.0, 0.12), rnorm(20, 1.9, 0.09))
# Factor
gen = gl(3,20, 60, paste0('g_', 1:3))
#Bloqueo
procedencia = gl(4, 5, 60, paste0('l_',1:4))
data = data.frame(gen, procedencia, diam_geom)
head(data)
## gen procedencia diam_geom
## 1 g_1 l_1 1.743952
## 2 g_1 l_1 1.776982
## 3 g_1 l_1 1.955871
## 4 g_1 l_1 1.807051
## 5 g_1 l_1 1.812929
## 6 g_1 l_2 1.971506
library(collapsibleTree)
collapsibleTreeSummary(data, c('procedencia', 'gen', 'diam_geom'), collapsed = FALSE)
library(lattice)
## Warning: package 'lattice' was built under R version 4.1.3
bwplot(diam_geom ~ gen | procedencia,
data)
mod3 = aov(diam_geom ~ procedencia * gen, data)
summary(mod3)
## Df Sum Sq Mean Sq F value Pr(>F)
## procedencia 3 0.0310 0.01034 1.135 0.344
## gen 2 0.3233 0.16164 17.732 1.71e-06 ***
## procedencia:gen 6 0.0407 0.00678 0.744 0.617
## Residuals 48 0.4376 0.00912
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1