1 Introducción

Imagine que desea comparar cuatro fertilizantes para determinar cuál produce el mayor rendimiento en un cultivo. Sin embargo, sabe que una parte del terreno recibe más sol y otra tiene un suelo ligeramente más fértil. Si simplemente aplica los fertilizantes al azar, podría concluir erróneamente que un fertilizante es mejor, cuando en realidad el efecto se debe a la posición en el terreno.

Este mismo problema aparece en muchos contextos:

El diseño de cuadrado latino surge precisamente para resolver este problema. Su idea es sencilla pero poderosa: comparar los tratamientos controlando simultáneamente dos fuentes de variación extraña.

De esta forma, cada tratamiento se evalúa una sola vez en cada fila y una sola vez en cada columna, logrando una comparación más justa, precisa y confiable.

En esta clase aprenderemos:

2 ¿Por qué no basta con un diseño completamente al azar?

Suponga que queremos comparar cuatro tratamientos: A, B, C y D.

Si los asignamos completamente al azar sobre el terreno, podríamos obtener algo como:

Columna 1 Columna 2 Columna 3 Columna 4
Fila 1 A A B C
Fila 2 D B C A
Fila 3 B C D D
Fila 4 C D A B

Observe que algunos tratamientos aparecen varias veces en la misma fila o en la misma columna.

Esto puede ser un problema porque:

Por ello necesitamos una forma de organizar el experimento para que todos los tratamientos tengan las mismas oportunidades.

3 Idea del diseño de cuadrado latino

El diseño de cuadrado latino organiza los tratamientos de manera que:

Por ejemplo:

Columna 1 Columna 2 Columna 3 Columna 4
Fila 1 A B C D
Fila 2 B C D A
Fila 3 C D A B
Fila 4 D A B C

Ahora todos los tratamientos están distribuidos de manera equilibrada.

4 Idea básica del diseño de cuadrado latino

Un diseño de cuadrado latino de orden \(t\):

4.1 Estructura general

Ejemplo de un cuadrado latino de orden 4:

Fila / Columna 1 2 3 4
1 A B C D
2 B C D A
3 C D A B
4 D A B C

Aquí:

  • Tratamientos: A, B, C, D
  • Bloque 1: filas
  • Bloque 2: columnas

5 ¿Cuándo usar este diseño?

Se usa cuando:

  1. hay un solo factor de interés: tratamiento;
  2. existen dos fuentes de variación que se desean controlar;
  3. el número de tratamientos es igual al número de niveles de cada factor de bloqueo;
  4. es razonable asumir que no hay interacción entre tratamiento y bloques.

6 Ventajas y limitaciones

6.1 Ventajas

  • controla dos fuentes de variación extraña;
  • mejora la precisión de la comparación entre tratamientos;
  • reduce el error experimental;
  • es eficiente cuando el número de tratamientos no es grande.

6.2 Limitaciones

  • requiere el mismo número de filas, columnas y tratamientos;
  • se vuelve poco práctico cuando el número de tratamientos es grande;
  • no permite estimar interacciones entre tratamientos y bloques;
  • si faltan observaciones, el análisis se complica.

7 Modelo estadístico

Sea \(Y_{ijk}\) la respuesta observada para el tratamiento \(k\), ubicado en la fila \(i\) y columna \(j\).

El modelo aditivo del diseño de cuadrado latino es:

\[ Y_{ijk} = \mu + \alpha_i + \beta_j + \tau_k + \varepsilon_{ijk} \]

donde:

Se asume que:

\[ \varepsilon_{ijk} \sim N(0,\sigma^2) \]

independientes.

7.1 Restricciones usuales

\[ \sum_i \alpha_i = 0, \qquad \sum_j \beta_j = 0, \qquad \sum_k \tau_k = 0 \]

8 Hipótesis de interés

La hipótesis principal es sobre tratamientos:

\[ H_0: \tau_1 = \tau_2 = \cdots = \tau_t = 0 \]

es decir, no hay diferencias entre tratamientos.

\[ H_1: \text{al menos uno de los efectos de tratamiento es diferente} \]

También pueden evaluarse los efectos de filas y columnas, aunque usualmente el interés central está en los tratamientos.

9 Grados de libertad y tabla ANOVA

Para un cuadrado latino de orden \(t\):

9.1 Tabla ANOVA general

Fuente de variación GL SC CM F
Filas \(t-1\) \(SC_F\) \(CM_F\) \(CM_F/CM_E\)
Columnas \(t-1\) \(SC_C\) \(CM_C\) \(CM_C/CM_E\)
Tratamientos \(t-1\) \(SC_T\) \(CM_T\) \(CM_T/CM_E\)
Error \((t-1)(t-2)\) \(SC_E\) \(CM_E\)
Total \(t^2-1\) \(SC_{Total}\)

