Diseño factorial simple en bloques al azar

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.2.3
collapsibleTreeSummary (data, 
                        c('procedencia', 
                          'gen', 
                          'diam_geom'), 
                        collapsed = FALSE)
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.2.2
ggplot(data)+ 
  aes (gen, diam_geom)+
  geom_point (size=7, 
              color = 'yellow')+
  facet_wrap(~procedencia)+
  theme_dark()

Análisis inferencial

\[H_0: \mu_{g_1}=\mu_{g2}= \mu_{g_3}\]

mod = aov (diam_geom ~ procedencia + gen, data)
summary(mod)
##             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

Advertencia No revisar el p valor de los bloques, es decir la procedencia. Solo interpretar el de los factores, es decir genes.

Como el p valor de los genes es mayor al 5% (0.007*100= 7%>5%), no rechaza la hipotesis nula, es decir que los fenotipos son iguales, el diametro es el mismo.

Eficiencia de bloqueo

¿Valió la pena preocuparse por la procedencia de las semillas? En el análsis inferencial H= a Fvalue H = 1.334 (F-value de los bloques), cuando H>1,sugiere que si valió la pena bloquear

En ocasiones cuando el H es menor a 1, si vale la pena bloquear (artículo Carlos y profe), porque el bloqueo suele resolver problemas que no están a la vista.

Pero quitar el bloqueo porque no “sirve” y dejar el análisis sin bloqueo está mal, se deben dejar así sean menores que 1 y la conclusión al final se mantenga.

Revisión de supuestos

#Paso 1. Extraer residuales del modelo.
res_mod1 = mod$residuals
#2. Probando normalidad

shapiro.test(res_mod1)
## 
##  Shapiro-Wilk normality test
## 
## data:  res_mod1
## W = 0.97765, p-value = 0.9725
# Paso 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.709, df = 2, p-value = 0.4255

Si los patrones de los residuales son exponenciales, en forma de cono, puede ser peligroso, lo ideal es que sean constantes

En ambos casos los p-value de las pruebas fueron >5%

plot(data$diam_geom, 
     res_mod1, 
     pch=16)

Cuando los resultados de los residuales tiene patrón, se vuelve cuestionable el análisis, en este caso, no existe ese patrón, se ve algo constante.

Conclusión

Valió la pena bloquear, no hay diferencias en el genotipo, se cumplieron todos los supuestos.

Diseño #3. Factorial simple, bloques generalizados y al azar.

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)

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
#Usando * es porque hay repeticiones, + cuando no#

Primero ver la interacción=(procedencia:gen),da mayor al 5%, es decir que no hay interacción, si se cumple esto se puede ver el genotipo, si no, se deben ir a los gráficos (el de las cajas).