Nota de integridad de datos. El segundo registro de Empeltre (4.9 tn/ha) se transcribió provisionalmente de una fotografía parcialmente ilegible. Debe contrastarse con el dato original antes de considerar definitivo el informe. La interpretación del DCA también presupone unidades experimentales independientes y asignación aleatoria efectiva, condiciones aún no documentadas en el enunciado disponible.

1 Planteamiento del problema

Se estudia el rendimiento agrícola de cuatro variedades de aceituna de mesa (Sevillana, Ascolana, Empeltre e Ilo), medido en toneladas por hectárea (tn/ha). Se dispone de 24 rendimientos, con tamaños muestrales desiguales por variedad. El propósito es determinar si el rendimiento medio depende de la variedad y precisar qué diferencias entre pares tienen respaldo estadístico.

Pregunta de investigación: ¿son iguales los rendimientos medios de las cuatro variedades en las condiciones del experimento?

Objetivo general: comparar los rendimientos medios mediante ANOVA de un factor para un diseño completamente aleatorizado (DCA) desbalanceado, examinar sus supuestos y comprobar la sensibilidad de las conclusiones ante varianzas desiguales.

Objetivos específicos: describir las distribuciones por variedad; contrastar la igualdad global de medias; evaluar residuos y dispersión; efectuar comparaciones de Tukey; contrastar los resultados con Welch y comparaciones por pares ajustadas mediante Holm.

2 Datos y diseño experimental

El factor fijo es variedad, con cuatro niveles. La variable respuesta es rendimiento (tn/ha). Se analiza como un DCA de un factor, sujeto a confirmar la unidad experimental y el proceso real de aleatorización. El desbalance (6, 7, 5 y 6 observaciones) no invalida por sí mismo un ANOVA de un factor.

sevillana <- c(4.9, 8.5, 9.4, 9.8, 7.4, 10.5)
ascolana  <- c(8.2, 4.7, 5.6, 6.3, 6.4, 5.2, 4.8)
empeltre  <- c(5.6, 4.9, 5.0, 4.7, 4.3)  # Segundo valor provisional
ilo       <- c(10.6, 9.6, 8.7, 10.3, 9.9, 8.7)

valores <- list(Sevillana=sevillana, Ascolana=ascolana,
                Empeltre=empeltre, Ilo=ilo)
datos <- data.frame(
  variedad=factor(rep(names(valores), lengths(valores)),
                  levels=names(valores)),
  rendimiento=unlist(valores, use.names=FALSE)
)
stopifnot(nrow(datos)==24, !anyNA(datos))
knitr::kable(data.frame(Variedad=names(valores),
                        n=lengths(valores),
                        Mediciones=vapply(valores,
                            function(x) paste(x, collapse=", "), character(1))),
             caption="Datos originales por variedad (tn/ha).")
Datos originales por variedad (tn/ha).
Variedad n Mediciones
Sevillana Sevillana 6 4.9, 8.5, 9.4, 9.8, 7.4, 10.5
Ascolana Ascolana 7 8.2, 4.7, 5.6, 6.3, 6.4, 5.2, 4.8
Empeltre Empeltre 5 5.6, 4.9, 5, 4.7, 4.3
Ilo Ilo 6 10.6, 9.6, 8.7, 10.3, 9.9, 8.7

2.1 Modelo e hipótesis

El modelo es:

\[Y_{ij}=\mu+\tau_i+\varepsilon_{ij},\quad i=1,\ldots,4;\quad j=1,\ldots,n_i.\]

En el ANOVA clásico se supone que los errores son independientes, normalmente distribuidos, con media cero y varianza común \(\sigma^2\). La hipótesis nula es \(H_0:\mu_1=\mu_2=\mu_3=\mu_4\), frente a que al menos una media difiere. Se adopta \(\alpha=0.05\).

3 Análisis exploratorio

resumen <- do.call(rbind, lapply(valores, function(x) {
  data.frame(n=length(x), Media=mean(x), Mediana=median(x),
             DE=sd(x), Varianza=var(x), Min=min(x), Max=max(x),
             CV_porcentaje=100*sd(x)/mean(x))
}))
knitr::kable(round(resumen, 3), caption="Estadísticos descriptivos de rendimiento.")
Estadísticos descriptivos de rendimiento.
n Media Mediana DE Varianza Min Max CV_porcentaje
Sevillana 6 8.417 8.95 2.031 4.126 4.9 10.5 24.133
Ascolana 7 5.886 5.60 1.220 1.488 4.7 8.2 20.726
Empeltre 5 4.900 4.90 0.474 0.225 4.3 5.6 9.680
Ilo 6 9.633 9.75 0.799 0.639 8.7 10.6 8.296

