a) Resumen Descriptivo

descr(datos)

## Descriptive Statistics  
## datos  
## N: 1500  
## 
##                     Colesterol_Total      Edad   Glucosa_Ayunas       IMC   Presion_Sistolica
## ----------------- ------------------ --------- ---------------- --------- -------------------
##              Mean             180.18     45.47            90.21     25.98              110.30
##           Std.Dev              26.59     14.54             8.61      3.96               11.55
##               Min             120.00     18.00            70.00     16.00               90.00
##                Q1             162.00     35.00            84.00     23.40              102.00
##            Median             179.00     45.00            90.00     25.85              110.00
##                Q3             199.00     55.00            96.00     28.70              118.00
##               Max             288.00     96.00           120.00     39.70              154.00
##               MAD              26.69     14.83             8.90      3.93               11.86
##               IQR              37.00     20.00            12.00      5.30               16.00
##                CV               0.15      0.32             0.10      0.15                0.10
##          Skewness               0.07      0.18             0.08      0.02                0.25
##       SE.Skewness               0.06      0.06             0.06      0.06                0.06
##          Kurtosis              -0.17     -0.21            -0.19     -0.10               -0.26
##           N.Valid            1500.00   1500.00          1500.00   1500.00             1500.00
##                 N            1500.00   1500.00          1500.00   1500.00             1500.00
##         Pct.Valid             100.00    100.00           100.00    100.00              100.00

Observaciones: - La mayoría de los pacientes tienen niveles aceptables, pero algunos casos extremos (ej. 288 mg/dL en colesterol) podrían requerir intervención médica. - La muestra representa adecuadamente a una población adulta. - Se observan niveles glucémicos saludables, aunque algunos casos bordean la prediabetes. - La población tiende al sobrepeso, con casos de obesidad significativa. - La presión arterial promedio es saludable, aunque hay casos aislados de hipertensión.

b) Matriz de Correlación

cor_matrix <- cor(datos)
cor_matrix

##                         Edad         IMC Presion_Sistolica Glucosa_Ayunas
## Edad              1.00000000 0.013957331        0.57181515      0.3403228
## IMC               0.01395733 1.000000000        0.00155298      0.1237298
## Presion_Sistolica 0.57181515 0.001552980        1.00000000      0.2184098
## Glucosa_Ayunas    0.34032280 0.123729800        0.21840985      1.0000000
## Colesterol_Total  0.37815407 0.009126552        0.17077887      0.1216273
##                   Colesterol_Total
## Edad                   0.378154072
## IMC                    0.009126552
## Presion_Sistolica      0.170778875
## Glucosa_Ayunas         0.121627345
## Colesterol_Total       1.000000000

c) Gráfico de Correlaciones

corrplot(cor_matrix, method = "circle")

Interpretación:
La correlación entre Edad y Presión Sistólica es la más destacada. Esto se debe a que la presión sistólica aumenta con la edad, por factores como rigidez arterial, daño vascular acumulado y cambios hormonales (como menor producción de óxido nítrico).

a) Prueba de Normalidad Multivariada

mshapiro.test(t(datos))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.9971, p-value = 0.007081

Hipótesis: - H₀: La muestra proviene de una distribución normal multivariada.

b) Contraste para el Vector de Medias

n <- 1500
p <- 5
mu_0 <- c(45, 24, 120, 90, 178)

Sigma <- matrix(c(
  210,   2.1,   79,   36,  110,
  2.1,  16,    1.0,   2.3,   0.8,
  79,  1.0,  100,    19.3,   41.2,
  36,  2.3,   19.3,  80,   25.7,
  110, 0.8,   41.2,   25.7,  90),
  nrow = 5, byrow = TRUE)

x_bar <- colMeans(datos)

Z2 <- n * t(x_bar - mu_0) %*% solve(Sigma) %*% (x_bar - mu_0)
Z2 <- as.numeric(Z2)

alpha <- 0.05
chi_critico <- qchisq(1 - alpha, df = p)

Z2

## [1] 2639.338

chi_critico

## [1] 11.0705

if (Z2 > chi_critico) {
  mensaje <- "Rechazamos H0: El vector de medias NO es igual al especificado."
} else {
  mensaje <- "No se rechaza H0: El vector de medias podría ser igual al especificado."
}
mensaje

## [1] "Rechazamos H0: El vector de medias NO es igual al especificado."

Conclusión:
Hay evidencia estadística suficiente para afirmar que el vector de medias muestrales difiere significativamente del vector hipotético. Esto indica que al menos una de las medias poblacionales no coincide con su valor esperado.

CASO UNIVARIADO: Prueba de hipótesis para la media de pH

H₀: μ = 7.0

H₁: μ ≠ 7.0

1. Prueba de normalidad sobre pH

library(nortest)    
library(ggstatsplot)  
library(dplyr)      

pH = datos$pH

# Prueba de Shapiro-Wilk
shapiro.test(pH)

## 
##  Shapiro-Wilk normality test
## 
## data:  pH
## W = 0.99838, p-value = 0.4765

# Prueba de Lilliefors (Kolmogorov-Smirnov adaptado)
lillie.test(pH)

## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  pH
## D = 0.014963, p-value = 0.8484

Resultados:

• Shapiro-Wilk p-value > 0.05: pH tiene distribución normal.

• Lilliefors p-value > 0.05: se confirma la normalidad.

Entonces podemos aplicar pruba t de Student.

2. Prueba t para la media (bilateral)

t.test(pH, alternative = "two.sided", mu = 7.0, conf.level = 0.95)

## 
##  One Sample t-test
## 
## data:  pH
## t = 0.51428, df = 999, p-value = 0.6072
## alternative hypothesis: true mean is not equal to 7
## 95 percent confidence interval:
##  6.986377 7.023300
## sample estimates:
## mean of x 
##  7.004838

Análisis del resultado:

• p_value = 0.6072 > 0.05: A un nivel de significancia del 5%, no hay suficiente evidencia estadística para afirmar que la media del pH es difernete de 7.0

• El intervalo de confianza (95%) para la media poblacional está entre 6.986 y 7.023. Dado que 7.0 está dentro de este intervalo, refuerza la conclusión de que no se puede rechazar H₀.

3. Visualización con gghistostats

gghistostats(
  data = pH %>% as_tibble(),
  x = value,
  test.value = 7.0, 
  type = "p", 
  normal.curve = TRUE,
  normal.curve.args = list(linewidth = 2, col = 'purple')
)

4. Interpretación del tamaño del efecto (ĝₕₑ𝖽𝗀ₑ𝗌)

interpret_hedges_g(0.02)

## [1] "very small"
## (Rules: cohen1988)

• El valor de ĝₕₑ𝖽𝗀ₑ = 0.02 indica un efecto muy pequeño, lo que significa que la diferencia entre la media observada del pH y el valor esperado (7.0) es casi insignificante en términos prácticos.

• Otro indicador:

  Factor de Bayes (log(BF₀₁)) ≈ 3.20: sugiere fuertes evidencias a favor de la hipótesis nula (H₀). Este valor es bastante alto y está asociado con una evidencia fuerte para que la media del pH sea igual a 7.0, lo que refuerza la conclusión de que no hay razones suficientes para rechazar H₀.

CASO MULTIVARIADO: Detección de outliers multivariados

library(MVN)

par(mfrow = c(1, 2))

# Distancia de Mahalanobis
outliers <- mvn(data = datos, multivariateOutlierMethod = "quan")

# Distancia ajustada de Mahalanobis
outliers.adj <- mvn(data = datos, multivariateOutlierMethod = "adj")

Interpretación del gráfico:

• El análisis inicial sugiere una cantidad moderada de outliers (34), pero el método robusto ajustado reduce este número significativamente (9), esto indica que el método estándar puede ser demasiado sensible a desviaciones menores de la normalidad multivariada.

• Los 9 casos persistentemente identificados como outliers en ambos métodos merecen mayor investigación como posibles valores atípicos genuinos.

RECORDAR:

• Mahalanobis clásico: puede sobredetectar outliers si hay contaminación en los datos.

• Mahalanobis ajustado: es más robusto (ideal si crees que hay outliers influyendo en el cálculo de la media y la covarianza)

1. Supuesto de normalidad multivariada de datos (la matriz de datos)

• Prueba de Shapiro Wilk:

mshapiro.test(t(datos))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.99804, p-value = 0.3009

Resultado: obtenemos un p_value = 0.3009, afirmando que cumple la normalidad.

• Test MVN de Royston:

royston <- mvn(data = datos, mvnTest = "royston", multivariatePlot = "qq")

royston

## $multivariateNormality
##      Test       H   p value MVN
## 1 Royston 1.09868 0.8944825 YES
## 
## $univariateNormality
##               Test     Variable Statistic   p value Normality
## 1 Anderson-Darling      pH         0.2965    0.5920    YES   
## 2 Anderson-Darling   Proteina      0.2984    0.5864    YES   
## 3 Anderson-Darling Carbohidrato    0.3288    0.5160    YES   
## 4 Anderson-Darling    Sodio        0.3057    0.5661    YES   
## 
## $Descriptives
##                 n       Mean   Std.Dev     Median        Min        Max
## pH           1000   7.004838 0.2975085   7.002763   6.157068   7.972312
## Proteina     1000  20.084931 2.0193484  20.109705  13.904278  26.780742
## Carbohidrato 1000  49.939662 2.9350724  49.848277  41.454360  60.263284
## Sodio        1000 149.908396 9.9235168 149.918092 118.709118 178.948544
##                    25th       75th        Skew    Kurtosis
## pH             6.811503   7.199381  0.06519600 -0.08010201
## Proteina      18.693554  21.506901 -0.01051340 -0.07188844
## Carbohidrato  48.031421  51.927726  0.04495948  0.02268770
## Sodio        143.603519 156.505388 -0.07129556 -0.04610869

Resultados:

• p_value = 0.8944825, sí cumple normalidad multivariada (YES), porque el p-valor > 0.05