10 Ejemplo conceptual

Supongamos que se quieren comparar 4 tratamientos (A, B, C y D) para evaluar rendimiento.

Se sospecha que la respuesta puede cambiar por:

Entonces, se organiza el experimento como un cuadrado latino de orden 4.

11 Base de datos de ejemplo en R

A continuación se construye una base de datos sencilla con un cuadrado latino de orden 4.

datos <- data.frame(
  fila = factor(rep(1:4, each = 4)),
  columna = factor(rep(1:4, times = 4)),
  tratamiento = factor(c("A","B","C","D",
                         "B","C","D","A",
                         "C","D","A","B",
                         "D","A","B","C")),
  respuesta = c(20, 24, 19, 18,
                22, 25, 23, 21,
                18, 20, 17, 19,
                16, 18, 17, 15)
)

datos
##    fila columna tratamiento respuesta
## 1     1       1           A        20
## 2     1       2           B        24
## 3     1       3           C        19
## 4     1       4           D        18
## 5     2       1           B        22
## 6     2       2           C        25
## 7     2       3           D        23
## 8     2       4           A        21
## 9     3       1           C        18
## 10    3       2           D        20
## 11    3       3           A        17
## 12    3       4           B        19
## 13    4       1           D        16
## 14    4       2           A        18
## 15    4       3           B        17
## 16    4       4           C        15

12 Análisis del diseño en R

12.1 Ajuste del modelo

modelo <- aov(respuesta ~ fila + columna + tratamiento, data = datos)
summary(modelo)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## fila         3   84.5  28.167   30.73 0.000488 ***
## columna      3   28.5   9.500   10.36 0.008679 ** 
## tratamiento  3    5.5   1.833    2.00 0.215553    
## Residuals    6    5.5   0.917                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

12.2 Interpretación

La tabla ANOVA permite responder:

  • ¿hay diferencias significativas entre tratamientos?
  • ¿las filas explican parte importante de la variabilidad?
  • ¿las columnas explican parte importante de la variabilidad?

La decisión se toma con el p-valor:

  • si \(p < 0.05\), se rechaza \(H_0\);
  • si \(p \ge 0.05\), no se rechaza \(H_0\).

13 Medias por tratamiento

aggregate(respuesta ~ tratamiento, data = datos, mean)
##   tratamiento respuesta
## 1           A     19.00
## 2           B     20.50
## 3           C     19.25
## 4           D     19.25

También se puede usar:

tapply(datos$respuesta, datos$tratamiento, mean)
##     A     B     C     D 
## 19.00 20.50 19.25 19.25

14 Comparaciones múltiples

Si el efecto de tratamiento es significativo, se pueden comparar medias.

TukeyHSD(modelo, which = "tratamiento")
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = respuesta ~ fila + columna + tratamiento, data = datos)
## 
## $tratamiento
##      diff        lwr     upr     p adj
## B-A  1.50 -0.8435897 3.84359 0.2211120
## C-A  0.25 -2.0935897 2.59359 0.9812251
## D-A  0.25 -2.0935897 2.59359 0.9812251
## C-B -1.25 -3.5935897 1.09359 0.3396748
## D-B -1.25 -3.5935897 1.09359 0.3396748
## D-C  0.00 -2.3435897 2.34359 1.0000000

15 Verificación de supuestos

15.1 Gráfico de residuos

par(mfrow = c(1,2))
plot(modelo, which = 1)
plot(modelo, which = 2)

par(mfrow = c(1,1))

15.2 Normalidad de residuos

shapiro.test(residuals(modelo))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(modelo)
## W = 0.97335, p-value = 0.8896

15.3 Homogeneidad de varianzas

library(car)

leveneTest(respuesta ~ tratamiento, data = datos)
## Levene's Test for Homogeneity of Variance (center = median)
##       Df F value Pr(>F)
## group  3  0.3836 0.7668
##       12

Nota: En el diseño de cuadrado latino, la revisión de supuestos se basa principalmente en los residuos del modelo ajustado.

16 Interpretación escrita de ejemplo

Suponga que al analizar la tabla ANOVA se obtiene un p-valor menor que 0.05 para tratamientos.

Una posible redacción sería:

Se encontró evidencia estadísticamente significativa para afirmar que existen diferencias en la respuesta media entre los tratamientos evaluados, después de controlar el efecto de las filas y las columnas del experimento.

