Con el propósito de solucionar un método de análisis químico (húmedo), se dispone de una gran cantidad de muestra que ha sido bien molida y mezclada. Se recogieron datos adecuados para poder realizar un ANOVA.

Para empezar, se disolvieron 10 g de muestra en 100 mL de agua.El siguiente paso fue tomar varias alícuotas diferentes y recoger múltiples lecturas de cada una de ellas. Se colocaron seis alícuotas en cada uno de los seis frascos y se repitieron seis mediciones en cada uno de ellos. Estas alícuotas dieron los resultados mostrados en la Tabla 9-2.

X= Número de Alícuota

Considerando los datos realice:

Análisis de Anova, evaluando todos los supuestos. Analice los resultados obtenidos Genere algunas ideas que usted considere puedan describir el comportamiento de los datos y formule algunas mejoras para el mencionado experimento.

Hipótesis

Ho= μ1=μ2=μ3=μ4=μ5=μ6 las medias de las alícuotas son iguales

HA= μ1≠μn Algunas de las medias de las alícuotas son diferentes

library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v ggplot2 3.3.5     v purrr   0.3.4
## v tibble  3.1.4     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
ID<-c(1,2,3,4,5,6)
X1<-c(7.25,7.68,7.76,8.10,7.50,7.58)
X2<-c(10.07,9.02,9.51,10.64,10.27,9.64)
X3<-c(5.96,6.66,5.87,6.95,6.54,6.29)
X4<-c(7.10,6.10,6.27,5.99,6.32,5.54)
X5<-c(5.74,6.90,6.29,6.37,5.99,6.58)
X6<-c(4.74,6.75,6.71,6.51,5.95,6.50)
df<-data.frame(ID=ID,X1=X1,X2=X2,X3=X3,X4=X4,X5=X5,X6=X6)
df
##   ID   X1    X2   X3   X4   X5   X6
## 1  1 7.25 10.07 5.96 7.10 5.74 4.74
## 2  2 7.68  9.02 6.66 6.10 6.90 6.75
## 3  3 7.76  9.51 5.87 6.27 6.29 6.71
## 4  4 8.10 10.64 6.95 5.99 6.37 6.51
## 5  5 7.50 10.27 6.54 6.32 5.99 5.95
## 6  6 7.58  9.64 6.29 5.54 6.58 6.50
df <- df %>% gather(key="alicuota",value="score",X1,X2,X3,X4,X5,X6) %>% convert_as_factor(ID,alicuota)
df
##    ID alicuota score
## 1   1       X1  7.25
## 2   2       X1  7.68
## 3   3       X1  7.76
## 4   4       X1  8.10
## 5   5       X1  7.50
## 6   6       X1  7.58
## 7   1       X2 10.07
## 8   2       X2  9.02
## 9   3       X2  9.51
## 10  4       X2 10.64
## 11  5       X2 10.27
## 12  6       X2  9.64
## 13  1       X3  5.96
## 14  2       X3  6.66
## 15  3       X3  5.87
## 16  4       X3  6.95
## 17  5       X3  6.54
## 18  6       X3  6.29
## 19  1       X4  7.10
## 20  2       X4  6.10
## 21  3       X4  6.27
## 22  4       X4  5.99
## 23  5       X4  6.32
## 24  6       X4  5.54
## 25  1       X5  5.74
## 26  2       X5  6.90
## 27  3       X5  6.29
## 28  4       X5  6.37
## 29  5       X5  5.99
## 30  6       X5  6.58
## 31  1       X6  4.74
## 32  2       X6  6.75
## 33  3       X6  6.71
## 34  4       X6  6.51
## 35  5       X6  5.95
## 36  6       X6  6.50

Resumen Estadístico

df %>% group_by(alicuota) %>% get_summary_stats(score,type="mean_sd")
## # A tibble: 6 x 5
##   alicuota variable     n  mean    sd
##   <fct>    <chr>    <dbl> <dbl> <dbl>
## 1 X1       score        6  7.64 0.284
## 2 X2       score        6  9.86 0.583
## 3 X3       score        6  6.38 0.418
## 4 X4       score        6  6.22 0.513
## 5 X5       score        6  6.31 0.413
## 6 X6       score        6  6.19 0.767

Visualización