Lectura: Ilo presenta el mayor promedio (9.633 tn/ha), seguida de Sevillana (8.417), Ascolana (5.886) y Empeltre (4.900). Sevillana presenta la mayor dispersión relativa (CV = 24.13%), mientras que Ilo presenta una dispersión observada menor (CV = 8.30%). La variabilidad en este conjunto no equivale a estabilidad varietal entre diferentes ambientes o años.

boxplot(rendimiento ~ variedad, data=datos, col="gray90", border="gray35",
        outline=FALSE, xlab="Variedad", ylab="Rendimiento (tn/ha)",
        main="Rendimiento por variedad")
set.seed(123)
stripchart(rendimiento ~ variedad, data=datos, vertical=TRUE,
           method="jitter", jitter=.1, add=TRUE,
           pch=19, col="navy")
Distribución del rendimiento por variedad, con cada observación individual.

Distribución del rendimiento por variedad, con cada observación individual.

4 ANOVA clásico

modelo <- aov(rendimiento ~ variedad, data=datos)
summary(modelo)
            Df Sum Sq Mean Sq F value   Pr(>F)    
variedad     3   82.9   27.63    16.4 0.000013 ***
Residuals   20   33.7    1.68                     
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
eta2 <- sum((fitted(modelo)-mean(datos$rendimiento))^2) /
        sum((datos$rendimiento-mean(datos$rendimiento))^2)
cat("Eta cuadrado =", round(eta2, 4), "\n")
Eta cuadrado = 0.7113

El ANOVA detecta diferencias entre medias de rendimiento: \(F(3,20)=16.42\), \(p=0.0000128\). Se rechaza la hipótesis de igualdad global. El valor \(\eta^2=0.7113\) indica que alrededor del 71.13% de la variación total observada se asocia con el factor variedad dentro de este modelo, sin constituir una atribución causal universal.

medias <- tapply(datos$rendimiento, datos$variedad, mean)
n <- table(datos$variedad)
CME <- deviance(modelo)/df.residual(modelo)
EE <- sqrt(CME/as.numeric(n))
margen <- qt(.975, df.residual(modelo))*EE
x <- seq_along(medias)
plot(x, medias, ylim=range(medias-margen, medias+margen), xaxt="n",
     pch=19, col="navy", cex=1.3, xlab="Variedad",
     ylab="Rendimiento medio (tn/ha)", main="Medias e IC individuales del 95%")
axis(1, at=x, labels=names(medias))
arrows(x, medias-margen, x, medias+margen,
       angle=90, code=3, length=.07, lwd=1.5)
abline(h=mean(datos$rendimiento), lty=2, col="gray50")
Medias por variedad e intervalos individuales del 95% calculados con el error común del ANOVA. No son intervalos simultáneos de comparaciones por pares.

Medias por variedad e intervalos individuales del 95% calculados con el error común del ANOVA. No son intervalos simultáneos de comparaciones por pares.

5 Diagnóstico de supuestos

shapiro.test(residuals(modelo))

    Shapiro-Wilk normality test

data:  residuals(modelo)
W = 0.947, p-value = 0.23
fligner.test(rendimiento ~ variedad, data=datos)

    Fligner-Killeen test of homogeneity of variances

data:  rendimiento by variedad
Fligner-Killeen:med chi-squared = 5.86, df = 3, p-value = 0.12

Shapiro–Wilk (\(W=0.94672\), \(p=0.2299\)) no encontró evidencia significativa contra la normalidad de los residuos. Fligner–Killeen (\(\chi^2=5.8551\), \(p=0.1189\)) no encontró evidencia significativa contra la homogeneidad de varianzas. No rechazar estas hipótesis no demuestra su cumplimiento. Dado el tamaño muestral y la diferencia de varianzas muestrales entre Sevillana y Empeltre, conviene analizar la sensibilidad.

