Se desea conocer el efecto de las cepas de inoculantes de Rhizobium, fijadoras de nitrógeno atmosférico,sobre el contenido del nitrógeno de plantas de trébol rojo. Para ello se dispone de 30 macetas de trébol rojo en un invernadero. Se asigna al azar 5 macetas para cada una de las cepas y se procede a inocularlas. Los resultados son los siguientes (en mg de nitrógeno/kg de materia seca):
library(agricolae)
trt<-c("cepa1","cepa2","cepa3","cepa4","cepa5","cepa6")
diseño<-design.crd(trt,r=5,serie=2,seed = 3,kinds = "Super-Duper",randomization=TRUE)
libro_campo<-diseño$book
libro_campo
factor de estudio: Cepas de inoculantes de Rhizobium
Tratamiento: Cepa 1 , Cepa 2 , Cepa 3 , Cepa 4 , Cepa 5 , Cepa 6
Unidad experimental: 1 maceta de trébol rojo
Variable respuesta: Mg de nitrógeno en las macetas de trébol
\(Y_{ij} = \mu + t_{i} + \epsilon_{ij }\) , con i : 1 , 2 , 3 , 4 , 5 , 6
\(Y_{ij}\) : Mg de Nitrógeno en la j-ésima maceta de trébol rojo que recibió el tratamiento i.
\(\mu\) : : Promedio general de Mg de Nitrógeno común a los seis tratamientos
\(t_{i}\) : Efecto del i-ésimo tratamiento
\(\epsilon_{ij}\): Efecto del error experimental al usar el i-ésimo tratamiento en la j-ésima repetición
library(readxl)
datos <- read_excel("diseños.xlsx")
RPTA <- as.numeric(datos$rpta)
TRATAMIENTOS <- factor(datos$trat)
modelo <- lm(RPTA~TRATAMIENTOS)
modelo
##
## Call:
## lm(formula = RPTA ~ TRATAMIENTOS)
##
## Coefficients:
## (Intercept) TRATAMIENTOScepa2 TRATAMIENTOScepa3 TRATAMIENTOScepa4
## 31.22 -5.24 -13.58 -11.30
## TRATAMIENTOScepa5 TRATAMIENTOScepa6
## -17.96 -12.52
\(H_{0}\) : Todos los tratamientos tienen efectos parecidos
\(H_{1}\) : Al menos un tratamiento difiere de los demás
\(H_{0}\) : \(\mu_{1}\) = \(\mu_{2}\) = \(\mu_{3}\) = \(\mu_{4}\) = \(\mu_{5}\) = \(\mu_{6}\)
\(H_{1}\) : \(\mu_{i}\) \(\neq\) \(\mu_{j}\) al menos para un par ( i , j )
Equivalentemente se puede plantear la hipótesis estadística
\(H_{0}\) : \(T_{1}\) = \(T_{1}\) = \(T_{2}\) = \(T_{3}\) = \(T_{4}\) = \(T_{5}\) = \(T_{6}\)
\(H_{1}\) : \(T_{i}\) \(\neq\) \(T_{j}\)
Para generar la tabla ANOVA empleamos la función aov
modelo <- aov(RPTA ~ TRATAMIENTOS)
summary(modelo)
## Df Sum Sq Mean Sq F value Pr(>F)
## TRATAMIENTOS 5 1034.1 206.82 84.63 1.92e-14 ***
## Residuals 24 58.6 2.44
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Puesto que p_value < 0.05 , se rechaza \(H_{0}\) por lo tanto se concluye que los Mg de nitrógeno promedios en los 6 tratamientos son diferentes ; es decir al menos uno de los 6 tratamientos ( Cepa1 , Cepa2 , Cepa3 , Cepa4 , Cepa5 , Cepa 6 ) afecta de manera significativa en la cantidad de Mg de nitrógeno producidos por las plantas de trébol rojo
A un nivel de significación de 5% existe suficiente evidencia estadística para rechazar \(H_{0}\) ; por lo que al menos un tratamiento difiere en los demás
Los supuestos que deben cumplir los residuales son los siguientes :
El supuesto del modelo aditivo se verifica mediante una gráfica
par(mfrow = c(2,2))
plot(modelo)
Se puede observar en la primera gráfica que los residuos no forman algun patron , lo cual podemos conluir que los errores \(\epsilon_{ij}\) tiene variancia constante \(\sigma^{2}\)
Algo similar podemos observar en la tercera gráfica que corresponde a los valores predichos versus la raíz cuadrada de los valores absolutos de los residuales estandarizados . En la gráfica se puede observar que no se presenta una tendencia marcada , lo que significa que si cumple el supuesto de homogeneidad de varianzas
mean(modelo$residuals)
## [1] 8.140913e-17
La media de \(\epsilon_{ij}\) = 8.140913e-17 que practicamente es igual a cero , lo cual si cumple con el supuesto
\(H_{0}\) : Los seis tratamientos tienen igual varianza
\(H_{1}\) : Al menos uno de los tratamientos tiene varianza distinta
bartlett.test(RPTA~TRATAMIENTOS)
##
## Bartlett test of homogeneity of variances
##
## data: RPTA by TRATAMIENTOS
## Bartlett's K-squared = 1.0124, df = 5, p-value = 0.9616
Dado que el p-value = 0.9616 > 0.05 , a un nivel de significancia de 5% no existe evidencia suficiente para rechazar \(H_{0}\) , por lo que los 6 tratamientos tienen igual varianza
library(car)
leveneTest(RPTA~TRATAMIENTOS)
\(H_{0}\) : Cumple con la homogeneidad de varianzas
\(H_{1}\) : No cumple con la homogeneidad de varianzas
Para un nivel de significancia de 5% ( p-value = 0.979 > 0.05 ) , no existe evidencia suficiente para rechazar \(H_{0}\) , por lo cual si cumple con el supuesto de homogeneidad de varianzas
\(H_{0}\) : La dsitribución de los errores es normal
\(H_{1}\) : La distribución de los errores no es normal
modelo |> rstandard() |> shapiro.test()
##
## Shapiro-Wilk normality test
##
## data: rstandard(modelo)
## W = 0.8946, p-value = 0.006211
Un zootecnista está interesado en comparar las medias de los tiempos en segundos de coagulación de animales sometidos a 4 dietas diferentes: D1, D2, D3 y D4. Se realizó un experimento con 24 animales de características similares. Las dietas fueron asignadas aleatoriamente entre los animales. Las muestras fueron tomadas en orden aleatorio. Los datos obtenidos de este experimento, se muestran a continuación:
(valores <- read_excel("diseño.xlsx"))
factor de estudio: 4 tipo de dietas
Tratamiento: Dieta 1 , Dieta 2 , Dieta 3 , Dieta 4
Unidad experimental: Un animal
Variable respuesta: Tiempo en segundos de coagulación
\(Y_{ij}\) : Tiempo en segundos de coagulación en el j-ésimo animal sometido a la dieta i
\(T_{i}\) : Efecto de la i-ésima dieta
Mediante el siguiente código vamos a poder concluir si los errores se distribuyen normalmente
valor <- read_excel("diseños_1.xlsx")
respuesta <- as.numeric(valor$valores)
tratamiento <- factor(valor$dietas)
mod <- lm(respuesta~tratamiento)
mod |> rstandard() |> shapiro.test()
##
## Shapiro-Wilk normality test
##
## data: rstandard(mod)
## W = 0.98217, p-value = 0.9322
\(H_{0}\) : La distribución de los errores es normal
\(H_{1}\) : La distribución de los errores no es normal
Con un nivel de significancia de 5% , p-value = 0.9322 > 0.05 , no existe evidencia científica para rechazar \(H_{0}\) por lo tanto los residuales se distribuyen normalmente
par(mfrow=c(2,2))
plot(mod)
En la primera gráfica se puede observarque los residuos no forman algun patron , lo cual podemos conluir que los errores \(\epsilon_{ij}\) tiene variancia constante \(\sigma^{2}\)
De la misma manera podemos concluir en la tercera gráfica que corresponde a los valores predichos versus la raíz cuadrada de los valores absolutos de los residuales estandarizados . En la gráfica se puede observar que no se presenta una tendencia marcada , lo que significa que si cumple el supuesto de homogeneidad de varianzas
library(car)
leveneTest(respuesta~tratamiento)
\(H_{0}\) : Cumple la homogeneidad de varianzas
\(H_{1}\) : No cumple la homogeneidad de varianzas
Con un nivel de significancia de 5% ( p-valor = 0.5926 > 0.05 ) , podemos concluir que no se rechaza \(H_{0}\) , es decir si cumple la homogeneidad de varianzas
modelo |> aov() |> summary()
## Df Sum Sq Mean Sq F value Pr(>F)
## TRATAMIENTOS 5 1034.1 206.82 84.63 1.92e-14 ***
## Residuals 24 58.6 2.44
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Cuadrado medio del tratamiento: 206.82
Suma de cuadrados de Error : 58.6
Cuadro medio del error : 2.44