[avehtari]

This PR was assisted by Claude

Summary

Current master: .ess_quantile() computes the ESS of the indicator x <= quantile(x, prob) for every prob. For discrete draws, or for draws with many ties, this makes the two tails behave differently. If more than 5% of the draws are at the maximum, the 95% quantile is the maximum, the indicator is TRUE for every draw, and ess_quantile(x, 0.95) and ess_tail(x) are NA. The lower tail is NA only if more than 95% of the draws are at the maximum.

New: This PR uses x >= quantile(x, prob) for prob > 0.5. Because type 7 quantiles satisfy quantile(-x, p) = -quantile(x, 1 - p), the upper-tail indicator of x is then the lower-tail indicator of -x, and ess_quantile(x, p) equals ess_quantile(-x, 1 - p).

For the tail-ESS, the upper-tail indicator changes from x <= q95 to x >= q95, as in Stan (stan::analyze::split_rank_normalized_ess(), used by CmdStan stansummary). With stan-dev/stan#3454, both then give the same tail-ESS for discrete x and -x, and a finite tail-ESS for discrete draws with more than 5% at the maximum.

Effect

• Continuous draws without ties at the quantile: no change. x >= q is the complement of x <= q, and an indicator and its complement have the same ESS.
• Draws with more than 5% at the maximum: ess_quantile(x, 0.95) and ess_tail(x) are finite.
• Draws with more than 95% at one value: ess_tail(x) is NA, at the minimum and at the maximum alike.
• mcse_quantile(), and pareto_smooth() without r_eff (it uses ess_tail()), change only for draws with ties at the quantile.

Tests

A new test in test-convergence.R uses 4 chains of 1000 Bernoulli(0.3) draws. It checks that ess_quantile(x, 0.95) and ess_tail(x) are not NA, that ess_quantile(x, p) equals ess_quantile(-x, 1 - p) for p = 0.05 and 0.95, and that ess_tail(x) equals ess_tail(-x). It also checks that for continuous draws the result equals the ESS of x <= q. With master, 4 of the 6 expectations fail.

Copyright and Licensing

Copyright holder: Aalto University

By submitting this pull request, the copyright holder is agreeing to license the submitted work under the following licenses:

• Code: BSD 3-clause (https://opensource.org/licenses/BSD-3-Clause)
• Documentation: CC-BY 4.0 (https://creativecommons.org/licenses/by/4.0/)

[github-actions]

This is how benchmark results would change (along with a 95% confidence interval in relative change) if 3dbe017 is merged into master:

• ✔️as_draws_array: 163ms -> 164ms [-0.67%, +2.1%]
• ✔️as_draws_df: 80.9ms -> 80.7ms [-1.35%, +0.78%]
• ✔️as_draws_list: 179ms -> 180ms [-0.45%, +1.16%]
• ✔️as_draws_matrix: 25.9ms -> 25.9ms [-0.93%, +0.91%]
• ✔️as_draws_rvars: 132ms -> 133ms [-1.28%, +1.48%]
• ✔️summarise_draws_100_variables: 773ms -> 777ms [-0.07%, +1.06%]
• ✔️summarise_draws_10_variables: 86.5ms -> 86.4ms [-0.76%, +0.45%]
  Further explanation regarding interpretation and methodology can be found in the documentation.

[florence-bockting]

Looks good to me. I just add the entry to the NEWS.md