Clase 02/06
Wilder Clavijo
2023-06-02
13/16 de Junio - Evaluacion del diseño del dia de hoy. Repeated Measures ANOVA in R. (Fuente de Datos). https://www.datanovia.com/en/lessons/repeated-measures-anova-in-r/
Diseño en medidas repetidas (One Way).
- One Way / Two Way / Three Way.
Se caracterizan por involucrar el tiempo, para evaluar estados, cambios y evolucion a nivel agronomico. (Evaluacion Temporal) - Factor Intrasujetos = Tiempo * Factores Entre-sujetos = FSCA/FSBA/FCCA/FCBA.
data("selfesteem", package = "datarium")
datos = selfesteem
head(datos,3)## id t1 t2 t3
## 1 1 4.005027 5.182286 7.107831
## 2 2 2.558124 6.912915 6.308434
## 3 3 3.244241 4.443434 9.778410One Way: Solo posee un factor (tiempo). Hipotesis Nula: \[H_0= \mu_0 = \mu_1 = \mu_3\] EJ= Aceite en limonaria en 3 tiempos equidistates. (30,60,90). - Se recomienda en tiempos equidistantes en todo momento.
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(tidyr)
# install.packages("datarium")
data("selfesteem", package = "datarium")
datos = selfesteem
datos = datos %>%
gather(key = "tiempo",
value = "rto",
t1, t2, t3) %>%
mutate_at(vars(id, tiempo), as.factor)
View(datos)Se convierte los datos en un sentido largo, empleando unicamente 2 columnas para agruparlos con la funcion “Gather”
- Se deben indicar los factores (Id de cada parcela y el tiempo) Con la funcion “mutate_at”
Resumen Estadistico:
datos %>%
group_by(tiempo) %>%
summarise(media = mean(rto),
desv = sd(rto),
n = n(),
cv = 100*desv/media)## # A tibble: 3 × 5
## tiempo media desv n cv
## <fct> <dbl> <dbl> <int> <dbl>
## 1 t1 3.14 0.552 10 17.6
## 2 t2 4.93 0.863 10 17.5
## 3 t3 7.64 1.14 10 15.0“cv” = Coeficiente de variacion - Suelen ser buenos en valores menores al 20% Indican una buena relacion Agronomica, en homogeneidad de datos
Valores superiores indican incongruencias, datos atipicos, etc…
boxplot(datos$rto ~ datos$tiempo)
Deteccion de Outliers
library(rstatix)##
## Attaching package: 'rstatix'## The following object is masked from 'package:stats':
##
## filterdatos %>%
group_by(tiempo) %>%
identify_outliers(rto)## # A tibble: 2 × 5
## tiempo id rto is.outlier is.extreme
## <fct> <fct> <dbl> <lgl> <lgl>
## 1 t1 6 2.05 TRUE FALSE
## 2 t2 2 6.91 TRUE FALSELa confirmacion la da la seccion “is.extreme” dado que a pesar de existir 2 Outlier, no son graves. - SE DEBE NOTIFICAR CUANDO SI EXISTEN OUTLIERS. Se tiene en cuenta el Coeficiente de variacoion > 20%.
#Supuesto de Normalidad: Se puede extraer de los residuales como forma alterna.
datos %>%
group_by(tiempo) %>%
shapiro_test(rto)## # A tibble: 3 × 4
## tiempo variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 t1 rto 0.967 0.859
## 2 t2 rto 0.876 0.117
## 3 t3 rto 0.923 0.380No existen problemas con el supuesto de Normalidad.
Supuesto de Esfericidad.
Prueba de Esfericidad (Mauchly’s Test)
Exclusivo de los diseños en medidas repetidas. Tiene que ver con la igualdad de varianzas
Compara la variabilidad de los datos entre tiempos sucesivos (t1 ~ t2 , t2 ~ t3, sucesivamente), y se espera que sean similares. (H_O).
res.aov <- anova_test(data = datos,
dv = rto,
wid = id,
within = tiempo)
res.aov$`Mauchly's Test for Sphericity`## Effect W p p<.05
## 1 tiempo 0.551 0.092get_anova_table(res.aov)## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 tiempo 2 18 55.469 2.01e-08 * 0.829El analisis funciona supiniendo las varianzas iguales (9.2% > 5%) es correcto. (De lo contratio se necesita otra metodologia).
“ges” = Tamaño del efecto. Para este caso el p_valor rechaza la H_0 - El Aceite producido en los 3 tiempos son diferentes.
La mayor produccion esta en el t3, sin Outliers.
Pruebas de Comparacion (Post - Analisis de Varianza)
datos %>%
pairwise_t_test(
rto ~ tiempo,
paired = TRUE,
p.adjust.method =
"bonferroni"
)## # A tibble: 3 × 10
## .y. group1 group2 n1 n2 statistic df p p.adj p.adj.signif
## * <chr> <chr> <chr> <int> <int> <dbl> <dbl> <dbl> <dbl> <chr>
## 1 rto t1 t2 10 10 -4.97 9 0.000772 2e-3 **
## 2 rto t1 t3 10 10 -13.2 9 0.000000334 1e-6 ****
## 3 rto t2 t3 10 10 -4.87 9 0.000886 3e-3 **- Se utiliza el valor del p_valor ajustado.
Se comprueba que existen las diferencias entre varianzas de los tiempos contiguas. Se comprueba que el rendimiento si aumenta en el tiempo.
- Se debe tener en cuenta los tiempos de la planta, ya que en base al ciclo de la planta pueden existir diferentes variaciones.
- El cultivo puede envejecer, o tener una germinacion del dato a medir tardia.
Medidas Repetidas en 2 Vias (Two Ways)
data("selfesteem2", package = "datarium")
datos2 = selfesteem2
datos2$treatment = gl(2,12,24, c('con fert', 'sin fert'))
View(datos2)Se añade un Factor adicional “trt / fertilizacion”
datos2 = datos2 %>%
gather(key='tiempo',
value = 'rto',
t1,t2,t3)datos2 %>%
group_by(treatment, tiempo) %>%
summarise(media = mean(rto),
desv = sd(rto),
n = n(),
cv = 100*desv/media)## `summarise()` has grouped output by 'treatment'. You can override using the
## `.groups` argument.## # A tibble: 6 × 6
## # Groups: treatment [2]
## treatment tiempo media desv n cv
## <fct> <chr> <dbl> <dbl> <int> <dbl>
## 1 con fert t1 88 8.08 12 9.18
## 2 con fert t2 83.8 10.2 12 12.2
## 3 con fert t3 78.7 10.5 12 13.4
## 4 sin fert t1 87.6 7.62 12 8.70
## 5 sin fert t2 87.8 7.42 12 8.45
## 6 sin fert t3 87.7 8.14 12 9.28Buenos cv, lo que indica una buena homogeneidad de datos.
Visualizacion:
library(ggplot2)
ggplot(datos2)+
aes(tiempo,rto, fill=treatment)+
geom_boxplot()
Se evidencia una alta variabildad al tener fertilizacion. - Pueden existir muchas razones por las cuales no hay diferencia (lixiviacion, evaporacion, existencia de arcillas 2:1, tiempos muy tempranos de medicion).
Comprobacion de Outliers.
datos2 %>%
group_by(treatment, tiempo) %>%
identify_outliers(rto)## [1] treatment tiempo id rto is.outlier is.extreme
## <0 rows> (or 0-length row.names)No existen datos atipicos.
Supuesto de Normalidad
datos2 %>%
group_by(treatment, tiempo) %>%
shapiro_test(rto)## # A tibble: 6 × 5
## treatment tiempo variable statistic p
## <fct> <chr> <chr> <dbl> <dbl>
## 1 con fert t1 rto 0.828 0.0200
## 2 con fert t2 rto 0.868 0.0618
## 3 con fert t3 rto 0.887 0.107
## 4 sin fert t1 rto 0.919 0.279
## 5 sin fert t2 rto 0.923 0.316
## 6 sin fert t3 rto 0.886 0.104Existen anormalidades en el primer tiempo con fertilizacion, sin embargo, se considera una prueba de normalidad de mala calidad. - La normalidad se aplica sobre los Residuales.
Analisis de Varianzas.
res.aov <- anova_test(
data = datos2,
dv = rto,
wid = id,
within = c(treatment,
tiempo)
)
get_anova_table(res.aov)## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 treatment 1.00 11.00 15.541 2.00e-03 * 0.059
## 2 tiempo 1.31 14.37 27.369 5.03e-05 * 0.049
## 3 treatment:tiempo 2.00 22.00 30.424 4.63e-07 * 0.050Si existe la interaccion, por ende se deben concentrar en la interaccion. - Cuando existe la interaccion, NO se deben hacer comparaciones posteriores (bontferroni).
Grafico de interaccion (Cuando existe)
interaction.plot(datos2$tiempo,
datos2$treatment,
datos2$rto)
* Forma Alternativa (reproducir)
datos2 %>%
group_by(tiempo, treatment) %>%
summarise(mean_rto = mean(rto)) %>%
ggplot()+
aes(tiempo, mean_rto,
color=treatment,
group=treatment)+
geom_point(size=5)+
geom_line(linewidth=3)## `summarise()` has grouped output by 'tiempo'. You can override using the
## `.groups` argument.
Las lineas se cruzan, por ende existe la intreraccion. - Se evidencia la
diferencia de rendimiento a traves del tiempo.
Test de Esfericidad
res.aov$`Mauchly's Test for Sphericity`## Effect W p p<.05
## 1 tiempo 0.469 0.023 *
## 2 treatment:tiempo 0.616 0.089- Las varianzas son iguales por interaccion, pero NO por tratamiento.