bayestestR icon indicating copy to clipboard operation
bayestestR copied to clipboard

Robust handling of priors (get_priors and describe_prior)

Open DominiqueMakowski opened this issue 3 years ago • 8 comments

Currently, we report the location, scale and df of priors. But as we know it doesn't well apply for all cases (for instance the shift parameter in shifted lognormal models etc.).

The problem is that all distribution types have different parameters, so it's quite tough to put them all in the same dataframe with a concise presentation... but we should probably try to rework that a bit to avoid strange outputs...

I'm not sure what would be the most logical implementation though...

Uniform priors example with brms

library(brms)

model <- brm(mpg ~ wt, data=mtcars, prior = set_prior("uniform(-10,10)", lb = -10, ub = 10), refresh= 0)

insight::get_priors(model)
#>     Parameter Distribution df        Location Scale
#> 1 (Intercept)    student_t  3            19.2   5.4
#> 2          wt      uniform NA uniform(-10,10)  10.0
bayestestR::describe_prior(model)
#>     Parameter Prior_Distribution Prior_df  Prior_Location Prior_Scale
#> 1 (Intercept)          student_t        3            19.2         5.4
#> 2          wt            uniform       NA uniform(-10,10)        10.0
bayestestR::describe_posterior(model, priors = TRUE)  # Messed up output
#> Warning: Note that the default rope range for linear models might change in
#> future versions (see https://github.com/easystats/bayestestR/issues/364).Please
#> set it explicitly to preserve current results.
#> Summary of Posterior Distribution
#> 
#> Parameter   | ... |                              Prior
#> --------------------------------------------------------------------------------------
#> (Intercept) |  ... |   Student_t (df=3) ( 19.2 +- 5.40)
#> wt          |  ... | Uniform (uniform(-10,10) +- 10.00)

Created on 2021-04-02 by the reprex package (v1.0.0)

Uniform priors example with rstanarm

library(rstanarm)

model <- stan_glm(mpg ~ wt, data=mtcars, prior = NULL, refresh= 0)

insight::get_priors(model)
#>     Parameter Distribution Location Scale Adjusted_Scale
#> 1 (Intercept)       normal 20.09062   2.5       15.06737
#> 2          wt      uniform       NA    NA             NA
bayestestR::describe_prior(model)
#>     Parameter Prior_Distribution Prior_Location Prior_Scale
#> 1 (Intercept)             normal       20.09062    15.06737
#> 2          wt            uniform             NA          NA
bayestestR::describe_posterior(model, priors = TRUE)  # Messed up output
#> Warning: Note that the default rope range for linear models might change in
#> future versions (see https://github.com/easystats/bayestestR/issues/364).Please
#> set it explicitly to preserve current results.
#> Summary of Posterior Distribution
#> 
#> Parameter   | ... |                   Prior
#> ---------------------------------------------------
#> (Intercept) |  ... | Normal (20.09 +- 15.07)
#> wt          |  ... |                 Uniform

Created on 2021-04-02 by the reprex package (v1.0.0)

DominiqueMakowski avatar Apr 02 '21 02:04 DominiqueMakowski

library(brms)
model <- brm(mpg ~ wt, data=mtcars, prior = set_prior("uniform(-10,10)", lb = -10, ub = 10), refresh= 0)

insight::get_priors(model)
#>     Parameter Distribution df Location Scale
#> 1 (Intercept)    student_t  3     19.2   5.4
#> 2          wt      uniform NA       NA    NA
bayestestR::describe_prior(model)
#>     Parameter Prior_Distribution Prior_df Prior_Location Prior_Scale
#> 1 (Intercept)          student_t        3           19.2         5.4
#> 2          wt            uniform       NA             NA          NA
bayestestR::describe_posterior(model, priors = TRUE)  # Messed up output
#> Summary of Posterior Distribution
#> 
#> Parameter   | Median |         95% CI |   pd |          ROPE | % in ROPE |  Rhat |     ESS |                            Prior
#> -----------------------------------------------------------------------------------------------------------------------------
#> (Intercept) |  37.31 | [33.30, 41.03] | 100% | [-0.60, 0.60] |        0% | 1.002 | 3397.00 | Student_t (df=3) (19.20 +- 5.40)
#> wt          |  -5.36 | [-6.47, -4.18] | 100% | [-0.60, 0.60] |        0% | 1.002 | 3389.00 |                          Uniform