Si además las filas resultan significativas:

El factor fila también explicó una proporción importante de la variación en la respuesta, lo cual justifica su inclusión como factor de bloqueo.

17 Cálculo manual de la lógica ANOVA

En un diseño de cuadrado latino, la variabilidad total se descompone en:

\[ SC_{Total} = SC_{Filas} + SC_{Columnas} + SC_{Tratamientos} + SC_{Error} \]

La idea es separar la parte de la variación que se debe a:

Esto permite una comparación más limpia entre tratamientos, porque se descuenta la variación atribuible a los dos bloques.

18 Ejemplo con datos simulados

set.seed(123)

fila <- factor(rep(1:5, each = 5))
columna <- factor(rep(1:5, times = 5))
tratamiento <- factor(c("A","B","C","D","E",
                        "B","C","D","E","A",
                        "C","D","E","A","B",
                        "D","E","A","B","C",
                        "E","A","B","C","D"))

ef_fila <- c(-2, -1, 0, 1, 2)
ef_col <- c(-1, 1, 0, 2, -2)
ef_trat <- c(0, 2, 4, 1, 3)

respuesta <- 20 +
  ef_fila[as.numeric(fila)] +
  ef_col[as.numeric(columna)] +
  ef_trat[as.numeric(tratamiento)] +
  rnorm(25, 0, 1.2)

datos2 <- data.frame(fila, columna, tratamiento, respuesta)

head(datos2)
##   fila columna tratamiento respuesta
## 1    1       1           A  16.32743
## 2    1       2           B  20.72379
## 3    1       3           C  23.87045
## 4    1       4           D  21.08461
## 5    1       5           E  19.15515
## 6    2       1           B  22.05808
modelo2 <- aov(respuesta ~ fila + columna + tratamiento, data = datos2)
summary(modelo2)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## fila         4  34.49   8.621   6.367 0.005482 ** 
## columna      4  57.21  14.303  10.563 0.000664 ***
## tratamiento  4  76.24  19.059  14.075 0.000175 ***
## Residuals   12  16.25   1.354                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
aggregate(respuesta ~ tratamiento, data = datos2, mean)
##   tratamiento respuesta
## 1           A  19.26079
## 2           B  22.14506
## 3           C  24.49009
## 4           D  21.07851
## 5           E  22.82558
TukeyHSD(modelo2, which = "tratamiento")
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = respuesta ~ fila + columna + tratamiento, data = datos2)
## 
## $tratamiento
##           diff           lwr        upr     p adj
## B-A  2.8842667  0.5384032209  5.2301303 0.0142678
## C-A  5.2292970  2.8834334284  7.5751605 0.0000992
## D-A  1.8177236 -0.5281398902  4.1635872 0.1623042
## E-A  3.5647872  1.2189236906  5.9106507 0.0030070
## C-B  2.3450302 -0.0008333207  4.6908937 0.0500963
## D-B -1.0665431 -3.4124066393  1.2793204 0.6105951
## E-B  0.6805205 -1.6653430585  3.0263840 0.8819775
## D-C -3.4115733 -5.7574368468 -1.0657098 0.0042430
## E-C -1.6645097 -4.0103732660  0.6813538 0.2225582
## E-D  1.7470636 -0.5987999474  4.0929271 0.1881257

18.1 Explicación del diseño

“Tenemos varios tratamientos, pero además sabemos que la respuesta puede cambiar por dos causas externas. En vez de ignorarlas, organizamos el experimento para que cada tratamiento pase una vez por cada fila y una vez por cada columna.”

19 Errores frecuentes

  1. confundir filas y columnas con tratamientos;
  2. creer que filas y columnas son factores de interés principal;
  3. olvidar que el número de tratamientos debe coincidir con el número de filas y columnas;
  4. interpretar los bloques como si fueran réplicas;
  5. no verificar que cada tratamiento aparece una sola vez en cada fila y columna;
  6. hacer comparaciones de medias sin revisar primero si el efecto de tratamiento es significativo.

20 Actividad de clase propuesta

20.1 Actividad 1. Identificación de elementos

A partir del siguiente contexto, identifique:

  • tratamiento;
  • bloque por filas;
  • bloque por columnas;
  • variable respuesta.

Contexto: Se evalúan 4 métodos de secado para semillas. El laboratorio sospecha que la humedad del ambiente cambia según el día y también según la posición del equipo en el cuarto.

20.2 Actividad 2. Construcción del diseño

Construya un cuadrado latino de orden 4 para los tratamientos A, B, C y D.

20.3 Actividad 3. Interpretación