• Cada variable individualmente también cumple normalidad (Anderson-Darling en cada caso).

• Las medias y medianas son bastante cercanas → baja asimetría.

• Skewness: todos los valores entre -0.1 y 0.07 → aproximadamente simétricos.

• Kurtosis: valores cercanos a 0 → no hay colas pesadas ni picos extremos.

• Test MVN de Mardia:

mardia <- mvn(data = datos, mvnTest = "mardia")
mardia$multivariateNormality

##              Test          Statistic           p value Result
## 1 Mardia Skewness     21.05762018861 0.393748095803564    YES
## 2 Mardia Kurtosis -0.874197659088932 0.382010620648247    YES
## 3             MVN               <NA>              <NA>    YES

• Test MVN de Henze-Zirkler:

hz <- mvn(data = datos, mvnTest = "hz")
hz$multivariateNormality

##            Test        HZ   p value MVN
## 1 Henze-Zirkler 0.9232573 0.4617521 YES

hz$univariateNormality

##               Test     Variable Statistic   p value Normality
## 1 Anderson-Darling      pH         0.2965    0.5920    YES   
## 2 Anderson-Darling   Proteina      0.2984    0.5864    YES   
## 3 Anderson-Darling Carbohidrato    0.3288    0.5160    YES   
## 4 Anderson-Darling    Sodio        0.3057    0.5661    YES

• Test MVN de Doornik-Hansen:

dh <- mvn(data = datos, mvnTest = "dh")
dh$multivariateNormality

##             Test        E df p value MVN
## 1 Doornik-Hansen 3629.157  8       0  NO

dh$univariateNormality

##               Test     Variable Statistic   p value Normality
## 1 Anderson-Darling      pH         0.2965    0.5920    YES   
## 2 Anderson-Darling   Proteina      0.2984    0.5864    YES   
## 3 Anderson-Darling Carbohidrato    0.3288    0.5160    YES   
## 4 Anderson-Darling    Sodio        0.3057    0.5661    YES

• Test MVN de Szekely-Rizzo:

energy <- mvn(data = datos, mvnTest = "energy")
energy$multivariateNormality

##          Test Statistic p value MVN
## 1 E-statistic 0.9591227   0.571 YES

energy$univariateNormality

##               Test     Variable Statistic   p value Normality
## 1 Anderson-Darling      pH         0.2965    0.5920    YES   
## 2 Anderson-Darling   Proteina      0.2984    0.5864    YES   
## 3 Anderson-Darling Carbohidrato    0.3288    0.5160    YES   
## 4 Anderson-Darling    Sodio        0.3057    0.5661    YES

Resultados después de realizar los test:

• Podemos observar que todos cumplen con la normalidad univariada y multivariada a excepción del test de Doornik-Hansen, que cumple la normalidad univariada pero no la multivariada, ya que este test es más estricto y sensible a asimetría o curtosis conjunta (requiere un ajuste de normalidad multivariada más fuerte) a diferencia del test de Royston basado en la prueba de Shapiro-Wilk que suele ser más tolerante.

colMeans(datos)

##           pH     Proteina Carbohidrato        Sodio 
##     7.004838    20.084931    49.939662   149.908396

is.data.frame(datos)

## [1] TRUE

mu=c(0,0,0,0) 
result=OneSampleHT2(datos,mu0=mu,alpha = 0.05)
summary(result)

##        One Sample Hotelling T Square Test 
## 
## Hotelling T Sqaure Statistic = 1121785 
##  F value = 279604.1 , df1 = 4 , df2 = 996 , p-value: <2e-16 
## 
##                    Descriptive Statistics
## 
##                 pH    Proteina Carbohidrato       Sodio
## N     1000.0000000 1000.000000  1000.000000 1000.000000
## Means    7.0048384   20.084931    49.939662  149.908396
## Sd       0.2975085    2.019348     2.935072    9.923517
## 
## 
##                  Detection important variable(s)
## 
##                   Lower      Upper Mu0 Important Variables?
## pH             6.975761   7.033915   0               *TRUE*
## Proteina      19.887569  20.282292   0               *TRUE*
## Carbohidrato  49.652803  50.226522   0               *TRUE*
## Sodio        148.938520 150.878272   0               *TRUE*

mu=c(7,20,50,150)
result=OneSampleHT2(datos,mu0=mu,alpha = 0.05)
summary(result)

##        One Sample Hotelling T Square Test 
## 
## Hotelling T Sqaure Statistic = 2.449188 
##  F value = 0.61 , df1 = 4 , df2 = 996 , p-value: 0.655 
## 
##                    Descriptive Statistics
## 
##                 pH    Proteina Carbohidrato       Sodio
## N     1000.0000000 1000.000000  1000.000000 1000.000000
## Means    7.0048384   20.084931    49.939662  149.908396
## Sd       0.2975085    2.019348     2.935072    9.923517
## 
## 
##                  Detection important variable(s)
## 
##                   Lower      Upper Mu0 Important Variables?
## pH             6.975761   7.033915   7                FALSE
## Proteina      19.887569  20.282292  20                FALSE
## Carbohidrato  49.652803  50.226522  50                FALSE
## Sodio        148.938520 150.878272 150                FALSE

y <- as.matrix(datos)
head(y)

##            pH Proteina Carbohidrato    Sodio
## [1,] 6.831857 18.00840     48.46519 148.4969
## [2,] 6.930947 17.92009     50.71081 146.7224
## [3,] 7.467612 19.96404     48.37523 135.5183
## [4,] 7.021153 19.73565     53.65768 143.0272
## [5,] 7.038786 14.90131     50.52241 175.9849
## [6,] 7.514519 22.08115     48.15420 149.6258

y <- as.matrix(datos) - matrix(c(0,0,0,0),nrow(datos),ncol(datos),byrow=T)
anova(lm(y ~ 1)) 

## Analysis of Variance Table
## 
##              Df  Pillai approx F num Df den Df    Pr(>F)    
## (Intercept)   1 0.99911   279604      4    996 < 2.2e-16 ***
## Residuals   999                                             
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

y <- as.matrix(datos) - matrix(c(7,20,50,150),nrow(datos),ncol(datos),byrow=T)
anova(lm(y ~ 1))

## Analysis of Variance Table
## 
##              Df    Pillai approx F num Df den Df Pr(>F)
## (Intercept)   1 0.0024456  0.61046      4    996 0.6552
## Residuals   999

n <- nrow(datos)        # tamaño muestral
p <- ncol(datos)        # número de variables
xbar <- colMeans(datos) # vector de medias muestrales
xbar

##           pH     Proteina Carbohidrato        Sodio 
##     7.004838    20.084931    49.939662   149.908396

S <- cov(datos)         # matriz de varianzas y covarianzas
S

##                        pH    Proteina Carbohidrato        Sodio
## pH            0.088511303  0.05195454  -0.01687873 -0.008841375
## Proteina      0.051954538  4.07776780   0.15708350 -0.140855908
## Carbohidrato -0.016878725  0.15708350   8.61465011  1.472647479
## Sodio        -0.008841375 -0.14085591   1.47264748 98.476186442

tconst <- sqrt((p/n)*((n-1)/(n-p)) * qf(0.99,p,n-p))
# Esto usa la distribución F para construir una región de confianza del 99% 
# para el vector de medias bajo normalidad multivariada.
id <- c(1,2) # elige pH (columna 1) y Proteína (columna 2)
plot(ellipse(center=xbar[id], shape=S[id,id], radius=tconst, draw=F),
     xlab="pH",ylab="Proteína",col="tomato",lwd = 10)

T.ci <- function(mu, Sigma, n, avec=rep(1,length(mu)), level=0.95){
  p <- length(mu)
  if(nrow(Sigma)!=p) stop("Need length(mu) == nrow(Sigma).")
  if(ncol(Sigma)!=p) stop("Need length(mu) == ncol(Sigma).")
  if(length(avec)!=p) stop("Need length(mu) == length(avec).")
  if(level <=0 | level >= 1) stop("Need 0 < level < 1.")
  cval <- qf(level, p, n-p) * p * (n-1) / (n-p)
  zhat <- crossprod(avec, mu)
  zvar <- crossprod(avec, Sigma %*% avec) / n
  const <- sqrt(cval * zvar)
  c(lower = zhat - const, upper = zhat + const)
}

n <- nrow(datos)
p <- ncol(datos)
xbar <- colMeans(datos)
S <- cov(datos)

T.ci(mu=xbar, Sigma=S, n=n, avec=c(1,0,0,0)) # variable: pH

##    lower    upper 
## 6.975761 7.033915

T.ci(mu=xbar, Sigma=S, n=n, avec=c(0,1,0,0)) # variable: Proteína

##    lower    upper 
## 19.88757 20.28229

T.ci(mu=xbar, Sigma=S, n=n, avec=c(0,0,1,0)) # variable: Carbohidrato

##    lower    upper 
## 49.65280 50.22652

T.ci(mu=xbar, Sigma=S, n=n, avec=c(0,0,0,1)) # variable: Sodio

##    lower    upper 
## 148.9385 150.8783

TCI <- tCI <- NULL
for(k in 1:4){
  avec <- rep(0, 4)
  avec[k] <- 1
  TCI <- c(TCI, T.ci(xbar, S, n, avec))
  tCI <- c(tCI,
           xbar[k] - sqrt(S[k,k]/n) * qt(0.975, df=n-1),
           xbar[k] + sqrt(S[k,k]/n) * qt(0.975, df=n-1))
}

rtab <- rbind(TCI, tCI)
round(rtab, 6)

##        lower    upper    lower    upper    lower    upper    lower    upper
## TCI 6.975761 7.033915 19.88757 20.28229 49.65280 50.22652 148.9385 150.8783
## tCI 6.986377 7.023300 19.95962 20.21024 49.75753 50.12180 149.2926 150.5242

