Se sabe que el dióxido de carbono tiene un efecto crítico en el crecimiento microbiológico. Cantidades pequeñas de \(CO_2\) estimulan el crecimiento de muchos microorganismos, mientras que altas concentraciones inhiben el crecimiento de la mayor parte de ellos. Este último efecto se utiliza comercialmente cuando se almacenan productos alimenticios perecederos. Se realizó un estudio para investigar el efecto de \(CO_2\) sobre la tasa de crecimiento de Pseudomonas fragi, un corruptor de alimentos. Se administró dióxido de carbono a cinco presiones atmósfericas diferentes. La respuesta anotada es el cambio porcentual en la masa celular después de un tiempo de crecimiento de una hora. Se utilizaron diez cultivos en cada nivel. Se obtuvieron los siguientes datos:
# install.packages("readxl")#solo una vez
library(readxl)
ejercicio1 <- read_excel("C:/Users/LENOVO/Desktop/Ejercicio1.xlsx")
ejercicio1 # para visualizar los datos
## # A tibble: 50 × 2
## Crecimiento Presion
## <dbl> <chr>
## 1 62.6 A
## 2 59.6 A
## 3 64.5 A
## 4 59.3 A
## 5 58.6 A
## 6 64.6 A
## 7 50.9 A
## 8 56.2 A
## 9 52.3 A
## 10 62.8 A
## # ℹ 40 more rows
conteo_valorespresion <- table(ejercicio1$Presion)
conteo_valorespresion
##
## A B C D E
## 10 10 10 10 10
Dado que el número de observaciones por tratamiento es el mismo, se puede concluir que es un diseño balanceado.
#install.packages("summarytools")# solo una vez
library(summarytools)
summarytools::descr(ejercicio1[,1])# todas las filas primera columna datos[,1]
## Descriptive Statistics
## ejercicio1$Crecimiento
## N: 50
##
## Crecimiento
## ----------------- -------------
## Mean 36.71
## Std.Dev 15.99
## Min 7.80
## Q1 22.80
## Median 36.75
## Q3 49.90
## Max 64.60
## MAD 19.87
## IQR 26.67
## CV 0.44
## Skewness 0.08
## SE.Skewness 0.34
## Kurtosis -1.13
## N.Valid 50.00
## Pct.Valid 100.00
A partir de los resultados que arroja el programa se puede concluir lo siguiente: el promedio de la tasa de crecimiento es de 36.71, con desviación estándar de 15.99, el valor mínimo de crecimiento es de 7.80, el valor máximo de crecimiento es de 64.60, el 50% de las observaciones(25) presentaron un crecimiento entre 7.80 y 36.75, mientras que el restante 50% (25) presentó un crecimiento entre 36.75 y 64.60, el coeficiente de asimetría fue de 0.08 presentando una asimetría positiva leve, además el coeficiente de curtosis fue de -0.13, indicando que la distribución de los datos es platicúrtica.
# Calcular estadísticas descriptivas por categoría
resultados_descriptivos <- aggregate(Crecimiento ~ Presion, data = ejercicio1, summary)
# Imprimir los resultados descriptivos
print(resultados_descriptivos)
## Presion Crecimiento.Min. Crecimiento.1st Qu. Crecimiento.Median
## 1 A 50.900 56.800 59.450
## 2 B 35.200 43.025 48.000
## 3 C 27.000 31.125 38.400
## 4 D 19.200 22.650 24.250
## 5 E 7.800 11.850 17.000
## Crecimiento.Mean Crecimiento.3rd Qu. Crecimiento.Max.
## 1 59.140 62.750 64.600
## 2 46.040 49.800 50.900
## 3 36.450 40.150 45.500
## 4 25.470 29.425 32.700
## 5 16.440 21.025 24.900
LLevando a cabo el análisis descriptivo por tratamiento se obtuvieron los siguientes resultados:
Para la presión A: el promedio de crecimiento es de 59.14, el valor mínimo de crecimiento es 50.9, el valor máximo de crecimiento es 64.6, el 50% de las observaciones(5) presentaron un crecimiento entre 50.9 y 59.450, mientras que el restante 50% (5) presento un desgaste entre 59.450 y 64.6.
Para la presión B: el promedio de crecimiento es de 46.04, el valor mínimo de crecimiento es 35.2, el valor máximo de crecimiento es 50.9, el 50% de las observaciones(5) presentaron un crecimiento entre 35.2 y 48, mientras que el restante 50% (5) presento un crecimiento entre 48 y 50.9.
Para la presión C: el promedio de crecimiento es de 36.45, el valor mínimo de crecimiento es 27, el valor máximo de crecimiento es 45.5, el 50% de las observaciones(5) presentaron un crecimiento entre 27 y 38.4, mientras que el restante 50% (5) presento un crecimiento entre 38.4 y 45.5.
Para la presión D: el promedio de crecimiento es de 25.47, el valor mínimo de crecimiento es 19.2, el valor máximo de crecimiento es 32.7, el 50% de las observaciones(5) presentaron un crecimiento entre 19.2 y 24.25, mientras que el restante 50% (5) presento un crecimiento entre 24.25 y 32.7.
Para la presión E: el promedio de crecimiento es de 16.44, el valor mínimo de crecimiento es 7.8, el valor máximo de crecimiento es 24.9, el 50% de las observaciones(5) presentaron un crecimiento entre 7.8 y 17, mientras que el restante 50% (5) presento un crecimiento entre 17 y 24.9.
A continuación se llevará a cabo el ANOVA de un factor o modelo factorial de un solo factor el cual nos permite estudiar si existen diferencias significativas entre las medias de crecimiento de los cinco niveles de presión.
Las hipótesis contrastadas en el ANOVA son: No hay diferencias entre las medias de las diferentes presiones de \(CO_2\) y existen por lo menos dos presiones de \(CO_2\) con diferencias significativas en el promedio de la tasa de crecimiento de Pseudomonas fragi. Que se plantearían de la siguiente forma:
\(H_0:μ_A=μ_B=μ_C=μ_D=μ_E=μ\)
\(H_a:μ_i=μ_j\) para \(i≠j\) donde \(i,j=A,B,C,D,E\) respectivamente
# Realizar el ANOVA
modelo_anova <- aov(Crecimiento ~ Presion, data = ejercicio1)
# Resumen del ANOVA
resumen_anova <- summary(modelo_anova)
# Imprimir el resumen del ANOVA
print(resumen_anova)
## Df Sum Sq Mean Sq F value Pr(>F)
## Presion 4 11274 2818.6 101.6 <2e-16 ***
## Residuals 45 1248 27.7
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
En este caso el ANOVA ha resultado significativo valor F = 101.6 y p-valor 2e-16, es menor que alpha (0.05), no se dispone de evidencia suficiente para considerar que los promedios son iguales, (es decir se rechaza \(H_0\)) por lo que se supone que existen por lo menos dos presiones de (\(CO_2\)) con diferencias significativas en el promedio de crecimiento de Pseudomonas fragi.Para estimar cuales son las presiones de (\(CO_2\)) cuyas medias de crecimiento son diferentes, hacemos uso de los diagramas de cajas y bigotes.
# Crear el diagrama de cajas por categorías
boxplot(ejercicio1$Crecimiento ~ ejercicio1$Presion, data = ejercicio1, col = c("green", "red", "blue","orange","pink"), ylab = "Tasa de crecimiento", xlab = "Presión")
Con base en la información proporcionada por la gráfica, se puede observar claramente que se registran diferencias significativas en los promedios de crecimiento en todas las condiciones de presión de dióxido de carbono (\(CO_2\)). Esto quiere decir, que el nivel de presión de \(CO_2\) parece ser un factor crítico que afecta significativamente al crecimiento de Pseudomonas fragi.
Como se ha rechazado la hipótesis de igualdad de medias con el test ANOVA, el interés radica en averiguar cuál o cuáles pares de medias son diferentes entre sí.
library(lsmeans)
## Loading required package: emmeans
## The 'lsmeans' package is now basically a front end for 'emmeans'.
## Users are encouraged to switch the rest of the way.
## See help('transition') for more information, including how to
## convert old 'lsmeans' objects and scripts to work with 'emmeans'.
# Calcular el LSD
lsd <- lsmeans(modelo_anova, pairwise ~ Presion)
# Ver los resultados
summary(lsd)
## $lsmeans
## Presion lsmean SE df lower.CL upper.CL
## A 59.1 1.67 45 55.8 62.5
## B 46.0 1.67 45 42.7 49.4
## C 36.5 1.67 45 33.1 39.8
## D 25.5 1.67 45 22.1 28.8
## E 16.4 1.67 45 13.1 19.8
##
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## A - B 13.10 2.36 45 5.562 <.0001
## A - C 22.69 2.36 45 9.634 <.0001
## A - D 33.67 2.36 45 14.296 <.0001
## A - E 42.70 2.36 45 18.130 <.0001
## B - C 9.59 2.36 45 4.072 0.0017
## B - D 20.57 2.36 45 8.734 <.0001
## B - E 29.60 2.36 45 12.568 <.0001
## C - D 10.98 2.36 45 4.662 0.0003
## C - E 20.01 2.36 45 8.496 <.0001
## D - E 9.03 2.36 45 3.834 0.0034
##
## P value adjustment: tukey method for comparing a family of 5 estimates
Los resultados de las comparaciones de contrastes (LSD) indican las diferencias específicas entre los pares de niveles y sugieren cuáles de ellos son significativamente diferentes entre sí.En este caso, las comparaciones de contrastes muestran las diferencias específicas entre los pares de niveles y sugieren cuáles de ellos son significativamente diferentes entre sí. Por ejemplo, las diferencias entre A-B, A-C, A-D, A-E son todas altamente significativas, lo que sugiere que estos pares son diferentes en términos de Presión. Por otro lado, las diferencias entre B-C, B-D, B-E, C-D, C-E, y D-E también son significativas, aunque en menor medida.En general, los niveles A, B, C, D y E presentan diferencias notables.
# Realizar un ANOVA
modelo_anova <- aov( Crecimiento~ Presion, data = ejercicio1)
# Realizar la prueba de Tukey
resultado_tukey <- TukeyHSD(modelo_anova)
# Ver los resultados
print(resultado_tukey)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Crecimiento ~ Presion, data = ejercicio1)
##
## $Presion
## diff lwr upr p adj
## B-A -13.10 -19.7921 -6.407896 0.0000133
## C-A -22.69 -29.3821 -15.997896 0.0000000
## D-A -33.67 -40.3621 -26.977896 0.0000000
## E-A -42.70 -49.3921 -36.007896 0.0000000
## C-B -9.59 -16.2821 -2.897896 0.0016698
## D-B -20.57 -27.2621 -13.877896 0.0000000
## E-B -29.60 -36.2921 -22.907896 0.0000000
## D-C -10.98 -17.6721 -4.287896 0.0002615
## E-C -20.01 -26.7021 -13.317896 0.0000000
## E-D -9.03 -15.7221 -2.337896 0.0034105
Luego de realizar la prueba de comparaciones múltiples de Tukey, los resultados en cada comparación indican que hay diferencias estadísticamente significativas en el crecimiento entre todos los pares de grupos (A, B, C, D y E). Los valores p ajustados cercanos a cero (menores a 0.05) en todos los casos respaldan la significancia estadística de estas diferencias.
plot(TukeyHSD(modelo_anova))
Lo mencionado anteriormente se confirma con la gráfica de la diferencia de medias(intervalos de confianza) en los distintos niveles de presión.
library(agricolae)
metodos.Duncan <-duncan.test(modelo_anova, trt = "Presion", group = T, console = T)
##
## Study: modelo_anova ~ "Presion"
##
## Duncan's new multiple range test
## for Crecimiento
##
## Mean Square Error: 27.73418
##
## Presion, means
##
## Crecimiento std r se Min Max Q25 Q50 Q75
## A 59.14 4.804674 10 1.665358 50.9 64.6 56.800 59.45 62.750
## B 46.04 5.052656 10 1.665358 35.2 50.9 43.025 48.00 49.800
## C 36.45 5.933942 10 1.665358 27.0 45.5 31.125 38.40 40.150
## D 25.47 4.483315 10 1.665358 19.2 32.7 22.650 24.25 29.425
## E 16.44 5.894480 10 1.665358 7.8 24.9 11.850 17.00 21.025
##
## Alpha: 0.05 ; DF Error: 45
##
## Critical Range
## 2 3 4 5
## 4.743560 4.988480 5.149154 5.265336
##
## Means with the same letter are not significantly different.
##
## Crecimiento groups
## A 59.14 a
## B 46.04 b
## C 36.45 c
## D 25.47 d
## E 16.44 e
De la salida anterior los tratamientos que comparten al menos una letra en la columna grupos se consideran no significativamente diferentes entre sí. En este caso los tratamientos(nivel de presión de \(CO_2\))no cumplen la anterior condición lo que indica que los grupos son estadísticamente diferentes entre sí en términos de crecimiento.
Lo cual se observa en la siguiente gráfica:
#out <- duncan.test(model, "virus", main = "yield of sweetpotato, Dealt with different virus")
plot(metodos.Duncan, variation="IQR" )
La validez de los resultados obtenidos en cualquier análisis de varianza queda condicionado a que los supuestos del modelo se cumplan. Estos supuestos son: normalidad, varianza constante (igual varianza de los tratamientos) e independencia.
library(car)
## Loading required package: carData
residuos<-residuals(modelo_anova) #Creando un objeto llamado residuos que contiene los residuos el modelo
par(mfrow=c(1,3)) #Para dividir el área del gráfico en dos partes (una fila y dos columnas)
dplot<-density(residuos) #Creando un objeto llamado dplot que recibe un Density_Plot de los residuos
plot(dplot, #Graficando el objeto dplot
main="Curva de densidad observada", #Título principal de la gráfica
xlab = "Residuos", #Etiqueta del eje x
ylab = "Densidad") #Etiqueta del eje y
polygon(dplot, #Añadiendo el poligono
col = "green", #Definiendo el color del poligono
border = "black") #Color del borde del poligono
qqPlot(residuos, #Un gráfico Cuantil-Cuantil de los residuos
pch =20, #Forma de los puntos
main="QQ-Plot de los residuos", #Título principal
xlab = "Cuantiles teoricos", #Etiqueta eje x
ylab="Cuantiles observados de los residuos") #Etiqueta eje y
## [1] 17 28
mtext(side=3, at=par("usr")[1], adj=0.7, cex=0.6, col="gray40", line=-21, #Posición del texto
text=paste("Derly Marin --", #Texto
format(Sys.time(),
"%d/%m/%Y %H:%M:%S --"), #Fecha y Hora
R.version.string)) #Versión de R
boxplot(residuos, col = c("pink"), ylab = "Residuos", main="Box-plot de los residuos")
Sin embargo, para confirmar de manera más sólida que los residuos siguen una distribución normal, se realiza la prueba de Shapiro-Wilk. Para la prueba Shapiro-Wilk para ratificar el cumplimiento del supuesto de normalidad de los residuos, evaluando las hipótesis:
\(H_0\): Los residuos de la variable crecimiento se distribuyen normalmente con media cero y varianza constante \(e i\) \(N(0,1)\).
\(H_a\): Los residuos de la variable crecimiento no siguen la distribución normal.
shapiro.test(residuals(modelo_anova)) #Prueba Shapiro-Wilk para los residuos de la variable
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_anova)
## W = 0.9627, p-value = 0.1153
Por lo tanto no hay evidencia estadística suficiente para rechazar \(H_0\), es decir se acepta la hipótesis nula, debido a que el valor de p (p-value = 0.1153) es mayor al valor del nivel de significancia (alfa=0.05), por lo que se concluye que los residuos de la variable tasa de crecimiento están normalmente distribuidos con media cero y varianza constante.
Este es el supuesto de mayor relevancia que los residuos deben cumplir para que el modelo empleado sea válido.
boxplot(residuos ~ ejercicio1$Presion,
main = "Boxplot de Residuos por presión",
xlab = "Presión",
col="lightgreen",
ylab = "Residuos")
En la siguiente gráfica, se representan los valores predichos por el modelo para la variable crecimiento en función de la raíz cuadrada de los residuos estandarizados. En esta gráfica, no se observa ninguna tendencia aparente en la distribución de los valores, lo que sugiere que no hay evidencia de incumplimiento del supuesto de homogeneidad de varianzas.
Ahora graficamos los residuos
library(car)
color_1 <-colorRampPalette(c("green", "blue", "green"))
plot(residuos, main = "Prueba de independencia", pch=20,cex = 2, col=color_1(120), ylab = "Residuos", xlab = " ")
En la grafica anterior se observan dispersos los puntos sin seguir un patron, esto es un indicio de homogeneidad de varianzas (entre más dispersos menos correlacionados)
Sin embargo, para validar de manera más sólida la homogeneidad de varianzas, se llevó a cabo la prueba de bartlett y la prueba de Levene.
Donde las hipótesis correspondientes son:
\(H_0\): Los residuos de la variable crecimiento son iguales para los ditintos niveles de presión de \(CO_2\).
\(H_a\):Existen por lo menos dos varianzas distintas para los ditintos niveles de presión de \(CO_2\). Es decir,
\(H_0: σ^2_A=σ^2_B=σ^2_C=σ^2_D=σ^2_E=σ2\)
\(H_a: σ^2_i=σ^2_j\) \(para i,j∈[A,B,C,D,E]\) \(e\) \(i≠j\)
library(stats)
bartlett.test(residuos ~ ejercicio1$Presion)
##
## Bartlett test of homogeneity of variances
##
## data: residuos by ejercicio1$Presion
## Bartlett's K-squared = 1.0701, df = 4, p-value = 0.899
De acuerdo al valor arrojado por la prueba de bartlett, valor de p (0.899) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, llegando a la misma conclusión anterior (grafica anterior).
ejercicio1$Presion <- as.factor(ejercicio1$Presion)
library(stats)
leveneTest(residuos ~ ejercicio1$Presion)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 4 0.1926 0.941
## 45
De acuerdo al valor arrojado por la prueba de levene, valor de p (0.941) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, confirmando las conclusiones obtenidas en los items anteriores.
Independencia de los residuos
Donde las hipótesis son:
\(H_0\): Los residuos entre los tratamientos son independientes.
\(H_a\):Los residuos entre los tratamientos no son independientes.
durbinWatsonTest(modelo_anova)
## lag Autocorrelation D-W Statistic p-value
## 1 -0.001714932 1.947068 0.426
## Alternative hypothesis: rho != 0
Los resultados indican que, en este análisis, los residuos están cerca de ser independientes. El valor del estadístico Durbin-Watson (DW) cercano a 2 sugiere que los residuos no muestran una autocorrelación significativa. Además, el p-value (0.456) es mayor que un nivel de significancia típico como α = 0.05, lo que sugiere que no hay suficiente evidencia para afirmar que existe una autocorrelación significativa en los residuos.
Comparación de cuatro Métodos de ensamble. Un equipo de mejora investiga el efecto de cuatro Métodos de ensamble A, B, C y D, sobre el tiempo de ensamble en minutos. En primera instancia, la estrategia experimental es aplicar cuatro veces los cuatro Metodos de ensamble en orden completamente aleatorio (las 16 pruebas en orden aleatorio). Los tiempos de ensamble obtenidos se muestran en la tabla 3.1. Si se usa el diseño completamente al azar (DCA), se supone que, además del Metodo de ensamble, no existe ningún otro factor que influya de manera significativa sobre la variable de respuesta (tiempo de ensamble).
# install.packages("readxl")#solo una vez
library(readxl)
ejercicio2 <- read_excel("C:/Users/LENOVO/Desktop/Ejercicio2.xlsx")
ejercicio2 # para visualizar los datos
## # A tibble: 16 × 2
## Desgaste Metodo
## <dbl> <chr>
## 1 6 A
## 2 8 A
## 3 7 A
## 4 8 A
## 5 7 B
## 6 9 B
## 7 10 B
## 8 8 B
## 9 11 C
## 10 16 C
## 11 11 C
## 12 13 C
## 13 10 D
## 14 12 D
## 15 11 D
## 16 9 D
conteo_valoresmetodo <- table(ejercicio2$Metodo)
conteo_valoresmetodo
##
## A B C D
## 4 4 4 4
Dado que el número de observaciones por tratamiento es el mismo, se puede concluir que es un diseño balanceado.
#install.packages("summarytools")# solo una vez
library(summarytools)
summarytools::descr(ejercicio2[,1])# todas las filas primera columna datos[,1]
## Descriptive Statistics
## ejercicio2$Desgaste
## N: 16
##
## Desgaste
## ----------------- ----------
## Mean 9.75
## Std.Dev 2.57
## Min 6.00
## Q1 8.00
## Median 9.50
## Q3 11.00
## Max 16.00
## MAD 2.22
## IQR 3.00
## CV 0.26
## Skewness 0.68
## SE.Skewness 0.56
## Kurtosis -0.11
## N.Valid 16.00
## Pct.Valid 100.00
A partir de los resultados que arroja el programa podemos concluir: el promedio de desgaste es de 9.75, con desviación estándar de 2.57, el valor mínimo de desgaste es de 6, el valor máximo de desgaste es de 16, el 50% de las observaciones(8) presentaron un desgaste entre 6 y 9.5, mientras que el restante 50% (8) presento un desgaste entre9.5 y 16, el coeficiente de asimetría fue de 0.68 presentando una asimetría positiva leve, además el coeficiente de curtosis fue de -0.11, indicando que la distribución es levemente platicurtica.
# Calcular estadísticas descriptivas por categoría
resultados_descriptivos <- aggregate(Desgaste ~ Metodo, data = ejercicio2, summary)
# Imprimir los resultados descriptivos
print(resultados_descriptivos)
## Metodo Desgaste.Min. Desgaste.1st Qu. Desgaste.Median Desgaste.Mean
## 1 A 6.00 6.75 7.50 7.25
## 2 B 7.00 7.75 8.50 8.50
## 3 C 11.00 11.00 12.00 12.75
## 4 D 9.00 9.75 10.50 10.50
## Desgaste.3rd Qu. Desgaste.Max.
## 1 8.00 8.00
## 2 9.25 10.00
## 3 13.75 16.00
## 4 11.25 12.00
LLevando a cabo el análisis descriptivo por tratamiento se obtuvieron los sigueintes resultados:
Para el método de ensamble A: el promedio de desgaste es de 7.25, el valor mínimo de desgaste es 6, el valor máximo de desgaste es 8, el 50% de las observaciones(2) presentaron un desgaste entre 6 y 7.5, mientras que el restante 50% (2) presento un desgaste entre 7.5 y 8.
Para el método de ensamble B: el promedio de desgaste es de 8.5, el valor mínimo de desgaste es 7, el valor máximo de desgaste es 10, el 50% de las observaciones(2) presentaron un desgaste entre 7 y 8.5, mientras que el restante 50% (2) presento un desgaste entre 8.5 y 10.
Para el método de ensamble C: el promedio de desgaste es de 12.75, el valor mínimo de desgaste es 11, el valor máximo de desgaste es 16, el 50% de las observaciones(2) presentaron un desgaste entre 11 y 12, mientras que el restante 50% (2) presento un desgaste entre 12 y 16.
Para el método de ensamble D: el promedio de desgaste es 10.5, el valor mínimo de desgaste es 9, el valor máximo de desgaste es 12, el 50% de las observaciones(2) presentaron un desgaste entre 9 y 10.5, mientras que el restante 50% (2) presento un desgaste entre 10.5 y 12.
A continuación se llevará a cabo el ANOVA de un factor o modelo factorial de un solo factor el cual nos permite estudiar si existen diferencias significativas entre las medias de desgaste de los cuatro métodos de ensamble.
Las hipótesis contrastadas en el ANOVA son: No hay diferencias entre las medias de los diferentes métodos de ensamble y existen por lo menos dos métodos de ensamble con diferencias significativas en el promedio de desgaste. Que se plantearían de la siguiente forma:
\(H_0:μ_A=μ_B=μ_C=μ_D=μ\)
\(H_a:μ_i=μ_j\) para \(i≠j\) donde \(i,j=A,B,C,D.\) respectivamente.
# Realizar el ANOVA
modelo_anova <- aov(Desgaste ~ Metodo, data = ejercicio2)
# Resumen del ANOVA
resumen_anova <- summary(modelo_anova)
# Imprimir el resumen del ANOVA
print(resumen_anova)
## Df Sum Sq Mean Sq F value Pr(>F)
## Metodo 3 69.5 23.167 9.424 0.00177 **
## Residuals 12 29.5 2.458
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Los resultados del análisis de varianza (ANOVA) ha resultado significativo valor F = 9.924 y p-valor 0.00177, es manor que alpha (0.05), no se dispone de evidencia suficiente para considerar que los promedios son iguales, (es decir se rechaza \(H_0\)) por lo que se supone que existen por lo menos dos métodos de ensamble con diferencias significativas en el promedio de desgaste. Para estimar cuales son los métodos de ensamble cuyos medias de desgaste son diferentes, hacemos uso de los diagramas de cajas y bigotes por cada método o tratamiento.
# Crear el diagrama de cajas por categorías
boxplot(ejercicio2$Desgaste ~ ejercicio2$Metodo, data = ejercicio2, col = c("green", "red", "blue","orange"), ylab = "Tiempo de ensamblaje (min)", xlab = "Método")
De la grafica anterior se puede evidenciar que existen diferencias en los promedios de desgaste entre el método de ensamble A y el método de ensamble C Y D, también se puede observar diferencias entre el método B y C (no se traslapan las gráficas).
Como se ha rechazado la hipótesis de igualdad de medias con el test ANOVA, el interés está en averiguar cuál o cuáles pares de medias son diferentes entre sí.
# Realizar la prueba de LSD
lsd <- lsmeans(modelo_anova, pairwise ~ Metodo)
# Mostrar los resultados de la prueba de LSD
print(lsd)
## $lsmeans
## Metodo lsmean SE df lower.CL upper.CL
## A 7.25 0.784 12 5.54 8.96
## B 8.50 0.784 12 6.79 10.21
## C 12.75 0.784 12 11.04 14.46
## D 10.50 0.784 12 8.79 12.21
##
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## A - B -1.25 1.11 12 -1.127 0.6805
## A - C -5.50 1.11 12 -4.961 0.0016
## A - D -3.25 1.11 12 -2.931 0.0533
## B - C -4.25 1.11 12 -3.833 0.0110
## B - D -2.00 1.11 12 -1.804 0.3181
## C - D 2.25 1.11 12 2.029 0.2309
##
## P value adjustment: tukey method for comparing a family of 4 estimates
Los resultados de las comparaciones de contrastes (LSD) indican las diferencias específicas entre los pares de niveles y sugieren cuáles de ellos son significativamente diferentes entre sí. En el caso de las comparaciones A - B, A - D, B-D y C - D, no se encontraron diferencias significativas, ya que sus valores p son mayores que el nivel de significancia de 0.05. Por otro lado, se identificaron diferencias significativas en las comparaciones A - C y B - C, ya que los valores p son menores al nivel de significancia preestablecido.
# Realizar un ANOVA
modelo_anova <- aov(Desgaste~ Metodo, data = ejercicio2)
# Realizar la prueba de Tukey
resultado_tukey <- TukeyHSD(modelo_anova)
# Ver los resultados
print(resultado_tukey)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Desgaste ~ Metodo, data = ejercicio2)
##
## $Metodo
## diff lwr upr p adj
## B-A 1.25 -2.04155503 4.541555 0.6804513
## C-A 5.50 2.20844497 8.791555 0.0016206
## D-A 3.25 -0.04155503 6.541555 0.0533380
## C-B 4.25 0.95844497 7.541555 0.0110423
## D-B 2.00 -1.29155503 5.291555 0.3181239
## D-C -2.25 -5.54155503 1.041555 0.2309373
Luego de realizar las prueba de comparaciones múltiples de Tukey, se concluye que existen diferencia significativas entre los métodos de ensamble A y C con una diferencia de 5.50 cuyo intervalo de confianza del 95% para la diferencia es (2.2, 8.8) y un p-valor de 0.002, lo que resulta significativo (menor de 0.05), es decir no se dispone de evidencia suficiente para considerar que los promedios de desgaste en los dos tratamiento son iguales.
También se estima que existen diferencia significativas entre los métodos de ensamble C y B con una diferencia de 4.25 cuyo intervalo de confianza del 95% para la diferencia es (0.95, 7.54) y un p-valor de 0.011, lo que resulta significativo (menor de 0.05), es decir no se dispone de evidencia suficiente para considerar que los promedios de desgaste en los dos tratamiento son iguales.
Para la comparación de los demás tratamientos(métodos de ensamble) no resulto significativo.
plot(TukeyHSD(modelo_anova))
Lo mencionado anteriormente se confirma con la grafica de la diferencia
de medias(intervalos de confianza) en los distintos niveles del
método.
library(agricolae)
metodos.Duncan <-duncan.test(modelo_anova, trt = "Metodo", group = T, console = T)
##
## Study: modelo_anova ~ "Metodo"
##
## Duncan's new multiple range test
## for Desgaste
##
## Mean Square Error: 2.458333
##
## Metodo, means
##
## Desgaste std r se Min Max Q25 Q50 Q75
## A 7.25 0.9574271 4 0.7839537 6 8 6.75 7.5 8.00
## B 8.50 1.2909944 4 0.7839537 7 10 7.75 8.5 9.25
## C 12.75 2.3629078 4 0.7839537 11 16 11.00 12.0 13.75
## D 10.50 1.2909944 4 0.7839537 9 12 9.75 10.5 11.25
##
## Alpha: 0.05 ; DF Error: 12
##
## Critical Range
## 2 3 4
## 2.415602 2.528441 2.596810
##
## Means with the same letter are not significantly different.
##
## Desgaste groups
## C 12.75 a
## D 10.50 ab
## B 8.50 bc
## A 7.25 c
De la salida anterior los tratamientos que comparten al menos una letra en la columna grupos se consideran no significativamente diferentes entre sí. Estos tratamientos forman grupos estadísticamente similares.
Mientras que los tratamientos que tienen letras diferentes en la columna grupos se consideran significativamente diferentes entre sí. Si un tratamiento tiene una letra diferente de otro, significa que hay una diferencia estadísticamente significativa en sus medias, para nuestro caso los tratamientos(métodos de ensamble) C con A, B con C y D con A, no comparten letra por lo tanto se consideran significativamente diferentes.
Lo cual se observa en la siguiente gráfica:
#out <- duncan.test(model, "virus", main = "yield of sweetpotato, Dealt with different virus")
plot(metodos.Duncan, variation="IQR" )
La validez de los resultados obtenidos en cualquier análisis de varianza queda condicionado a que los supuestos del modelo se cumplan. Estos supuestos son: normalidad, varianza constante (igual varianza de los tratamientos) e independencia.
library(car)
residuos<-residuals(modelo_anova) #Creando un objeto llamado residuos que contiene los residuos el modelo
par(mfrow=c(1,3)) #Para dividir el área del gráfico en dos partes (una fila y dos columnas)
dplot<-density(residuos) #Creando un objeto llamado dplot que recibe un Density_Plot de los residuos
plot(dplot, #Graficando el objeto dplot
main="Curva de densidad observada", #Título principal de la gráfica
xlab = "Residuos", #Etiqueta del eje x
ylab = "Densidad") #Etiqueta del eje y
polygon(dplot, #Añadiendo el poligono
col = "blue", #Definiendo el color del poligono
border = "black") #Color del borde del poligono
qqPlot(residuos, #Un grafico Cuantil-Cuantil de los residuos
pch =20, #Forma de los puntos
main="QQ-Plot de los residuos", #Titulo principal
xlab = "Cuantiles teóricos", #Etiqueta eje x
ylab="Cuantiles observados de los residuos") #Etiqueta eje y
## [1] 10 9
mtext(side=3, at=par("usr")[1], adj=0.7, cex=0.6, col="gray40", line=-21, #Posición del texto
text=paste("Derly Marin --", #Texto
format(Sys.time(),
"%d/%m/%Y %H:%M:%S --"), #Fecha y Hora
R.version.string)) #Versión de R
boxplot(residuos, col = c("orange"), ylab = "residuos", main="Box-plot de los residuos")
Sin embargo, para confirmar de manera más sólida que los residuos siguen una distribución normal, se realiza la prueba de Shapiro-Wilk.
Para la prueba Shapiro-Wilk para ratificar el cumplimiento del supuesto de normalidad de los residuos, evaluando las hipótesis:
\(H_0\): Los residuos de la variable desgaste se distribuyen normalmente con media cero y varianza constante \(e i\) \(N(0,1)\).
\(H_a\): los residuos de la variable desgaste no siguen la distribución normal.
shapiro.test(residuals(modelo_anova)) #Prueba Shapiro-Wilk para los residuos de la variable
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_anova)
## W = 0.93389, p-value = 0.2808
Por lo tanto no hay evidencia estadística suficiente para rechazar \(H_0\),es decir, se acepta la hipótesis nula, debido a que el valor de p (p-value = 0.2808) es mayor al valor del nivel de significancia (α =0.05), por lo que se concluye que los residuos de la variable porcentaje de parasitismo están normalmente distribuidos con media cero y varianza constante.
Este es el supuesto de mayor relevancia que los residuos deben cumplir para que el modelo empleado sea válido.
boxplot(residuos ~ ejercicio2$Metodo,
main = "Boxplot de Residuos por método de ensamble",
xlab = "Método de ensamble",
col="lightgreen",
ylab = "Residuos")
En la siguiente gráfica, se representan los valores predichos por el modelo para la variable desgaste del material en función de la raíz cuadrada de los residuos estandarizados. En esta gráfica, no se observa ninguna tendencia aparente en la distribución de los valores, lo que sugiere que no hay evidencia de incumplimiento del supuesto de homogeneidad de varianzas.
Ahora graficamos los residuos
library(car)
color_1 <-colorRampPalette(c("violet", "blue", "violet"))
plot(residuos, main = "Prueba de independencia", pch=20,cex = 2, col=color_1(120), ylab = "Residuos", xlab = " ")
En la grafica anterior se observan dispersos los puntos sin seguir un
patrón, esto es un indicio de homogeneidad de varianzas (entre más
dispersos menos correlacionados).
Sin embargo, para validar de manera más sólida la homogeneidad de varianzas, se llevó a cabo la prueba de bartlett y la prueba de Levene.
Donde las hipótesis correspondientes son: \(H_0\): Los residuos de la variable crecimiento son iguales para los ditintos niveles de presión de \(CO_2\).
\(H_a\):Existen por lo menos dos varianzas distintas para los ditintos niveles de presión de \(CO_2\). Es decir,
\(H_0: σ^2_A=σ^2_B=σ^2_C=σ^2_D=σ2\)
\(H_a: σ^2_i=σ^2_j\) \(para i,j∈[A,B,C,D]\) \(e\) \(i≠j\)
# Verifica la longitud de las variables
length(residuos)
## [1] 16
length(ejercicio2$Metodo)
## [1] 16
# Verifica la cantidad de valores faltantes en cada variable
sum(is.na(residuos))
## [1] 0
sum(is.na(ejercicio2$Metodo))
## [1] 0
# Calcular los residuos del modelo ANOVA
residuos <- residuals(modelo_anova)
# Realizar la prueba de Bartlett para evaluar la homogeneidad de varianza
bartlett.test(residuos ~ ejercicio2$Metodo)
##
## Bartlett test of homogeneity of variances
##
## data: residuos by ejercicio2$Metodo
## Bartlett's K-squared = 2.4853, df = 3, p-value = 0.4779
De acuerdo al valor arrojado por la prueba de bartlett, valor de p (0.4779) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, llegando a la misma conclusión anterior (grafica anterior).
library(stats)
leveneTest(residuos ~ ejercicio2$Metodo)
## Warning in leveneTest.default(y = y, group = group, ...): group coerced to
## factor.
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 3 0.9474 0.4485
## 12
De acuerdo al valor arrojado por la prueba de levene, valor de p (0.4485) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, confirmando las conclusiones obtenidas en los items anteriores.
Independencia de los residuos
Prueba de Durbin Watson
Donde las hipótesis son:
\(H_0\): Los residuos entre los tratamientos son independientes.
\(H_a\):Los residuos entre los tratamientos no son independientes.
durbinWatsonTest(modelo_anova)
## lag Autocorrelation D-W Statistic p-value
## 1 -0.5084746 2.887712 0.272
## Alternative hypothesis: rho != 0
Puesto que el valor de DW es aproximadamente igual a 2 (1-Rho) donde Rho es la autocorrelación de la muestra de los residuos.
Se debe tener en cuenta que si el valor del estadístico Durbin Watson (DW) está próximo a 2 entonces los residuos no están autocorrelacionados. Si su valor es 0 hay autocorrelación perfecta positiva. Si tiene un valor de 4 existe autocorrelación perfecta negativa.
Al realizar la prueba de independencia de residuos para la variable desgaste se determinó que los residuos no están correlacionados, debido a que el DW está próximo a 2 y el valor de p (p-value = 0.252) es superior al nivel de significancia de 5% (α=0.05) por lo que se concluye que existe independencia de los residuos.
Se hace un estudio sobre la efectividad de tres marcas de spray para matar moscas.Para ello, cada producto se aplica a un grupo de 100 moscas, y se cuenta el número de moscas muertas expresado en porcentajes. Se hacen seis réplicas y los resultados obtenidos se muestran a continuación.
\(y_{ij}\)=\(\mu\)+\(\tau_i\)+\(\epsilon_j\)
# install.packages("readxl")#solo una vez
library(readxl)
ejercicio3 <- read_excel("C:/Users/LENOVO/Desktop/Ejercicio3.xlsx")
ejercicio3 # para visualizar los datos
## # A tibble: 18 × 2
## Porcentaje Spray
## <dbl> <chr>
## 1 72 A
## 2 65 A
## 3 67 A
## 4 75 A
## 5 62 A
## 6 73 A
## 7 55 B
## 8 59 B
## 9 68 B
## 10 70 B
## 11 53 B
## 12 50 B
## 13 64 C
## 14 74 C
## 15 61 C
## 16 58 C
## 17 51 C
## 18 69 C
conteo_valoresspray <- table(ejercicio3$Spray)
conteo_valoresspray
##
## A B C
## 6 6 6
Dado que el número de observaciones por tratamiento es el mismo, se puede concluir que es un diseño balanceado.
#install.packages("summarytools")# solo una vez
library(summarytools)
summarytools::descr(ejercicio3[,1])# todas las filas primera columna datos[,1]
## Descriptive Statistics
## ejercicio3$Porcentaje
## N: 18
##
## Porcentaje
## ----------------- ------------
## Mean 63.67
## Std.Dev 8.01
## Min 50.00
## Q1 58.00
## Median 64.50
## Q3 70.00
## Max 75.00
## MAD 8.90
## IQR 11.50
## CV 0.13
## Skewness -0.25
## SE.Skewness 0.54
## Kurtosis -1.33
## N.Valid 18.00
## Pct.Valid 100.00
A partir de los resultados que arroja el programa podemos concluir: el promedio de porcentaje de moscas muertas es de 63.67, con desviación estándar de 8.01, el valor mínimo de porcentaje de moscas muertas es de 50, el valor máximo de porcentaje de moscas muertas es de 75, el 50% de las observaciones(9) presentaron un porcentaje entre 50 y 64.5, mientras que el restante 50% (9) presento un desgaste entre 64.5 y 75, el coeficiente de asimetría fue de -0.25 presentando una asimetría negativa leve, además el coeficiente de curtosis fue de -1.33, indicando que la distribución es levemente platicúrtica (hay una menor concentración de datos en torno a la media).
# Calcular estadísticas descriptivas por categoría
resultados_descriptivos <- aggregate(Porcentaje ~ Spray, data = ejercicio3, summary)
# Imprimir los resultados descriptivos
print(resultados_descriptivos)
## Spray Porcentaje.Min. Porcentaje.1st Qu. Porcentaje.Median Porcentaje.Mean
## 1 A 62.00000 65.50000 69.50000 69.00000
## 2 B 50.00000 53.50000 57.00000 59.16667
## 3 C 51.00000 58.75000 62.50000 62.83333
## Porcentaje.3rd Qu. Porcentaje.Max.
## 1 72.75000 75.00000
## 2 65.75000 70.00000
## 3 67.75000 74.00000
LLevando a cabo el análisis descriptivo por tratamiento se obtuvieron los siguentes resultados:
Para el spray A: el promedio de porcentaje de moscas muertas es de 69, el valor mínimo de porcentaje es 62, el valor máximo de porcentaje es 75, el 50% de las observaciones(3) presentaron un porcentaje entre 62 y 69.5, mientras que el restante 50% (3) presento un porcentaje entre 69.5 y 75.
Para el spray B: el promedio de porcentaje de moscas muertas es de 59.17, el valor mínimo de porcentaje es 50, el valor máximo de porcentaje es 70, el 50% de las observaciones(3) presentaron un porcentaje entre 50 y 57, mientras que el restante 50% (3) presento un porcentaje entre 57 y 70.
Para el spray C: el promedio de porcentaje de moscas muertas es de 62.8333, el valor mínimo de porcentaje es 51, el valor máximo de porcentaje es 74, el 50% de las observaciones(3) presentaron un porcentaje entre 51 y 62.5, mientras que el restante 50% (3) presento un porcentaje entre 62.5 y 74.
A continuación se llevará a cabo el ANOVA de un factor o modelo factorial de un solo factor el cual nos permite estudiar si existen diferencias significativas entre las medias de porcentaje de moscas muertas de las tres marcas de spray.
Las hipótesis contrastadas en el ANOVA son: No hay diferencias entre las medias de los diferentes spray y existen por lo menos dos spray con diferencias significativas en el promedio de porcentje de moscas muertas. Que se plantearían de la siguiente forma:
\(H_0:μ_A=μ_B=μ_C=μ\)
\(H_a:μ_i=μ_j\) para \(i≠j\) donde \(i,j=A,B,C\) respectivamente.
# Realizar el ANOVA
modelo_anova <- aov(Porcentaje ~ Spray, data = ejercicio3)
# Resumen del ANOVA
resumen_anova <- summary(modelo_anova)
# Imprimir el resumen del ANOVA
print(resumen_anova)
## Df Sum Sq Mean Sq F value Pr(>F)
## Spray 2 296.3 148.17 2.793 0.0931 .
## Residuals 15 795.7 53.04
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
En este caso el ANOVA, el p-valor 0.0931 es mayor que el nivel de significancia ( (α= 0.05), es decir, que se acepta la hipótesis nula (\(H_0\)).Esto significa que según los datos y el nivel de significancia seleccionado, no se tiene suficiente evidencia para rechazar la hipótesis nula y, por lo tanto, no se puede concluir que hay una diferencia o efecto significativo entre las marcas de spray.
# Realizar el ANOVA
modelo_anova <- aov(Porcentaje ~ Spray, data = ejercicio3)
library(lsmeans)
# Realizar la prueba de LSD
lsd <- lsmeans(modelo_anova, pairwise ~ Spray)
# Mostrar los resultados de la prueba de LSD
print(lsd)
## $lsmeans
## Spray lsmean SE df lower.CL upper.CL
## A 69.0 2.97 15 62.7 75.3
## B 59.2 2.97 15 52.8 65.5
## C 62.8 2.97 15 56.5 69.2
##
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## A - B 9.83 4.2 15 2.339 0.0808
## A - C 6.17 4.2 15 1.467 0.3341
## B - C -3.67 4.2 15 -0.872 0.6655
##
## P value adjustment: tukey method for comparing a family of 3 estimates
No es posible identificar un spray que sea estadísticamente superior en términos de efectividad. Esto significa que, según los resultados del análisis, todas las marcas de spray podrían ser igualmente efectivas en promedio, al menos dentro del margen de error permitido por el nivel de significancia.
# Calcular intervalos de confianza al 95% para cada marca (A, B y C)
intervalos_confianza <- confint(modelo_anova, level = 0.95)
nombres_marcas <- c("A", "B", "C")
rownames(intervalos_confianza) <- nombres_marcas
print(intervalos_confianza)
## 2.5 % 97.5 %
## A 62.66248 75.3375206
## B -18.79594 -0.8707257
## C -15.12927 2.7959409
library(gplots)
##
## Attaching package: 'gplots'
## The following object is masked from 'package:stats':
##
## lowess
# Crear un gráfico de medias con intervalos de confianza (95%)
plotmeans(Porcentaje ~ Spray, data = ejercicio3,
xlab = "Spray", ylab = "Porcentaje",
main = "Gráfico de Medias con 95% IC",
cex.axis = 0.6)
# Crear el diagrama de cajas por categorías
boxplot(ejercicio3$Porcentaje ~ ejercicio3$Spray, data = ejercicio3, col = c("pink", "lightblue", "violet"), ylab = "Porcentaje", xlab = "Marca de spray")
En el diagrama de cajas se observa que los tres tratamientos(A,B y C) se superponen entre sí, esto sugiere que no hay evidencia clara de diferencias significativas en la distribución de los datos en las marcas de los spray.
En un centro de investigación se realiza un estudio para comparar varios tratamientos que, al aplicarse previamente a los frijoles crudos, reducen su tiempo de cocción. Estos tratamientos son a base de bicarbonato de sodio (NaHCO3) y cloruro de sodio o sal común (NaCl). El primer tratamiento es el de control, que consiste en no aplicar ningún tratamiento. El tratamiento \(T_{B}\) es el remojo en agua con bicarbonato de sodio, el \(T_{C}\) es remojar en agua con sal común y el \(T_{D}\) es remojar en agua con una combinación de ambos ingredientes en proporciones iguales. La variable de respuesta es el tiempo de cocción en minutos. Los datos se muestran en la siguiente tabla:
El experimentador debe utilizar un Diseño completamente al azar (DCA). En lo que respecta a aleatorizacion de experimento se debe tener en cuenta que el experimentador asigne de forma aleatoria de los tratamientos a los grupos de frijoles, además de aleatorizar el lote de frijoles para que las diferencias en la calidad de los frijoles no afecte el resultado.En cuanto al material experimental, es esencial que los recipientes y utensilios utilizados sean consistentes en términos de material y estado, de manera que todos sean idénticos.
Deben estar fijos factores como la variedad de los frijoles utilizados, la temperatura ambiente, el nivel de temperatura, el tiempo de remojo de los frijoles, la cantidad y calidad de agua empleada en el proceso y el tamaño de los grupos experimentales.
Las hipótesis contrastadas en el ANOVA son:Los tratamientos no reducen el tiempo de cocción a los frijoles y existen por lo menos dos tratamientos con diferencias significativas en la reducción del tiempo de cocción de los frijoles. Que plantearían de la siguiente forma:
\(H_0:μ_A=μ_B=μ_C=μ_D=μ\).
\(H_a:μ_i=μ_j\) para \(i≠j\) donde
\(i,j=A,B,C,D\) respectivamente.
# install.packages("readxl")#solo una vez
library(readxl)
ejercicio4 <- read_excel("C:/Users/LENOVO/Desktop/Ejercicio4.xlsx")
ejercicio4 # para visualizar los datos
## # A tibble: 28 × 2
## Tiempo Tratamiento
## <dbl> <chr>
## 1 213 A
## 2 214 A
## 3 204 A
## 4 208 A
## 5 212 A
## 6 200 A
## 7 207 A
## 8 76 B
## 9 85 B
## 10 74 B
## # ℹ 18 more rows
conteo_valorestratamiento <- table(ejercicio4$Tratamiento)
conteo_valorestratamiento
##
## A B C D
## 7 7 7 7
Dado que el número de observaciones por tratamiento es el mismo, se puede concluir que es un diseño balanceado.
#install.packages("summarytools")# solo una vez
library(summarytools)
summarytools::descr(ejercicio4[,1])# todas las filas primera columna datos[,1]
## Descriptive Statistics
## ejercicio4$Tiempo
## N: 28
##
## Tiempo
## ----------------- --------
## Mean 109.18
## Std.Dev 59.00
## Min 55.00
## Q1 74.50
## Median 82.00
## Q3 146.00
## Max 214.00
## MAD 13.34
## IQR 44.25
## CV 0.54
## Skewness 1.02
## SE.Skewness 0.44
## Kurtosis -0.86
## N.Valid 28.00
## Pct.Valid 100.00
A partir de los resultados que arroja el programa podemos concluir: el promedio de tiempo de cocción es de 109.18, con desviación estándar de 59, el valor mínimo de tiempo de cocción es de 55, el valor máximo de tiempo es de 214, el 50% de las observaciones(14) presentaron un tiempo entre 55 y 82, mientras que el restante 50% (14) presento un desgaste entre 82 y 214, el coeficiente de asimetría fue de 1.02 presentando una asimetría positiva, además el coeficiente de curtosis fue de -0.86, indicando que la distribución es platicúrtica.
# Calcular estadísticas descriptivas por categoría
resultados_descriptivos <- aggregate(Tiempo ~ Tratamiento, data = ejercicio4, summary)
# Imprimir los resultados descriptivos
print(resultados_descriptivos)
## Tratamiento Tiempo.Min. Tiempo.1st Qu. Tiempo.Median Tiempo.Mean
## 1 A 200.00000 205.50000 208.00000 208.28571
## 2 B 74.00000 75.50000 78.00000 78.85714
## 3 C 55.00000 62.00000 63.00000 64.00000
## 4 D 79.00000 83.00000 85.00000 85.57143
## Tiempo.3rd Qu. Tiempo.Max.
## 1 212.50000 214.00000
## 2 82.00000 85.00000
## 3 65.50000 75.00000
## 4 88.50000 92.00000
LLevando a cabo el análisis descriptivo por tratamiento se obtuvieron los siguientes resultados:
Para el tratamiento A: el promedio de tiempo de cocción es de 208.28571, el valor mínimo de tiempo de cocción es 200, el valor máximo de tiempo de cocción es 214, el 50% de las observaciones presentaron un tiempo de cocción entre 200 y 208, mientras que el restante 50% presento un tiempo de cocción entre 208 y 214.
Para el tratamiento B: el promedio de tiempo de cocción es de 78.85714, el valor mínimo de tiempo de cocción es 74, el valor máximo de tiempo de cocción es 85, el 50% de las observaciones presentaron un tiempo de cocción entre 74 y 78, mientras que el restante 50% presento un tiempo de cocción entre 78 y 85.
Para el tratamiento C: el promedio de tiempo de cocción es de 64, el valor mínimo de tiempo de cocción es 55, el valor máximo de tiempo de cocción es 75, el 50% de las observaciones presentaron un tiempo de cocción entre 55 y 63, mientras que el restante 50% presento un tiempo de cocción entre 63 y 75.
Para el tratamiento D: el promedio de tiempo de cocción es de 85.57143, el valor mínimo de tiempo de cocción es 79, el valor máximo de tiempo de cocción es 92, el 50% de las observaciones presentaron un tiempo de cocción entre 79 y 85, mientras que el restante 50% presento un tiempo de cocción entre 85 y 92.
# Realizar el ANOVA
modelo_anova <- aov(Tiempo ~ Tratamiento, data = ejercicio4)
# Resumen del ANOVA
resumen_anova <- summary(modelo_anova)
# Imprimir el resumen del ANOVA
print(resumen_anova)
## Df Sum Sq Mean Sq F value Pr(>F)
## Tratamiento 3 93380 31127 1233 <2e-16 ***
## Residuals 24 606 25
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Los resultados del análisis de varianza (ANOVA) arrojan un p-valor 2e-16, es menor al nivel de significancia (α=0.05). En consecuencia, no se dispone de evidencia suficiente para considerar que los promedios son iguales, es decir, se rechaza \(H_0\), por lo que se supone que existen por lo menos dos tratamientos con diferencias significativas en el promedio tiempo de cocción. Para estimar cuales son los tratamientos cuyos medias de tiempo de cocción son diferentes, hacemos uso de los diagramas de cajas y bigotes por cada método o tratamiento.
library(gplots)
# Crear un gráfico de medias con intervalos de confianza (95%)
plotmeans(Tiempo ~ Tratamiento, data = ejercicio4,
xlab = "Tratamiento", ylab = "Tiempo",
main = "Grafico de Medias con 95% IC", cex.axis = 0.6)
De la gráfica anterior, se puede apreciar que la media del tratamiento A es la más alta, ya que se trata del grupo de control. En contraste, el tratamiento C (agua con sal común) muestra la media más baja. El tratamiento B (agua con bicarbonato de sodio) se ubica en un nivel intermedio, mientras que el tratamiento D (combinación de ambos compuestos) se encuentra entre los tratamientos A y C en términos de su efectividad para reducir el tiempo de cocción.
# Crear el diagrama de cajas por categorías
boxplot(ejercicio4$Tiempo ~ ejercicio4$Tratamiento, data = ejercicio4, col = c("pink", "lightblue", "violet","lightgreen"), ylab = "Tiempo de cocción", xlab = "Tratamiento")
El diagrama de cajas revela notables diferencias significativas en los
promedios de los tratamientos en todas los tiempos de cocción (no se
traslapan las cajas). Esto quiere decir, que el cada tratamiento parece
ser un factor crítico que afecta significativamente el tiempo de cocción
de los frijoles.
Como se ha rechazado la hipótesis de igualdad de medias con el test ANOVA, el interés radica en averiguar cuál o cuáles pares de medias son diferentes entre sí.
# Realizar la prueba de LSD
lsd <- lsmeans(modelo_anova, pairwise ~ Tratamiento)
# Mostrar los resultados de la prueba de LSD
print(lsd)
## $lsmeans
## Tratamiento lsmean SE df lower.CL upper.CL
## A 208.3 1.9 24 204.4 212.2
## B 78.9 1.9 24 74.9 82.8
## C 64.0 1.9 24 60.1 67.9
## D 85.6 1.9 24 81.7 89.5
##
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## A - B 129.43 2.69 24 48.187 <.0001
## A - C 144.29 2.69 24 53.719 <.0001
## A - D 122.71 2.69 24 45.688 <.0001
## B - C 14.86 2.69 24 5.531 0.0001
## B - D -6.71 2.69 24 -2.500 0.0855
## C - D -21.57 2.69 24 -8.031 <.0001
##
## P value adjustment: tukey method for comparing a family of 4 estimates
Los resultados de las comparaciones de contrastes (LSD) indican las diferencias específicas entre los pares de niveles y sugieren cuáles de ellos son significativamente diferentes entre sí.En todos los casos de comparación entre A y los demás niveles (B, C y D), se observaron diferencias altamente significativas. Además, se identificaron diferencias significativas tanto entre B y C como entre C y D (p-valor es menor a 0.05).Sin embargo, la comparación entre B-D no alcanzó significancia estadística, ya que el valor p (0.0855) es ligeramente mayor que 0.05.
# Realizar un ANOVA
modelo_anova <- aov(Tiempo ~ Tratamiento, data = ejercicio4)
# Realizar la prueba de Tukey
resultado_tukey <- TukeyHSD(modelo_anova)
# Ver los resultados
print(resultado_tukey)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Tiempo ~ Tratamiento, data = ejercicio4)
##
## $Tratamiento
## diff lwr upr p adj
## B-A -129.428571 -136.8380353 -122.019108 0.0000000
## C-A -144.285714 -151.6951781 -136.876250 0.0000000
## D-A -122.714286 -130.1237496 -115.304822 0.0000000
## C-B -14.857143 -22.2666067 -7.447679 0.0000610
## D-B 6.714286 -0.6951781 14.123750 0.0854646
## D-C 21.571429 14.1619647 28.980892 0.0000002
Luego de realizar las prueba de comparaciones múltiples de Tukey, se concluye que existen diferencias significativas de los tratamientos B,C y D con el tratamiento A, ya que se trata del control y esto se evidencia en cada p-valor de estas comparaciones es menor (α=0.05),es decir,no se dispone de evidencia suficiente para considerar que los promedios de tiempo en los dos tratamientos son iguales para cada caso.
También se estima que existen diferencia significativas entre los tratamientos C y B con una diferencia de -14.9 cuyo intervalo de confianza del 95% para la diferencia es (-22.27,-7.45) y un p-valor de 0.00006, lo que resulta significativo (menor de 0.05), es decir no se dispone de evidencia suficiente para considerar que los promedios de desgaste en los dos tratamiento son iguales.
También se estima que existen diferencia significativas entre los tratamientos D y C con una diferencia de 21.6 cuyo intervalo de confianza del 95% para la diferencia es (14.16,28.9) y un p-valor de 0.0000002, lo que resulta significativo (menor de 0.05), es decir, no se dispone de evidencia suficiente para considerar que los promedios de desgaste en los dos tratamiento son iguales.
Para la comparación de los tratamientos(D-B) no resulto significativo.
plot(resultado_tukey)
Lo mencionado anteriormente se confirma con la gráfica de la diferencia de medias(intervalos de confianza) en los distintos niveles del tratamiento.
library(agricolae)
metodos.Duncan <-duncan.test(modelo_anova, trt = "Tratamiento", group = T, console = T)
##
## Study: modelo_anova ~ "Tratamiento"
##
## Duncan's new multiple range test
## for Tiempo
##
## Mean Square Error: 25.25
##
## Tratamiento, means
##
## Tiempo std r se Min Max Q25 Q50 Q75
## A 208.28571 5.122313 7 1.899248 200 214 205.5 208 212.5
## B 78.85714 4.180453 7 1.899248 74 85 75.5 78 82.0
## C 64.00000 6.082763 7 1.899248 55 75 62.0 63 65.5
## D 85.57143 4.503967 7 1.899248 79 92 83.0 85 88.5
##
## Alpha: 0.05 ; DF Error: 24
##
## Critical Range
## 2 3 4
## 5.543512 5.822354 6.001384
##
## Means with the same letter are not significantly different.
##
## Tiempo groups
## A 208.28571 a
## D 85.57143 b
## B 78.85714 c
## C 64.00000 d
De la salida anterior los tratamientos que comparten al menos una letra en la columna grupos se consideran no significativamente diferentes entre sí. En este caso los tratamientos no cumplen la anterior condición lo que indica que los grupos son estadísticamente diferentes entre sí en términos de tiempo de cocción.
Lo cual se observa en la siguiente gráfica:
#out <- duncan.test(model, "virus", main = "yield of sweetpotato, Dealt with different virus")
plot(metodos.Duncan, variation="IQR" )
La validez de los resultados obtenidos en cualquier análisis de varianza queda condicionado a que los supuestos del modelo se cumplan. Estos supuestos son: normalidad, varianza constante (igual varianza de los tratamientos) e independencia.
library(car)
residuos<-residuals(modelo_anova) #Creando un objeto llamado residuos que contiene los residuos el modelo
par(mfrow=c(1,3)) #Para dividir el área del gráfico en dos partes (una fila y dos columnas)
dplot<-density(residuos) #Creando un objeto llamado dplot que recibe un Density_Plot de los residuos
plot(dplot, #Graficando el objeto dplot
main="Curva de densidad observada", #Título principal de la gráfica
xlab = "Residuos", #Etiqueta del eje x
ylab = "Densidad") #Etiqueta del eje y
polygon(dplot, #Añadiendo el poligono
col = "yellow", #Definiendo el color del poligono
border = "black") #Color del borde del poligono
qqPlot(residuos, #Un gráfico Cuantil-Cuantil de los residuos
pch =20, #Forma de los puntos
main="QQ-Plot de los residuos", #Título principal
xlab = "Cuantiles teoricos", #Etiqueta eje x
ylab="Cuantiles observados de los residuos") #Etiqueta eje y
## [1] 15 17
mtext(side=3, at=par("usr")[1], adj=0.7, cex=0.6, col="gray40", line=-21, #Posición del texto
text=paste("Derly Marin --", #Texto
format(Sys.time(),
"%d/%m/%Y %H:%M:%S --"), #Fecha y Hora
R.version.string)) #Versión de R
boxplot(residuos, col = c("purple"), ylab = "Residuos", main="Box-plot de los residuos")
Sin embargo, para confirmar de manera más sólida que los residuos siguen una distribución normal, se realiza la prueba de Shapiro-Wilk. Para la prueba Shapiro-Wilk para ratificar el cumplimiento del supuesto de normalidad de los residuos, evaluando las hipótesis:
\(H_0\): Los residuos de la variable tiempo de cocción se distribuyen normalmente con media cero y varianza constante \(e i\) \(N(0,1)\).
\(H_a\): Los residuos de la variable tiempo de cocción no siguen la distribución normal.
shapiro.test(residuals(modelo_anova))
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_anova)
## W = 0.98274, p-value = 0.9101
Por lo tanto no hay evidencia estadística suficiente para rechazar \(H_0\) ,es decir, se acepta la hipótesis nula, debido a que el valor de p (p-value = 0.9101) es mayor al valor del nivel de significancia (α =0.05), por lo que se concluye que los residuos de la variable porcentaje de parasitismo están normalmente distribuidos con media cero y varianza constante.
Este es el supuesto de mayor relevancia que los residuos deben cumplir para que el modelo empleado sea válido.
boxplot(residuos ~ ejercicio4$Tratamiento,
main = "Boxplot de Residuos por tratamiento",
xlab = "Tratamiento",
col="violet",
ylab = "Residuos")
En la siguiente gráfica, se representan los valores predichos por el
modelo para la variable tiempo de cocción en función de la raíz cuadrada
de los residuos estandarizados. En esta gráfica, no se observa ninguna
tendencia aparente en la distribución de los valores, lo que sugiere que
no hay evidencia de incumplimiento del supuesto de homogeneidad de
varianzas.
Ahora graficamos los residuos
library(car)
color_1 <-colorRampPalette(c("red", "purple", "red"))
plot(residuos, main = "Prueba de independencia", pch=20,cex = 2, col=color_1(120), ylab = "Residuos", xlab = " ")
En la grafica anterior se observan dispersos los puntos sin seguir un
patron, esto es un indicio de homogeneidad de varianzas (entre más
dispersos menos correlacionados).
Sin embargo, para validar de manera más sólida la homogeneidad de varianzas, se llevó a cabo la prueba de bartlett y la prueba de Levene.
Donde las hipótesis correspondientes son:
\(H_0\): Los residuos de la variable tiempo son iguales para los ditintos tratamientos.
\(H_a\): Existen por lo menos dos varianzas distintas para los ditintos tratamientos. Es decir,
\(H_0: σ^2_A=σ^2_B=σ^2_C=σ^2_D=σ2\)
\(H_a: σ^2_i=σ^2_j\) \(para i,j∈[A,B,C,D]\) \(e\) \(i≠j\)
library(stats)
bartlett.test(residuos ~ ejercicio4$Tratamiento)
##
## Bartlett test of homogeneity of variances
##
## data: residuos by ejercicio4$Tratamiento
## Bartlett's K-squared = 0.93367, df = 3, p-value = 0.8173
De acuerdo al valor arrojado por la prueba de bartlett, valor de p (0.8173) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, llegando a la misma conclusión anterior (grafica anterior).
ejercicio4$Tratamiento <- as.factor(ejercicio4$Tratamiento)
library(stats)
leveneTest(residuos ~ ejercicio4$Tratamiento)
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 3 0.0606 0.98
## 24
De acuerdo al valor arrojado por la prueba de levene, valor de p (0.98) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, confirmando las conclusiones obtenidas en los items anteriores.
Independencia de los residuos
Prueba de Durbin Watson
Donde las hipótesis son:
\(H_0\):Los residuos entre los tratamientos son independientes.
\(H_a\):Los residuos entre los tratamientos no son independientes
durbinWatsonTest(modelo_anova)
## lag Autocorrelation D-W Statistic p-value
## 1 -0.1042635 2.139489 0.846
## Alternative hypothesis: rho != 0
Los resultados indican que, en este análisis, los residuos están cerca de ser independientes. El valor del estadístico Durbin-Watson (DW) cercano a 2 sugiere que los residuos no muestran una autocorrelación significativa. Además, el p-value (0.822) es mayor que un nivel de significancia típico como α = 0.05, lo que sugiere que no hay suficiente evidencia para afirmar que existe una autocorrelación significativa en los residuos.
Una compañía farmacéutica desea evaluar el efecto que tiene la cantidad de almidón en la dureza de las tabletas. Se decidió producir lotes con una cantidad determinada de almidón, y que las cantidades de almidón a aprobar fueran 2%, 5% y 10%. La variable de respuesta sería el promedio de la dureza de 20 tabletas de cada lote. Se hicieron 4 réplicas por tratamiento y se obtuvieron los siguientes resultados:
# install.packages("readxl")#solo una vez
library(readxl)
ejercicio5 <- read_excel("C:/Users/LENOVO/Desktop/Ejercicio5.xlsx")
ejercicio5 # para visualizar los datos
## # A tibble: 12 × 2
## Dureza Almidon
## <dbl> <chr>
## 1 4.3 A
## 2 5.2 A
## 3 4.8 A
## 4 4.5 A
## 5 6.5 B
## 6 7.3 B
## 7 6.9 B
## 8 6.1 B
## 9 9 C
## 10 7.8 C
## 11 8.5 C
## 12 8.1 C
conteo_valoresalmidon <- table(ejercicio5$Almidon)
conteo_valoresalmidon
##
## A B C
## 4 4 4
Dado que el número de observaciones por tratamiento es el mismo, se puede concluir que es un diseño balanceado.
#install.packages("summarytools")# solo una vez
library(summarytools)
summarytools::descr(ejercicio5[,1])# todas las filas primera columna datos[,1]
## Descriptive Statistics
## ejercicio5$Dureza
## N: 12
##
## Dureza
## ----------------- --------
## Mean 6.58
## Std.Dev 1.62
## Min 4.30
## Q1 5.00
## Median 6.70
## Q3 7.95
## Max 9.00
## MAD 2.15
## IQR 2.77
## CV 0.25
## Skewness -0.05
## SE.Skewness 0.64
## Kurtosis -1.60
## N.Valid 12.00
## Pct.Valid 100.00
A partir de los resultados que arroja el programa podemos concluir: el promedio de dureza es de 6.58, con desviación estándar de 1.62, el valor mínimo de dureza es de 4.3, el valor máximo de desgaste es de 9, el 50% de las observaciones(6) presentaron una dureza entre 4.3 y 6.7, mientras que el restante 50% (6) presento un desgaste entre 6.7 y 9, el coeficiente de asimetría fue de -0.05 presentando una asimetría negativa leve, además el coeficiente de curtosis fue de ,-1.6 indicando que la distribución es levemente platicúrtica.
# Calcular estadísticas descriptivas por categoría
resultados_descriptivos <- aggregate(Dureza ~ Almidon, data = ejercicio5, summary)
# Imprimir los resultados descriptivos
print(resultados_descriptivos)
## Almidon Dureza.Min. Dureza.1st Qu. Dureza.Median Dureza.Mean Dureza.3rd Qu.
## 1 A 4.300 4.450 4.650 4.700 4.900
## 2 B 6.100 6.400 6.700 6.700 7.000
## 3 C 7.800 8.025 8.300 8.350 8.625
## Dureza.Max.
## 1 5.200
## 2 7.300
## 3 9.000
LLevando a cabo el análisis descriptivo por tratamiento se obtuvieron los siguientes resultados:
Para la cantidad de almidón A: el promedio de dureza es de 4.7, el valor mínimo de dureza es 4.3, el valor máximo de dureza es 5.2, el 50% de las observaciones(2) presentaron una dureza entre 4.3 y 4.65, mientras que el restante 50% (2) presento un desgaste entre 4.65 y 5.2.
Para la cantidad de almidón B: el promedio de dureza es de 6.7, el valor mínimo de dureza es 6.1, el valor máximo de dureza es 7.3, el 50% de las observaciones(2) presentaron una dureza entre 6.1 y 6.7, mientras que el restante 50% (2) presento un desgaste entre 6.7 y 7.3.
Para la cantidad de almidón C: el promedio de dureza es de 8.35, el valor mínimo de dureza es 7.8, el valor máximo de dureza es 9, el 50% de las observaciones(2) presentaron una dureza entre 7.8 y 8.3, mientras que el restante 50% (2) presento un desgaste entre 8.3 y 9.
A continuación se llevará a cabo el ANOVA de un factor o modelo factorial de un solo factor el cual nos permite estudiar si existen diferencias significativas entre las medias de dureza de las tabletas de las tres cantidades de almidón.
Las hipótesis contrastadas en el ANOVA son: No hay diferencias entre las medias de las diferentes cantidades de almidón y existen por lo menos dos cantidades de almidón con diferencias significativas en el promedio de dureza de las tabletas.
Que se plantearían de la siguiente forma:
\(H_0:μ_A=μ_B=μ_C=μ\)
\(H_a:μ_i=μ_j\) para \(i≠j\) donde \(i,j=A,B,C\) respectivamente.
# Realizar el ANOVA
modelo_anova <- aov(Dureza ~ Almidon, data = ejercicio5)
# Resumen del ANOVA
resumen_anova <- summary(modelo_anova)
# Imprimir el resumen del ANOVA
print(resumen_anova)
## Df Sum Sq Mean Sq F value Pr(>F)
## Almidon 2 26.73 13.36 58.1 7.16e-06 ***
## Residuals 9 2.07 0.23
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
¿Hay evidencia suficiente de que el almidón influye en la dureza en las tabletas?
Los resultados del análisis de varianza (ANOVA) arrojan un p-valor 7.16e-06, es menor al nivel de significancia (α=0.05). En consecuencia, no se dispone de evidencia suficiente para considerar que los promedios son iguales, es decir, se rechaza \(H_0\) por lo que se supone que existen por lo menos dos % de almidón con diferencias significativas en el promedio de la dureza de las tabletas. Para estimar cuales son los tratamientos cuyos medias de tiempo de cocción son diferentes, hacemos uso de los diagramas de cajas y bigotes por cada método o tratamiento.
boxplot(ejercicio5$Dureza ~ ejercicio5$Almidon, data = ejercicio5, col = c("pink", "lightblue", "blue"), ylab = "Dureza", xlab = "Porcentaje de almidón")
De la grafica anterior se puede evidenciar que existen diferencias en los promedios de dureza entre todos los porcentajes de almidón (no se traslapan las gráficas).De forma que la dureza de las tabletas varía entre los diferentes porcentajes de almidón, con mayor grado de dureza al incrementar el % de almidón.
Como se ha rechazado la hipótesis de igualdad de medias con el test ANOVA, el interés está en averiguar cuál o cuáles pares de medias son diferentes entre sí.
# Realizar la prueba de LSD
lsd <- lsmeans(modelo_anova, pairwise ~ Almidon)
# Mostrar los resultados de la prueba de LSD
print(lsd)
## $lsmeans
## Almidon lsmean SE df lower.CL upper.CL
## A 4.70 0.24 9 4.16 5.24
## B 6.70 0.24 9 6.16 7.24
## C 8.35 0.24 9 7.81 8.89
##
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## A - B -2.00 0.339 9 -5.898 0.0006
## A - C -3.65 0.339 9 -10.763 <.0001
## B - C -1.65 0.339 9 -4.866 0.0023
##
## P value adjustment: tukey method for comparing a family of 3 estimates
Los resultados del análisis de comparaciones de contrastes (LSD) destacan diferencias estadísticamente significativas entre todos los pares de niveles evaluados (A - B, A - C y B - C), con valores de p menores al preestablecido (0.05).
TukeyHSD(modelo_anova)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Dureza ~ Almidon, data = ejercicio5)
##
## $Almidon
## diff lwr upr p adj
## B-A 2.00 1.0531848 2.946815 0.0006016
## C-A 3.65 2.7031848 4.596815 0.0000052
## C-B 1.65 0.7031848 2.596815 0.0022940
Luego de realizar la prueba de comparaciones múltiples de Tukey, los resultados en cada comparación indican que hay diferencias estadísticamente significativas en la dureza entre todos los tratamientos (A, B y C). Los valores p ajustados cercanos a cero (menores a 0.05) en todos los casos respaldan la significancia estadística de estas diferencias.
plot(TukeyHSD(modelo_anova))
Lo mencionado anteriormente se confirma con la gráfica de la diferencia de medias(intervalos de confianza) en los distintos niveles de almidón.
library(agricolae)
metodos.Duncan <-duncan.test(modelo_anova, trt = "Almidon", group = T, console = T)
##
## Study: modelo_anova ~ "Almidon"
##
## Duncan's new multiple range test
## for Dureza
##
## Mean Square Error: 0.23
##
## Almidon, means
##
## Dureza std r se Min Max Q25 Q50 Q75
## A 4.70 0.3915780 4 0.2397916 4.3 5.2 4.450 4.65 4.900
## B 6.70 0.5163978 4 0.2397916 6.1 7.3 6.400 6.70 7.000
## C 8.35 0.5196152 4 0.2397916 7.8 9.0 8.025 8.30 8.625
##
## Alpha: 0.05 ; DF Error: 9
##
## Critical Range
## 2 3
## 0.7671348 0.8006971
##
## Means with the same letter are not significantly different.
##
## Dureza groups
## C 8.35 a
## B 6.70 b
## A 4.70 c
De la salida anterior los tratamientos que comparten al menos una letra en la columna grupos se consideran no significativamente diferentes entre sí. En este caso los tratamientos (% de almidón) no cumplen la anterior condición lo que indica que los grupos son estadísticamente diferentes entre sí en términos de dureza.
Lo cual se observa en la siguiente gráfica:
#out <- duncan.test(model, "virus", main = "yield of sweetpotato, Dealt with different virus")
plot(metodos.Duncan, variation="IQR" )
La validez de los resultados obtenidos en cualquier análisis de varianza queda condicionado a que los supuestos del modelo se cumplan. Estos supuestos son: normalidad, varianza constante (igual varianza de los tratamientos) e independencia.
library(car)
residuos<-residuals(modelo_anova) #Creando un objeto llamado residuos que contiene los residuos el modelo
par(mfrow=c(1,3)) #Para dividir el área del gráfico en dos partes (una fila y dos columnas)
dplot<-density(residuos) #Creando un objeto llamado dplot que recibe un Density_Plot de los residuos
plot(dplot, #Graficando el objeto dplot
main="Curva de densidad observada", #Título principal de la gráfica
xlab = "Residuos", #Etiqueta del eje x
ylab = "Densidad") #Etiqueta del eje y
polygon(dplot, #Añadiendo el poligono
col = "orange", #Definiendo el color del poligono
border = "black") #Color del borde del poligono
qqPlot(residuos, #Un gráfico Cuantil-Cuantil de los residuos
pch =20, #Forma de los puntos
main="QQ-Plot de los residuos", #Título principal
xlab = "Cuantiles teoricos", #Etiqueta eje x
ylab="Cuantiles observados de los residuos") #Etiqueta eje y
## [1] 9 8
mtext(side=3, at=par("usr")[1], adj=0.7, cex=0.6, col="gray40", line=-21, #Posición del texto
text=paste("Derly Marin --", #Texto
format(Sys.time(),
"%d/%m/%Y %H:%M:%S --"), #Fecha y Hora
R.version.string)) #Versión de R
boxplot(residuos, col = c("violet"), ylab = "Residuos", main="Box-plot de los residuos")
Sin embargo, para confirmar de manera más sólida que los residuos siguen una distribución normal, se realiza la prueba de Shapiro-Wilk. Para la prueba Shapiro-Wilk para ratificar el cumplimiento del supuesto de normalidad de los residuos, evaluando las hipótesis:
\(H_0\):Los residuos de la variable dureza se distribuyen normalmente con media cero y varianza constante \(ei N(0,1)\).
\(H_a\): Los residuos de la variable dureza no siguen la distribución normal.
shapiro.test(residuals(modelo_anova)) #Prueba Shapiro-Wilk para los residuos de la variable
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_anova)
## W = 0.93444, p-value = 0.4295
Por lo tanto no hay evidencia estadística suficiente para rechazar \(H_0\),es decir se acepta la hipótesis nula, debido a que el valor de p (p-value = 0.4295) es mayor al valor del nivel de significancia (alfa=0.05), por lo que se concluye que los residuos de la variable dureza están normalmente distribuidos con media cero y varianza constante.
Este es el supuesto de mayor relevancia que los residuos deben cumplir para que el modelo empleado sea válido.
boxplot(residuos ~ ejercicio5$Almidon,
main = "Boxplot de Residuos por metodo de ensamble",
xlab = "Porcentaje de almidón",
col="blue",
ylab = "Residuos")
En la siguiente gráfica, se representan los valores predichos por el modelo para la variable dureza en función de la raíz cuadrada de los residuos estandarizados. En esta gráfica, no se observa ninguna tendencia aparente en la distribución de los valores, lo que sugiere que no hay evidencia de incumplimiento del supuesto de homogeneidad de varianzas.
Ahora graficamos los residuos
library(car)
color_1 <-colorRampPalette(c("yellow ", "black", "purple"))
plot(residuos, main = "Prueba de independencia", pch=20,cex = 2, col=color_1(120), ylab = "Residuos", xlab = " ")
En la grafica anterior se observan dispersos los puntos sin seguir un patrón, esto es un indicio de homogeneidad de varianzas (entre más dispersos menos correlacionados)
Sin embargo, para validar de manera más sólida la homogeneidad de varianzas, se llevó a cabo la prueba de bartlett y la prueba de Levene.
Donde las hipótesis correspondientes son:
\(H_0\): Los residuos de la variable dureza son iguales para las distintas cantidades de almidón.
\(H_a\): Existen por lo menos dos varianzas distintas para las distintas cantidades de almidón.
Es decir,
\(H_0: σ^2_A=σ^2_B=σ^2_C=σ2\)
\(H_a: σ^2_i=σ^2_j\) \(para i,j∈[A,B,C]\) \(e\) \(i≠j\)
library(stats)
bartlett.test(residuos ~ ejercicio5$Almidon)
##
## Bartlett test of homogeneity of variances
##
## data: residuos by ejercicio5$Almidon
## Bartlett's K-squared = 0.25398, df = 2, p-value = 0.8807
De acuerdo al valor arrojado por la prueba de bartlett, valor de p (0.8807) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, llegando a la misma conclusión anterior (gráfica anterior).
library(stats)
leveneTest(residuos ~ ejercicio5$Almidon)
## Warning in leveneTest.default(y = y, group = group, ...): group coerced to
## factor.
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.2667 0.7718
## 9
De acuerdo al valor arrojado por la prueba de levene, valor de p (0.7718) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, confirmando las conclusiones obtenidas en los items anteriores.
Independencia de los residuos
Prueba de Durbin Watson
Donde las hipótesis son:
\(H_0\): Los residuos entre los tratamientos son independientes.
\(H_a\):Los residuos entre los tratamientos no son independientes.
durbinWatsonTest(modelo_anova)
## lag Autocorrelation D-W Statistic p-value
## 1 -0.5398551 2.972222 0.168
## Alternative hypothesis: rho != 0
Puesto que el valor de DW es aproximadamente igual a 2 (1-Rho) donde Rho es la autocorrelación de la muestra de los residuos.
Al realizar la prueba de independencia de residuos para la variable desgaste se determinó que los residuos no están correlacionados, debido a que el DW está próximo a 2 y el valor de p (p-value = 0.468) es superior al nivel de significancia de 5% (α=0.05) por lo que se concluye que existe independencia de los residuos.
Los datos que se presentan enseguida son rendimientos en toneladas por hectárea de un pasto con tres niveles de fertilización nitrogenada. El diseño fue completamente aleatorizado, con cinco repeticiones por tratamiento.
# install.packages("readxl")#solo una vez
library(readxl)
ejercicio6 <- read_excel("C:/Users/LENOVO/Desktop/Ejercicio6.xlsx")
ejercicio6 # para visualizar los datos
## # A tibble: 15 × 2
## Rendimiento Nitrogeno
## <dbl> <chr>
## 1 14.8 A
## 2 14.7 A
## 3 14.7 A
## 4 14.5 A
## 5 15.1 A
## 6 25.2 B
## 7 25.4 B
## 8 25.1 B
## 9 25.0 B
## 10 25.3 B
## 11 32.6 C
## 12 32.5 C
## 13 32.3 C
## 14 32.7 C
## 15 32.1 C
conteo_valoresnitrogeno <- table(ejercicio6$Nitrogeno)
conteo_valoresnitrogeno
##
## A B C
## 5 5 5
Dado que el número de observaciones por tratamiento es el mismo, se puede concluir que es un diseño balanceado.
#install.packages("summarytools")# solo una vez
library(summarytools)
summarytools::descr(ejercicio6[,1])# todas las filas primera columna datos[,1]
## Descriptive Statistics
## ejercicio6$Rendimiento
## N: 15
##
## Rendimiento
## ----------------- -------------
## Mean 24.13
## Std.Dev 7.51
## Min 14.51
## Q1 14.82
## Median 25.15
## Q3 32.26
## Max 32.67
## MAD 10.84
## IQR 17.24
## CV 0.31
## Skewness -0.20
## SE.Skewness 0.58
## Kurtosis -1.69
## N.Valid 15.00
## Pct.Valid 100.00
A partir de los resultados que arroja el programa podemos concluir: el promedio de rendimiento es de 24.13, con desviación estándar de 7.51, el valor mínimo de rendimiento es de 14.51, el valor máximo de desgaste es de 32.67, el 50% de las observaciones presentaron un rendimiento entre 14.51 y 25.15, mientras que el restante 50% presentó un rendimiento entre 25.15 y 32.67, el coeficiente de asimetría fue de -0.20 presentando una asimetría negativa, además el coeficiente de curtosis fue de -1.69, indicando que la distribución es platicúrtica.
# Calcular estadísticas descriptivas por categoría
resultados_descriptivos <- aggregate(Rendimiento ~ Nitrogeno, data = ejercicio6, summary)
# Imprimir los resultados descriptivos
print(resultados_descriptivos)
## Nitrogeno Rendimiento.Min. Rendimiento.1st Qu. Rendimiento.Median
## 1 A 14.51410 14.67600 14.72000
## 2 B 25.03100 25.13100 25.15100
## 3 C 32.11100 32.25600 32.46000
## Rendimiento.Mean Rendimiento.3rd Qu. Rendimiento.Max.
## 1 14.75962 14.82300 15.06500
## 2 25.19620 25.26700 25.40100
## 3 32.42020 32.60500 32.66900
LLevando a cabo el análisis descriptivo por tratamiento se obtuvieron os siguientes resultados:
Para nivel de nitrógeno A: el promedio de rendimiento es de 14.8, el valor mínimo de rendimiento es 14.5141, el valor máximo de rendimiento es 15.07, el 50% de las observaciones presentaron un rendimiento entre 14.5141 y 14.72, mientras que el restante 50% presentó un rendimiento entre 14.72 y 15.07.
Para nivel de nitrógeno B: el promedio de rendimiento es de 25.2, el valor mínimo de rendimiento es 25.031, el valor máximo de rendimiento es 25.401, el 50% de las observaciones presentaron un rendimiento entre 25.031 y 25.151, mientras que el restante 50% presentó un rendimiento entre 25.151 y 25.401.
Para nivel de nitrógeno C: el promedio de rendimiento es de 32.4202, el valor mínimo de rendimiento es 32.111, el valor máximo de rendimiento es 32.669, el 50% de las observaciones presentaron un rendimiento entre 32.111 y 32.46, mientras que el restante 50% presentó un rendimiento entre 32.46 y 32.669.
A continuación se llevará a cabo el ANOVA de un factor o modelo factorial de un solo factor el cual nos permite estudiar si existen diferencias significativas entre las medias rendimiento de los tres niveles de nitrógeno.
Las hipótesis contrastadas en el ANOVA son: No hay diferencias entre las medias de los diferentes niveles de nitrógeno y existen por lo menos dos niveles de nitrógeno con diferencias significativas en el promedio de erndimiento. Que se plantearían de la siguiente forma:
\(H_0:μ_A=μ_B=μ_C=μ\)
\(H_a:μ_i=μ_j\) para \(i≠j\) donde \(i,j=A,B,C\) respectivamente
# Realizar el ANOVA
modelo_anova <- aov(Rendimiento ~ Nitrogeno, data = ejercicio6)
# Resumen del ANOVA
resumen_anova <- summary(modelo_anova)
# Imprimir el resumen del ANOVA
print(resumen_anova)
## Df Sum Sq Mean Sq F value Pr(>F)
## Nitrogeno 2 788.3 394.2 10132 <2e-16 ***
## Residuals 12 0.5 0.0
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
En este caso el ANOVA ha resultado significativo valor F = 10132 y p-valor 2e-16 ,es menor que alpha (0.05), no se dispone de evidencia suficiente para considerar que los promedios son iguales,es decir, se rechaza \(H_0\) por lo que se supone que existen por lo menos dos niveles de nitrógeno con diferencias significativas en el promedio de rendimiento.
Para estimar cuales son los niveles de nitrógeno cuyos medias de rendimiento son diferentes, hacemos uso de los diagramas de cajas y bigotes por cada método o tratamiento.
Las diferencias muestrales son estadísticamente significativas, lo que sugiere que hay variabilidad en los datos, pero se requiere realizar pruebas adicionales, como pruebas post hoc o análisis de comparaciones múltiples, para identificar cuáles niveles de nitrógeno específicos difieren entre sí y cuantificar la magnitud de esas diferencias. Estas pruebas adicionales pueden ayudarnos a determinar si las diferencias observadas en las muestras son representativas de diferencias poblacionales significativas.
library(gplots)
# Crear un gráfico de medias con intervalos de confianza (95%)
plotmeans(Rendimiento ~ Nitrogeno, data = ejercicio6,
xlab = "Nitrógeno", ylab = "Rendimiento",
main = "Grafico de Medias con 95% IC", cex.axis = 0.6)
De la gráfica anterior, se puede apreciar que la media del nivel C es la
más alta. En contraste, el nivel A que presenta la media más baja y en
con una media intermedia esta en nivel B. Esas diferencias en los
niveles de fertilización nitrogenada parecen tener un impacto notable en
el rendimiento del pasto.
# Crear el diagrama de cajas por categorías
boxplot(ejercicio6$Rendimiento ~ ejercicio6$Nitrogeno, data = ejercicio6, col = c("pink", "lightblue", "orange"), ylab = "Rendimiento", xlab = "Niveles de Nitrógeno")
Con base en la información proporcionada por la gráfica, se puede
observar claramente que se registran diferencias significativas en los
promedios de rendimiento en todos los niveles de nitrógeno (no se
traslapan las gráficas).
Como se ha rechazado la hipótesis de igualdad de medias con el test ANOVA, el interés está en averiguar cuál o cuáles pares de medias son diferentes entre sí.
# Realizar la prueba de LSD
lsd <- lsmeans(modelo_anova, pairwise ~ Nitrogeno)
# Mostrar los resultados de la prueba de LSD
print(lsd)
## $lsmeans
## Nitrogeno lsmean SE df lower.CL upper.CL
## A 14.8 0.0882 12 14.6 15.0
## B 25.2 0.0882 12 25.0 25.4
## C 32.4 0.0882 12 32.2 32.6
##
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## A - B -10.44 0.125 12 -83.662 <.0001
## A - C -17.66 0.125 12 -141.570 <.0001
## B - C -7.22 0.125 12 -57.909 <.0001
##
## P value adjustment: tukey method for comparing a family of 3 estimates
Los resultados del análisis de comparaciones de contrastes (LSD) indican diferencias estadísticamente significativas en todos los pares de niveles evaluados (A - B, A - C y B - C). Estas diferencias son consistentes con un valor de p menor a 0.05.
TukeyHSD(modelo_anova)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Rendimiento ~ Nitrogeno, data = ejercicio6)
##
## $Nitrogeno
## diff lwr upr p adj
## B-A 10.43658 10.10377 10.76939 0
## C-A 17.66058 17.32777 17.99339 0
## C-B 7.22400 6.89119 7.55681 0
Luego de realizar la prueba de comparaciones múltiples de Tukey, los resultados en cada comparación indican que hay diferencias estadísticamente significativas en el rendimiento entre todos niveles (A, B y C). Los valores p ajustados dan cero (menores a 0.05), es decir, no se dispone de evidencia suficiente para considerar que los promedios de rendimiento en los tratamientos son iguales.
plot(resultado_tukey)
Lo mencionado anteriormente se confirma con la grafica de la diferencia
de medias(intervalos de confianza) en los distintos niveles de
nitrógeno.
library(agricolae)
metodos.Duncan <-duncan.test(modelo_anova, trt = "Nitrogeno", group = T, console = T)
##
## Study: modelo_anova ~ "Nitrogeno"
##
## Duncan's new multiple range test
## for Rendimiento
##
## Mean Square Error: 0.03890497
##
## Nitrogeno, means
##
## Rendimiento std r se Min Max Q25 Q50 Q75
## A 14.75962 0.2037867 5 0.08820995 14.5141 15.065 14.676 14.720 14.823
## B 25.19620 0.1418986 5 0.08820995 25.0310 25.401 25.131 25.151 25.267
## C 32.42020 0.2346289 5 0.08820995 32.1110 32.669 32.256 32.460 32.605
##
## Alpha: 0.05 ; DF Error: 12
##
## Critical Range
## 2 3
## 0.2718019 0.2844986
##
## Means with the same letter are not significantly different.
##
## Rendimiento groups
## C 32.42020 a
## B 25.19620 b
## A 14.75962 c
De la salida anterior los tratamientos que comparten al menos una letra en la columna grupos se consideran no significativamente diferentes entre sí. En este caso los tratamientos(nivel de nittrógeno) no cumplen la anterior condición lo que indica que los grupos son estadísticamente diferentes entre sí en términos de rendimiento.
Lo cual se observa en la siguiente gráfica:
plot(metodos.Duncan, variation="IQR" )
La validez de los resultados obtenidos en cualquier análisis de varianza queda condicionado a que los supuestos del modelo se cumplan. Estos supuestos son: normalidad, varianza constante (igual varianza de los tratamientos) e independencia.
library(car)
residuos<-residuals(modelo_anova) #Creando un objeto llamado residuos que contiene los residuos el modelo
par(mfrow=c(1,3)) #Para dividir el área del gráfico en dos partes (una fila y dos columnas)
dplot<-density(residuos) #Creando un objeto llamado dplot que recibe un Density_Plot de los residuos
plot(dplot, #Graficando el objeto dplot
main="Curva de densidad observada", #Título principal de la gráfica
xlab = "Residuos", #Etiqueta del eje x
ylab = "Densidad") #Etiqueta del eje y
polygon(dplot, #Añadiendo el poligono
col = "pink", #Definiendo el color del poligono
border = "black") #Color del borde del poligono
qqPlot(residuos, #Un grafico Cuantil-Cuantil de los residuos
pch =20, #Forma de los puntos
main="QQ-Plot de los residuos", #Titulo principal
xlab = "Cuantiles teóricos", #Etiqueta eje x
ylab="Cuantiles observados de los residuos") #Etiqueta eje y
## [1] 15 5
mtext(side=3, at=par("usr")[1], adj=0.7, cex=0.6, col="gray40", line=-21, #Posición del texto
text=paste("Derly Marin --", #Texto
format(Sys.time(),
"%d/%m/%Y %H:%M:%S --"), #Fecha y Hora
R.version.string)) #Versión de R
boxplot(residuos, col = c("orange"), ylab = "Residuos", main="Box-plot de los residuos")
Sin embargo, para confirmar de manera más sólida que los residuos siguen una distribución normal, se realiza la prueba de Shapiro-Wilk. Para la prueba Shapiro-Wilk para ratificar el cumplimiento del supuesto de normalidad de los residuos, evaluando las hipótesis:
\(H_0\): Los residuos de la variable rendimiento se distribuyen normalmente con media cero y varianza constante \(e i\) \(N(0,1)\).
\(H_a\): Los residuos de la variable rendimiento no siguen la distribución normal.
shapiro.test(residuals(modelo_anova)) #Prueba Shapiro-Wilk para los residuos de la variable
##
## Shapiro-Wilk normality test
##
## data: residuals(modelo_anova)
## W = 0.97219, p-value = 0.8891
Por lo tanto no hay evidencia estadística suficiente para rechazar \(H_0\), es decir, se acepta la hipótesis nula, debido a que el valor de p (p-value = 0.8891) es mayor al valor del nivel de significancia (alfa=0.05), por lo que se concluye que los residuos de la variable rendimiento están normalmente distribuidos con media cero y varianza constante.
Este es el supuesto de mayor relevancia que los residuos deben cumplir para que el modelo empleado sea válido.
boxplot(residuos ~ ejercicio6$Nitrogeno,
main = "Boxplot de Residuos por niveles de nitrógeno",
xlab = "Nivel de nitrógeno",
col="grey",
ylab = "Residuos")
En la siguiente gráfica, se representan los valores predichos por el modelo para la variable rendimiento en función de la raíz cuadrada de los residuos estandarizados. En esta gráfica, no se observa ninguna tendencia aparente en la distribución de los valores, lo que sugiere que no hay evidencia de incumplimiento del supuesto de homogeneidad de varianzas.
Ahora graficamos los residuos
library(car)
color_1 <-colorRampPalette(c("blue", "green", "blue"))
plot(residuos, main = "Prueba de independencia", pch=20,cex = 2, col=color_1(120), ylab = "Residuos", xlab = " ")
En la grafica anterior se observan dispersos los puntos sin seguir un patrón, esto es un indicio de homogeneidad de varianzas (entre más dispersos menos correlacionados).
Sin embargo, para validar de manera más sólida la homogeneidad de varianzas, se llevó a cabo la prueba de bartlett y la prueba de Levene.
Donde las hipótesis correspondientes son:
\(H_0\): Los residuos de la variable rendimiento son iguales para los ditintos niveles de nitrógeno.
\(H_a\): Existen por lo menos dos varianzas distintas para los ditintos niveles de nitrógeno. Es decir,
\(H_0: σ^2_A=σ^2_B=σ^2_C=σ2\)
\(H_a: σ^2_i=σ^2_j\) \(para i,j∈[A,B,C]\) \(e\) \(i≠j\)
library(stats)
bartlett.test(residuos ~ ejercicio6$Nitrogeno)
##
## Bartlett test of homogeneity of variances
##
## data: residuos by ejercicio6$Nitrogeno
## Bartlett's K-squared = 0.8865, df = 2, p-value = 0.6419
De acuerdo al valor arrojado por la prueba de bartlett, valor de p (0.6419) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, llegando a la misma conclusión anterior (gráfica anterior).
library(stats)
leveneTest(residuos ~ ejercicio6$Nitrogeno)
## Warning in leveneTest.default(y = y, group = group, ...): group coerced to
## factor.
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.5372 0.5978
## 12
De acuerdo al valor arrojado por la prueba de levene, valor de p (0.5978) mayor a 0.05 se acepta la hipótesis nula, por lo que existe homogeneidad de varianzas, confirmando las conclusiones obtenidas en los items anteriores
Independencia de los residuos
Prueba de Durbin Watson
Donde las hipótesis son:
\(H_0\): Los residuos entre los tratamientos son independientes.
\(H_a\): Los residuos entre los tratamientos no son independientes.
durbinWatsonTest(modelo_anova)
## lag Autocorrelation D-W Statistic p-value
## 1 -0.4464995 2.679612 0.36
## Alternative hypothesis: rho != 0
Los resultados indican que, en este análisis, los residuos están cerca de ser independientes. El valor del estadístico Durbin-Watson (DW) cercano a 2 sugiere que los residuos no muestran una autocorrelación significativa. Además, el p-value (0.394) es mayor que un nivel de significancia típico como α = 0.05, lo que sugiere que no hay suficiente evidencia para afirmar que existe una autocorrelación significativa en los residuos.