This report documents the results of a simulation based calibration (SBC) run for RBesT. The calibration data will be generated whenever relevant changes to the gMAP function were made. The calibration runs are performed for typical use cases of gMAP. These include the three likelihoods (binomial, gaussian & Poisson), a sparse (\(2\) trials) and dense (\(10\) trials) data situation and finally a run with a very/less conservative prior choice for between-trial heterogeniety parameter.

The calibration data presented here has been generated at and with the RBesT git version as:

## Created:  2024-11-21 10:30:15 UTC
## git hash: 5f5afb101f98533e4229f9dc267b3540b172c573
## MD5:      c8f566f591c59887653bac53ab03bb45

The MD5 hash of the calibration data file presented here must match the above listed MD5:

## /Users/weberse2/rwork/BesT/inst/sbc/calibration.rds 
##                  "c8f566f591c59887653bac53ab03bb45"

Introduction

Simulation based calibration (SBC) is a necessary condition which must be met for any Bayesian analysis with proper priors. The details are presented in Talts, et. al (see https://arxiv.org/abs/1804.06788).

Self-consistency of any Bayesian analysis with a proper prior:

\[ p(\theta) = \iint \mbox{d}\tilde{y} \, \mbox{d}\tilde{\theta} \, p(\theta|\tilde{y}) \, p(\tilde{y}|\tilde{\theta}) \, p(\tilde{\theta}) \] \[ \Leftrightarrow p(\theta) = \iint \mbox{d}\tilde{y} \, \mbox{d}\tilde{\theta} \, p(\theta,\tilde{y},\tilde{\theta}) \]

SBC procedure:

Repeat \(s=1, ..., S\) times:

  1. Sample from the prior \[\tilde{\theta} \sim p(\theta)\]

  2. Sample fake data \[\tilde{y} \sim p(y|\tilde{\theta})\]

  3. Obtain \(L\) posterior samples \[\{\theta_1, ..., \theta_L\} \sim p(\tilde{\theta}|\tilde{y})\]

  4. Calculate the rank \(r_s\) of the prior draw \(\tilde{\theta}\) wrt to the posterior sample \(\{\theta_1, ..., \theta_L\} \sim p(\tilde{\theta}|\tilde{y})\) which falls into the range \([0,L]\) out of the possible \(L+1\) ranks. The rank is calculated as \[r_s = \sum_{l=1}^L \mathbb{I}[ \theta_l < \tilde{\theta}]\]

The \(S\) ranks then form a uniform \(0-1\) density and the count in each bin has a binomial distribution with probability of \[p(r \in \mbox{Any Bin}) =\frac{(L+1)}{S}.\]

Hierarchical intercept only (random-effects intercept) model for binomial, gaussian and Poisson likelihood

Likelihood:

  • Binary \[y_i|\theta_{j} \sim \mbox{Bernoulli}(\theta_j), \] \[g(\theta) = \mbox{logit}(\theta)\]
  • Normal \[y_i|\theta_{j} \sim \mbox{Normal}(\theta_j, \sigma^2), \] \[g(\theta) = \theta\]
  • Poisson \[y_i|\theta_{j} \sim \mbox{Poisson}(\theta_j), \] \[g(\theta) = \log(\theta)\]

Hierarchical prior:

\[ g(\theta_j)|\mu,\tau \sim \mbox{Normal}(\mu, \tau^2)\]

\[\mu \sim \mbox{Normal}(m_\mu, s^2_\mu)\] \[\tau \sim \mbox{Normal}^+(0, s^2_\tau)\]

The fake data simulation function returns for binomial and Poisson data the sum of the responses while for normal the mean summary is used. Please refer to the sbc_tools.R and make_reference_rankhist.R R programs for the implementation details.

The reference runs are created with \(L=1023\) posterior draws for each replication and a total of \(S=10^4\) replications are run per case. For the evaluation here the results are reduced to \(B=L'+1=64\) bins to ensure a sufficiently large sample size per bin.

SBC results

Sampler Diagnostics Overview

family data_scenario sd_tau N total_divergent total_divergent_sim_fraction min_ess max_Rhat total_large_Rhat min_lp_ess_bulk min_lp_ess_tail
binomial dense 0.5 10000 322 0.028 109 1.010 0 234 142
binomial dense 1.0 10000 353 0.029 58 1.032 0 136 129
binomial sparse 0.5 10000 1229 0.058 9 1.217 1 19 9
binomial sparse 1.0 10000 2373 0.113 12 1.279 1 109 58
gaussian dense 0.5 10000 521 0.028 30 1.071 0 227 164
gaussian dense 1.0 10000 339 0.023 42 1.039 0 145 72
gaussian sparse 0.5 10000 21734 0.467 3 1.840 1 12 21
gaussian sparse 1.0 10000 41989 0.402 2 3.895 1 9 17
poisson dense 0.5 10000 742 0.057 64 1.035 0 147 50
poisson dense 1.0 10000 1888 0.110 12 1.233 1 123 185
poisson sparse 0.5 10000 26234 0.396 2 2.722 1 6 48
poisson sparse 1.0 10000 39505 0.411 2 11.114 1 4 16

Note: Large Rhat is defined as exceeding 1.2.

Summary Statistics

\(\chi^2\) Statistic, \(\mu\)

data_scenario likelihood sd_tau parameter statistic df p.value
1 dense binomial 0.5 mu 53.901 63 0.786
13 dense binomial 1 mu 49.830 63 0.886
25 dense gaussian 0.5 mu 64.499 63 0.424
37 dense gaussian 1 mu 49.498 63 0.893
49 dense poisson 0.5 mu 49.395 63 0.895
61 dense poisson 1 mu 68.493 63 0.296
73 sparse binomial 0.5 mu 61.363 63 0.535
78 sparse binomial 1 mu 61.222 63 0.540
83 sparse gaussian 0.5 mu 50.317 63 0.876
88 sparse gaussian 1 mu 49.126 63 0.900
93 sparse poisson 0.5 mu 58.509 63 0.637
98 sparse poisson 1 mu 57.498 63 0.672

\(\chi^2\) Statistic, \(\tau\)

data_scenario likelihood sd_tau parameter statistic df p.value
2 dense binomial 0.5 tau 63.091 63 0.473
14 dense binomial 1 tau 58.048 63 0.653
26 dense gaussian 0.5 tau 34.266 63 0.999
38 dense gaussian 1 tau 88.768 63 0.018
50 dense poisson 0.5 tau 62.938 63 0.479
62 dense poisson 1 tau 61.850 63 0.517
74 sparse binomial 0.5 tau 78.451 63 0.091
79 sparse binomial 1 tau 82.598 63 0.049
84 sparse gaussian 0.5 tau 71.245 63 0.223
89 sparse gaussian 1 tau 60.966 63 0.549
94 sparse poisson 0.5 tau 52.416 63 0.827
99 sparse poisson 1 tau 50.355 63 0.875

\(\chi^2\) Statistic, group estimates \(\theta\)

data_scenario likelihood sd_tau parameter statistic df p.value
3 dense binomial 0.5 theta\[10\] 48.013 63 0.919
4 dense binomial 0.5 theta\[1\] 74.291 63 0.156
5 dense binomial 0.5 theta\[2\] 60.979 63 0.549
6 dense binomial 0.5 theta\[3\] 62.451 63 0.496
7 dense binomial 0.5 theta\[4\] 55.142 63 0.749
8 dense binomial 0.5 theta\[5\] 54.362 63 0.773
9 dense binomial 0.5 theta\[6\] 68.800 63 0.288
10 dense binomial 0.5 theta\[7\] 67.597 63 0.323
11 dense binomial 0.5 theta\[8\] 104.563 63 0.001
12 dense binomial 0.5 theta\[9\] 80.819 63 0.065
15 dense binomial 1 theta\[10\] 65.562 63 0.388
16 dense binomial 1 theta\[1\] 63.437 63 0.461
17 dense binomial 1 theta\[2\] 53.325 63 0.802
18 dense binomial 1 theta\[3\] 61.670 63 0.524
19 dense binomial 1 theta\[4\] 59.840 63 0.590
20 dense binomial 1 theta\[5\] 58.944 63 0.622
21 dense binomial 1 theta\[6\] 79.194 63 0.082
22 dense binomial 1 theta\[7\] 45.120 63 0.957
23 dense binomial 1 theta\[8\] 55.667 63 0.733
24 dense binomial 1 theta\[9\] 79.411 63 0.079
27 dense gaussian 0.5 theta\[10\] 66.355 63 0.362
28 dense gaussian 0.5 theta\[1\] 74.342 63 0.155
29 dense gaussian 0.5 theta\[2\] 59.840 63 0.590
30 dense gaussian 0.5 theta\[3\] 86.797 63 0.025
31 dense gaussian 0.5 theta\[4\] 96.845 63 0.004
32 dense gaussian 0.5 theta\[5\] 56.000 63 0.722
33 dense gaussian 0.5 theta\[6\] 81.382 63 0.059
34 dense gaussian 0.5 theta\[7\] 58.342 63 0.643
35 dense gaussian 0.5 theta\[8\] 77.542 63 0.103
36 dense gaussian 0.5 theta\[9\] 58.419 63 0.640
39 dense gaussian 1 theta\[10\] 78.784 63 0.087
40 dense gaussian 1 theta\[1\] 70.234 63 0.248
41 dense gaussian 1 theta\[2\] 64.486 63 0.424
42 dense gaussian 1 theta\[3\] 55.398 63 0.741
43 dense gaussian 1 theta\[4\] 56.474 63 0.706
44 dense gaussian 1 theta\[5\] 61.939 63 0.514
45 dense gaussian 1 theta\[6\] 45.504 63 0.953
46 dense gaussian 1 theta\[7\] 44.352 63 0.964
47 dense gaussian 1 theta\[8\] 59.878 63 0.588
48 dense gaussian 1 theta\[9\] 66.189 63 0.367
51 dense poisson 0.5 theta\[10\] 48.166 63 0.916
52 dense poisson 0.5 theta\[1\] 64.038 63 0.440
53 dense poisson 0.5 theta\[2\] 64.397 63 0.427
54 dense poisson 0.5 theta\[3\] 62.528 63 0.493
55 dense poisson 0.5 theta\[4\] 66.432 63 0.360
56 dense poisson 0.5 theta\[5\] 68.557 63 0.295
57 dense poisson 0.5 theta\[6\] 74.291 63 0.156
58 dense poisson 0.5 theta\[7\] 59.366 63 0.607
59 dense poisson 0.5 theta\[8\] 63.987 63 0.442
60 dense poisson 0.5 theta\[9\] 62.464 63 0.495
63 dense poisson 1 theta\[10\] 78.861 63 0.086
64 dense poisson 1 theta\[1\] 53.837 63 0.788
65 dense poisson 1 theta\[2\] 73.510 63 0.172
66 dense poisson 1 theta\[3\] 63.181 63 0.470
67 dense poisson 1 theta\[4\] 62.528 63 0.493
68 dense poisson 1 theta\[5\] 87.270 63 0.023
69 dense poisson 1 theta\[6\] 60.890 63 0.552
70 dense poisson 1 theta\[7\] 78.643 63 0.088
71 dense poisson 1 theta\[8\] 60.019 63 0.583
72 dense poisson 1 theta\[9\] 92.864 63 0.009
75 sparse binomial 0.5 theta\[1\] 61.837 63 0.518
76 sparse binomial 0.5 theta\[2\] 55.142 63 0.749
77 sparse binomial 0.5 theta\[3\] 50.918 63 0.863
80 sparse binomial 1 theta\[1\] 59.226 63 0.612
81 sparse binomial 1 theta\[2\] 66.189 63 0.367
82 sparse binomial 1 theta\[3\] 76.493 63 0.118
85 sparse gaussian 0.5 theta\[1\] 51.674 63 0.845
86 sparse gaussian 0.5 theta\[2\] 64.678 63 0.418
87 sparse gaussian 0.5 theta\[3\] 68.442 63 0.298
90 sparse gaussian 1 theta\[1\] 56.768 63 0.697
91 sparse gaussian 1 theta\[2\] 82.701 63 0.049
92 sparse gaussian 1 theta\[3\] 63.155 63 0.471
95 sparse poisson 0.5 theta\[1\] 71.552 63 0.215
96 sparse poisson 0.5 theta\[2\] 65.882 63 0.377
97 sparse poisson 0.5 theta\[3\] 72.179 63 0.200
100 sparse poisson 1 theta\[1\] 67.162 63 0.337
101 sparse poisson 1 theta\[2\] 54.413 63 0.771
102 sparse poisson 1 theta\[3\] 56.294 63 0.712

Session Info

## R version 4.4.2 (2024-10-31)
## Platform: aarch64-apple-darwin20
## Running under: macOS Sonoma 14.6.1
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRblas.0.dylib 
## LAPACK: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.0
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## time zone: Europe/Zurich
## tzcode source: internal
## 
## attached base packages:
## [1] tools     stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] purrr_1.0.2        rstan_2.32.6       StanHeaders_2.32.9 here_1.0.1        
##  [5] ggplot2_3.5.1      broom_1.0.6        tidyr_1.3.1        dplyr_1.1.4       
##  [9] assertthat_0.2.1   knitr_1.47        
## 
## loaded via a namespace (and not attached):
##  [1] utf8_1.2.4         generics_0.1.3     digest_0.6.35      magrittr_2.0.3    
##  [5] evaluate_0.24.0    grid_4.4.2         fastmap_1.2.0      rprojroot_2.0.4   
##  [9] jsonlite_1.8.8     pkgbuild_1.4.4     backports_1.5.0    gridExtra_2.3     
## [13] fansi_1.0.6        QuickJSR_1.2.2     scales_1.3.0       codetools_0.2-20  
## [17] cli_3.6.2          rlang_1.1.4        munsell_0.5.1      withr_3.0.0       
## [21] yaml_2.3.8         inline_0.3.19      parallel_4.4.2     colorspace_2.1-0  
## [25] curl_5.2.1         vctrs_0.6.5        R6_2.5.1           matrixStats_1.3.0 
## [29] stats4_4.4.2       lifecycle_1.0.4    V8_4.4.2           pkgconfig_2.0.3   
## [33] RcppParallel_5.1.7 pillar_1.9.0       gtable_0.3.5       loo_2.7.0         
## [37] glue_1.7.0         Rcpp_1.0.12        xfun_0.45          tibble_3.2.1      
## [41] tidyselect_1.2.1   htmltools_0.5.8.1  rmarkdown_2.27     compiler_4.4.2