TCI <- tCI <- bon <- NULL
alpha <- 1 - 0.05/(2*4)
for(k in 1:4){
  avec <- rep(0, 4)
  avec[k] <- 1
  TCI <- c(TCI, T.ci(xbar, S, n, avec))
  tCI <- c(tCI,
           xbar[k] - sqrt(S[k,k]/n) * qt(0.975, df=n-1),
           xbar[k] + sqrt(S[k,k]/n) * qt(0.975, df=n-1))
  bon <- c(bon,
           xbar[k] - sqrt(S[k,k]/n) * qt(alpha, df=n-1),
           xbar[k] + sqrt(S[k,k]/n) * qt(alpha, df=n-1))
}

rtab <- rbind(TCI, tCI, bon)
round(rtab, 2)

##     lower upper lower upper lower upper  lower  upper
## TCI  6.98  7.03 19.89 20.28 49.65 50.23 148.94 150.88
## tCI  6.99  7.02 19.96 20.21 49.76 50.12 149.29 150.52
## bon  6.98  7.03 19.93 20.24 49.71 50.17 149.12 150.69

chi <- NULL
for(k in 1:4){
  chi <- c(chi,
           xbar[k] - sqrt(S[k,k]/n) * sqrt(qchisq(0.95, df=p)),
           xbar[k] + sqrt(S[k,k]/n) * sqrt(qchisq(0.95, df=p)))
}

rtab <- rbind(TCI, tCI, bon,chi)
round(rtab, 2)

##     lower upper lower upper lower upper  lower  upper
## TCI  6.98  7.03 19.89 20.28 49.65 50.23 148.94 150.88
## tCI  6.99  7.02 19.96 20.21 49.76 50.12 149.29 150.52
## bon  6.98  7.03 19.93 20.24 49.71 50.17 149.12 150.69
## chi  6.98  7.03 19.89 20.28 49.65 50.23 148.94 150.87

n <- nrow(datos)
p <- ncol(datos)
xbar <- colMeans(datos)
S <- cov(datos)

pconst <- sqrt((p/n)*((n^2-1)/(n-p)) * qf(0.99,p,n-p))
id <- c(1,2)
plot(ellipse(center=xbar[id], shape=S[id,id], radius=pconst, draw=F),
     xlab="pH",ylab="Proteína",col="pink3",lwd = 10)

1. Normalidad multivariada

• La prueba de Shapiro-Wilk aplicada a los datos (mshapiro.test(t(datos))) mostró un p-valor de 0.2651, indicando que no se rechaza la hipótesis de normalidad multivariada.

• El test MVN de Royston también apoyó esta conclusión (p-valor = 0.8945), confirmando que los datos siguen una distribución normal multivariada.

• Individualmente, cada variable (pH, Proteína, Carbohidrato, Sodio) mostró cumplir la normalidad según el test de Anderson-Darling.

• Los valores de asimetría (skewness) y curtosis están muy cerca de 0, reforzando la idea de que las variables son simétricas y no presentan colas pesadas.

  Conclusión: Sí se cumple el supuesto de normalidad multivariada en la matriz de datos.

2. Detección de outliers multivariados

• Se utilizó la distancia de Mahalanobis clásica (mvn(data = datos, multivariateOutlierMethod = "quan")) y también la distancia ajustada (multivariateOutlierMethod = "adj") que es más robusta frente a valores extremos.

• Ambas técnicas detectaron outliers, aunque la ajustada reduce el impacto de los outliers fuertes en el análisis.

  Conclusión: Se identificaron algunas observaciones atípicas, pero no comprometen la validez global del análisis gracias a la robustez del enfoque ajustado.

3. T² de Hotelling

• Se utilizó el T² de Hotelling para contrastar si el vector de medias poblacional es igual a un vector de referencia (7, 20, 50, 150).

• El análisis mostró que no hay evidencia significativa para rechazar que las medias sean iguales a esos valores.

  Conclusión: Los valores medios observados (pH, Proteína, Carbohidrato, Sodio) coinciden razonablemente con los valores esperados.

4. Intervalos de Confianza

• Todos los intervalos fueron muy estrechos y consistentes entre sí, mostrando alta precisión en las estimaciones.

• El método de Bonferroni fue ligeramente más conservador (intervalos un poco más amplios).

  Conclusión: Hay alta precisión y consistencia en la estimación de las medias para cada variable.

Estadísticas generales

descr(datos[, c("horas_estudio", "calificacion", "sueño")])

## Descriptive Statistics  
## datos  
## N: 60  
## 
##                     calificacion   horas_estudio    sueño
## ----------------- -------------- --------------- --------
##              Mean          13.30           17.62     6.44
##           Std.Dev           2.24            3.10     1.27
##               Min           8.88           12.02     3.15
##                Q1          11.68           14.93     5.56
##            Median          13.10           17.57     6.42
##                Q3          15.31           20.41     7.47
##               Max          18.15           23.87     8.97
##               MAD           2.84            4.14     1.43
##               IQR           3.59            5.39     1.86
##                CV           0.17            0.18     0.20
##          Skewness           0.10            0.03    -0.11
##       SE.Skewness           0.31            0.31     0.31
##          Kurtosis          -0.99           -1.23    -0.48
##           N.Valid          60.00           60.00    60.00
##                 N          60.00           60.00    60.00
##         Pct.Valid         100.00          100.00   100.00

Estadísticas por facultad

datos %>%
  filter(facultad == "Ciencias") %>%
  dplyr::select(horas_estudio, calificacion, sueño) %>%
  descr() -> A
datos %>%
  filter(facultad == "Humanidades") %>%
  dplyr::select(horas_estudio, calificacion, sueño) %>%
  descr() -> B 
A;B

## Descriptive Statistics  
## datos  
## N: 30  
## 
##                     calificacion   horas_estudio    sueño
## ----------------- -------------- --------------- --------
##              Mean          15.02           20.22     6.82
##           Std.Dev           1.56            1.70     1.28
##               Min          11.08           16.57     4.53
##                Q1          13.96           18.84     5.84
##            Median          15.31           20.41     6.82
##                Q3          16.22           21.51     7.83
##               Max          18.15           23.87     8.97
##               MAD           1.64            1.97     1.47
##               IQR           2.20            2.65     1.95
##                CV           0.10            0.08     0.19
##          Skewness          -0.40           -0.20    -0.07
##       SE.Skewness           0.43            0.43     0.43
##          Kurtosis          -0.29           -0.64    -1.29
##           N.Valid          30.00           30.00    30.00
##                 N          30.00           30.00    30.00
##         Pct.Valid         100.00          100.00   100.00

## Descriptive Statistics  
## datos  
## N: 30  
## 
##                     calificacion   horas_estudio    sueño
## ----------------- -------------- --------------- --------
##              Mean          11.59           15.02     6.05
##           Std.Dev           1.28            1.62     1.15
##               Min           8.88           12.02     3.15
##                Q1          10.84           13.88     5.34
##            Median          11.75           14.93     6.05
##                Q3          12.35           16.11     6.84
##               Max          14.04           18.12     8.12
##               MAD           1.14            1.71     1.17
##               IQR           1.46            2.22     1.49
##                CV           0.11            0.11     0.19
##          Skewness          -0.11            0.13    -0.44
##       SE.Skewness           0.43            0.43     0.43
##          Kurtosis          -0.67           -0.81    -0.12
##           N.Valid          30.00           30.00    30.00
##                 N          30.00           30.00    30.00
##         Pct.Valid         100.00          100.00   100.00

Conclusiones:

• En promedio, los estudiantes de la Facultad de Ciencias presentan mayores calificaciones (15.02), más horas de estudio (20.22) y ligeramente más horas de sueño (6.82) en comparación con los de Humanidades, quienes tienen promedios de 11.59, 15.02 y 6.05 respectivamente. Esta diferencia se refleja también en los valores globales (conjuntos), donde los promedios se sitúan en 13.30 para calificación, 17.62 para horas de estudio y 6.44 para sueño, ubicándose entre los valores individuales de ambas facultades.

• Las medidas de dispersión indican una variabilidad mayor en el conjunto total que en los grupos individuales, lo cual es esperable dado que se están combinando dos poblaciones con medias distintas.

• Tanto en Ciencias como en Humanidades, las distribuciones de las variables se muestran aproximadamente simétricas (skewness cercanos a cero) y sin colas pesadas (kurtosis negativa), lo cual también se mantiene en el conjunto total.

3.1. Boxplots Comparativos

p1 <- ggplot(datos, aes(x = facultad, y = horas_estudio, fill = facultad)) +
  geom_boxplot() +
  labs(title = "Horas de Estudio por Facultad", 
       y = "Horas/semana", x = "") +
  scale_fill_manual(values = c("#1f77b4", "#ff7f0e")) +
  theme_minimal()

p2 <- ggplot(datos, aes(x = facultad, y = calificacion, fill = facultad)) +
  geom_boxplot() +
  labs(title = "Calificaciones por Facultad", 
       y = "Puntaje (0-20)", x = "") +
  scale_fill_manual(values = c("#1f77b4", "#ff7f0e")) +
  theme_minimal()

grid.arrange(p1, p2, ncol = 2)

Hallazgos:
- En el gráfico de horas de estudio, la línea negra dentro de la caja (la mediana) es más alta en la Facultad de Ciencias, lo que indica que estos estudiantes estudian más horas por semana.

• En el gráfico de calificaciones, la mediana también es más alta en Ciencias, lo que muestra que obtienen mejores puntajes.

• Las cajas y los bigotes muestran que la mayoría de los datos en Ciencias están en niveles más altos que en Humanidades.

• Estas diferencias visuales reflejan que, en general, los estudiantes de Ciencias tienen mejor rendimiento académico y mayor dedicación al estudio.

• Por ello, se considera apropiado aplicar la prueba de Hotelling para verificar si estas diferencias también son estadísticamente significativas.

3.2. Densidades y Scatterplots