Created on 2021-04-03 by the reprex package (v2.0.0)

strengejacke avatar Apr 03 '21 13:04 strengejacke

library(rstanarm)

model <- stan_glm(mpg ~ wt, data=mtcars, prior = NULL, refresh= 0)

insight::get_priors(model)
#>     Parameter Distribution Location Scale Adjusted_Scale
#> 1 (Intercept)       normal 20.09062   2.5       15.06737
#> 2          wt      uniform       NA    NA             NA
bayestestR::describe_prior(model)
#>     Parameter Prior_Distribution Prior_Location Prior_Scale
#> 1 (Intercept)             normal       20.09062    15.06737
#> 2          wt            uniform             NA          NA
bayestestR::describe_posterior(model, priors = TRUE)  # Messed up output
#> Summary of Posterior Distribution
#> 
#> Parameter   | Median |         95% CI |   pd |          ROPE | % in ROPE |  Rhat |     ESS |                   Prior
#> --------------------------------------------------------------------------------------------------------------------
#> (Intercept) |  37.25 | [33.51, 41.20] | 100% | [-0.60, 0.60] |        0% | 1.000 | 3511.00 | Normal (20.09 +- 15.07)
#> wt          |  -5.33 | [-6.47, -4.23] | 100% | [-0.60, 0.60] |        0% | 1.000 | 3601.00 |                 Uniform

Created on 2021-04-03 by the reprex package (v2.0.0)

strengejacke avatar Apr 03 '21 13:04 strengejacke

For Student's, I have now an additional "df" column. Actually, we would need scale/shape for Gamma/Beta priors or so...

strengejacke avatar Apr 03 '21 13:04 strengejacke

See my comment here. We always had a basic support that worked well for rstanarm and simple brms models, but I think for a more comprehensive support we need to rewrite get_priors() for brms. We should take some more time for this, and for now, run tests for describe_prior() and check_prior() only on GitHub CI, not CRAN.

strengejacke avatar Apr 03 '21 19:04 strengejacke

@strengejacke I think it'd be good to have an effects and components argument for describe_priors, following the standards of describe_posteriors (so that we can match them easily). However, I wonder how to best implement it, since either I use clean_parameters() (which doesn't allow for smart filtering) or find_parameters() which doesn't returned clean names ^^

I'm not sure how to use these arguments:

https://github.com/easystats/bayestestR/blob/fbc883b7bc61bade683bf79ad81ddd3582c27b7d/R/describe_prior.R#L42-L48

DominiqueMakowski avatar Apr 06 '21 01:04 DominiqueMakowski

This is not yet possible until https://github.com/easystats/insight/pull/335 is resolved. The idea is to revise get_priors() to match all priors with the output of clean_names (see https://github.com/easystats/insight/issues/332#issuecomment-812887151). I would address this for the next release cycle.

strengejacke avatar Apr 06 '21 06:04 strengejacke

Right I see, I was trying to fix this bug in the checks before CRAN submission:

── Failure (test-describe_prior.R:85:5): describe_prior ────────────────────────
describe_prior(mod_brms) (`actual`) not equal to structure(...) (`expected`).

  `attr(actual, 'row.names')`: 1 2 3 4
`attr(expected, 'row.names')`: 1 2 3  

`actual$Parameter`:   "(Intercept)" "cyl" "wt" "sigma"
`expected$Parameter`: "(Intercept)" "cyl" "wt"        

`actual$Prior_Distribution`:   "student_t" "uniform" "uniform" "student_t"
`expected$Prior_Distribution`: "student_t" "uniform" "uniform"            

  `actual$Prior_df`: 3   3
`expected$Prior_df`: 3    

  `actual$Prior_Location`: 19.2   0  
`expected$Prior_Location`: 19.2   0 0

  `actual$Prior_Scale`: 5.4   5.4
`expected$Prior_Scale`: 5.4      

Any workaround?

DominiqueMakowski avatar Apr 06 '21 06:04 DominiqueMakowski

I suggest skipping these tests for now (I added packageVersion("insight") > "0.13.2", and those tests do not run on CRAN anyway, I think), and once we resolved the PR, we can rework describe_prior() here. It should work now, better than the current CRAN version, but not fully for complex prior settings.

strengejacke avatar Apr 06 '21 06:04 strengejacke