qqnorm(residuals(modelo), pch=19, col="navy",
       main="Q-Q de residuos")
qqline(residuals(modelo), col="red3", lwd=2)
Gráfico Q-Q de los residuos del ANOVA.

Gráfico Q-Q de los residuos del ANOVA. 

plot(fitted(modelo), residuals(modelo), pch=19, col="navy",
     xlab="Valores ajustados (tn/ha)", ylab="Residuos",
     main="Residuos frente a ajustados")
abline(h=0, col="red3", lty=2)
Residuos frente a valores ajustados del ANOVA.

Residuos frente a valores ajustados del ANOVA. 

boxplot(residuals(modelo) ~ datos$variedad, col="gray90",
        border="gray35", outline=FALSE,
        xlab="Variedad", ylab="Residuos", main="Residuos por variedad")
set.seed(123)
stripchart(residuals(modelo) ~ datos$variedad, vertical=TRUE,
           method="jitter", jitter=.1, add=TRUE, pch=19, col="navy")
abline(h=0, col="red3", lty=2)
Residuos por variedad: inspección visual de su dispersión.

Residuos por variedad: inspección visual de su dispersión.

La independencia debe justificarse mediante el protocolo experimental (asignación aleatoria y unidades independientes); ningún gráfico o test de los anteriores puede verificarla por sí solo.

6 Comparaciones múltiples de Tukey

tukey <- TukeyHSD(modelo, "variedad", conf.level=.95)
knitr::kable(round(as.data.frame(tukey$variedad), 5),
             caption="Diferencias entre medias, intervalos simultáneos del 95% y valores p ajustados (Tukey).")
Diferencias entre medias, intervalos simultáneos del 95% y valores p ajustados (Tukey).
diff lwr upr p adj
Ascolana-Sevillana -2.53095 -4.55080 -0.5111 0.01096
Empeltre-Sevillana -3.51667 -5.71507 -1.3183 0.00122
Ilo-Sevillana 1.21667 -0.87943 3.3128 0.38804
Empeltre-Ascolana -0.98571 -3.11154 1.1401 0.57481
Ilo-Ascolana 3.74762 1.72777 5.7675 0.00024
Ilo-Empeltre 4.73333 2.53493 6.9317 0.00004

Tukey identifica diferencias significativas entre Sevillana y Ascolana; Sevillana y Empeltre; Ilo y Ascolana; e Ilo y Empeltre. No identifica diferencias significativas entre Ilo y Sevillana, ni entre Ascolana y Empeltre.

t <- tukey$variedad
dif <- t[,"diff"]; li <- t[,"lwr"]; ls <- t[,"upr"]
y <- seq_along(dif)
color <- ifelse(t[,"p adj"] < .05, "navy", "gray50")
op <- par(mar=c(5,13,4,2))
plot(dif, y, xlim=range(c(li,ls,0)), ylim=c(.5,length(y)+.5),
     yaxt="n", pch=19, col=color, xlab="Diferencia de medias (tn/ha)",
     ylab="", main="Tukey: intervalos simultáneos del 95%")
axis(2, at=y, labels=rownames(t), las=1, cex.axis=.85)
abline(v=0, lty=2, col="red3")
arrows(li,y,ls,y,angle=90,code=3,length=.05,col=color,lwd=1.6)
Diferencias de medias e intervalos simultáneos del 95% de Tukey. Los intervalos que no cruzan cero representan comparaciones significativas al 5%.

Diferencias de medias e intervalos simultáneos del 95% de Tukey. Los intervalos que no cruzan cero representan comparaciones significativas al 5%. 

par(op)

7 Sensibilidad: Welch y Welch–Holm

El ANOVA de Welch no presupone varianzas iguales entre grupos; los contrastes por pares de Welch con ajuste de Holm tampoco requieren una varianza agrupada común. Welch–Holm no es lo mismo que Games–Howell.

welch <- oneway.test(rendimiento ~ variedad, data=datos,
                     var.equal=FALSE)
holm <- pairwise.t.test(datos$rendimiento, datos$variedad,
                        pool.sd=FALSE, p.adjust.method="holm")
print(welch)

    One-way analysis of means (not assuming equal variances)

data:  rendimiento and variedad
F = 46, num df = 3.0, denom df = 10.7, p-value = 0.0000021
print(holm)

    Pairwise comparisons using t tests with non-pooled SD 

