Medidas repetidas
2024-10-01
\[\text{Valeria Riveros Parra y Juan Felipe Paz Paiba}\] \[\text{Analisis de varianza para medidas repetidas de una vía}\]
library(ggpubr)## Warning: package 'ggpubr' was built under R version 4.2.3## Loading required package: ggplot2## Warning: package 'ggplot2' was built under R version 4.2.3library(rstatix)## Warning: package 'rstatix' was built under R version 4.2.3##
## Attaching package: 'rstatix'## The following object is masked from 'package:stats':
##
## filterlibrary(readxl)## Warning: package 'readxl' was built under R version 4.2.3library(multcomp)## Loading required package: mvtnorm## Warning: package 'mvtnorm' was built under R version 4.2.3## Loading required package: survival## Loading required package: TH.data## Warning: package 'TH.data' was built under R version 4.2.3## Loading required package: MASS## Warning: package 'MASS' was built under R version 4.2.3##
## Attaching package: 'MASS'## The following object is masked from 'package:rstatix':
##
## select##
## Attaching package: 'TH.data'## The following object is masked from 'package:MASS':
##
## geyserlibrary(lme4)## Warning: package 'lme4' was built under R version 4.2.3## Loading required package: Matrix## Warning: package 'Matrix' was built under R version 4.2.3library(dplyr)## Warning: package 'dplyr' was built under R version 4.2.3##
## Attaching package: 'dplyr'## The following object is masked from 'package:MASS':
##
## select## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union# Datos en formato ancho
data= read_xlsx ("C:/Users/juanf/downloads/m_repetidas.xlsx", sheet= "biomasa_1")
data## # A tibble: 30 × 13
## `0_días` `05_días` `10_días` `15_días` `20_días` `25_días` `30_días`
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0.00001 2.95 3.31 3.12 3.27 3.86 3.65
## 2 0.00001 2.53 3.02 3.35 3.41 3.73 3.93
## 3 0.00001 2.68 3.07 3.39 2.72 3.52 4.00
## 4 0.00001 2.77 2.91 3.11 3.42 2.99 3.60
## 5 0.00001 2.83 3.15 3.46 3.15 3.08 3.31
## 6 0.00001 2.88 3.18 3.27 3.75 3.24 3.56
## 7 0.00001 2.92 3.21 2.67 3.37 3.81 3.75
## 8 0.00001 2.86 2.96 3.24 2.9 3.45 3.86
## 9 0.00001 2.99 3.27 3.04 3.51 3.03 3.22
## 10 0.00001 3.02 3.14 3.19 3.44 3.57 3.96
## # ℹ 20 more rows
## # ℹ 6 more variables: `35_días` <dbl>, `40_días` <dbl>, `45_días` <dbl>,
## # `50_días` <dbl>, `55_días` <dbl>, `60_días` <dbl>#Datos en formato largo
# Gather columns t1, t2 and t3 into long format
# Convert id and time into factor variables
data_larga <- data %>%
gather(key = "time", value = "Biomasa", `0_días`, `05_días`, `10_días`, `15_días`, `20_días`, `25_días`, `30_días`, `35_días`, `40_días`, `45_días`, `50_días`, `55_días`, `60_días`) %>%
convert_as_factor(time)
data_larga$id= rep(seq(1,30),13)
head(data_larga, 5) ## # A tibble: 5 × 3
## time Biomasa id
## <fct> <dbl> <int>
## 1 0_días 0.00001 1
## 2 0_días 0.00001 2
## 3 0_días 0.00001 3
## 4 0_días 0.00001 4
## 5 0_días 0.00001 5# Cargar la librería necesaria
library(collapsibleTree)## Warning: package 'collapsibleTree' was built under R version 4.2.3# Crear un árbol colapsible
tree <- collapsibleTree(data_larga, c("time", "Biomasa", "id"), collapsed = TRUE)
tree#Resumen estadistico
data_larga %>%
group_by(time) %>%
get_summary_stats(Biomasa, type = "mean_sd")## # A tibble: 13 × 5
## time variable n mean sd
## <fct> <fct> <dbl> <dbl> <dbl>
## 1 0_días Biomasa 30 0 0
## 2 05_días Biomasa 30 2.86 0.175
## 3 10_días Biomasa 30 3.11 0.144
## 4 15_días Biomasa 30 3.20 0.238
## 5 20_días Biomasa 30 3.31 0.257
## 6 25_días Biomasa 30 3.56 0.233
## 7 30_días Biomasa 30 3.72 0.242
## 8 35_días Biomasa 30 3.95 0.224
## 9 40_días Biomasa 30 4.32 0.2
## 10 45_días Biomasa 30 4.57 0.305
## 11 50_días Biomasa 30 4.83 0.34
## 12 55_días Biomasa 30 5.18 0.301
## 13 60_días Biomasa 30 5.83 0.304# Gráfico de líneas
ggplot(data_larga, aes(x = time, y = Biomasa)) +
geom_line() +
geom_point() + scale_y_continuous(breaks = seq(0, max(data_larga$Biomasa, na.rm = TRUE), by = 5)) +
labs(title = "Biomasa a lo largo del tiempo",
x = "Tiempo",
y = "Biomasa") +
theme_minimal()
# Gráfico de cajas (boxplot)
ggplot(data_larga, aes(x = time, y = Biomasa)) +
geom_boxplot() + scale_fill_manual(values = rainbow(n_distinct(data_larga$time))) +
labs(title = "Distribución de Biomasa a lo largo del tiempo",
x = "Tiempo",
y = "Biomasa") +
theme_minimal()
#Supuestos
#outliers
data_larga %>%
group_by(time) %>%
identify_outliers(Biomasa)## # A tibble: 12 × 5
## time Biomasa id is.outlier is.extreme
## <fct> <dbl> <int> <lgl> <lgl>
## 1 10_días 3.45 14 TRUE FALSE
## 2 10_días 2.78 29 TRUE FALSE
## 3 20_días 2.72 3 TRUE FALSE
## 4 20_días 2.75 16 TRUE FALSE
## 5 25_días 2.99 4 TRUE FALSE
## 6 25_días 3.08 5 TRUE FALSE
## 7 25_días 3.03 9 TRUE FALSE
## 8 30_días 3.13 13 TRUE FALSE
## 9 35_días 3.40 1 TRUE FALSE
## 10 35_días 3.54 10 TRUE FALSE
## 11 35_días 4.41 25 TRUE FALSE
## 12 45_días 5.35 16 TRUE FALSE#supuestos de la normalidad (sin embargo solo para determinar la normalidad de la respuesta, luego se realizar al ralizar el análisis de varianza)
hist(data_larga$Biomasa)
#Supuesto de esfericidad (una manera de comparar la variabilidad entre los pares de tiempo
#Comparacion de varianzas entre pares de tiempos
res.aov <- anova_test(data = data_larga, dv = Biomasa, wid = id, within = time) #`Mauchly's Test for Sphericity`
res.aov #Los datos proveen evidencia en contra de la hipotesis nula, lo cual rechaza la hipotesis de esfericidad## ANOVA Table (type III tests)
##
## $ANOVA
## Effect DFn DFd F p p<.05 ges
## 1 time 12 348 1083.509 3.19e-267 * 0.97
##
## $`Mauchly's Test for Sphericity`
## Effect W p p<.05
## 1 time 0.000513 6.82e-11 *
##
## $`Sphericity Corrections`
## Effect GGe DF[GG] p[GG] p[GG]<.05 HFe DF[HF] p[HF]
## 1 time 0.518 6.22, 180.25 1.05e-139 * 0.675 8.1, 234.78 3.59e-181
## p[HF]<.05
## 1 *get_anova_table(res.aov) #se rechaza la hipotesis de esfericidad## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 time 6.22 180.25 1083.509 1.05e-139 * 0.97#transformación de los datos para cumplir supuestos
# Aplicar la transformación de Box-Cox
boxcox_result <- boxcox(Biomasa ~ time + (1 | id), data = data_larga, lambda = seq(-2, 2, by = 0.1))
best_lambda <- boxcox_result$x[which.max(boxcox_result$y)]
data_larga <- data_larga %>%
mutate(Biomasa_boxcox = ifelse(best_lambda == 0, log(Biomasa), (Biomasa^best_lambda - 1) / best_lambda))
head(data_larga, 5)## # A tibble: 5 × 4
## time Biomasa id Biomasa_boxcox
## <fct> <dbl> <int> <dbl>
## 1 0_días 0.00001 1 -4999999999.
## 2 0_días 0.00001 2 -4999999999.
## 3 0_días 0.00001 3 -4999999999.
## 4 0_días 0.00001 4 -4999999999.
## 5 0_días 0.00001 5 -4999999999.#boxplot datos tranformados con raíz cuadrada
bxp2 = ggboxplot(data_larga, x= "time", y= "Biomasa_boxcox", add= "point")
bxp2
#Comparacion de varianzas entre pares de tiempos con boxcox
res.aov2 <- anova_test(data = data_larga, dv = Biomasa_boxcox, wid = id, within = time)
res.aov2## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 time 12 348 NaN NaN <NA> 0get_anova_table(res.aov2)## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 time 12 348 NaN NaN <NA> 0#transformación de los datos para cumplir supuestos 2
# Aplicar funcion logaritmica a los datos
data_larga$log_Biomasa= log(data_larga$Biomasa)
#boxplot datos tranformados con logaritmo
bxp3 = ggboxplot(data_larga, x= "time", y= "log_Biomasa", add= "point")
bxp3
#Comparacion de varianzas entre pares de tiempos con logaritmo de la Biomasa
res.aov3 <- anova_test(data = data_larga, dv = log_Biomasa, wid = id, within = time)
res.aov3## ANOVA Table (type III tests)
##
## $ANOVA
## Effect DFn DFd F p p<.05 ges
## 1 time 12 348 106162.5 0 * 1
##
## $`Mauchly's Test for Sphericity`
## Effect W p p<.05
## 1 time 0.000699 6.46e-10 *
##
## $`Sphericity Corrections`
## Effect GGe DF[GG] p[GG] p[GG]<.05 HFe DF[HF] p[HF]
## 1 time 0.527 6.33, 183.48 9.881313e-324 * 0.69 8.28, 240.22 0
## p[HF]<.05
## 1 *get_anova_table(res.aov3)## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 time 6.33 183.48 106162.5 9.881313e-324 * 1#transformación de los datos para cumplir supuestos 3
# Aplicar funcion logaritmica a los datos
data_larga$sqrt_Biomasa= sqrt(data_larga$Biomasa)
#boxplot datos tranformados con raíz cuadrada
bxp4 = ggboxplot(data_larga, x= "time", y= "sqrt_Biomasa", add= "point")
bxp4
res.aov4 <- anova_test(data = data_larga, dv = sqrt_Biomasa, wid = id, within = time)
res.aov4## ANOVA Table (type III tests)
##
## $ANOVA
## Effect DFn DFd F p p<.05 ges
## 1 time 12 348 3038.58 0 * 0.989
##
## $`Mauchly's Test for Sphericity`
## Effect W p p<.05
## 1 time 0.000799 1.67e-09 *
##
## $`Sphericity Corrections`
## Effect GGe DF[GG] p[GG] p[GG]<.05 HFe DF[HF] p[HF]
## 1 time 0.534 6.41, 185.89 6.35e-185 * 0.702 8.42, 244.3 3.45e-242
## p[HF]<.05
## 1 *get_anova_table(res.aov4)## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 time 6.41 185.89 3038.58 6.35e-185 * 0.989# Datos en formato ancho con razon de cumplir supuestos
supuestos_1= read_xlsx ("C:/Users/juanf/downloads/m_repetidas.xlsx", sheet= "supuestos_1")
supuestos_1## # A tibble: 30 × 12
## `05_días` `10_días` `15_días` `20_días` `25_días` `30_días` `35_días`
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 2.95 3.31 3.12 3.27 3.86 3.65 3.40
## 2 2.53 3.02 3.35 3.41 3.73 3.93 3.68
## 3 2.68 3.07 3.39 2.72 3.52 4.00 3.89
## 4 2.77 2.91 3.11 3.42 2.99 3.60 4.07
## 5 2.83 3.15 3.46 3.15 3.08 3.31 3.68
## 6 2.88 3.18 3.27 3.75 3.24 3.56 3.64
## 7 2.92 3.21 2.67 3.37 3.81 3.75 3.98
## 8 2.86 2.96 3.24 2.9 3.45 3.86 3.90
## 9 2.99 3.27 3.04 3.51 3.03 3.22 4.02
## 10 3.02 3.14 3.19 3.44 3.57 3.96 3.54
## # ℹ 20 more rows
## # ℹ 5 more variables: `40_días` <dbl>, `45_días` <dbl>, `50_días` <dbl>,
## # `55_días` <dbl>, `60_días` <dbl>#Datos para cumplir supuestos en formato largo
# Gather columns t1, t2 and t3 into long format
# Convert id and time into factor variables
data_supuestos_larga <- supuestos_1 %>%
gather(key = "time", value = "Biomasa", `05_días`, `10_días`, `15_días`, `20_días`, `25_días`, `30_días`, `35_días`, `40_días`, `45_días`, `50_días`, `55_días`, `60_días`) %>%
convert_as_factor(time)
data_supuestos_larga$id= rep(seq(1,30),12)
head(data_supuestos_larga, 5) ## # A tibble: 5 × 3
## time Biomasa id
## <fct> <dbl> <int>
## 1 05_días 2.95 1
## 2 05_días 2.53 2
## 3 05_días 2.68 3
## 4 05_días 2.77 4
## 5 05_días 2.83 5# Gráfico de cajas (boxplot)
ggplot(data_supuestos_larga, aes(x = time, y = Biomasa)) +
geom_boxplot() + scale_fill_manual(values = rainbow(n_distinct(data_larga$time))) +
labs(title = "Distribución de Biomasa a lo largo del tiempo",
x = "Tiempo",
y = "Biomasa") +
theme_minimal()
#Supuesto de esfericidad (una manera de comparar la variabilidad entre los pares de tiempo
#Comparacion de varianzas entre pares de tiempos
res.aov5 <- anova_test(data = data_supuestos_larga, dv = Biomasa, wid = id, within = time) #`Mauchly's Test for Sphericity`
res.aov5 #Los datos proveen evidencia en contra de la hipotesis nula, lo cual rechaza la hipotesis de esfericidad## ANOVA Table (type III tests)
##
## $ANOVA
## Effect DFn DFd F p p<.05 ges
## 1 time 11 319 419.484 2.49e-182 * 0.926
##
## $`Mauchly's Test for Sphericity`
## Effect W p p<.05
## 1 time 0.005 1.08e-06 *
##
## $`Sphericity Corrections`
## Effect GGe DF[GG] p[GG] p[GG]<.05 HFe DF[HF] p[HF]
## 1 time 0.56 6.16, 178.65 2.78e-103 * 0.728 8, 232.11 2.11e-133
## p[HF]<.05
## 1 *get_anova_table(res.aov5) #se rechaza la hipotesis de esfericidad## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 time 6.16 178.65 419.484 2.78e-103 * 0.926\[H_0= \mu_{t_0}=\mu_{t_5}=\mu_{t_{10}}=\mu_{t_{15}}=\mu_{t_{20}}=\mu_{t_{25}}=\mu_{t_{30}}=\mu_{t_{35}}=\mu_{t_{40}}=\mu_{t_{45}}=\mu_{t_{50}}=\mu_{t_{55}}=\mu_{t_{60}}\]
# Obtener la tabla de ANOVA con los datos de la Comparacion de varianzas entre pares de tiempos de los datos sin transformación
anova_table <- get_anova_table(res.aov)
print(anova_table)## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 time 6.22 180.25 1083.509 1.05e-139 * 0.97# pairwise comparisons
pwc <- data_larga %>%
pairwise_t_test(
Biomasa ~ time, paired = TRUE,
p.adjust.method = "bonferroni"
)
pwc## # A tibble: 78 × 10
## .y. group1 group2 n1 n2 statistic df p p.adj
## * <chr> <chr> <chr> <int> <int> <dbl> <dbl> <dbl> <dbl>
## 1 Biomasa 0_días 05_días 30 30 -89.9 29 4.94e-37 3.85e-35
## 2 Biomasa 0_días 10_días 30 30 -118. 29 1.82e-40 1.42e-38
## 3 Biomasa 0_días 15_días 30 30 -73.5 29 1.64e-34 1.28e-32
## 4 Biomasa 0_días 20_días 30 30 -70.6 29 5.33e-34 4.16e-32
## 5 Biomasa 0_días 25_días 30 30 -83.8 29 3.79e-36 2.96e-34
## 6 Biomasa 0_días 30_días 30 30 -84.1 29 3.33e-36 2.60e-34
## 7 Biomasa 0_días 35_días 30 30 -96.7 29 5.90e-38 4.6 e-36
## 8 Biomasa 0_días 40_días 30 30 -118. 29 1.69e-40 1.32e-38
## 9 Biomasa 0_días 45_días 30 30 -82.1 29 6.85e-36 5.34e-34
## 10 Biomasa 0_días 50_días 30 30 -77.8 29 3.21e-35 2.5 e-33
## # ℹ 68 more rows
## # ℹ 1 more variable: p.adj.signif <chr>\[\text{Los datos proveen evidencia en contra de H_0 lo cual indica que la biomasa no se mantuvo constante con el paso del tiempo}\]
# Modelo Lineal Generalizado Mixto
model <- lmer(sqrt_Biomasa ~ time + (1 | id), data = data_larga)
# Comparaciones por pares con letras de agrupamiento (Tukey)
comp <- glht(model, linfct = mcp(time = "Tukey"))
summary(comp)## Warning in RET$pfunction("adjusted", ...): Completion with error > abseps
## Warning in RET$pfunction("adjusted", ...): Completion with error > abseps
## Warning in RET$pfunction("adjusted", ...): Completion with error > abseps
## Warning in RET$pfunction("adjusted", ...): Completion with error > abseps##
## Simultaneous Tests for General Linear Hypotheses
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: lmer(formula = sqrt_Biomasa ~ time + (1 | id), data = data_larga)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## 05_días - 0_días == 0 1.68839 0.01522 110.901 <0.01 ***
## 10_días - 0_días == 0 1.75960 0.01522 115.579 <0.01 ***
## 15_días - 0_días == 0 1.78315 0.01522 117.125 <0.01 ***
## 20_días - 0_días == 0 1.81362 0.01522 119.127 <0.01 ***
## 25_días - 0_días == 0 1.88254 0.01522 123.654 <0.01 ***
## 30_días - 0_días == 0 1.92437 0.01522 126.402 <0.01 ***
## 35_días - 0_días == 0 1.98328 0.01522 130.271 <0.01 ***
## 40_días - 0_días == 0 2.07517 0.01522 136.307 <0.01 ***
## 45_días - 0_días == 0 2.13406 0.01522 140.175 <0.01 ***
## 50_días - 0_días == 0 2.19237 0.01522 144.005 <0.01 ***
## 55_días - 0_días == 0 2.27128 0.01522 149.188 <0.01 ***
## 60_días - 0_días == 0 2.41033 0.01522 158.322 <0.01 ***
## 10_días - 05_días == 0 0.07122 0.01522 4.678 <0.01 ***
## 15_días - 05_días == 0 0.09476 0.01522 6.224 <0.01 ***
## 20_días - 05_días == 0 0.12523 0.01522 8.226 <0.01 ***
## 25_días - 05_días == 0 0.19415 0.01522 12.753 <0.01 ***
## 30_días - 05_días == 0 0.23599 0.01522 15.501 <0.01 ***
## 35_días - 05_días == 0 0.29490 0.01522 19.370 <0.01 ***
## 40_días - 05_días == 0 0.38678 0.01522 25.405 <0.01 ***
## 45_días - 05_días == 0 0.44567 0.01522 29.274 <0.01 ***
## 50_días - 05_días == 0 0.50399 0.01522 33.104 <0.01 ***
## 55_días - 05_días == 0 0.58289 0.01522 38.287 <0.01 ***
## 60_días - 05_días == 0 0.72194 0.01522 47.420 <0.01 ***
## 15_días - 10_días == 0 0.02354 0.01522 1.546 0.9457
## 20_días - 10_días == 0 0.05401 0.01522 3.548 0.0226 *
## 25_días - 10_días == 0 0.12293 0.01522 8.075 <0.01 ***
## 30_días - 10_días == 0 0.16477 0.01522 10.823 <0.01 ***
## 35_días - 10_días == 0 0.22368 0.01522 14.692 <0.01 ***
## 40_días - 10_días == 0 0.31556 0.01522 20.728 <0.01 ***
## 45_días - 10_días == 0 0.37446 0.01522 24.596 <0.01 ***
## 50_días - 10_días == 0 0.43277 0.01522 28.426 <0.01 ***
## 55_días - 10_días == 0 0.51168 0.01522 33.609 <0.01 ***
## 60_días - 10_días == 0 0.65072 0.01522 42.743 <0.01 ***
## 20_días - 15_días == 0 0.03047 0.01522 2.002 0.7326
## 25_días - 15_días == 0 0.09939 0.01522 6.528 <0.01 ***
## 30_días - 15_días == 0 0.14123 0.01522 9.277 <0.01 ***
## 35_días - 15_días == 0 0.20014 0.01522 13.146 <0.01 ***
## 40_días - 15_días == 0 0.29202 0.01522 19.181 <0.01 ***
## 45_días - 15_días == 0 0.35091 0.01522 23.050 <0.01 ***
## 50_días - 15_días == 0 0.40923 0.01522 26.880 <0.01 ***
## 55_días - 15_días == 0 0.48813 0.01522 32.063 <0.01 ***
## 60_días - 15_días == 0 0.62718 0.01522 41.196 <0.01 ***
## 25_días - 20_días == 0 0.06892 0.01522 4.527 <0.01 ***
## 30_días - 20_días == 0 0.11075 0.01522 7.275 <0.01 ***
## 35_días - 20_días == 0 0.16966 0.01522 11.144 <0.01 ***
## 40_días - 20_días == 0 0.26155 0.01522 17.180 <0.01 ***
## 45_días - 20_días == 0 0.32044 0.01522 21.048 <0.01 ***
## 50_días - 20_días == 0 0.37875 0.01522 24.878 <0.01 ***
## 55_días - 20_días == 0 0.45766 0.01522 30.061 <0.01 ***
## 60_días - 20_días == 0 0.59671 0.01522 39.195 <0.01 ***
## 30_días - 25_días == 0 0.04184 0.01522 2.748 0.2267
## 35_días - 25_días == 0 0.10075 0.01522 6.618 <0.01 ***
## 40_días - 25_días == 0 0.19263 0.01522 12.653 <0.01 ***
## 45_días - 25_días == 0 0.25152 0.01522 16.521 <0.01 ***
## 50_días - 25_días == 0 0.30984 0.01522 20.352 <0.01 ***
## 55_días - 25_días == 0 0.38874 0.01522 25.535 <0.01 ***
## 60_días - 25_días == 0 0.52779 0.01522 34.668 <0.01 ***
## 35_días - 30_días == 0 0.05891 0.01522 3.870 <0.01 **
## 40_días - 30_días == 0 0.15079 0.01522 9.905 <0.01 ***
## 45_días - 30_días == 0 0.20969 0.01522 13.773 <0.01 ***
## 50_días - 30_días == 0 0.26800 0.01522 17.604 <0.01 ***
## 55_días - 30_días == 0 0.34691 0.01522 22.786 <0.01 ***
## 60_días - 30_días == 0 0.48596 0.01522 31.920 <0.01 ***
## 40_días - 35_días == 0 0.09188 0.01522 6.035 <0.01 ***
## 45_días - 35_días == 0 0.15078 0.01522 9.904 <0.01 ***
## 50_días - 35_días == 0 0.20909 0.01522 13.734 <0.01 ***
## 55_días - 35_días == 0 0.28800 0.01522 18.917 <0.01 ***
## 60_días - 35_días == 0 0.42704 0.01522 28.050 <0.01 ***
## 45_días - 40_días == 0 0.05889 0.01522 3.868 <0.01 **
## 50_días - 40_días == 0 0.11721 0.01522 7.699 <0.01 ***
## 55_días - 40_días == 0 0.19611 0.01522 12.882 <0.01 ***
## 60_días - 40_días == 0 0.33516 0.01522 22.015 <0.01 ***
## 50_días - 45_días == 0 0.05831 0.01522 3.830 <0.01 **
## 55_días - 45_días == 0 0.13722 0.01522 9.013 <0.01 ***
## 60_días - 45_días == 0 0.27627 0.01522 18.147 <0.01 ***
## 55_días - 50_días == 0 0.07891 0.01522 5.183 <0.01 ***
## 60_días - 50_días == 0 0.21796 0.01522 14.316 <0.01 ***
## 60_días - 55_días == 0 0.13905 0.01522 9.133 <0.01 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)# Obtener las letras de agrupamiento
cld_result <- cld(comp, level = 0.05)## Warning in RET$pfunction("adjusted", ...): Completion with error > abseps
## Warning in RET$pfunction("adjusted", ...): Completion with error > abseps
## Warning in RET$pfunction("adjusted", ...): Completion with error > abseps
## Warning in RET$pfunction("adjusted", ...): Completion with error > abseps
## Warning in RET$pfunction("adjusted", ...): Completion with error > absepsprint(cld_result)## 0_días 05_días 10_días 15_días 20_días 25_días 30_días 35_días 40_días 45_días
## "a" "b" "c" "cd" "d" "e" "e" "f" "g" "h"
## 50_días 55_días 60_días
## "i" "j" "k"# A mayor numero de * la diferencia es mas significativa#Las letras de agrupamiento, indican los grupos de tiempos que no difieren significativamente entre sí en términos de la media de la biomasa. En este caso, los tiempos se agrupan en 8 grupos, cada uno identificado por una letra diferente:
#Grupo "a": 0 días;
#Grupo "b": 5 días;
#Grupo "c": 10 días;
#Grupo "cd": 15 días (No hay diferencias significativas entre 10 - 15 y 15 - 20, pero si 10 - 20 )
#Grupo "d": 20 días
#Grupo "e": 25 días y 30 significa que no hay diferencia significativa en la biomasa entre estos dos tiempos
#Grupo "f": 35 días
#Grupo "g": 40 días
#Grupo "h": 45 días
#Grupo "i": 50 días
#Grupo "j": 55 días
#Grupo "k": 60 días# Crea el gráfico de cajas
bxp <- ggplot(data_larga, aes(x = time, y = biomasa)) +
geom_boxplot()
# Realiza la prueba de comparaciones múltiples
pwc <- data_larga %>%
pairwise_t_test(Biomasa ~ time, paired = TRUE, p.adjust.method = "bonferroni")
# Asegúrate de que el objeto "pwc" sea numérico
pwc$p.adj <- as.numeric(pwc$p.adj)
# Por ejemplo, si quieres utilizar la variable "biomasa" como "y.position"
data_larga$y.position <- data_larga$Biomasa
# Luego, puedes utilizar la variable "y.position" en tu función
# Por ejemplo, en un gráfico de cajas
bxp <- ggplot(data_larga, aes(x = time, y = y.position)) +
geom_boxplot()
# Realiza la prueba de comparaciones múltiples
pwc <- data_larga %>%
pairwise_t_test(Biomasa ~ time, paired = TRUE, p.adjust.method = "bonferroni")
# Agrega las posiciones x e y a los resultados de la prueba de comparaciones múltiples
pwc <- pwc %>% add_xy_position(x = "time")
# Agrega las etiquetas de valor p al gráfico
bxp_final <- bxp +
stat_pvalue_manual(pwc) +
labs(
subtitle = get_test_label(res.aov, detailed = TRUE),
caption = get_pwc_label(pwc)
)
# Muestra el gráfico
print(bxp_final)