1 Load packages
library(tidyverse) # data wrangling
2 Motivation
Let’s try to see how much the results of Jamovi (2.2.5) and rstanarm (2.21.1) converge. It’s probably difficult to say because the defaults are different, and it may not be straight forward to translate back and forth.
Anyhow, let’s see.
3 data
The ubiqious mtcars
.
data(mtcars)
write_csv(mtcars, file = "mtcars.csv")
4 Model 1
Let’s try: mpg ~ am
.
4.1 rstanarm
library(rstanarm)
m1 <- stan_glm(mpg ~ am, data = mtcars, refresh = 0)
posterior_interval(m1)
#> 5% 95%
#> (Intercept) 15.166342 19.115190
#> am 4.331449 10.211882
#> sigma 4.080859 6.218033
I did not touch the defaults, which can be accessed like this:
prior_summary(m1)
#> Priors for model 'm1'
#> ------
#> Intercept (after predictors centered)
#> Specified prior:
#> ~ normal(location = 20, scale = 2.5)
#> Adjusted prior:
#> ~ normal(location = 20, scale = 15)
#>
#> Coefficients
#> Specified prior:
#> ~ normal(location = 0, scale = 2.5)
#> Adjusted prior:
#> ~ normal(location = 0, scale = 30)
#>
#> Auxiliary (sigma)
#> Specified prior:
#> ~ exponential(rate = 1)
#> Adjusted prior:
#> ~ exponential(rate = 0.17)
#> ------
#> See help('prior_summary.stanreg') for more details
4.2 Jamovi
I did not touch the defaults, which is a JZS prior with scale r = 0.354.
Jamovi report a 95% CI of [3.07, 9.88] which is quite close to what rstanarm found.
5 Model 2
mpg ~ am + hp
5.1 rstanarm
m2 <- stan_glm(mpg ~ am + hp, data = mtcars, refresh = 0)
posterior_interval(m2)
#> 5% 95%
#> (Intercept) 24.04358663 29.18365031
#> am 3.39356991 7.20358097
#> hp -0.07265686 -0.04540518
#> sigma 2.41869351 3.74714944
5.2 Jamovi
Jamovi provides R code for Bayes regresison, but it not yet kind of released.
jsq::blinReg(
data = data,
dep = mpg,
covs = vars(am, hp),
postSumm = TRUE,
modelTerms = list(
list(
components="am",
isNuisance=FALSE),
list(
components="hp",
isNuisance=FALSE)))
The source code can be accessed here though.
Anyhow, it produces (in Jamovi) this output:
#
# Model Comparison
# ───────────────────────────────────────────────────────────────────────────
# Models P(M) P(M|data) BF-M BF-10 R²
# ───────────────────────────────────────────────────────────────────────────
# Null model 0.2500 2.919e-8 8.756e-8 1.000 0.0000
# am + hp 0.2500 0.998333 1796.888017 3.421e+7 0.7820
# hp 0.2500 0.001664 0.005001 57020.790 0.6024
# am 0.2500 2.549e-6 7.648e-6 87.348 0.3598
# ───────────────────────────────────────────────────────────────────────────
#
#
#
#
# Posterior Summaries of Coefficients
# ──────────────────────────────────────────────────────────────────────────────────────────────────────────
# Coefficient Mean SD P(incl) P(incl|data) BF-inclusion Lower Upper
# ──────────────────────────────────────────────────────────────────────────────────────────────────────────
# Intercept 20.09062 0.514278 1.0000 1.0000 1.000 19.04175 21.13950
# am 5.11035 1.062349 0.5000 0.9983 599.882 2.94368 7.27703
# hp -0.05703 0.007732 0.5000 1.0000 387822.393 -0.07280 -0.04126
# ──────────────────────────────────────────────────────────────────────────────────────────────────────────
#
#
#
Which is in turn quite comparable to rstanarm.
5.3 Interim conclusion
This little exercise shows that the two software packages may yield similar results, at least in some cases.
6 Reproducibility
#> ─ Session info ───────────────────────────────────────────────────────────────────────────────────────────────────────
#> setting value
#> version R version 4.1.3 (2022-03-10)
#> os macOS Big Sur/Monterey 10.16
#> system x86_64, darwin17.0
#> ui X11
#> language (EN)
#> collate en_US.UTF-8
#> ctype en_US.UTF-8
#> tz Europe/Berlin
#> date 2022-05-09
#>
#> ─ Packages ───────────────────────────────────────────────────────────────────────────────────────────────────────────
#> package * version date lib source
#> assertthat 0.2.1 2019-03-21 [1] CRAN (R 4.1.0)
#> backports 1.4.1 2021-12-13 [1] CRAN (R 4.1.0)
#> blogdown 1.8 2022-02-16 [2] CRAN (R 4.1.2)
#> bookdown 0.26.2 2022-05-02 [1] Github (rstudio/bookdown@6adacc3)
#> brio 1.1.3 2021-11-30 [1] CRAN (R 4.1.0)
#> broom 0.8.0 2022-04-13 [1] CRAN (R 4.1.2)
#> bslib 0.3.1 2021-10-06 [1] CRAN (R 4.1.0)
#> cachem 1.0.6 2021-08-19 [1] CRAN (R 4.1.0)
#> callr 3.7.0 2021-04-20 [1] CRAN (R 4.1.0)
#> cellranger 1.1.0 2016-07-27 [1] CRAN (R 4.1.0)
#> cli 3.3.0 2022-04-25 [1] CRAN (R 4.1.2)
#> codetools 0.2-18 2020-11-04 [2] CRAN (R 4.1.3)
#> colorout * 1.2-2 2022-01-04 [1] Github (jalvesaq/colorout@79931fd)
#> colorspace 2.0-3 2022-02-21 [1] CRAN (R 4.1.2)
#> crayon 1.5.1 2022-03-26 [1] CRAN (R 4.1.2)
#> DBI 1.1.2 2021-12-20 [1] CRAN (R 4.1.0)
#> dbplyr 2.1.1 2021-04-06 [1] CRAN (R 4.1.0)
#> desc 1.4.1 2022-03-06 [1] CRAN (R 4.1.2)
#> devtools 2.4.3 2021-11-30 [1] CRAN (R 4.1.0)
#> digest 0.6.29 2021-12-01 [1] CRAN (R 4.1.0)
#> dplyr * 1.0.9 2022-04-28 [1] CRAN (R 4.1.2)
#> ellipsis 0.3.2 2021-04-29 [1] CRAN (R 4.1.0)
#> evaluate 0.15 2022-02-18 [1] CRAN (R 4.1.2)
#> fansi 1.0.3 2022-03-24 [1] CRAN (R 4.1.2)
#> fastmap 1.1.0 2021-01-25 [2] CRAN (R 4.1.0)
#> forcats * 0.5.1 2021-01-27 [1] CRAN (R 4.1.0)
#> fs 1.5.2 2021-12-08 [1] CRAN (R 4.1.0)
#> generics 0.1.2 2022-01-31 [1] CRAN (R 4.1.2)
#> ggplot2 * 3.3.6 2022-05-03 [1] CRAN (R 4.1.2)
#> glue 1.6.2 2022-02-24 [1] CRAN (R 4.1.2)
#> gtable 0.3.0 2019-03-25 [1] CRAN (R 4.1.0)
#> haven 2.4.3 2021-08-04 [1] CRAN (R 4.1.0)
#> hms 1.1.1 2021-09-26 [1] CRAN (R 4.1.0)
#> htmltools 0.5.2 2021-08-25 [1] CRAN (R 4.1.0)
#> httr 1.4.3 2022-05-04 [1] CRAN (R 4.1.2)
#> jquerylib 0.1.4 2021-04-26 [1] CRAN (R 4.1.0)
#> jsonlite 1.8.0 2022-02-22 [1] CRAN (R 4.1.2)
#> knitr 1.39 2022-04-26 [1] CRAN (R 4.1.2)
#> lifecycle 1.0.1 2021-09-24 [1] CRAN (R 4.1.0)
#> lubridate 1.8.0 2021-10-07 [1] CRAN (R 4.1.0)
#> magrittr 2.0.3 2022-03-30 [1] CRAN (R 4.1.2)
#> memoise 2.0.0 2021-01-26 [2] CRAN (R 4.1.0)
#> modelr 0.1.8 2020-05-19 [1] CRAN (R 4.1.0)
#> munsell 0.5.0 2018-06-12 [1] CRAN (R 4.1.0)
#> pillar 1.7.0 2022-02-01 [1] CRAN (R 4.1.2)
#> pkgbuild 1.2.0 2020-12-15 [2] CRAN (R 4.1.0)
#> pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 4.1.0)
#> pkgload 1.2.4 2021-11-30 [1] CRAN (R 4.1.0)
#> prettyunits 1.1.1 2020-01-24 [1] CRAN (R 4.1.0)
#> processx 3.5.3 2022-03-25 [1] CRAN (R 4.1.2)
#> ps 1.7.0 2022-04-23 [1] CRAN (R 4.1.2)
#> purrr * 0.3.4 2020-04-17 [1] CRAN (R 4.1.0)
#> R6 2.5.1 2021-08-19 [1] CRAN (R 4.1.0)
#> Rcpp 1.0.8.3 2022-03-17 [1] CRAN (R 4.1.2)
#> readr * 2.1.2 2022-01-30 [1] CRAN (R 4.1.2)
#> readxl 1.3.1 2019-03-13 [1] CRAN (R 4.1.0)
#> remotes 2.4.0 2021-06-02 [2] CRAN (R 4.1.0)
#> reprex 2.0.1 2021-08-05 [1] CRAN (R 4.1.0)
#> rlang 1.0.2 2022-03-04 [1] CRAN (R 4.1.2)
#> rmarkdown 2.14 2022-04-25 [1] CRAN (R 4.1.2)
#> rprojroot 2.0.3 2022-04-02 [1] CRAN (R 4.1.2)
#> rstudioapi 0.13 2020-11-12 [1] CRAN (R 4.1.0)
#> rvest 1.0.2 2021-10-16 [1] CRAN (R 4.1.0)
#> sass 0.4.1 2022-03-23 [1] CRAN (R 4.1.2)
#> scales 1.2.0 2022-04-13 [1] CRAN (R 4.1.3)
#> sessioninfo 1.1.1 2018-11-05 [2] CRAN (R 4.1.0)
#> stringi 1.7.6 2021-11-29 [1] CRAN (R 4.1.0)
#> stringr * 1.4.0 2019-02-10 [1] CRAN (R 4.1.0)
#> testthat 3.1.4 2022-04-26 [1] CRAN (R 4.1.2)
#> tibble * 3.1.7 2022-05-03 [1] CRAN (R 4.1.2)
#> tidyr * 1.2.0 2022-02-01 [1] CRAN (R 4.1.2)
#> tidyselect 1.1.2 2022-02-21 [1] CRAN (R 4.1.2)
#> tidyverse * 1.3.1 2021-04-15 [1] CRAN (R 4.1.0)
#> tzdb 0.1.2 2021-07-20 [2] CRAN (R 4.1.0)
#> usethis 2.0.1 2021-02-10 [2] CRAN (R 4.1.0)
#> utf8 1.2.2 2021-07-24 [1] CRAN (R 4.1.0)
#> vctrs 0.4.1 2022-04-13 [1] CRAN (R 4.1.2)
#> withr 2.5.0 2022-03-03 [1] CRAN (R 4.1.2)
#> xfun 0.30 2022-03-02 [1] CRAN (R 4.1.2)
#> xml2 1.3.3 2021-11-30 [1] CRAN (R 4.1.0)
#> yaml 2.3.5 2022-02-21 [1] CRAN (R 4.1.2)
#>
#> [1] /Users/sebastiansaueruser/Library/R/x86_64/4.1/library
#> [2] /Library/Frameworks/R.framework/Versions/4.1/Resources/library