data:  datos$rendimiento and datos$variedad 

         Sevillana Ascolana Empeltre
Ascolana 0.09      -        -       
Empeltre 0.03      0.17     -       
Ilo      0.22      0.000243 0.000009

P value adjustment method: holm

Welch confirma la diferencia global: \(F(3,10.686)=45.999\), \(p=0.000002099\). Los F de los dos enfoques no deben compararse como tamaños de efecto.

pares <- data.frame(
  Par=c("Ilo - Sevillana", "Ilo - Ascolana", "Ilo - Empeltre",
        "Sevillana - Ascolana", "Sevillana - Empeltre",
        "Ascolana - Empeltre"),
  Tukey_p=c(.3880384,.0002396,.0000378,.0109566,.0012178,.5748123),
  Welch_Holm_p=c(.21743,.00024,.0000087,.08600,.02868,.17386)
)
pares$Tukey <- ifelse(pares$Tukey_p<.05,"Significativa","No significativa")
pares$Welch_Holm <- ifelse(pares$Welch_Holm_p<.05,"Significativa","No significativa")
knitr::kable(pares, digits=6,
             caption="Contraste de decisiones de Tukey y Welch–Holm al 5%.")
Contraste de decisiones de Tukey y Welch–Holm al 5%.
Par Tukey_p Welch_Holm_p Tukey Welch_Holm
Ilo - Sevillana 0.388038 0.217430 No significativa No significativa
Ilo - Ascolana 0.000240 0.000240 Significativa Significativa
Ilo - Empeltre 0.000038 0.000009 Significativa Significativa
Sevillana - Ascolana 0.010957 0.086000 Significativa No significativa
Sevillana - Empeltre 0.001218 0.028680 Significativa Significativa
Ascolana - Empeltre 0.574812 0.173860 No significativa No significativa

Resultado clave: cinco de las seis decisiones coinciden. La diferencia Sevillana–Ascolana es significativa con Tukey (\(p=0.01096\)), pero no con Welch–Holm (\(p=0.08600\)). Esta comparación es sensible al método, de modo que no debería presentarse como un hallazgo robusto.

8 Discusión

Ilo obtuvo el mayor rendimiento medio observado, con variación relativa baja, pero no se demuestra una diferencia estadísticamente significativa entre Ilo y Sevillana con ninguno de los dos procedimientos de comparaciones múltiples. Tampoco se demuestra equivalencia entre ambas. La evidencia para Ilo frente a Ascolana y Empeltre, y para Sevillana frente a Empeltre, se mantiene bajo los dos enfoques.

El rendimiento de 4.9 tn/ha registrado para Sevillana contribuye a su mayor dispersión; no procede eliminar ninguna observación sin evidencia externa de error de medición o registro. Para recomendar una variedad de manera generalizable se requerirían datos del protocolo, replicación en diferentes condiciones y, potencialmente, varias campañas.

9 Conclusión final

  1. Se identificaron diferencias globales significativas en rendimiento entre las cuatro variedades (ANOVA clásico y Welch, ambos \(p<0.001\)).
  2. El orden de medias observado fue Ilo (9.633) > Sevillana (8.417) > Ascolana (5.886) > Empeltre (4.900) tn/ha; este orden descriptivo no implica que todos los pares difieran estadísticamente.
  3. Las comparaciones Ilo–Ascolana, Ilo–Empeltre y Sevillana–Empeltre son significativas con ambos procedimientos. Ilo–Sevillana y Ascolana–Empeltre no muestran diferencias significativas.
  4. La comparación Sevillana–Ascolana requiere cautela, pues su decisión cambia entre Tukey y Welch–Holm.
  5. La interpretación definitiva exige verificar el segundo dato de Empeltre (4.9) y documentar la unidad experimental, independencia y asignación aleatoria.

10 Reproducibilidad y publicación

Este documento contiene todos los datos, el código R y las gráficas. No necesita archivos de imágenes externos. Los resultados se producen en el momento de compilar el .Rmd a HTML con rmarkdown::render() o el botón Knit en RStudio.

Programas empleados: R, rmarkdown y knitr. Los análisis y gráficos utilizan funciones incluidas en R base. Para RPubs, publique el HTML generado, no directamente el archivo .Rmd.