setwd("G:/Mi unidad/Agrosavia/SOCODEVI/Informe final")
datos<-read.table("Barril_grano4.csv", header=T, sep=';')
datos$Vereda<-as.factor(datos$Vereda)
datos$id<-as.factor(datos$id)
datos$Rep<-as.factor(datos$Rep)
datos$Estado<-as.factor(datos$Estado)
library(ggplot2)
library(Rmisc)
## Cargando paquete requerido: lattice
## Cargando paquete requerido: plyr
library(dplyr)
##
## Adjuntando el paquete: 'dplyr'
## The following objects are masked from 'package:plyr':
##
## arrange, count, desc, failwith, id, mutate, rename, summarise,
## summarize
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ lubridate 1.9.4 ✔ tibble 3.2.1
## ✔ purrr 1.0.4 ✔ tidyr 1.3.1
## ✔ readr 2.1.5
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::arrange() masks plyr::arrange()
## ✖ purrr::compact() masks plyr::compact()
## ✖ dplyr::count() masks plyr::count()
## ✖ dplyr::desc() masks plyr::desc()
## ✖ dplyr::failwith() masks plyr::failwith()
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::id() masks plyr::id()
## ✖ dplyr::lag() masks stats::lag()
## ✖ dplyr::mutate() masks plyr::mutate()
## ✖ dplyr::rename() masks plyr::rename()
## ✖ dplyr::summarise() masks plyr::summarise()
## ✖ dplyr::summarize() masks plyr::summarize()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(ggpubr)
##
## Adjuntando el paquete: 'ggpubr'
##
## The following object is masked from 'package:plyr':
##
## mutate
library(rstatix)
##
## Adjuntando el paquete: 'rstatix'
##
## The following objects are masked from 'package:plyr':
##
## desc, mutate
##
## The following object is masked from 'package:stats':
##
## filter
library(emmeans)
## Welcome to emmeans.
## Caution: You lose important information if you filter this package's results.
## See '? untidy'
##Summary statistics
summ<-datos %>%
group_by(Vereda, Estado) %>%
get_summary_stats(Cd_dry, type = "mean_sd")
summ %>% as_tibble() %>% print(n=Inf)
## # A tibble: 6 × 6
## Vereda Estado variable n mean sd
## <fct> <fct> <fct> <dbl> <dbl> <dbl>
## 1 Caño Baul Baba Cd_dry 3 12.3 2.82
## 2 Caño Baul Fermentado Cd_dry 3 9.93 2.65
## 3 Mataredonda Baba Cd_dry 3 17.4 0.179
## 4 Mataredonda Fermentado Cd_dry 3 15.6 0.761
## 5 Mira Valle Baba Cd_dry 3 4.07 0.368
## 6 Mira Valle Fermentado Cd_dry 3 3.07 0.242
##Visualization
bxp <- ggboxplot(
datos, x = "Vereda", y = "Cd_dry",
color = "Estado", palette = "jco"
)
bxp
##Check assumptions
##Outliers
datos %>%
group_by(Vereda, Estado) %>%
identify_outliers(Cd_dry)
## [1] Vereda Estado Hora Tejido Rep
## [6] id Ensayo Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm<-datos %>%
group_by(Vereda, Estado) %>%
shapiro_test(Cd_dry)
norm %>% as_tibble() %>% print(n=Inf)
## # A tibble: 6 × 5
## Vereda Estado variable statistic p
## <fct> <fct> <chr> <dbl> <dbl>
## 1 Caño Baul Baba Cd_dry 0.987 0.780
## 2 Caño Baul Fermentado Cd_dry 0.895 0.370
## 3 Mataredonda Baba Cd_dry 0.986 0.775
## 4 Mataredonda Fermentado Cd_dry 1.00 0.967
## 5 Mira Valle Baba Cd_dry 0.940 0.527
## 6 Mira Valle Fermentado Cd_dry 0.884 0.337
##Create QQ plot for each cell of design:
ggqqplot(datos, "Cd_dry", ggtheme = theme_bw()) +
facet_grid(Estado~ Vereda, labeller = "label_both")
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev<-datos %>%
group_by(Estado) %>%
levene_test(Cd_dry ~ Vereda)
lev %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 5
## Estado df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 Baba 2 6 2.86 0.134
## 2 Fermentado 2 6 1.25 0.351
##Computation
res.aov <- anova_test(
data = datos, dv = Cd_dry, wid = id,
within = Estado, between = c(Vereda)
)
res.aov
## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Vereda 2 6 51.368 0.000168 * 0.941
## 2 Estado 1 6 38.606 0.000802 * 0.292
## 3 Vereda:Estado 2 6 2.085 0.205000 0.043
get_anova_table(res.aov)
## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Vereda 2 6 51.368 0.000168 * 0.941
## 2 Estado 1 6 38.606 0.000802 * 0.292
## 3 Vereda:Estado 2 6 2.085 0.205000 0.043
#Table by error
res.aov.error <- aov(Cd_dry ~ Estado*Vereda + Error(id/Estado), datos)
res.aov.error
##
## Call:
## aov(formula = Cd_dry ~ Estado * Vereda + Error(id/Estado), data = datos)
##
## Grand Mean: 10.39023
##
## Stratum 1: id
##
## Terms:
## Vereda Residuals
## Sum of Squares 505.3332 29.5126
## Deg. of Freedom 2 6
##
## Residual standard error: 2.21783
## 2 out of 4 effects not estimable
## Estimated effects may be unbalanced
##
## Stratum 2: id:Estado
##
## Terms:
## Estado Estado:Vereda Residuals
## Sum of Squares 13.032389 1.407580 2.025469
## Deg. of Freedom 1 2 6
##
## Residual standard error: 0.5810148
## Estimated effects may be unbalanced
summary(res.aov.error)
##
## Error: id
## Df Sum Sq Mean Sq F value Pr(>F)
## Vereda 2 505.3 252.67 51.37 0.000168 ***
## Residuals 6 29.5 4.92
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: id:Estado
## Df Sum Sq Mean Sq F value Pr(>F)
## Estado 1 13.032 13.032 38.606 0.000802 ***
## Estado:Vereda 2 1.408 0.704 2.085 0.205370
## Residuals 6 2.025 0.338
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Emmeans
emmip(res.aov.error, Vereda~Estado, CIs = T, xlab="Estado", ylab="Almendra Cd [ppm]")
## Note: re-fitting model with sum-to-zero contrasts
emm_Vereda <- emmeans(res.aov.error, pairwise ~ Vereda)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
emm_Vereda
## $emmeans
## Vereda emmean SE df lower.CL upper.CL
## Caño Baul 11.12 0.905 6 8.90 13.33
## Mataredonda 16.49 0.905 6 14.27 18.70
## Mira Valle 3.57 0.905 6 1.35 5.78
##
## Results are averaged over the levels of: Estado
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## Caño Baul - Mataredonda -5.37 1.28 6 -4.193 0.0135
## Caño Baul - Mira Valle 7.55 1.28 6 5.895 0.0026
## Mataredonda - Mira Valle 12.92 1.28 6 10.088 0.0001
##
## Results are averaged over the levels of: Estado
## P value adjustment: tukey method for comparing a family of 3 estimates
emm_gen_Estado <- emmeans(res.aov.error, pairwise ~ Estado | Vereda)
## Note: re-fitting model with sum-to-zero contrasts
emm_gen_Estado
## $emmeans
## Vereda = Caño Baul:
## Estado emmean SE df lower.CL upper.CL
## Baba 12.30 0.936 6.82 10.077 14.53
## Fermentado 9.93 0.936 6.82 7.706 12.16
##
## Vereda = Mataredonda:
## Estado emmean SE df lower.CL upper.CL
## Baba 17.35 0.936 6.82 15.127 19.58
## Fermentado 15.62 0.936 6.82 13.395 17.84
##
## Vereda = Mira Valle:
## Estado emmean SE df lower.CL upper.CL
## Baba 4.07 0.936 6.82 1.844 6.29
## Fermentado 3.07 0.936 6.82 0.842 5.29
##
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## Vereda = Caño Baul:
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado 2.37 0.474 6 4.998 0.0025
##
## Vereda = Mataredonda:
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado 1.73 0.474 6 3.651 0.0107
##
## Vereda = Mira Valle:
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado 1.00 0.474 6 2.112 0.0791
emm_Estado <- emmeans(res.aov.error, pairwise ~ Estado)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
emm_Estado
## $emmeans
## Estado emmean SE df lower.CL upper.CL
## Baba 11.24 0.54 6.82 9.96 12.5
## Fermentado 9.54 0.54 6.82 8.25 10.8
##
## Results are averaged over the levels of: Vereda
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado 1.7 0.274 6 6.213 0.0008
##
## Results are averaged over the levels of: Vereda
emm_Estado_Vereda <- emmeans(res.aov.error, pairwise ~ Vereda | Estado)
## Note: re-fitting model with sum-to-zero contrasts
emm_Estado_Vereda
## $emmeans
## Estado = Baba:
## Vereda emmean SE df lower.CL upper.CL
## Caño Baul 12.30 0.936 6.82 10.077 14.53
## Mataredonda 17.35 0.936 6.82 15.127 19.58
## Mira Valle 4.07 0.936 6.82 1.844 6.29
##
## Estado = Fermentado:
## Vereda emmean SE df lower.CL upper.CL
## Caño Baul 9.93 0.936 6.82 7.706 12.16
## Mataredonda 15.62 0.936 6.82 13.395 17.84
## Mira Valle 3.07 0.936 6.82 0.842 5.29
##
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## Estado = Baba:
## contrast estimate SE df t.ratio p.value
## Caño Baul - Mataredonda -5.05 1.32 6.82 -3.815 0.0166
## Caño Baul - Mira Valle 8.23 1.32 6.82 6.219 0.0012
## Mataredonda - Mira Valle 13.28 1.32 6.82 10.035 0.0001
##
## Estado = Fermentado:
## contrast estimate SE df t.ratio p.value
## Caño Baul - Mataredonda -5.69 1.32 6.82 -4.298 0.0092
## Caño Baul - Mira Valle 6.86 1.32 6.82 5.185 0.0034
## Mataredonda - Mira Valle 12.55 1.32 6.82 9.483 0.0001
##
## P value adjustment: tukey method for comparing a family of 3 estimates
##Splitting dataframe by Vereda
## Vereda=="Caño Baul"
datos.vereda1<-filter(datos, Vereda=="Caño Baul")
##Check assumptions
##Outliers
datos.vereda1 %>%
group_by(Estado) %>%
identify_outliers(Cd_dry)
## [1] Estado Hora Tejido Rep id
## [6] Ensayo Vereda Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm1<-datos.vereda1 %>%
group_by(Estado) %>%
shapiro_test(Cd_dry)
norm1 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 4
## Estado variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 Baba Cd_dry 0.987 0.780
## 2 Fermentado Cd_dry 0.895 0.370
##Create QQ plot for each cell of design:
ggqqplot(datos.vereda1, "Cd_dry", ggtheme = theme_bw())
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev1<-datos.vereda1 %>%
levene_test(Cd_dry~Estado)
lev1 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 1 × 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 1 4 0.0166 0.904
##Computation
res.aov1 <- anova_test(
data = datos.vereda1, dv = Cd_dry, wid = id,
within = Estado
)
get_anova_table(res.aov1)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Estado 1 2 11.366 0.078 0.22
## Vereda=="Mataredonda"
datos.vereda2<-filter(datos, Vereda=="Mataredonda")
##Check assumptions
##Outliers
datos.vereda2 %>%
group_by(Estado) %>%
identify_outliers(Cd_dry)
## [1] Estado Hora Tejido Rep id
## [6] Ensayo Vereda Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm2<-datos.vereda2 %>%
group_by(Estado) %>%
shapiro_test(Cd_dry)
norm2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 4
## Estado variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 Baba Cd_dry 0.986 0.775
## 2 Fermentado Cd_dry 1.00 0.967
##Create QQ plot for each cell of design:
ggqqplot(datos.vereda2, "Cd_dry", ggtheme = theme_bw())
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev2<-datos.vereda2 %>%
levene_test(Cd_dry ~ Estado)
lev2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 1 × 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 1 4 2.20 0.212
##Computation
res.aov2 <- anova_test(
data = datos.vereda2, dv = Cd_dry, wid = id,
within = Estado
)
get_anova_table(res.aov2)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Estado 1 2 19.935 0.047 * 0.786
## vereda==Mira Valle
datos.vereda3<-filter(datos, Vereda=="Mira Valle")
##Check assumptions
##Outliers
datos.vereda3 %>%
group_by(Estado) %>%
identify_outliers(Cd_dry)
## [1] Estado Hora Tejido Rep id
## [6] Ensayo Vereda Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm2<-datos.vereda3 %>%
group_by(Estado) %>%
shapiro_test(Cd_dry)
norm2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 4
## Estado variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 Baba Cd_dry 0.940 0.527
## 2 Fermentado Cd_dry 0.884 0.337
##Create QQ plot for each cell of design:
ggqqplot(datos.vereda3, "Cd_dry", ggtheme = theme_bw())
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev2<-datos.vereda3 %>%
levene_test(Cd_dry ~ Estado)
lev2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 1 × 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 1 4 0.214 0.668
##Computation
res.aov2 <- anova_test(
data = datos.vereda3, dv = Cd_dry, wid = id,
within = Estado
)
get_anova_table(res.aov2)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Estado 1 2 33.444 0.029 * 0.795
###Análisis Testa
##Summary statistics
summ<-datos %>%
group_by(Vereda, Estado) %>%
get_summary_stats(Cd_dry_testa, type = "mean_sd")
summ %>% as_tibble() %>% print(n=Inf)
## # A tibble: 6 × 6
## Vereda Estado variable n mean sd
## <fct> <fct> <fct> <dbl> <dbl> <dbl>
## 1 Caño Baul Baba Cd_dry_testa 3 9.54 3.06
## 2 Caño Baul Fermentado Cd_dry_testa 3 23.5 6.34
## 3 Mataredonda Baba Cd_dry_testa 3 14.9 1.63
## 4 Mataredonda Fermentado Cd_dry_testa 3 25.4 3.27
## 5 Mira Valle Baba Cd_dry_testa 3 1.67 0.353
## 6 Mira Valle Fermentado Cd_dry_testa 3 3.88 0.338
##Visualization
bxp <- ggboxplot(
datos, x = "Vereda", y = "Cd_dry_testa",
color = "Estado", palette = "jco"
)
bxp
##Check assumptions
##Outliers
datos %>%
group_by(Vereda, Estado) %>%
identify_outliers(Cd_dry_testa)
## [1] Vereda Estado Hora Tejido Rep
## [6] id Ensayo Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm<-datos %>%
group_by(Vereda, Estado) %>%
shapiro_test(Cd_dry_testa)
norm %>% as_tibble() %>% print(n=Inf)
## # A tibble: 6 × 5
## Vereda Estado variable statistic p
## <fct> <fct> <chr> <dbl> <dbl>
## 1 Caño Baul Baba Cd_dry_testa 0.791 0.0934
## 2 Caño Baul Fermentado Cd_dry_testa 0.996 0.887
## 3 Mataredonda Baba Cd_dry_testa 1.00 1.00
## 4 Mataredonda Fermentado Cd_dry_testa 0.845 0.227
## 5 Mira Valle Baba Cd_dry_testa 0.765 0.0333
## 6 Mira Valle Fermentado Cd_dry_testa 0.762 0.0261
##Create QQ plot for each cell of design:
ggqqplot(datos, "Cd_dry_testa", ggtheme = theme_bw()) +
facet_grid(Estado~ Vereda, labeller = "label_both")
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev<-datos %>%
group_by(Estado) %>%
levene_test(Cd_dry_testa ~ Vereda)
lev %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 5
## Estado df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 Baba 2 6 0.625 0.567
## 2 Fermentado 2 6 1.67 0.265
##Computation
res.aov <- anova_test(
data = datos, dv = Cd_dry_testa, wid = id,
within = Estado, between = c(Vereda)
)
res.aov
## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Vereda 2 6 30.987 0.000688 * 0.889
## 2 Estado 1 6 75.157 0.000130 * 0.738
## 3 Vereda:Estado 2 6 11.542 0.009000 * 0.464
get_anova_table(res.aov)
## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Vereda 2 6 30.987 0.000688 * 0.889
## 2 Estado 1 6 75.157 0.000130 * 0.738
## 3 Vereda:Estado 2 6 11.542 0.009000 * 0.464
#Table by error
res.aov.error <- aov(Cd_dry_testa ~ Estado*Vereda + Error(id/Estado), datos)
res.aov.error
##
## Call:
## aov(formula = Cd_dry_testa ~ Estado * Vereda + Error(id/Estado),
## data = datos)
##
## Grand Mean: 13.1604
##
## Stratum 1: id
##
## Terms:
## Vereda Residuals
## Sum of Squares 1010.7522 97.8554
## Deg. of Freedom 2 6
##
## Residual standard error: 4.038469
## 2 out of 4 effects not estimable
## Estimated effects may be unbalanced
##
## Stratum 2: id:Estado
##
## Terms:
## Estado Estado:Vereda Residuals
## Sum of Squares 355.6148 109.2273 28.3898
## Deg. of Freedom 1 2 6
##
## Residual standard error: 2.175234
## Estimated effects may be unbalanced
summary(res.aov.error)
##
## Error: id
## Df Sum Sq Mean Sq F value Pr(>F)
## Vereda 2 1010.8 505.4 30.99 0.000688 ***
## Residuals 6 97.9 16.3
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: id:Estado
## Df Sum Sq Mean Sq F value Pr(>F)
## Estado 1 355.6 355.6 75.16 0.00013 ***
## Estado:Vereda 2 109.2 54.6 11.54 0.00878 **
## Residuals 6 28.4 4.7
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Emmeans
emmip(res.aov.error, Vereda~Estado, CIs = T, xlab="Estado", ylab="Testa Cd [ppm]")
## Note: re-fitting model with sum-to-zero contrasts
emm_Vereda <- emmeans(res.aov.error, pairwise ~ Vereda)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
emm_Vereda
## $emmeans
## Vereda emmean SE df lower.CL upper.CL
## Caño Baul 16.52 1.65 6 12.49 20.56
## Mataredonda 20.18 1.65 6 16.15 24.22
## Mira Valle 2.78 1.65 6 -1.26 6.81
##
## Results are averaged over the levels of: Estado
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## Caño Baul - Mataredonda -3.66 2.33 6 -1.569 0.3278
## Caño Baul - Mira Valle 13.75 2.33 6 5.896 0.0026
## Mataredonda - Mira Valle 17.41 2.33 6 7.465 0.0007
##
## Results are averaged over the levels of: Estado
## P value adjustment: tukey method for comparing a family of 3 estimates
emm_gen_Estado <- emmeans(res.aov.error, pairwise ~ Estado | Vereda)
## Note: re-fitting model with sum-to-zero contrasts
emm_gen_Estado
## $emmeans
## Vereda = Caño Baul:
## Estado emmean SE df lower.CL upper.CL
## Baba 9.54 1.87 9.21 5.319 13.76
## Fermentado 23.51 1.87 9.21 19.286 27.73
##
## Vereda = Mataredonda:
## Estado emmean SE df lower.CL upper.CL
## Baba 14.94 1.87 9.21 10.719 19.16
## Fermentado 25.42 1.87 9.21 21.202 29.64
##
## Vereda = Mira Valle:
## Estado emmean SE df lower.CL upper.CL
## Baba 1.67 1.87 9.21 -2.555 5.89
## Fermentado 3.88 1.87 9.21 -0.337 8.11
##
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## Vereda = Caño Baul:
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado -13.97 1.78 6 -7.864 0.0002
##
## Vereda = Mataredonda:
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado -10.48 1.78 6 -5.902 0.0011
##
## Vereda = Mira Valle:
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado -2.22 1.78 6 -1.249 0.2581
emm_Estado <- emmeans(res.aov.error, pairwise ~ Estado)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
emm_Estado
## $emmeans
## Estado emmean SE df lower.CL upper.CL
## Baba 8.72 1.08 9.21 6.28 11.2
## Fermentado 17.61 1.08 9.21 15.17 20.0
##
## Results are averaged over the levels of: Vereda
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado -8.89 1.03 6 -8.669 0.0001
##
## Results are averaged over the levels of: Vereda
emm_Estado_Vereda <- emmeans(res.aov.error, pairwise ~ Vereda | Estado)
## Note: re-fitting model with sum-to-zero contrasts
emm_Estado_Vereda
## $emmeans
## Estado = Baba:
## Vereda emmean SE df lower.CL upper.CL
## Caño Baul 9.54 1.87 9.21 5.319 13.76
## Mataredonda 14.94 1.87 9.21 10.719 19.16
## Mira Valle 1.67 1.87 9.21 -2.555 5.89
##
## Estado = Fermentado:
## Vereda emmean SE df lower.CL upper.CL
## Caño Baul 23.51 1.87 9.21 19.286 27.73
## Mataredonda 25.42 1.87 9.21 21.202 29.64
## Mira Valle 3.88 1.87 9.21 -0.337 8.11
##
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## Estado = Baba:
## contrast estimate SE df t.ratio p.value
## Caño Baul - Mataredonda -5.40 2.65 9.21 -2.039 0.1573
## Caño Baul - Mira Valle 7.87 2.65 9.21 2.973 0.0369
## Mataredonda - Mira Valle 13.27 2.65 9.21 5.012 0.0018
##
## Estado = Fermentado:
## contrast estimate SE df t.ratio p.value
## Caño Baul - Mataredonda -1.92 2.65 9.21 -0.724 0.7561
## Caño Baul - Mira Valle 19.62 2.65 9.21 7.409 0.0001
## Mataredonda - Mira Valle 21.54 2.65 9.21 8.133 <.0001
##
## P value adjustment: tukey method for comparing a family of 3 estimates
##Splitting dataframe by Vereda
## Vereda=="Caño Baul"
datos.vereda1<-filter(datos, Vereda=="Caño Baul")
##Check assumptions
##Outliers
datos.vereda1 %>%
group_by(Estado) %>%
identify_outliers(Cd_dry_testa)
## [1] Estado Hora Tejido Rep id
## [6] Ensayo Vereda Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm1<-datos.vereda1 %>%
group_by(Estado) %>%
shapiro_test(Cd_dry_testa)
norm1 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 4
## Estado variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 Baba Cd_dry_testa 0.791 0.0934
## 2 Fermentado Cd_dry_testa 0.996 0.887
##Create QQ plot for each cell of design:
ggqqplot(datos.vereda1, "Cd_dry_testa", ggtheme = theme_bw())
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev1<-datos.vereda1 %>%
levene_test(Cd_dry_testa~Estado)
lev1 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 1 × 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 1 4 0.787 0.425
##Computation
res.aov1 <- anova_test(
data = datos.vereda1, dv = Cd_dry_testa, wid = id,
within = Estado
)
get_anova_table(res.aov1)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Estado 1 2 36.148 0.027 * 0.747
## Vereda=="Mataredonda"
datos.vereda2<-filter(datos, Vereda=="Mataredonda")
##Check assumptions
##Outliers
datos.vereda2 %>%
group_by(Estado) %>%
identify_outliers(Cd_dry_testa)
## [1] Estado Hora Tejido Rep id
## [6] Ensayo Vereda Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm2<-datos.vereda2 %>%
group_by(Estado) %>%
shapiro_test(Cd_dry_testa)
norm2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 4
## Estado variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 Baba Cd_dry_testa 1.00 1.00
## 2 Fermentado Cd_dry_testa 0.845 0.227
##Create QQ plot for each cell of design:
ggqqplot(datos.vereda2, "Cd_dry_testa", ggtheme = theme_bw())
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev2<-datos.vereda2 %>%
levene_test(Cd_dry_testa ~ Estado)
lev2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 1 × 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 1 4 0.286 0.621
##Computation
res.aov2 <- anova_test(
data = datos.vereda2, dv = Cd_dry_testa, wid = id,
within = Estado
)
get_anova_table(res.aov2)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Estado 1 2 27.293 0.035 * 0.861
## vereda==Mira Valle
datos.vereda3<-filter(datos, Vereda=="Mira Valle")
##Check assumptions
##Outliers
datos.vereda3 %>%
group_by(Estado) %>%
identify_outliers(Cd_dry_testa)
## [1] Estado Hora Tejido Rep id
## [6] Ensayo Vereda Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm2<-datos.vereda3 %>%
group_by(Estado) %>%
shapiro_test(Cd_dry_testa)
norm2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 4
## Estado variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 Baba Cd_dry_testa 0.765 0.0333
## 2 Fermentado Cd_dry_testa 0.762 0.0261
##Create QQ plot for each cell of design:
ggqqplot(datos.vereda3, "Cd_dry_testa", ggtheme = theme_bw())
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev2<-datos.vereda3 %>%
levene_test(Cd_dry_testa ~ Estado)
lev2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 1 × 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 1 4 0.00112 0.975
##Computation
res.aov2 <- anova_test(
data = datos.vereda3, dv = Cd_dry_testa, wid = id,
within = Estado
)
get_anova_table(res.aov2)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Estado 1 2 122.98 0.008 * 0.939
###Análisis Testa
##Summary statistics
summ<-datos %>%
group_by(Vereda, Estado) %>%
get_summary_stats(ITF, type = "mean_sd")
summ %>% as_tibble() %>% print(n=Inf)
## # A tibble: 6 × 6
## Vereda Estado variable n mean sd
## <fct> <fct> <fct> <dbl> <dbl> <dbl>
## 1 Caño Baul Baba ITF 3 0.765 0.117
## 2 Caño Baul Fermentado ITF 3 2.37 0.196
## 3 Mataredonda Baba ITF 3 0.86 0.085
## 4 Mataredonda Fermentado ITF 3 1.62 0.136
## 5 Mira Valle Baba ITF 3 0.417 0.13
## 6 Mira Valle Fermentado ITF 3 1.27 0.092
##Visualization
bxp <- ggboxplot(
datos, x = "Vereda", y = "ITF",
color = "Estado", palette = "jco"
)
bxp
##Check assumptions
##Outliers
datos %>%
group_by(Vereda, Estado) %>%
identify_outliers(ITF)
## [1] Vereda Estado Hora Tejido Rep
## [6] id Ensayo Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm<-datos %>%
group_by(Vereda, Estado) %>%
shapiro_test(ITF)
norm %>% as_tibble() %>% print(n=Inf)
## # A tibble: 6 × 5
## Vereda Estado variable statistic p
## <fct> <fct> <chr> <dbl> <dbl>
## 1 Caño Baul Baba ITF 0.999 0.956
## 2 Caño Baul Fermentado ITF 0.965 0.641
## 3 Mataredonda Baba ITF 1.00 0.988
## 4 Mataredonda Fermentado ITF 0.823 0.170
## 5 Mira Valle Baba ITF 0.810 0.140
## 6 Mira Valle Fermentado ITF 0.790 0.0907
##Create QQ plot for each cell of design:
ggqqplot(datos, "ITF", ggtheme = theme_bw()) +
facet_grid(Estado~ Vereda, labeller = "label_both")
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev<-datos %>%
group_by(Estado) %>%
levene_test(ITF ~ Vereda)
lev %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 5
## Estado df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 Baba 2 6 0.0641 0.939
## 2 Fermentado 2 6 0.321 0.737
##Computation
res.aov <- anova_test(
data = datos, dv = ITF, wid = id,
within = Estado, between = c(Vereda)
)
res.aov
## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Vereda 2 6 27.921 9.13e-04 * 0.884
## 2 Estado 1 6 846.885 1.09e-07 * 0.962
## 3 Vereda:Estado 2 6 52.210 1.60e-04 * 0.756
get_anova_table(res.aov)
## ANOVA Table (type II tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Vereda 2 6 27.921 9.13e-04 * 0.884
## 2 Estado 1 6 846.885 1.09e-07 * 0.962
## 3 Vereda:Estado 2 6 52.210 1.60e-04 * 0.756
#Table by error
res.aov.error <- aov(ITF ~ Estado*Vereda + Error(id/Estado), datos)
res.aov.error
##
## Call:
## aov(formula = ITF ~ Estado * Vereda + Error(id/Estado), data = datos)
##
## Grand Mean: 1.216818
##
## Stratum 1: id
##
## Terms:
## Vereda Residuals
## Sum of Squares 1.5747178 0.1691963
## Deg. of Freedom 2 6
##
## Residual standard error: 0.1679267
## 2 out of 4 effects not estimable
## Estimated effects may be unbalanced
##
## Stratum 2: id:Estado
##
## Terms:
## Estado Estado:Vereda Residuals
## Sum of Squares 5.173956 0.637941 0.036656
## Deg. of Freedom 1 2 6
##
## Residual standard error: 0.07816264
## Estimated effects may be unbalanced
summary(res.aov.error)
##
## Error: id
## Df Sum Sq Mean Sq F value Pr(>F)
## Vereda 2 1.5747 0.7874 27.92 0.000913 ***
## Residuals 6 0.1692 0.0282
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: id:Estado
## Df Sum Sq Mean Sq F value Pr(>F)
## Estado 1 5.174 5.174 846.88 1.09e-07 ***
## Estado:Vereda 2 0.638 0.319 52.21 0.00016 ***
## Residuals 6 0.037 0.006
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Emmeans
emmip(res.aov.error, Vereda~Estado, CIs = T, xlab="Estado", ylab="ITF")
## Note: re-fitting model with sum-to-zero contrasts
emm_Vereda <- emmeans(res.aov.error, pairwise ~ Vereda)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
emm_Vereda
## $emmeans
## Vereda emmean SE df lower.CL upper.CL
## Caño Baul 1.566 0.0686 6 1.398 1.73
## Mataredonda 1.242 0.0686 6 1.074 1.41
## Mira Valle 0.843 0.0686 6 0.675 1.01
##
## Results are averaged over the levels of: Estado
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## Caño Baul - Mataredonda 0.323 0.097 6 3.337 0.0361
## Caño Baul - Mira Valle 0.723 0.097 6 7.459 0.0007
## Mataredonda - Mira Valle 0.400 0.097 6 4.122 0.0146
##
## Results are averaged over the levels of: Estado
## P value adjustment: tukey method for comparing a family of 3 estimates
emm_gen_Estado <- emmeans(res.aov.error, pairwise ~ Estado | Vereda)
## Note: re-fitting model with sum-to-zero contrasts
emm_gen_Estado
## $emmeans
## Vereda = Caño Baul:
## Estado emmean SE df lower.CL upper.CL
## Baba 0.765 0.0756 8.48 0.592 0.937
## Fermentado 2.367 0.0756 8.48 2.194 2.540
##
## Vereda = Mataredonda:
## Estado emmean SE df lower.CL upper.CL
## Baba 0.860 0.0756 8.48 0.688 1.033
## Fermentado 1.624 0.0756 8.48 1.451 1.797
##
## Vereda = Mira Valle:
## Estado emmean SE df lower.CL upper.CL
## Baba 0.417 0.0756 8.48 0.244 0.590
## Fermentado 1.268 0.0756 8.48 1.095 1.441
##
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## Vereda = Caño Baul:
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado -1.602 0.0638 6 -25.108 <.0001
##
## Vereda = Mataredonda:
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado -0.764 0.0638 6 -11.964 <.0001
##
## Vereda = Mira Valle:
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado -0.851 0.0638 6 -13.333 <.0001
emm_Estado <- emmeans(res.aov.error, pairwise ~ Estado)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
emm_Estado
## $emmeans
## Estado emmean SE df lower.CL upper.CL
## Baba 0.681 0.0437 8.48 0.581 0.78
## Fermentado 1.753 0.0437 8.48 1.653 1.85
##
## Results are averaged over the levels of: Vereda
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## Baba - Fermentado -1.07 0.0368 6 -29.101 <.0001
##
## Results are averaged over the levels of: Vereda
emm_Estado_Vereda <- emmeans(res.aov.error, pairwise ~ Vereda | Estado)
## Note: re-fitting model with sum-to-zero contrasts
emm_Estado_Vereda
## $emmeans
## Estado = Baba:
## Vereda emmean SE df lower.CL upper.CL
## Caño Baul 0.765 0.0756 8.48 0.592 0.937
## Mataredonda 0.860 0.0756 8.48 0.688 1.033
## Mira Valle 0.417 0.0756 8.48 0.244 0.590
##
## Estado = Fermentado:
## Vereda emmean SE df lower.CL upper.CL
## Caño Baul 2.367 0.0756 8.48 2.194 2.540
## Mataredonda 1.624 0.0756 8.48 1.451 1.797
## Mira Valle 1.268 0.0756 8.48 1.095 1.441
##
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
##
## $contrasts
## Estado = Baba:
## contrast estimate SE df t.ratio p.value
## Caño Baul - Mataredonda -0.0959 0.107 8.48 -0.897 0.6564
## Caño Baul - Mira Valle 0.3474 0.107 8.48 3.249 0.0263
## Mataredonda - Mira Valle 0.4433 0.107 8.48 4.146 0.0071
##
## Estado = Fermentado:
## contrast estimate SE df t.ratio p.value
## Caño Baul - Mataredonda 0.7429 0.107 8.48 6.947 0.0002
## Caño Baul - Mira Valle 1.0989 0.107 8.48 10.276 <.0001
## Mataredonda - Mira Valle 0.3560 0.107 8.48 3.329 0.0234
##
## P value adjustment: tukey method for comparing a family of 3 estimates
##Splitting dataframe by Vereda
## Vereda=="Caño Baul"
datos.vereda1<-filter(datos, Vereda=="Caño Baul")
##Check assumptions
##Outliers
datos.vereda1 %>%
group_by(Estado) %>%
identify_outliers(ITF)
## [1] Estado Hora Tejido Rep id
## [6] Ensayo Vereda Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm1<-datos.vereda1 %>%
group_by(Estado) %>%
shapiro_test(ITF)
norm1 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 4
## Estado variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 Baba ITF 0.999 0.956
## 2 Fermentado ITF 0.965 0.641
##Create QQ plot for each cell of design:
ggqqplot(datos.vereda1, "ITF", ggtheme = theme_bw())
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev1<-datos.vereda1 %>%
levene_test(ITF~Estado)
lev1 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 1 × 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 1 4 0.362 0.580
##Computation
res.aov1 <- anova_test(
data = datos.vereda1, dv = ITF, wid = id,
within = Estado
)
get_anova_table(res.aov1)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Estado 1 2 1064.91 0.000938 * 0.974
## Vereda=="Mataredonda"
datos.vereda2<-filter(datos, Vereda=="Mataredonda")
##Check assumptions
##Outliers
datos.vereda2 %>%
group_by(Estado) %>%
identify_outliers(ITF)
## [1] Estado Hora Tejido Rep id
## [6] Ensayo Vereda Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm2<-datos.vereda2 %>%
group_by(Estado) %>%
shapiro_test(ITF)
norm2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 4
## Estado variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 Baba ITF 1.00 0.988
## 2 Fermentado ITF 0.823 0.170
##Create QQ plot for each cell of design:
ggqqplot(datos.vereda2, "ITF", ggtheme = theme_bw())
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev2<-datos.vereda2 %>%
levene_test(ITF ~ Estado)
lev2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 1 × 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 1 4 0.112 0.754
##Computation
res.aov2 <- anova_test(
data = datos.vereda2, dv = ITF, wid = id,
within = Estado
)
get_anova_table(res.aov2)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Estado 1 2 62.835 0.016 * 0.945
## vereda==Mira Valle
datos.vereda3<-filter(datos, Vereda=="Mira Valle")
##Check assumptions
##Outliers
datos.vereda3 %>%
group_by(Estado) %>%
identify_outliers(ITF)
## [1] Estado Hora Tejido Rep id
## [6] Ensayo Vereda Finca Conc_grano Matriz
## [11] Matter_dry Cd_wet Cd_dry id2 Cd_dry_testa
## [16] Estado2 ITF is.outlier is.extreme
## <0 rows> (o 0- extensión row.names)
##Normality assumption
##Compute Shapiro-Wilk test for each combinations of factor levels:
norm2<-datos.vereda3 %>%
group_by(Estado) %>%
shapiro_test(ITF)
norm2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 2 × 4
## Estado variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 Baba ITF 0.810 0.140
## 2 Fermentado ITF 0.790 0.0907
##Create QQ plot for each cell of design:
ggqqplot(datos.vereda3, "ITF", ggtheme = theme_bw())
##Homogneity of variance assumption
##Compute the Levene’s test at each level of the within-subjects factor, here time variable:
lev2<-datos.vereda3 %>%
levene_test(ITF ~ Estado)
lev2 %>% as_tibble() %>% print(n=Inf)
## # A tibble: 1 × 4
## df1 df2 statistic p
## <int> <int> <dbl> <dbl>
## 1 1 4 0.0749 0.798
##Computation
res.aov2 <- anova_test(
data = datos.vereda3, dv = ITF, wid = id,
within = Estado
)
get_anova_table(res.aov2)
## ANOVA Table (type III tests)
##
## Effect DFn DFd F p p<.05 ges
## 1 Estado 1 2 1368.573 0.00073 * 0.955