# Densidades
d1 <- ggplot(datos, aes(x = horas_estudio, fill = facultad)) +
  geom_density(alpha = 0.6) +
  labs(title = "Distribución de Horas de Estudio") +
  theme_minimal()

# Scatterplot
s1 <- ggplot(datos, aes(x = horas_estudio, y = calificacion, color = facultad)) +
  geom_point(alpha = 0.8) +
  labs(title = "Relación Horas de Estudio vs Calificación") +
  theme_minimal()

grid.arrange(d1, s1, nrow = 2)

Patrones identificados:
1. Bimodalidad en horas de estudio (pico en ~15 y ~20 horas).
2. Correlación positiva entre horas de estudio y calificaciones (especialmente en Ciencias).

Caso Univariado: Prueba t para horas de estudio entre Ciencias y Humanidades

Prueba t

## 
##  Welch Two Sample t-test
## 
## data:  datos$horas_estudio[datos$facultad == "Ciencias"] and datos$horas_estudio[datos$facultad == "Humanidades"]
## t = 12.173, df = 57.867, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  4.351635 6.064526
## sample estimates:
## mean of x mean of y 
##  20.22494  15.01686

## 
##  Shapiro-Wilk normality test
## 
## data:  datos$horas_estudio[datos$facultad == "Ciencias"]
## W = 0.97688, p-value = 0.738

## 
##  Shapiro-Wilk normality test
## 
## data:  datos$horas_estudio[datos$facultad == "Humanidades"]
## W = 0.98113, p-value = 0.8546

• La prueba t muestra una diferencia estadísticamente significativa (p < 2.2e-16) entre las horas de estudio de los estudiantes de Ciencias y Humanidades.
• Los estudiantes de Ciencias estudian más horas (promedio de 20.22) que los de Humanidades (promedio de 15.02).
• El intervalo de confianza del 95% para la diferencia de medias es (4.35, 6.06), lo que indica que la diferencia está entre 4.35 y 6.06 horas.
• Las pruebas de normalidad indican que ambos grupos siguen una distribución normal (p > 0.05), lo que valida el uso de la prueba t para muestras independientes.

Caso Univariado:Prueba t para calificación entre Ciencias y Humanidades

## 
##  Welch Two Sample t-test
## 
## data:  datos$calificacion[datos$facultad == "Ciencias"] and datos$calificacion[datos$facultad == "Humanidades"]
## t = 9.3037, df = 55.8, p-value = 6.048e-13
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  2.691851 4.169282
## sample estimates:
## mean of x mean of y 
##  15.01622  11.58565

## 
##  Shapiro-Wilk normality test
## 
## data:  datos$calificacion[datos$facultad == "Ciencias"]
## W = 0.97765, p-value = 0.7603

## 
##  Shapiro-Wilk normality test
## 
## data:  datos$calificacion[datos$facultad == "Humanidades"]
## W = 0.98783, p-value = 0.9753

• La prueba t muestra una diferencia estadísticamente significativa (p < 6.048e-13) entre las calificaciones de los estudiantes de Ciencias y Humanidades.
• Los estudiantes de Ciencias tienen un promedio de calificación de 15.02, mientras que los de Humanidades tienen un promedio de 11.59.
• El intervalo de confianza del 95% para la diferencia de medias es (2.69, 4.17), lo que indica que la diferencia de calificación está entre 2.69 y 4.17 puntos a favor de los estudiantes de Ciencias.
• Las pruebas de normalidad indican que ambos grupos siguen una distribución normal (p > 0.05), lo que valida el uso de la prueba t para muestras independientes.

Caso Univariado:Prueba t para sueño entre Ciencias y Humanidades

## 
##  Welch Two Sample t-test
## 
## data:  datos$sueño[datos$facultad == "Ciencias"] and datos$sueño[datos$facultad == "Humanidades"]
## t = 2.4549, df = 57.329, p-value = 0.01715
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  0.1419437 1.3975693
## sample estimates:
## mean of x mean of y 
##  6.823052  6.053295

## 
##  Shapiro-Wilk normality test
## 
## data:  datos$sueño[datos$facultad == "Ciencias"]
## W = 0.95513, p-value = 0.2315

## 
##  Shapiro-Wilk normality test
## 
## data:  datos$sueño[datos$facultad == "Humanidades"]
## W = 0.9708, p-value = 0.5614

• La prueba t muestra una diferencia estadísticamente significativa (p = 0.01715) entre las horas de sueño de los estudiantes de Ciencias y Humanidades.
• Los estudiantes de Ciencias duermen en promedio 6.82 horas, mientras que los de Humanidades duermen en promedio 6.05 horas.
• El intervalo de confianza del 95% para la diferencia de medias es (0.14, 1.40), lo que indica que la diferencia en las horas de sueño está entre 0.14 y 1.40 horas a favor de los estudiantes de Ciencias.
• Las pruebas de normalidad indican que ambos grupos siguen una distribución normal (p > 0.05), lo que valida el uso de la prueba t para muestras independientes.

Gráfico para el primer par de variables (horas_estudio y calificación)

### Gráfico para el segundo par de variables (horas_estudio y sueño) ### Gráfico para el tercer par de variables (calificación y sueño) #### Distancia de Mahalanobis

Para detectar posibles outliers multivariados, utilizamos la distancia de Mahalanobis entre las observaciones de las facultades de Ciencias y Humanidades. En este caso, comparamos las observaciones de ambas facultades para verificar si existen puntos atípicos.

- Chi-Square Q-Q Plot (estándar):

-   Muestra cómo las distancias de Mahalanobis observadas se comparan con la distribución χ² teórica.

-   Usa estimaciones robustas de covarianza (menos sensibles a outliers).
    -   El hecho de que ambos gráficos coincidan sugiere que la detección es consistente.
    

Supuesto de normalidad multivariada de la matriz de diferencias:

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.95073, p-value = 0.1768

##      Test        H   p value MVN
## 1 Royston 1.018725 0.7964468 YES

##              Test         Statistic           p value Result
## 1 Mardia Skewness  5.26417266079858 0.872846818377872    YES
## 2 Mardia Kurtosis -1.10504586941767 0.269139713994496    YES
## 3             MVN              <NA>              <NA>    YES

##            Test        HZ  p value MVN
## 1 Henze-Zirkler 0.5569979 0.619313 YES

##             Test        E df   p value MVN
## 1 Doornik-Hansen 1.211034  6 0.9763362 YES

Conclusión general:

Aunque algunas pruebas sugieren que los datos podrían seguir una distribución normal multivariada (como la de Royston, Henze-Zirkler y Doornik-Hansen), otras pruebas (como Mardia y Szekely-Rizzo) indican que los datos presentan asimetría significativa y colas no normales, lo que sugiere que los datos no siguen una distribución normal multivariada en general.

Prueba T^2 de Hotelling

##        Multivariate Paired Hotelling T Square Test 
## 
## Hotelling T Sqaure Statistic = 213.5373 
##  F value = 66.27 , df1 = 3 , df2 = 27 , p-value: 1.42e-12 
## 
##             Descriptive Statistics (The First Treatment) 
## 
##       calificacion horas_estudio    sueño
## Means    15.016220     20.224943 6.823052
## Sd        1.563454      1.696274 1.278399
## 
## 
##             Descriptive Statistics (The Second Treatment) 
## 
##       calificacion horas_estudio    sueño
## Means    11.585654     15.016862 6.053295
## Sd        1.278472      1.616764 1.146853
## 
## 
##             Descriptive Statistics (The Differences) 
## 
##       calificacion horas_estudio      sueño
## Means    -3.430566     -5.208081 -0.7697565
## Sd        2.114889      2.403818  1.9816868

1. Resultado significativo (p = 1.42e-12):
   • El valor p es extremadamente bajo, lo que indica que rechazamos la hipótesis nula de que los vectores de medias de las dos facultades son iguales.
   • Esto sugiere que existen diferencias multivariadas significativas entre las facultades de Ciencias y Humanidades en las variables consideradas: calificación, horas de estudio y sueño. la prueba T² de Hotelling muestra que hay diferencias multivariadas significativas entre las facultades de Ciencias y Humanidades en las tres variables analizadas.
2. Diferencias observadas entre las facultades:

• Calificación: Las diferencias en la media de calificación entre las dos facultades son de -3.43, con una desviación estándar de 2.11. Esto indica que los estudiantes de Ciencias tienen, en promedio, calificaciones más altas que los de Humanidades.

• Horas de estudio: La diferencia en las horas de estudio es de -5.21, con una desviación estándar de 2.40. Los estudiantes de Ciencias estudian más horas, en promedio, que los de Humanidades.

• Sueño: La diferencia en el promedio de horas de sueño es de -0.77, con una desviación estándar de 1.98. Los estudiantes de Humanidades tienen, en promedio, más horas de sueño que los de Ciencias.

# Usando la librería ICSNP
xbar <- colMeans(datos_ciencias[, c("calificacion", "horas_estudio", "sueño")])
xbar  # Promedio de las diferencias entre las categorías de ambas variables

##  calificacion horas_estudio         sueño 
##     15.016220     20.224943      6.823052

# Realizar la prueba T² de Hotelling respecto a un vector de medias conocido
HotellingsT2(datos_ciencias[, c("calificacion", "horas_estudio", "sueño")], mu=c(15, 18, 7))

## 
##  Hotelling's one sample T2-test
## 
## data:  datos_ciencias[, c("calificacion", "horas_estudio", "sueño")]
## T.2 = 17.022, df1 = 3, df2 = 27, p-value = 2.089e-06
## alternative hypothesis: true location is not equal to c(15,18,7)

a. Test de normalidad multivariada

mshapiro.test(t(tratA))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.99485, p-value = 0.2451

• Fertilizante A:Las papas tratadas con A muestran una distribución normal en sus características combinadas de tamaño y peso.

mshapiro.test(t(tratB))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.99635, p-value = 0.5491

• Fertilizante B:Las papas tratadas con B muestran una distribución normal en sus características combinadas de tamaño y peso