bxp<-ggboxplot(df,x="alicuota",y="score", add="point")
bxp

Outliers

df %>% group_by(alicuota) %>% identify_outliers(score)
## # A tibble: 4 x 5
##   alicuota ID    score is.outlier is.extreme
##   <fct>    <fct> <dbl> <lgl>      <lgl>     
## 1 X1       4      8.1  TRUE       FALSE     
## 2 X4       1      7.1  TRUE       FALSE     
## 3 X4       6      5.54 TRUE       FALSE     
## 4 X6       1      4.74 TRUE       FALSE

Verificación de Normalidad

df %>% group_by(alicuota) %>% shapiro_test(score)
## # A tibble: 6 x 4
##   alicuota variable statistic      p
##   <fct>    <chr>        <dbl>  <dbl>
## 1 X1       score        0.980 0.953 
## 2 X2       score        0.985 0.973 
## 3 X3       score        0.954 0.774 
## 4 X4       score        0.937 0.635 
## 5 X5       score        0.990 0.989 
## 6 X6       score        0.771 0.0320
ggqqplot(df,"score",facet.by="alicuota")

ANOVA

res.aov<-anova_test(data=df,dv=score,wid=ID,within=alicuota)
get_anova_table(res.aov)
## ANOVA Table (type III tests)
## 
##     Effect DFn DFd      F        p p<.05   ges
## 1 alicuota   5  25 46.592 7.12e-12     * 0.888

Comparaciones de medias

pwc<-df %>% pairwise_t_test(score~alicuota,paired=TRUE,p.adjust.method="bonferroni")
pwc
## # A tibble: 15 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 X1     X2         6     6   -9.07       5 0.000272  0.004    **          
##  2 score X1     X3         6     6    9.29       5 0.000244  0.004    **          
##  3 score X1     X4         6     6    4.88       5 0.005     0.069    ns          
##  4 score X1     X5         6     6    9.02       5 0.00028   0.004    **          
##  5 score X1     X6         6     6    6.06       5 0.002     0.026    *           
##  6 score X2     X3         6     6   14.2        5 0.0000311 0.000466 ***         
##  7 score X2     X4         6     6   12.7        5 0.0000534 0.000801 ***         
##  8 score X2     X5         6     6    9.65       5 0.000203  0.003    **          
##  9 score X2     X6         6     6    7.94       5 0.00051   0.008    **          
## 10 score X3     X4         6     6    0.488      5 0.646     1        ns          
## 11 score X3     X5         6     6    0.369      5 0.727     1        ns          
## 12 score X3     X6         6     6    0.631      5 0.556     1        ns          
## 13 score X4     X5         6     6   -0.258      5 0.806     1        ns          
## 14 score X4     X6         6     6    0.0533     5 0.96      1        ns          
## 15 score X5     X6         6     6    0.606      5 0.571     1        ns

Reporte Final

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

Conclusión

Inicialmente vemos en el diagrama de cajas y bigotes las diferencias en los datos obtenidos para las alícuotas 1 y 2, respecto a las demás. A partir de allí se evalúan los outlier donde presentan 4, sin embargo los mismos no son significativos, por lo que no es necesario descartarlos de la muestra. Los datos presentan un comportamiento normal según la prueba de Shapiro, con un valor de p >0.05, al 95% de confianza, en cambio la alícuota 6, presenta un p<0.05, por lo que sus datos no se ajustan a la normalidad, quizás por la presencia del dato anómalo. Con el ANOVA se demuestra un efecto significativo de la alícuota 1 con p<0.05 al 95% de confianza. En la prueba de Bonferroni, se observan diferencias significativas de algunos de los niveles críticos asociados a cada comparación con un valor de p<0.05. En los diagramas de cajas y bigotes se pueden apreciar mejor las diferencias significativas de X1 y X2 con relación a las demás X’s. Con todo lo expuesto anteriormente podemos rechazar la hipótesis nula de igualdad de medias y concluir que Las dos primeras muestras son diferentes. Con los resultados estadísticos obtenidos, se puede afirmar que los valores de X2 se encuentran fuera de rango, lo recomendable para mejorar el experimento podría ser repetir el análisis con esa alícuota. Algunas de las posibles causas de ese error podrían ser contaminación de material o presencias de errores sistemáticos y aleatorios.