Just curious if there’s any reason it was decided to have rstan::summary() show quantiles rather than highest posterior density intervals?
3 Likes
My guess is that it’s 1) easy 2) purely descriptive and 3) HPDI does not always yield one set of intervals (e.g., when bimodality is present), which complicates output.
Using HDInterval on stan samples is relatively simple.
3 Likes
Yeah I think what @Stephen_Martin said covers it.
Thanks @Stephen_Martin, that makes sense. In the case of a multimodal posterior, I wonder if that’s not particularly pertinent since, as I understand it, such posteriors usually indicate something has gone wrong with the model specification/sampling, no?
Not necessarily. I’ve had some funny models where it converged, but the posterior has one long tail with a small hump. Consequently, I get two intervals - One near the big mode and one near the small mode.
Like:
/\
___/ \----/`\___
By all diagnostics, it seemed fine - Just a funky posterior shape with one small hump that caused 2 intervals.
1 Like
. 4) HPDI is not invariant to transformations
3 Likes
@Stephen_Martin Nice art
Who needs bayesplot anyway? We need a bayesascii
7 Likes
And transformations of unimodal distributions can be multimodal. (Take the inverse logit of a diffuse Gaussian)
- They require density estimates, which are much much noisier than CDF estimates especially in the tails (which is the only region we care about).
Basically HPDIs are a computationally burdensome, minimally interpretable, and highly contextual bad solution to a non-existent problem and should be left in the darkest ditch of history.
5 Likes
@anon75146577 Do they require density estimates? I thought they were constructed by computing the # of samples that would be contained in the interval, then finding the narrowest interval that contains that number of samples?
1 Like
That’s a density estimator.
4 Likes
Ok, wondered if that might be the case. I’m used to thinking of density estimation as the result of a complicated thing like stats::density() computes, but makes sense that any time you try to relate a count of samples to an interval, this is also density estimation.
And if will still be unreliable in the tails. I mean if you want to commit 50% HPDIs it’s still a bad idea, but at least it’s not a bad estimator of a bad idea
2 Likes
devout :)
# install.packages("remotes")
# remotes::install_github("coolbutuseless/devout")
library(bayesplot)
library(devout)
x <- example_mcmc_draws()
p <- mcmc_dens(x, c("alpha", "beta[1]"))
ascii(width = 80)
p
dev.off()
#> alpha beta[1]
#> # XXXX # XXX
#> # XX XX # X XXXXX
#> # XX XX # X XX
#> # X XX # X XX
#> # XX XX # X X
#> # X X # XX XX
#> # X XX # X XX
#> # X X # XX X
#> # X X # XX XX
#> # XX XX # X XX
#> # X X # XX X
#> # XX XX # XX XX
#> # XX X # X XX
#> # XX XX # XX XX
#> # XX XX # XX XX
#> # XX XXX # XX XX
#> # XXX XXX # XXXXX XXX
#> #########O#########O#########O#########O###O#########O########O#########O#######
#> -50 -25 0 25 -0.5 0.0 0.5 1.0
11 Likes