Con base en una tabla ANOVA dada por el docente, responda:

  • ¿hay diferencias entre tratamientos?
  • ¿el bloqueo fue útil?
  • ¿qué conclusión práctica se puede comunicar?

21 Ejercicio propuesto con solución en R

21.1 Enunciado

Se comparan 4 tratamientos para evaluar el crecimiento de plántulas. Se controla el efecto de la fila del invernadero y de la columna. Los datos son:

ejercicio <- data.frame(
  fila = factor(rep(1:4, each = 4)),
  columna = factor(rep(1:4, times = 4)),
  tratamiento = factor(c("A","B","C","D",
                         "B","C","D","A",
                         "C","D","A","B",
                         "D","A","B","C")),
  respuesta = c(14, 18, 16, 15,
                17, 20, 19, 16,
                13, 15, 14, 16,
                12, 13, 15, 14)
)

ejercicio
##    fila columna tratamiento respuesta
## 1     1       1           A        14
## 2     1       2           B        18
## 3     1       3           C        16
## 4     1       4           D        15
## 5     2       1           B        17
## 6     2       2           C        20
## 7     2       3           D        19
## 8     2       4           A        16
## 9     3       1           C        13
## 10    3       2           D        15
## 11    3       3           A        14
## 12    3       4           B        16
## 13    4       1           D        12
## 14    4       2           A        13
## 15    4       3           B        15
## 16    4       4           C        14

21.2 Solución

modelo_ej <- aov(respuesta ~ fila + columna + tratamiento, data = ejercicio)
summary(modelo_ej)
##             Df Sum Sq Mean Sq F value  Pr(>F)   
## fila         3  45.19  15.062  23.323 0.00105 **
## columna      3  14.19   4.729   7.323 0.01978 * 
## tratamiento  3  10.69   3.563   5.516 0.03686 * 
## Residuals    6   3.88   0.646                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
aggregate(respuesta ~ tratamiento, data = ejercicio, mean)
##   tratamiento respuesta
## 1           A     14.25
## 2           B     16.50
## 3           C     15.75
## 4           D     15.25
TukeyHSD(modelo_ej, which = "tratamiento")
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = respuesta ~ fila + columna + tratamiento, data = ejercicio)
## 
## $tratamiento
##      diff        lwr       upr     p adj
## B-A  2.25  0.2828563 4.2171437 0.0285463
## C-A  1.50 -0.4671437 3.4671437 0.1327332
## D-A  1.00 -0.9671437 2.9671437 0.3740521
## C-B -0.75 -2.7171437 1.2171437 0.5844031
## D-B -1.25 -3.2171437 0.7171437 0.2253452
## D-C -0.50 -2.4671437 1.4671437 0.8154819

22 Recomendaciones para reportar resultados

Una estructura clara para escribir resultados es:

  1. describir el diseño usado;
  2. indicar la variable respuesta;
  3. reportar la significancia de tratamientos;
  4. mencionar si los bloques fueron relevantes;
  5. presentar medias ajustadas o comparaciones múltiples si aplica;
  6. cerrar con una conclusión sustantiva.

22.1 Ejemplo de redacción

Se empleó un diseño de cuadrado latino de orden 4 para comparar cuatro tratamientos, controlando simultáneamente el efecto de las filas y columnas. El análisis de varianza mostró diferencias estadísticamente significativas entre tratamientos (\(p < 0.05\)). Las comparaciones múltiples indicaron que el tratamiento C presentó la mayor media de respuesta.

23 Código final compacto

modelo_final <- aov(respuesta ~ fila + columna + tratamiento, data = datos)
summary(modelo_final)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## fila         3   84.5  28.167   30.73 0.000488 ***
## columna      3   28.5   9.500   10.36 0.008679 ** 
## tratamiento  3    5.5   1.833    2.00 0.215553    
## Residuals    6    5.5   0.917                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
TukeyHSD(modelo_final, which = "tratamiento")
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = respuesta ~ fila + columna + tratamiento, data = datos)
## 
## $tratamiento
##      diff        lwr     upr     p adj
## B-A  1.50 -0.8435897 3.84359 0.2211120
## C-A  0.25 -2.0935897 2.59359 0.9812251
## D-A  0.25 -2.0935897 2.59359 0.9812251
## C-B -1.25 -3.5935897 1.09359 0.3396748
## D-B -1.25 -3.5935897 1.09359 0.3396748
## D-C  0.00 -2.3435897 2.34359 1.0000000
shapiro.test(residuals(modelo_final))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(modelo_final)
## W = 0.97335, p-value = 0.8896