mshapiro.test(t(tratC))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.99072, p-value = 0.01835

• Fertilizante C:Las papas tratadas con C no muestran una distribución normal en sus características combinadas de tamaño y peso

mshapiro.test(t(tratD))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.99663, p-value = 0.6246

• Fertilizante D:Las papas tratadas con D muestran una distribución normal en sus características combinadas de tamaño y peso

Conclusión: 3 de 4 fertilizantes cumplen con el supuesto de normalidad, la unica que no cumple es el fertilizante C, eso nos indica que hay mayor variabilidad en cuanto al peso de las papas

b. Homogeneidad de matrices varianza-covarianza (Box’s M)

heplots::boxM(cbind(Largo, Ancho, Peso) ~ Tratamiento, data = datos)

## 
##  Box's M-test for Homogeneity of Covariance Matrices
## 
## data:  Y
## Chi-Sq (approx.) = 325.73, df = 18, p-value < 2.2e-16

Conclusión: Cada fertilizante produce no solo diferencias en los promedios (Largo, Ancho, Peso), sino también en la variabilidad individual de cada característica.Las relaciones entre variables morfológicas cambian según el fertilizante. Por ejemplo: en un fertilizante el Largo y Peso pueden estar fuertemente asociados, mientras que en otro no

c. Esfericidad de Bartlett

res1 <- Bsper(datos[, -1])
summary(res1)

##        Bartlett's Sphericity Test 
## 
## Chi-Squared Value = 213.2051 , df = 3  and p-value: <2e-16 
## 
##        Correlation Matrix 
## 
##           Largo     Ancho      Peso
## Largo 1.0000000 0.1607042 0.1259117
## Ancho 0.1607042 1.0000000 0.3226589
## Peso  0.1259117 0.3226589 1.0000000

Conclusión: Existen correlaciones significativas entre al menos algunas de las variables respuesta (Largo, Ancho, Peso), esto nos indica que los atributos físicos de las papas están interrelacionados.
- Analisis de relaciones: - Ancho-Peso: Las papas más anchas tienden a ser más pesadas (relación más fuerte) - Largo-Ancho: El largo no predice fuertemente el ancho (papas largas no necesariamente anchas) - Largo-Peso: El peso depende más del ancho que del largo

a) SUPUESTOS DE NORMALIDAD MULTIVARIADA

• Prueba utilizada: Prueba de Shapiro-Wilk Multivariada

  Hipótesis:

  • H₀:Los datos siguen una distribución normal multivariada.

  • H₁:: Los datos no siguen una distribución normal multivariada.

library(tidyverse)

parc1 = datos %>% filter(Parcela== "P1") %>%
  dplyr::select(Rendimiento, Altura_planta, N_macollos, Humedad_grano)
parc2 = datos %>% filter(Parcela== "P2") %>%
  dplyr::select( Rendimiento, Altura_planta, N_macollos, Humedad_grano)
parc3 = datos %>% filter(Parcela== "P3") %>%
  dplyr::select( Rendimiento, Altura_planta, N_macollos, Humedad_grano)
parc4 = datos %>% filter(Parcela== "P4") %>%
  dplyr::select(Rendimiento, Altura_planta, N_macollos, Humedad_grano)
parc5 = datos %>% filter(Parcela== "P5") %>%
  dplyr::select(Rendimiento, Altura_planta, N_macollos, Humedad_grano)
parc6 = datos %>% filter(Parcela== "P6") %>%
  dplyr::select(Rendimiento, Altura_planta, N_macollos, Humedad_grano)
parc7 = datos %>% filter(Parcela== "P7") %>%
  dplyr::select(Rendimiento, Altura_planta, N_macollos, Humedad_grano)
parc8 = datos %>% filter(Parcela== "P8") %>%
  dplyr::select(Rendimiento, Altura_planta, N_macollos, Humedad_grano)
parc9 = datos %>% filter(Parcela== "P9") %>%
  dplyr::select(Rendimiento, Altura_planta, N_macollos, Humedad_grano)
parc10 = datos %>% filter(Parcela== "P10") %>%
  dplyr::select(Rendimiento, Altura_planta, N_macollos, Humedad_grano)

#Ejecutamos el test

library(mvnormtest)
mshapiro.test(t(parc1))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.98598, p-value = 0.3727

mshapiro.test(t(parc2))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.97124, p-value = 0.0275

mshapiro.test(t(parc3))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.97823, p-value = 0.09675

mshapiro.test(t(parc4))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.98069, p-value = 0.1503

mshapiro.test(t(parc5))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.97796, p-value = 0.09208

mshapiro.test(t(parc6))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.97663, p-value = 0.07245

mshapiro.test(t(parc7))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.98803, p-value = 0.5101

mshapiro.test(t(parc8))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.97599, p-value = 0.06445

mshapiro.test(t(parc9))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.98004, p-value = 0.1338

mshapiro.test(t(parc10))

## 
##  Shapiro-Wilk normality test
## 
## data:  Z
## W = 0.98376, p-value = 0.2578

Conclusión general: Solo el grupo P2 presentó evidencia significativa (p < 0.05) para rechazar la normalidad multivariada. El resto de grupos no presentaron evidencia suficiente para rechazar H₀, por lo que se asume normalidad multivariada en general, aunque con una ligera violación en P2.

b) SUPUESTOS DE HOMOGENEIDAD DE MATRICES VARIANZA-COVARIANZAS

Prueba utilizada:

• Box’s M test:

  Hipótesis:

  • H₀:Las matrices varianza-covarianza son iguales entre los grupos.

  • H₁:Al menos una matriz de variaanza-covarianza difiere.

library(heplots)

res <- boxM(datos[, 3:6], datos[, "Parcela"])
res

## 
##  Box's M-test for Homogeneity of Covariance Matrices
## 
## data:  datos[, 3:6]
## Chi-Sq (approx.) = 96.973, df = 90, p-value = 0.289

summary(res)

## Summary for Box's M-test of Equality of Covariance Matrices
## 
## Chi-Sq:   96.97255 
## df:   90 
## p-value: 0.289 
## 
## log of Covariance determinants:
##       P1      P10       P2       P3       P4       P5       P6       P7 
## 12.19054 12.43337 13.07821 12.20480 12.81953 12.41135 12.58394 12.38739 
##       P8       P9 
## 12.97447 12.32000 
## 
## Eigenvalues:
##             P1          P10          P2           P3           P4           P5
## 1 1553.4461600 1940.0028127 1832.236290 1529.0320797 1933.3313377 2064.7884392
## 2   32.3004261   41.1206659   39.557862   31.8598596   40.5728053   29.9121505
## 3    4.4566978    4.1773853    5.038287    5.0984854    5.4225543    4.1302928
## 4    0.8805773    0.7533134    1.310077    0.8042227    0.8683702    0.9626653
##            P6           P7          P8          P9      pooled
## 1 1449.399640 1485.6962628 1882.276721 1335.152950 1700.122171
## 2   29.138199   30.3532268   31.817286   32.671696   34.225869
## 3    5.894682    6.2239226    6.584333    4.871831    5.279203
## 4    1.172263    0.8542285    1.093656    1.054659    1.004694
## 
## Statistics based on eigenvalues:
##                     P1          P10           P2           P3           P4
## product   1.969179e+05 2.510399e+05 4.784032e+05 1.997459e+05 3.693602e+05
## sum       1.591084e+03 1.986054e+03 1.878143e+03 1.566795e+03 1.980195e+03
## precision 7.185957e-01 6.282642e-01 1.012536e+00 6.795256e-01 7.346665e-01
## max       1.553446e+03 1.940003e+03 1.832236e+03 1.529032e+03 1.933331e+03
##                     P5           P6           P7           P8           P9
## product   2.455723e+05 2.918342e+05 2.397580e+05 4.312599e+05 2.241336e+05
## sum       2.099794e+03 1.485605e+03 1.523128e+03 1.921772e+03 1.373751e+03
## precision 7.605651e-01 9.454434e-01 7.326351e-01 9.105803e-01 8.440302e-01
## max       2.064788e+03 1.449400e+03 1.485696e+03 1.882277e+03 1.335153e+03
##                 pooled
## product   3.086292e+05
## sum       1.740632e+03
## precision 8.233462e-01
## max       1.700122e+03

heplots::boxM(cbind(Rendimiento, Altura_planta, N_macollos, Humedad_grano) 
              ~ Parcela, data = datos)

## 
##  Box's M-test for Homogeneity of Covariance Matrices
## 
## data:  Y
## Chi-Sq (approx.) = 96.973, df = 90, p-value = 0.289

library(biotools)

boxM(datos[3:6], datos[, 1])

## 
##  Box's M-test for Homogeneity of Covariance Matrices
## 
## data:  datos[3:6]
## Chi-Sq (approx.) = 44.222, df = 40, p-value = 0.2979

Conclusión: Dado que los p-valores son mayores a 0.05, no se rechaza H₀.
Se cumple el supuesto de homogeneidad de matrices de varianza-covarianza.

• Prueba de Ahmad et al

  Hipótesis:

  • H₀:Todas las matrices de covarianza son iguales

  • H₁:Al menos una matriz es diferente.

library(covTestR)
maiz<- unique(datos$Parcela)
maiz1 <- lapply(maiz,
    function(x){as.matrix(datos[datos$Parcela == x, 3:6])}
)

names(maiz1) <- maiz
Ahmad2017(maiz1)

## 
##  Ahmad 2017 Homogeneity of Covariance Matrices Test
## 
## data:  P1, P2, P3, P4, P5, P6, P7, P8, P9 and P10
## Standard Normal = 260653384, Mean = 0, Variance = 1, p-value < 2.2e-16
## alternative hypothesis: true difference in covariance matrices is not equal to 0

Conclusión: Dado que el p-valor es muy pequeño (menor a 0.05), se rechaza H₀.
No se cumple el supuesto de homogeneidad de matrices de covarianza, según el test de Ahmad (2017).

