Anova-Medidas Repetidas

Anova de 1 via Medidas Repetidas

library(tidyverse)

## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --

## v ggplot2 3.3.5     v purrr   0.3.4
## v tibble  3.1.5     v dplyr   1.0.7
## v tidyr   1.1.4     v stringr 1.4.0
## v readr   2.0.2     v forcats 0.5.1

## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()

library(ggpubr)
library(rstatix)

## 
## Attaching package: 'rstatix'

## The following object is masked from 'package:stats':
## 
##     filter

library(datarium)
data("selfesteem",package="datarium")
selfesteem

## # A tibble: 10 x 4
##       id    t1    t2    t3
##    <int> <dbl> <dbl> <dbl>
##  1     1  4.01  5.18  7.11
##  2     2  2.56  6.91  6.31
##  3     3  3.24  4.44  9.78
##  4     4  3.42  4.71  8.35
##  5     5  2.87  3.91  6.46
##  6     6  2.05  5.34  6.65
##  7     7  3.53  5.58  6.84
##  8     8  3.18  4.37  7.82
##  9     9  3.51  4.40  8.47
## 10    10  3.04  4.49  8.58

selfesteem <- selfesteem %>% gather(key="time",value="score",t1,t2,t3) %>% convert_as_factor(id,time)
selfesteem

## # A tibble: 30 x 3
##    id    time  score
##    <fct> <fct> <dbl>
##  1 1     t1     4.01
##  2 2     t1     2.56
##  3 3     t1     3.24
##  4 4     t1     3.42
##  5 5     t1     2.87
##  6 6     t1     2.05
##  7 7     t1     3.53
##  8 8     t1     3.18
##  9 9     t1     3.51
## 10 10    t1     3.04
## # ... with 20 more rows

RESUMEN ESTADÍSTICO

selfesteem %>% group_by(time) %>% get_summary_stats(score,type="mean_sd")

## # A tibble: 3 x 5
##   time  variable     n  mean    sd
##   <fct> <chr>    <dbl> <dbl> <dbl>
## 1 t1    score       10  3.14 0.552
## 2 t2    score       10  4.93 0.863
## 3 t3    score       10  7.64 1.14

Visualización

bxp<-ggboxplot(selfesteem,x="time",y="score", add="point")
bxp

SUPUESTOS DEL MODELO

Outliers

selfesteem %>% group_by(time) %>% identify_outliers(score)

## # A tibble: 2 x 5
##   time  id    score is.outlier is.extreme
##   <fct> <fct> <dbl> <lgl>      <lgl>     
## 1 t1    6      2.05 TRUE       FALSE     
## 2 t2    2      6.91 TRUE       FALSE

Normalidad

selfesteem %>% group_by(time) %>% shapiro_test(score)

## # A tibble: 3 x 4
##   time  variable statistic     p
##   <fct> <chr>        <dbl> <dbl>
## 1 t1    score        0.967 0.859
## 2 t2    score        0.876 0.117
## 3 t3    score        0.923 0.380

ggqqplot(selfesteem,"score",facet.by="time")

Cálculo del ANOVA

res.aov<-anova_test(data=selfesteem,dv=score,wid=id,within=time)
get_anova_table(res.aov)

## ANOVA Table (type III tests)
## 
##   Effect DFn DFd      F        p p<.05   ges
## 1   time   2  18 55.469 2.01e-08     * 0.829

Comparaciones de medias

pwc<-selfesteem %>% pairwise_t_test(score~time,paired=TRUE,p.adjust.method="bonferroni")
pwc

## # A tibble: 3 x 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 score t1     t2        10    10     -4.97     9 0.000772     2e-3 **          
## 2 score t1     t3        10    10    -13.2      9 0.000000334  1e-6 ****        
## 3 score t2     t3        10    10     -4.87     9 0.000886     3e-3 **

Reporte

pwc <- pwc %>% add_xy_position(x="time")
  bxp+
  stat_pvalue_manual(pwc) +
    labs(
         subtitle = get_test_label(res.aov,detailed=TRUE),
         caption=get_pwc_label(pwc)
         )