The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.
rstatix provides a simple, pipe-friendly framework for common statistical tests, coherent with the tidyverse design philosophy. Each function takes a data frame and a formula, and returns the result as a tidy tibble — so a test flows straight into dplyr pipelines, into reporting tables (gtsummary, gt, broom), and into ggpubr plots, without reshaping.
This vignette walks one analysis from start to finish: describe the data, check the assumptions, run the test, quantify the effect with a confidence interval, and report it. The same handful of verbs covers t-tests, Wilcoxon tests, ANOVA, Kruskal-Wallis, correlation, and their effect sizes.
We use ToothGrowth, a base R dataset: the length of odontoblasts (len) in 60 guinea pigs, each given one of three vitamin C doses (dose) by one of two delivery methods (supp).
df <- ToothGrowth
df$dose <- factor(df$dose)
head(df)
#> len supp dose
#> 1 4.2 VC 0.5
#> 2 11.5 VC 0.5
#> 3 7.3 VC 0.5
#> 4 5.8 VC 0.5
#> 5 6.4 VC 0.5
#> 6 10.0 VC 0.5get_summary_stats() returns the summaries you actually report, for one or many variables, grouped or not.
Before comparing the three doses, check_test_assumptions() tests normality (Shapiro-Wilk, per group) and homogeneity of variance (Levene), and returns a one-row recommendation: which omnibus test and which post-hoc test the data call for.
df %>% check_test_assumptions(len ~ dose)
#> # A tibble: 1 × 8
#> .y. normality.p homogeneity.p normal equal.variance significance omnibus
#> <chr> <dbl> <dbl> <lgl> <lgl> <dbl> <chr>
#> 1 len 0.164 0.528 TRUE TRUE 0.05 anova_test
#> # ℹ 1 more variable: posthoc <chr>Here both assumptions hold, so the recommendation is a one-way ANOVA followed by Tukey HSD. When normality fails it points to Kruskal-Wallis with Dunn’s test; when only the variances differ, to Welch ANOVA with Games-Howell. (Choosing a test by first testing its assumptions on the same data is a known trade-off, noted in ?check_test_assumptions — pre-specifying the test is always defensible.)
Follow the recommendation. anova_test() returns the ANOVA table; asking for the partial eta-squared effect size with a confidence interval adds pes and its bounds.
df %>% anova_test(len ~ dose, effect.size = "pes", ci = 0.95)
#> ANOVA Table (type II tests)
#>
#> Effect DFn DFd F p p<.05 pes conf.low conf.high
#> 1 dose 2 57 67.4 9.53e-16 * 0.703 0.567 0.785Dose has a large, highly significant effect on tooth length.
posthoc_test() runs the post-hoc appropriate to the same decision tree and shows which test it picked and why, so the routing is transparent.
df %>% posthoc_test(len ~ dose)
#> Post-hoc test chosen: Tukey HSD
#> Normality (Shapiro-Wilk, min across groups): p = 0.164 -> normal
#> Homogeneity of variance (Levene): p = 0.528 -> equal variances
#>
#> # A tibble: 3 × 9
#> term group1 group2 null.value estimate conf.low conf.high p.adj
#> * <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 dose 0.5 1 0 9.13 5.90 12.4 2.00e- 8
#> 2 dose 0.5 2 0 15.5 12.3 18.7 1.12e-11
#> 3 dose 1 2 0 6.37 3.14 9.59 4.25e- 5
#> # ℹ 1 more variable: p.adj.signif <chr>Every pairwise difference is significant after adjustment. If you prefer to drive the pairwise step yourself, the individual verbs compose the same way — tukey_hsd(), or dunn_test() / games_howell_test() — and any pairwise test can be piped through adjust_pvalue() and add_significance().
A p-value tells you an effect exists; an effect size tells you how large. For a pairwise mean comparison that is Cohen’s d:
df %>% cohens_d(len ~ dose)
#> # A tibble: 3 × 7
#> .y. group1 group2 effsize n1 n2 magnitude
#> * <chr> <chr> <chr> <dbl> <int> <int> <ord>
#> 1 len 0.5 1 -2.05 20 20 large
#> 2 len 0.5 2 -3.73 20 20 large
#> 3 len 1 2 -1.55 20 20 largeEach metric has a matching verb: eta_squared() for ANOVA (on the fitted model), cramer_v() for association, cliff_delta() and wilcox_effsize() for the non-parametric case, kruskal_effsize() for Kruskal-Wallis. Each can also return a confidence interval — ci = TRUE for the formula-based verbs, a confidence level such as ci = 0.95 for eta_squared().
For a manuscript, get_test_label(style = "apa") writes the APA-7 in-text sentence, including the effect-size confidence interval when the result carries one:
model <- df %>% anova_test(len ~ dose, effect.size = "pes", ci = 0.95)
get_test_label(model, type = "text", style = "apa")
#> [1] "F(2, 57) = 67.42, p < .001, eta2[p] = .70, 95% CI [.57, .79]"The test results are rstatix_test objects with tidy() and glance() methods, which hand them as plain tibbles to the reporting-table ecosystem (broom, gtsummary, gt):
The workflow is group_by()-aware: run the same test within every level of a grouping variable, then adjust and flag, in a single pipe. Here we compare the two delivery methods within each dose.
df %>%
group_by(dose) %>%
t_test(len ~ supp) %>%
adjust_pvalue(method = "holm") %>%
add_significance()
#> # A tibble: 3 × 11
#> dose .y. group1 group2 n1 n2 statistic df p p.adj
#> <fct> <chr> <chr> <chr> <int> <int> <dbl> <dbl> <dbl> <dbl>
#> 1 0.5 len OJ VC 10 10 3.17 15.0 0.00636 0.0127
#> 2 1 len OJ VC 10 10 4.03 15.4 0.00104 0.00312
#> 3 2 len OJ VC 10 10 -0.0461 14.0 0.964 0.964
#> # ℹ 1 more variable: p.adj.signif <chr>The delivery method matters at the two lower doses but not at the highest.
rstatix returns the data; ggpubr draws it. The tidy result of any test pipes straight into a plot: add_xy_position() places the significance brackets and stat_pvalue_manual() adds them, so the p-values on the figure are the ones you computed. See the ggpubr documentation.
Extended tutorials — comparing means, ANOVA designs, correlation analysis, effect sizes — are on the package website: https://rpkgs.datanovia.com/rstatix/.
sessionInfo()
#> R version 4.5.1 (2025-06-13)
#> Platform: aarch64-apple-darwin20
#> Running under: macOS Sequoia 15.6.1
#>
#> Matrix products: default
#> BLAS: /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRblas.0.dylib
#> LAPACK: /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.12.1
#>
#> locale:
#> [1] C/fr_FR.UTF-8/fr_FR.UTF-8/C/fr_FR.UTF-8/fr_FR.UTF-8
#>
#> time zone: Europe/Paris
#> tzcode source: internal
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] dplyr_1.2.1 rstatix_1.1.0
#>
#> loaded via a namespace (and not attached):
#> [1] vctrs_0.7.2 cli_3.6.6 knitr_1.51 rlang_1.2.0
#> [5] xfun_0.57 Formula_1.2-5 otel_0.2.0 purrr_1.2.2
#> [9] car_3.1-3 generics_0.1.4 jsonlite_2.0.0 glue_1.8.0
#> [13] backports_1.5.0 htmltools_0.5.9 sass_0.4.10 rmarkdown_2.31
#> [17] abind_1.4-8 carData_3.0-5 evaluate_1.0.5 jquerylib_0.1.4
#> [21] tibble_3.3.1 fastmap_1.2.0 yaml_2.3.12 lifecycle_1.0.5
#> [25] compiler_4.5.1 pkgconfig_2.0.3 tidyr_1.3.2 digest_0.6.39
#> [29] R6_2.6.1 utf8_1.2.6 tidyselect_1.2.1 pillar_1.11.1
#> [33] magrittr_2.0.4 bslib_0.10.0 withr_3.0.2 tools_4.5.1
#> [37] broom_1.0.13 cachem_1.1.0These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.