• Prueba Wrapper

  Hipótesis:

  • H₀:Todas las matrices de covarianza son iguales

  • H₁:Al menos una difiere.

## Prueba Wrapper
homogeneityCovariances(datos[, -1], group = Parcela, covTest = BoxesM)

## 
##  Boxes' M Homogeneity of Covariance Matrices Test
## 
## data:  P1, P2, P3, P4, P5, P6, P7, P8, P9 and P10
## Chi-Squared = 98.542, df = 45450, p-value = 1
## alternative hypothesis: true difference in covariance matrices is not equal to 0

library(DFA.CANCOR)
HOMOGENEITY(data = datos[, -1], groups = 'Parcela', 
            variables = c('Rendimiento', 'Altura_planta',
                          'N_macollos', 'Humedad_grano'))   

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       1540.17        136.46      35.27         12.45
## Altura_planta      136.46         43.02       9.50          2.86
## N_macollos          35.27          9.50       6.79          0.88
## Humedad_grano       12.45          2.86       0.88          1.10

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       1826.37         94.03      36.32         18.55
## Altura_planta       94.03         43.98       5.62          3.13
## N_macollos          36.32          5.62       6.16          0.87
## Humedad_grano       18.55          3.13       0.87          1.65

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       1518.89        114.46      42.78         12.78
## Altura_planta      114.46         40.33       6.04          2.27
## N_macollos          42.78          6.04       6.60          0.33
## Humedad_grano       12.78          2.27       0.33          0.97

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       1921.30        139.49      53.33         19.03
## Altura_planta      139.49         50.77       5.66          3.04
## N_macollos          53.33          5.66       6.99          0.44
## Humedad_grano       19.03          3.04       0.44          1.13

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       2057.74        106.74      52.11         14.40
## Altura_planta      106.74         35.01       6.20          1.89
## N_macollos          52.11          6.20       5.91          0.77
## Humedad_grano       14.40          1.89       0.77          1.13

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       1437.70        117.23      51.86          8.30
## Altura_planta      117.23         38.34       7.86          0.80
## N_macollos          51.86          7.86       8.34          0.35
## Humedad_grano        8.30          0.80       0.35          1.22

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       1471.62        134.55      46.46          6.79
## Altura_planta      134.55         42.80       6.00          0.89
## N_macollos          46.46          6.00       7.77          0.75
## Humedad_grano        6.79          0.89       0.75          0.94

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       1869.86        136.32      62.67         18.27
## Altura_planta      136.32         41.14       8.86          2.68
## N_macollos          62.67          8.86       9.42          1.13
## Humedad_grano       18.27          2.68       1.13          1.35

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       1323.85        115.48      34.29          9.28
## Altura_planta      115.48         42.28       7.35          1.88
## N_macollos          34.29          7.35       6.38          0.95
## Humedad_grano        9.28          1.88       0.95          1.24

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento       1932.99        105.68      44.69         11.20
## Altura_planta      105.68         46.43       6.82          2.41
## N_macollos          44.69          6.82       5.70          0.81
## Humedad_grano       11.20          2.41       0.81          0.94

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento      1690.049       120.045     45.978        13.104
## Altura_planta     120.045        42.411      6.989         2.184
## N_macollos         45.978         6.989      7.006         0.727
## Humedad_grano      13.104         2.184      0.727         1.166

##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento         1.000         0.448      0.423         0.295
## Altura_planta       0.448         1.000      0.405         0.311
## N_macollos          0.423         0.405      1.000         0.254
## Humedad_grano       0.295         0.311      0.254         1.000

##        Log Determinant
## P1              12.191
## P2              13.078
## P3              12.205
## P4              12.820
## P5              12.411
## P6              12.584
## P7              12.387
## P8              12.974
## P9              12.320
## P10             12.433
## Pooled          12.640

Conclusión: Se cumple el supuesto de homogeneidad de las matrices de varianzas-covarianzas entre los grupos (parcelas)

b) SUPUESTOS DE VARIABLES DEPENDIENTES CORRELACIONADAS

• Prueba de esfericidad de Bartlett

  Hipótesis:

  • H₀:No estan correlacionas det(Rp)=1 (La matriz de correlación es la matriz identidad)

  • H₁:Al menos una difiere.(La matriz de correlación no es la matriz identidad)

library(psych)
options(scipen = 0)
cortest.bartlett(cor(datos[,c(-1, -2)]), n = nrow(datos[, c(-1, -2)]))

## $chisq
## [1] 627.3808
## 
## $p.value
## [1] 2.888833e-132
## 
## $df
## [1] 6

library(MVTests)
res1 = Bsper(datos[, c(-1, -2)])
summary(res1)

##        Bartlett's Sphericity Test 
## 
## Chi-Squared Value = 627.3808 , df = 6  and p-value: <2e-16 
## 
##        Correlation Matrix 
## 
##               Rendimiento Altura_planta N_macollos Humedad_grano
## Rendimiento     1.0000000     0.4477554  0.4168608     0.2906873
## Altura_planta   0.4477554     1.0000000  0.4012487     0.3104042
## N_macollos      0.4168608     0.4012487  1.0000000     0.2462650
## Humedad_grano   0.2906873     0.3104042  0.2462650     1.0000000

Conclusión: Dado que el valor-p es mucho menor que 0.05, se rechaza la hipótesis nula, concluyéndose que existe una correlación significativa entre las variables. Por lo tanto, es adecuado aplicar análisis MANOVA .

a.COMPARACIÓN POR PARES

Este análisis busca evaluar si realmente Parcela tiene un efecto importante en un subconjunto de los datos (P1 y P3).

• Prueba utilizada: Pillai, Wilks, Hotelling-Lawley, Roy.

  Hipótesis Para el factor Parcela (P1 vs P3):

  • H₀:Los vectores de medias para las parcelas P1 y P3 son iguales.

  • H₁:Al menos un vector de medias difiere entre las parcelas P1 y P3.

  Hipótesis Para el factor Parcela (bloqueo:):

  • H₀:Los vectores de medias para las distintas variedades son iguales.

  • H₁:Al menos un vector de medias difiere entre variedades.

modelo1 = manova(cbind(Rendimiento, Altura_planta,
                          N_macollos, Humedad_grano) ~ Variedad + Parcela, data = datos,
                 subset = Parcela %in% c("P1", "P3"))

summary(modelo1, test = "Pillai")

##            Df  Pillai approx F num Df den Df Pr(>F)    
## Variedad    4 0.77046  11.5704     16    776 <2e-16 ***
## Parcela     1 0.02629   1.2894      4    191 0.2756    
## Residuals 194                                          
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary(modelo1, test = "Wilks")

##            Df   Wilks approx F num Df den Df Pr(>F)    
## Variedad    4 0.27270  19.3530     16 584.15 <2e-16 ***
## Parcela     1 0.97371   1.2894      4 191.00 0.2756    
## Residuals 194                                          
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary(modelo1, test = "Hotelling-Lawley")

##            Df Hotelling-Lawley approx F num Df den Df Pr(>F)    
## Variedad    4           2.5098  29.7253     16    758 <2e-16 ***
## Parcela     1           0.0270   1.2894      4    191 0.2756    
## Residuals 194                                                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary(modelo1, test = "Roy")

##            Df    Roy approx F num Df den Df Pr(>F)    
## Variedad    4 2.4461  118.636      4    194 <2e-16 ***
## Parcela     1 0.0270    1.289      4    191 0.2756    
## Residuals 194                                         
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary.aov(modelo1)

##  Response Rendimiento :
##              Df Sum Sq Mean Sq F value Pr(>F)    
## Variedad      4 142412   35603 43.0516 <2e-16 ***
## Parcela       1    126     126  0.1529 0.6962    
## Residuals   194 160435     827                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##  Response Altura_planta :
##              Df Sum Sq Mean Sq F value Pr(>F)    
## Variedad      4 3794.4  948.60 41.2802 <2e-16 ***
## Parcela       1   22.1   22.11  0.9624 0.3278    
## Residuals   194 4458.0   22.98                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##  Response N_macollos :
##              Df Sum Sq Mean Sq F value Pr(>F)    
## Variedad      4 513.60 128.400 30.6767 <2e-16 ***
## Parcela       1   6.07   6.069  1.4499   0.23    
## Residuals   194 812.00   4.186                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##  Response Humedad_grano :
##              Df  Sum Sq Mean Sq F value    Pr(>F)    
## Variedad      4  48.450 12.1126 15.0259 1.023e-10 ***
## Parcela       1   2.548  2.5479  3.1606     0.077 .  
## Residuals   194 156.387  0.8061                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Decisión sobre las hipótesis

Para Parcela: dado que p-valor > 0.05 en todos los tests, no se rechaza H0. Es decir, no hay diferencias significativas entre las parcelas P1 y P3 en el vector de respuestas.

Para Variedad: se rechaza H0 con evidencia muy fuerte, por lo tanto, sí hay diferencias entre las variedades.

Conclusión específica del análisis restringido

El bloqueo por Parcela no es necesario ni útil al evaluar únicamente P1 y P3, ya que no hay diferencias significativas entre ellas.

Variedad sigue mostrando un efecto significativo fuerte, lo que reafirma su importancia como factor de estudio.

b.COMPARACIÓN GENERAL POR PARES

