Clase
Stefania Hernández Azuero
2023-06-02
#Diseño en medidas repetidas (Repeated Measures ANOVA in R) https://www.datanovia.com/en/lessons/repeated-measures-anova-in-r/
Esta vez vamos a involucrar un diseño que tenga en cuenta el tiempo (Factor intrasujetos = Tiempo). Factores entresujetos = FSCA, FSBA, FCCA, FCBA.
Diseño en medidas repetidas: (i) Una vía, (ii) Dos vías y (iii) Tres vías.
#Datos de una vía
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)## Warning: package 'tidyr' was built under R version 4.2.3# 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)Por ejempo: -Una vía Resp: Aceite(Limonaria) Tiempo 1: Corte 1 (30dds) Tiempo 2: Corte 2 (60dds) Tiempo 3: Corte 3 (90dds) Ho: Miu1 = Miu2 = Miu3
La función gather, reune todos los datos en una sola matríz –> Convertir formato ancho en formato largo
#Resumen estadístico
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.0Coeficiente de variación > 20% significa que hay problemas en los datos –> Casi siempre asociados a datos que no son normales.
boxplot(datos$rto ~ datos$tiempo)
#Detección de outliers
library(rstatix)## Warning: package 'rstatix' was built under R version 4.2.3##
## 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 FALSE#Supuesto de normalidad
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.380En t1 el id si es outlier, pero nos indica que no es extremo. No me preocupo por eliminarlos ya que el coeficiente de varición es menor al 20%.
#Supuesto de esfericidad (relacionado a igualdad de varianzas) Compara si la variabilidad en los tiempos sucesivos es similar.
Primero, corremos ANOVA
res.aov <- anova_test(data = datos,
dv = rto,
wid = id,
within = tiempo)
#Esfericidad
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.829Si las varianzas no son iguales, sólo se coloca que de acuerdo con la prueba de Maucly no son iguales y se necesitaría otra metodología.
Se rechaza hipótesis nula.
#Método Bonferroni
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.s…¹
## * <chr> <chr> <chr> <int> <int> <dbl> <dbl> <dbl> <dbl> <chr>
## 1 rto t1 t2 10 10 -4.97 9 0.000772 0.002 **
## 2 rto t1 t3 10 10 -13.2 9 0.000000334 0.000001 ****
## 3 rto t2 t3 10 10 -4.87 9 0.000886 0.003 **
## # … with abbreviated variable name ¹p.adj.signif#Dos vías
data("selfesteem2", package = "datarium")
datos2 = selfesteem2
datos2$treatment = gl(2, 12, 24, c('con fert', 'sin fert'))
View(datos2)
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.28library(ggplot2)## Warning: package 'ggplot2' was built under R version 4.2.3ggplot(datos2)+
aes(tiempo, rto, fill=treatment)+
geom_boxplot()
Datos atípicos
datos2 %>%
group_by(treatment, tiempo) %>%
identify_outliers(rto)## [1] treatment tiempo id rto is.outlier is.extreme
## <0 rows> (or 0-length row.names)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.104La normalidad se le hace a los residuales, no al rendimiento.
Cuando hay interacción, no se mira para arriba. Solo se analiza la interacción.
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.050Gráfico de interacción
interaction.plot(datos2$tiempo,
datos2$treatment,
datos2$rto)
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.
res.aov$`Mauchly's Test for Sphericity`## Effect W p p<.05
## 1 tiempo 0.469 0.023 *
## 2 treatment:tiempo 0.616 0.089Se evidencia la interacción, hay un cre¿uce en el tiempo 1 debido a que las medias no varian significativamente.