• Prueba utilizada: Pillai, Wilks, Hotelling-Lawley, Roy.

  Hipótesis para el factor variedad:

  • H₀: Las dos parcelas en comparación no presentan diferencias multivariadas significativas.

  • H₁:Las dos parcelas en comparación presentan diferencias multivariadas significativas.

  Parc <- c("P1", "P2", "P3", "P4","P5", "P6", "P7", "P8","P9","P10")
  comb <- t(combn(length(Parc), 2))
  
  for(i in 1:nrow(comb)){
    modelo.comp = manova(cbind(Rendimiento, Altura_planta, N_macollos, Humedad_grano) ~
                      Variedad + Parcela, data = datos,
                     subset = Parcela %in% Parc[comb[i,]])
    print(paste("Trat: ",Parc[comb[i,]][1], "y", Parc[comb[i,]][2]))
    print(summary(modelo.comp, test = "Pillai"))
    cat("\n")
  }

  ## [1] "Trat:  P1 y P2"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.78501  11.8423     16    776 <2e-16 ***
  ## Parcela     1 0.02115   1.0318      4    191 0.3921    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P1 y P3"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.77046  11.5704     16    776 <2e-16 ***
  ## Parcela     1 0.02629   1.2894      4    191 0.2756    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P1 y P4"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.74119  11.0309     16    776 <2e-16 ***
  ## Parcela     1 0.03821   1.8971      4    191 0.1126    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P1 y P5"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.76736  11.5128     16    776 <2e-16 ***
  ## Parcela     1 0.02437   1.1927      4    191 0.3154    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P1 y P6"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.80918   12.299     16    776 <2e-16 ***
  ## Parcela     1 0.00981    0.473      4    191 0.7555    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P1 y P7"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.78455  11.8337     16    776 <2e-16 ***
  ## Parcela     1 0.02612   1.2805      4    191 0.2791    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P1 y P8"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.75323  11.2518     16    776 <2e-16 ***
  ## Parcela     1 0.02815   1.3832      4    191 0.2413    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P1 y P9"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.76537  11.4760     16    776 <2e-16 ***
  ## Parcela     1 0.03507   1.7354      4    191 0.1438    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P1 y P10"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.83646  12.8238     16    776 <2e-16 ***
  ## Parcela     1 0.01133   0.5474      4    191 0.7011    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P2 y P3"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.75917  11.3612     16    776 <2e-16 ***
  ## Parcela     1 0.02504   1.2263      4    191 0.3011    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P2 y P4"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.72754  10.7826     16    776 <2e-16 ***
  ## Parcela     1 0.01740   0.8457      4    191 0.4978    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P2 y P5"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.75013  11.1948     16    776 <2e-16 ***
  ## Parcela     1 0.03775   1.8732      4    191 0.1168    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P2 y P6"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.79209  11.9755     16    776 <2e-16 ***
  ## Parcela     1 0.01165   0.5628      4    191   0.69    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P2 y P7"
  ##            Df  Pillai approx F num Df den Df    Pr(>F)    
  ## Variedad    4 0.77985  11.7456     16    776 < 2.2e-16 ***
  ## Parcela     1 0.06858   3.5158      4    191  0.008537 ** 
  ## Residuals 194                                             
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P2 y P8"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.77536  11.6618     16    776 <2e-16 ***
  ## Parcela     1 0.00774   0.3727      4    191 0.8279    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P2 y P9"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.79693  12.0669     16    776 <2e-16 ***
  ## Parcela     1 0.03328   1.6437      4    191 0.1649    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P2 y P10"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.80098  12.1436     16    776 <2e-16 ***
  ## Parcela     1 0.05002   2.5142      4    191  0.043 *  
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P3 y P4"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.71839  10.6173     16    776 <2e-16 ***
  ## Parcela     1 0.01664   0.8081      4    191 0.5214    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P3 y P5"
  ##            Df  Pillai approx F num Df den Df    Pr(>F)    
  ## Variedad    4 0.73820  10.9764     16    776 < 2.2e-16 ***
  ## Parcela     1 0.07718   3.9937      4    191  0.003903 ** 
  ## Residuals 194                                             
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P3 y P6"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.79779  12.0831     16    776 < 2e-16 ***
  ## Parcela     1 0.04173   2.0796      4    191 0.08505 .  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P3 y P7"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.78097  11.7665     16    776 < 2e-16 ***
  ## Parcela     1 0.05456   2.7554      4    191 0.02926 *  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P3 y P8"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.75800  11.3395     16    776 <2e-16 ***
  ## Parcela     1 0.03872   1.9235      4    191 0.1081    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P3 y P9"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.73783  10.9696     16    776 <2e-16 ***
  ## Parcela     1 0.01244   0.6017      4    191 0.6619    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P3 y P10"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.77570  11.6681     16    776 < 2e-16 ***
  ## Parcela     1 0.05367   2.7082      4    191 0.03156 *  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P4 y P5"
  ##            Df  Pillai approx F num Df den Df    Pr(>F)    
  ## Variedad    4 0.79971  12.1195     16    776 < 2.2e-16 ***
  ## Parcela     1 0.08004   4.1543      4    191  0.002997 ** 
  ## Residuals 194                                             
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P4 y P6"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.75032  11.1982     16    776 <2e-16 ***
  ## Parcela     1 0.03495   1.7294      4    191 0.1451    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P4 y P7"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.75537   11.291     16    776 < 2e-16 ***
  ## Parcela     1 0.04736    2.374      4    191 0.05368 .  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P4 y P8"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.74562  11.1120     16    776 <2e-16 ***
  ## Parcela     1 0.01185   0.5725      4    191 0.6829    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P4 y P9"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.71801  10.6105     16    776 <2e-16 ***
  ## Parcela     1 0.00880   0.4241      4    191 0.7911    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P4 y P10"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.77177  11.5948     16    776 < 2e-16 ***
  ## Parcela     1 0.05866   2.9756      4    191 0.02053 *  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P5 y P6"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.77810  11.7129     16    776 <2e-16 ***
  ## Parcela     1 0.01136   0.5488      4    191 0.7001    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P5 y P7"
  ##            Df  Pillai approx F num Df den Df    Pr(>F)    
  ## Variedad    4 0.78193  11.7847     16    776 < 2.2e-16 ***
  ## Parcela     1 0.07529   3.8879      4    191  0.004642 ** 
  ## Residuals 194                                             
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P5 y P8"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.77050  11.5713     16    776 < 2e-16 ***
  ## Parcela     1 0.04389   2.1918      4    191 0.07143 .  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P5 y P9"
  ##            Df  Pillai approx F num Df den Df    Pr(>F)    
  ## Variedad    4 0.73789  10.9707     16    776 < 2.2e-16 ***
  ## Parcela     1 0.10189   5.4171      4    191 0.0003747 ***
  ## Residuals 194                                             
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P5 y P10"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.80519  12.2234     16    776 < 2e-16 ***
  ## Parcela     1 0.04119   2.0516      4    191 0.08881 .  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P6 y P7"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.80696  12.2571     16    776 < 2e-16 ***
  ## Parcela     1 0.05762   2.9197      4    191 0.02247 *  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P6 y P8"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.79644  12.0576     16    776 <2e-16 ***
  ## Parcela     1 0.01432   0.6937      4    191 0.5971    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P6 y P9"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.81416  12.3944     16    776 <2e-16 ***
  ## Parcela     1 0.04252   2.1206      4    191 0.0798 .  
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P6 y P10"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.79398  12.0112     16    776 <2e-16 ***
  ## Parcela     1 0.03006   1.4797      4    191 0.2099    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P7 y P8"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.76720  11.5100     16    776 < 2e-16 ***
  ## Parcela     1 0.04538   2.2698      4    191 0.06323 .  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P7 y P9"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.76977  11.5577     16    776 <2e-16 ***
  ## Parcela     1 0.03186   1.5716      4    191 0.1835    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P7 y P10"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.82999   12.699     16    776 <2e-16 ***
  ## Parcela     1 0.01024    0.494      4    191 0.7401    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P8 y P9"
  ##            Df  Pillai approx F num Df den Df Pr(>F)    
  ## Variedad    4 0.77083  11.5774     16    776 <2e-16 ***
  ## Parcela     1 0.02733   1.3417      4    191  0.256    
  ## Residuals 194                                          
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P8 y P10"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.76078  11.3909     16    776 < 2e-16 ***
  ## Parcela     1 0.04300   2.1458      4    191 0.07674 .  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  ## 
  ## [1] "Trat:  P9 y P10"
  ##            Df  Pillai approx F num Df den Df  Pr(>F)    
  ## Variedad    4 0.76825  11.5294     16    776 < 2e-16 ***
  ## Parcela     1 0.04856   2.4373      4    191 0.04857 *  
  ## Residuals 194                                           
  ## ---
  ## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Interpretación:

El bloqueo por parcela controla parte de la variabilidad experimental, pero no garantiza homogeneidad total entre todas las parcelas. Hay subconjuntos (por ejemplo, P2–P7 o P5–P9) donde el bloque ejerce un efecto significativo que habría que investigar o ajustar en el modelo. Sin embargo, la Variedad mantiene un efecto fuertemente significativo en todas las comparaciones.

Interacciones de Fertilizante con Riego

linearHypothesis(modelo, "FertilizanteMix:RiegoGoteo = 0")

## 
## Sum of squares and products for the hypothesis:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg     1677117.12   14558.500451 -137.27998949
## Altura_plantas       14558.50     126.377540   -1.19168230
## Humedad_grano         -137.28      -1.191682    0.01123702
## 
## Sum of squares and products for error:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   125198211.15     26679.1023    11266.1071
## Altura_plantas       26679.10    142573.5247      219.1194
## Humedad_grano        11266.11       219.1194     1465.1196
## 
## Multivariate Tests: 
##                  Df test stat approx F num Df den Df     Pr(>F)    
## Pillai            1 0.0140772 7.086766      3   1489 9.9546e-05 ***
## Wilks             1 0.9859228 7.086766      3   1489 9.9546e-05 ***
## Hotelling-Lawley  1 0.0142782 7.086766      3   1489 9.9546e-05 ***
## Roy               1 0.0142782 7.086766      3   1489 9.9546e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

• Existe una interacción significativa entre Mix y Goteo, lo que indica que su efecto combinado tiene un comportamiento único que no puede predecirse sumando sus efectos individuales.

linearHypothesis(modelo, "FertilizanteNitro:RiegoGoteo = 0")

## 
## Sum of squares and products for the hypothesis:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   258341.88079  -4660.9572517  40.329712717
## Altura_plantas    -4660.95725     84.0921435  -0.727621346
## Humedad_grano        40.32971     -0.7276213   0.006295865
## 
## Sum of squares and products for error:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   125198211.15     26679.1023    11266.1071
## Altura_plantas       26679.10    142573.5247      219.1194
## Humedad_grano        11266.11       219.1194     1465.1196
## 
## Multivariate Tests: 
##                  Df test stat approx F num Df den Df  Pr(>F)
## Pillai            1 0.0026618 1.324651      3   1489 0.26472
## Wilks             1 0.9973382 1.324651      3   1489 0.26472
## Hotelling-Lawley  1 0.0026689 1.324651      3   1489 0.26472
## Roy               1 0.0026689 1.324651      3   1489 0.26472

• No hay evidencia suficiente para concluir que la interacción entre FertilizanteNitro y RiegoGoteo tenga un efecto significativo. Esto implica que su efecto conjunto no se desvía significativamente de la suma de sus efectos individuales. Puedes interpretar los efectos de FertilizanteNitro de manera independiente del tipo de riego

linearHypothesis(modelo, "FertilizanteMix:RiegoGoteo - FertilizanteNitro:RiegoGoteo = 0")

## 
## Sum of squares and products for the hypothesis:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg    3280107.281   3622.5957628 -4.719400e+01
## Altura_plantas       3622.596      4.0008448 -5.212170e-02
## Humedad_grano         -47.194     -0.0521217  6.790245e-04
## 
## Sum of squares and products for error:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   125198211.15     26679.1023    11266.1071
## Altura_plantas       26679.10    142573.5247      219.1194
## Humedad_grano        11266.11       219.1194     1465.1196
## 
## Multivariate Tests: 
##                  Df test stat approx F num Df den Df     Pr(>F)    
## Pillai            1 0.0255716  13.0251      3   1489 2.1258e-08 ***
## Wilks             1 0.9744284  13.0251      3   1489 2.1258e-08 ***
## Hotelling-Lawley  1 0.0262426  13.0251      3   1489 2.1258e-08 ***
## Roy               1 0.0262426  13.0251      3   1489 2.1258e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

• El efecto de Fertilizante Mix es significativamente diferente al efecto de Fertilizante Nitro en el riego Goteo. Esto implica que los efectos de Fertilizante Mix y Fertilizante Nitro no son iguales cuando se usa el riego Goteo.

Efectos simples

# Comparar Mix vs Fosfo cuando el riego es por Goteo
linearHypothesis(modelo, "FertilizanteMix + FertilizanteMix:RiegoGoteo = 0")

## 
## Sum of squares and products for the hypothesis:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg    46000953.34    1590294.162   -71623.9316
## Altura_plantas     1590294.16      54977.894    -2476.1035
## Humedad_grano       -71623.93      -2476.103      111.5192
## 
## Sum of squares and products for error:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   125198211.15     26679.1023    11266.1071
## Altura_plantas       26679.10    142573.5247      219.1194
## Humedad_grano        11266.11       219.1194     1465.1196
## 
## Multivariate Tests: 
##                  Df test stat approx F num Df den Df     Pr(>F)    
## Pillai            1 0.4562412 416.4488      3   1489 < 2.22e-16 ***
## Wilks             1 0.5437588 416.4488      3   1489 < 2.22e-16 ***
## Hotelling-Lawley  1 0.8390507 416.4488      3   1489 < 2.22e-16 ***
## Roy               1 0.8390507 416.4488      3   1489 < 2.22e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

• Existe un efecto simple significativo al comparar el Fertilizante Mix vs Fosfo cuando el Riego es por Goteo.

# Comparar Nitro vs Fosfo cuando el riego es por Goteo
linearHypothesis(modelo, "FertilizanteNitro + FertilizanteNitro:RiegoGoteo = 0")

## 
## Sum of squares and products for the hypothesis:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg     3859602.58     230569.089   -22240.4525
## Altura_plantas      230569.09      13773.984    -1328.6241
## Humedad_grano       -22240.45      -1328.624      128.1577
## 
## Sum of squares and products for error:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   125198211.15     26679.1023    11266.1071
## Altura_plantas       26679.10    142573.5247      219.1194
## Humedad_grano        11266.11       219.1194     1465.1196
## 
## Multivariate Tests: 
##                  Df test stat approx F num Df den Df     Pr(>F)    
## Pillai            1 0.1802447 109.1319      3   1489 < 2.22e-16 ***
## Wilks             1 0.8197553 109.1319      3   1489 < 2.22e-16 ***
## Hotelling-Lawley  1 0.2198763 109.1319      3   1489 < 2.22e-16 ***
## Roy               1 0.2198763 109.1319      3   1489 < 2.22e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

• El efecto del Fertilizante Nitro es diferente al del Fosfo cuando se aplica riego por Goteo.

# Comparar Mix vs Fosfo cuando el riego es por Aspersión
linearHypothesis(modelo, "FertilizanteMix = 0")

## 
## Sum of squares and products for the hypothesis:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg    23297621.10    1031935.056   -49180.8584
## Altura_plantas     1031935.06      45708.099    -2178.3963
## Humedad_grano       -49180.86      -2178.396      103.8199
## 
## Sum of squares and products for error:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   125198211.15     26679.1023    11266.1071
## Altura_plantas       26679.10    142573.5247      219.1194
## Humedad_grano        11266.11       219.1194     1465.1196
## 
## Multivariate Tests: 
##                  Df test stat approx F num Df den Df     Pr(>F)    
## Pillai            1 0.3692771 290.5944      3   1489 < 2.22e-16 ***
## Wilks             1 0.6307229 290.5944      3   1489 < 2.22e-16 ***
## Hotelling-Lawley  1 0.5854823 290.5944      3   1489 < 2.22e-16 ***
## Roy               1 0.5854823 290.5944      3   1489 < 2.22e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

• Existe una una diferencia significativa entre el fertilizante Mix y Fosfo cuando el riego es por Aspersión.

# Comparar Nitro vs Fosfo cuando el riego es por Aspersión
linearHypothesis(modelo, "FertilizanteNitro = 0")

## 
## Sum of squares and products for the hypothesis:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg     6941238.33     268925.032   -28908.5154
## Altura_plantas      268925.03      10418.987    -1120.0053
## Humedad_grano       -28908.52      -1120.005      120.3967
## 
## Sum of squares and products for error:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   125198211.15     26679.1023    11266.1071
## Altura_plantas       26679.10    142573.5247      219.1194
## Humedad_grano        11266.11       219.1194     1465.1196
## 
## Multivariate Tests: 
##                  Df test stat approx F num Df den Df     Pr(>F)    
## Pillai            1 0.1775939 107.1804      3   1489 < 2.22e-16 ***
## Wilks             1 0.8224061 107.1804      3   1489 < 2.22e-16 ***
## Hotelling-Lawley  1 0.2159443 107.1804      3   1489 < 2.22e-16 ***
## Roy               1 0.2159443 107.1804      3   1489 < 2.22e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

• Nitro es significativamente diferente a Fosfo bajo riego por aspersión

ggplot(cultivos, aes(x = interaction(Fertilizante, Riego), y = Rendimiento_kg, fill = Fertilizante)) +
  geom_boxplot() +
  labs(title = "Boxplot del Rendimiento por combinación de Fertilizante y Riego", x = "Fertilizante:Riego") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

# Comparar Goteo vs Aspersión cuando el fertilizante es Mix
linearHypothesis(modelo,"RiegoGoteo + FertilizanteMix:RiegoGoteo = 0")

## 
## Sum of squares and products for the hypothesis:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg     47016693.9    1244669.498  -155036.5658
## Altura_plantas      1244669.5      32950.045    -4104.2717
## Humedad_grano       -155036.6      -4104.272      511.2298
## 
## Sum of squares and products for error:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   125198211.15     26679.1023    11266.1071
## Altura_plantas       26679.10    142573.5247      219.1194
## Humedad_grano        11266.11       219.1194     1465.1196
## 
## Multivariate Tests: 
##                  Df test stat approx F num Df den Df     Pr(>F)    
## Pillai            1 0.4950169 486.5379      3   1489 < 2.22e-16 ***
## Wilks             1 0.5049831 486.5379      3   1489 < 2.22e-16 ***
## Hotelling-Lawley  1 0.9802644 486.5379      3   1489 < 2.22e-16 ***
## Roy               1 0.9802644 486.5379      3   1489 < 2.22e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# Comparar Goteo vs Aspersión cuando el fertilizante es Nitro
linearHypothesis(modelo,"RiegoGoteo + FertilizanteNitro:RiegoGoteo = 0")

## 
## Sum of squares and products for the hypothesis:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg    17693324.30     739328.323    -93406.590
## Altura_plantas      739328.32      30893.367     -3903.062
## Humedad_grano       -93406.59      -3903.062       493.112
## 
## Sum of squares and products for error:
##                Rendimiento_kg Altura_plantas Humedad_grano
## Rendimiento_kg   125198211.15     26679.1023    11266.1071
## Altura_plantas       26679.10    142573.5247      219.1194
## Humedad_grano        11266.11       219.1194     1465.1196
## 
## Multivariate Tests: 
##                  Df test stat approx F num Df den Df     Pr(>F)    
## Pillai            1 0.4160653 353.6475      3   1489 < 2.22e-16 ***
## Wilks             1 0.5839347 353.6475      3   1489 < 2.22e-16 ***
## Hotelling-Lawley  1 0.7125202 353.6475      3   1489 < 2.22e-16 ***
## Roy               1 0.7125202 353.6475      3   1489 < 2.22e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

ggplot(cultivos, aes(x = interaction(Fertilizante, Riego), y = Rendimiento_kg, fill = Riego)) +
  geom_boxplot() +
  labs(title = "Boxplot por combinación Fertilizante × Riego",
       x = "Tratamientos", y = "Rendimiento (kg)") +
  theme_minimal()