Every fit so far has held the variance components \(\tau^2_k\) fixed at the lmer
estimates. That is a real assumption: it says you know exactly how much
subjects vary, and it makes the posterior for everything else narrower
than it should be. This chapter relaxes it.
The choice is made per random coefficient, through the prior family.
| Family | \(\tau^2_k\) | Consequence |
|---|---|---|
dNormal |
Fixed at the reference estimate | Joint posterior stays Gaussian; exact iid draws |
dIndependent_Normal_Gamma |
Sampled, Gamma prior on \(1/\tau^2_k\) | Target no longer Gaussian; Gibbs sweeps and a pilot stage |
Calibrate as usual, adding a prior weight for the dispersion:
ps <- Prior_Setup_GLMM(form, data = dat, pop.pwt = 0.01,
pop.dispersion.pwt = 0.2)
ps$pop.ing_prior
#> $`(Intercept)`
#> $`(Intercept)`$shape
#> [1] 3.25
#>
#> $`(Intercept)`$rate
#> [1] 3028.828
#>
#> $`(Intercept)`$disp_lower
#> [1] 632.3273
#>
#> $`(Intercept)`$disp_upper
#> [1] 2998.203
#>
#> $`(Intercept)`$max_disp_perc
#> [1] 0.99
#>
#>
#> $Days_c
#> $Days_c$shape
#> [1] 3.25
#>
#> $Days_c$rate
#> [1] 66.08751
#>
#> $Days_c$disp_lower
#> [1] 13.79707
#>
#> $Days_c$disp_upper
#> [1] 65.41931
#>
#> $Days_c$max_disp_perc
#> [1] 0.99pop.dispersion.pwt = 0.2 gives the dispersion prior the
weight of about 20 percent of the groups. With only 18 subjects there is
genuinely little information about a variance component, so this is not
a small choice — compare it against 0.05 and
0.5 before trusting the width of your intervals.
Each ING component carries disp_lower and
disp_upper, and the sampled \(\tau^2_k\) is confined to that
interval.
The bounds are quantiles of the limiting posterior of the variance component, not of the prior. That distinction matters: a Gamma prior on a precision has a fat tail toward zero variance, and a \(\tau^2\) driven to zero collapses the hierarchy — every group’s coefficient is pinned to the population mean, and the sampler’s coupling becomes total. Truncating away that region is what keeps the chain geometrically ergodic.
The lower bound does double duty. The sampler evaluates its
convergence rate at disp_lower, because smaller \(\tau^2\) means tighter coupling and a
slower chain, so a bound computed there is valid across the entire
window. That is conservative rather than correct-on-average, which is
the right direction for a certificate. See
vignette("Chapter-B05", package = "lmebayesCore") and
inst/ING_TRUNCATION_WINDOW.md in the package sources.
pf_ing <- pfamily_list(ps, ptypes = "dIndependent_Normal_Gamma")
class(pf_ing[["Days_c"]])
#> [1] "pfamily"ptypes also accepts a named vector when you want to
sample some components and fix others — for example
ptypes = c("(Intercept)" = "dNormal", "Days_c" = "dIndependent_Normal_Gamma").
set.seed(1)
fit <- lmerb(form, data = dat,
pfamily_list = pf_ing,
dispersion_ranef = ps$group.dispersion,
n = 300L, progbar = FALSE)fit
#> Call:
#> lmerb(formula = form, data = dat, pfamily_list = pf_ing, dispersion_ranef = ps$group.dispersion,
#> n = 300L, progbar = FALSE)
#>
#> Bayesian linear mixed model [300 draws, two-block Gibbs]
#> Formula: Reaction ~ Days_c + (1 + Days_c || Subject)
#>
#> Random effects (lmer reference; tau^2 sampled for non-dNormal components):
#> Groups Name Std.Dev.
#> Subject (Intercept) 37.533
#> Subject.1 Days_c 5.922
#> Residual 25.592
#> Number of obs: 180, groups: Subject, 18
#>
#> Posterior mean tau^2_k: (Intercept) = 1391, Days_c = 31.45
#>
#> --- Block 2 hyperparameters (gamma at lmer tau^2 plug-in; MCMC means when simulated) ---
#>
#> RE component parameter gamma @ lmer tau2 fixef.means draws SD
#> ------------ ----------- ------------ ------------ ----------
#> (Intercept) (Intercept) 298.5079 298.3910 8.2584
#> Days_c (Intercept) 10.3748 10.3266 1.5139Two things changed in the header. The engine is now
two-block Gibbs rather than exact iid, and the
random-effects line says the variance components are sampled rather than
fixed.
head(fit$popef.dispersion)
#> (Intercept) Days_c
#> [1,] 1600.298 46.59567
#> [2,] 1145.382 24.62881
#> [3,] 1420.746 33.92399
#> [4,] 1893.599 24.74103
#> [5,] 1399.663 30.99423
#> [6,] 1348.311 41.06287
round(fit$popef.dispersion.mean, 3)
#> (Intercept) Days_c
#> 1391.40 31.45fit$popef.dispersion holds the \(\tau^2_k\) draws — one column per random
coefficient. Now they have a distribution:
round(t(apply(fit$popef.dispersion, 2, quantile,
probs = c(0.025, 0.5, 0.975))), 2)
#> 2.5% 50% 97.5%
#> (Intercept) 757.35 1328.48 2506.01
#> Days_c 15.49 29.57 56.05On the standard-deviation scale, which is easier to read against
lmer:
round(t(apply(sqrt(fit$popef.dispersion), 2, quantile,
probs = c(0.025, 0.5, 0.975))), 3)
#> 2.5% 50% 97.5%
#> (Intercept) 27.519 36.448 50.055
#> Days_c 3.935 5.437 7.486Compare against the lmer point estimates:
print(lme4::VarCorr(fit$lmer), comp = "Std.Dev.")
#> Groups Name Std.Dev.
#> Subject (Intercept) 37.5334
#> Subject.1 Days_c 5.9221
#> Residual 25.5918The lmer numbers should sit inside the credible
intervals. They are point estimates with no uncertainty attached; the
intervals here are wide, and honestly so — 18 subjects is not much
information about how subjects vary.
set.seed(1)
fit_fixed <- lmerb(form, data = dat,
pfamily_list = pfamily_list(ps),
dispersion_ranef = ps$group.dispersion,
n = 300L, progbar = FALSE)sd_fixed <- vapply(fit_fixed$popef, function(m) sd(m[, 1]), numeric(1))
sd_ing <- vapply(fit$popef, function(m) sd(m[, 1]), numeric(1))
round(data.frame(
fixed_tau2 = sd_fixed,
sampled_tau2 = sd_ing,
ratio = sd_ing / sd_fixed
), 4)
#> fixed_tau2 sampled_tau2 ratio
#> (Intercept) 9.4427 8.2584 0.8746
#> Days_c 1.3637 1.5139 1.1101Posterior standard deviations on the population coefficients are wider once the variance components are uncertain. That widening is the point: fixing \(\tau^2\) at a point estimate and reporting the resulting intervals understates what you know.
The convergence theory of Chapter 06 assumed chains start at the exact posterior mean, making the mean-error term of the bound exactly zero. With sampled variance components the joint target is no longer Gaussian, so there is no exact mean to start from.
The engine’s answer is a pilot stage: a batch of short chains whose cross-chain mean becomes the start point for the main chains. The pilot draws are discarded. Their job is not to estimate anything but to produce a starting point whose distance from the truth can be bounded.
fit$pilot
#> $n
#> [1] 28
#>
#> $m_convergence
#> [1] 3
#>
#> $chisq
#> $chisq$Q
#> [1] 1.313762
#>
#> $chisq$df
#> [1] 2
#>
#> $chisq$p_value
#> [1] 0.518466
#>
#> $chisq$n_pilot
#> [1] 28
#>
#>
#> $draws
#> $draws$sweep_history
#>
#> --- two-block [pilot stage summary: fixef by sweep (3 sweeps)] ---
#> Block 2 fixed effects (mode and chain stats after each sweep):
#> Random effect Covariate mode/sweep mean sd
#> --------------------------------------------------------------------------------------------
#> (Intercept) (Intercept) mode 298.5079
#> (Intercept) (Intercept) sweep 1 297.8788 10.0647
#> (Intercept) (Intercept) sweep 2 298.4804 7.7030
#> (Intercept) (Intercept) sweep 3 298.1352 7.5755
#> Days_c (Intercept) mode 10.3748
#> Days_c (Intercept) sweep 1 10.7468 1.4401
#> Days_c (Intercept) sweep 2 10.2859 1.2920
#> Days_c (Intercept) sweep 3 10.6915 1.6056
#>
#> Block 2 RE precision (1/tau^2; chain mean and SD after each sweep):
#> Random effect sweep mean sd
#> ------------------------------------------------------------
#> (Intercept) sweep 1 0.0008 0.0003
#> (Intercept) sweep 2 0.0008 0.0002
#> (Intercept) sweep 3 0.0009 0.0002
#> Days_c sweep 1 0.0392 0.0129
#> Days_c sweep 2 0.0372 0.0110
#> Days_c sweep 3 0.0378 0.0122
#>
#>
#> $draws$n
#> [1] 28
#>
#> $draws$popef
#> $draws$popef$`(Intercept)`
#> (Intercept)
#> [1,] 287.5407
#> [2,] 293.7575
#> [3,] 296.7314
#> [4,] 306.5957
#> [5,] 298.1009
#> [6,] 296.9867
#> [7,] 300.9297
#> [8,] 297.9107
#> [9,] 301.2009
#> [10,] 297.0981
#> [11,] 288.0879
#> [12,] 292.5252
#> [13,] 294.8369
#> [14,] 311.5246
#> [15,] 299.7641
#> [16,] 310.8793
#> [17,] 302.9973
#> [18,] 281.8691
#> [19,] 311.0191
#> [20,] 292.4624
#> [21,] 298.8445
#> [22,] 286.6510
#> [23,] 306.5081
#> [24,] 290.0598
#> [25,] 303.5048
#> [26,] 304.3594
#> [27,] 294.8012
#> [28,] 300.2390
#>
#> $draws$popef$Days_c
#> (Intercept)
#> [1,] 11.222950
#> [2,] 9.469730
#> [3,] 12.518366
#> [4,] 10.985740
#> [5,] 9.645349
#> [6,] 12.074105
#> [7,] 9.177914
#> [8,] 11.955890
#> [9,] 9.035576
#> [10,] 11.793516
#> [11,] 12.708542
#> [12,] 11.710641
#> [13,] 11.549556
#> [14,] 13.425279
#> [15,] 11.027747
#> [16,] 12.867358
#> [17,] 11.961647
#> [18,] 8.672453
#> [19,] 11.032201
#> [20,] 10.974276
#> [21,] 10.265890
#> [22,] 7.943870
#> [23,] 6.874187
#> [24,] 9.209482
#> [25,] 10.843135
#> [26,] 9.708206
#> [27,] 11.730409
#> [28,] 8.978489
#>
#>
#> $draws$groupef
#> draw Subject (Intercept) Days_c
#> 1 1 308 342.1144 21.54528019
#> 2 1 309 218.6502 7.33574088
#> 3 1 310 245.5235 12.80660354
#> 4 1 330 303.2802 6.01433155
#> 5 1 331 319.4513 8.09794782
#> 6 1 332 309.1553 13.72657430
#> 7 1 333 322.9175 10.12226301
#> 8 1 334 282.8377 11.98467400
#> 9 1 335 263.0113 0.52480150
#> 10 1 337 352.6487 17.57947961
#> 11 1 349 272.8557 12.81494101
#> 12 1 350 312.3657 14.34972776
#> 13 1 351 289.8535 6.93684108
#> 14 1 352 333.7984 11.78890560
#> 15 1 369 302.6338 10.21994965
#> 16 1 370 307.6069 17.69706107
#> 17 1 371 298.3370 14.30061900
#> 18 1 372 310.0293 12.96791787
#> 19 2 308 332.3437 18.95022608
#> 20 2 309 223.0653 2.46005962
#> 21 2 310 230.0001 11.12481220
#> 22 2 330 316.4974 2.62483180
#> 23 2 331 308.9870 9.12659443
#> 24 2 332 291.8458 8.88688311
#> 25 2 333 323.4866 10.73135954
#> 26 2 334 293.2396 9.61014401
#> 27 2 335 258.9817 3.57923446
#> 28 2 337 354.8635 16.86939591
#> 29 2 349 293.0589 14.63719746
#> 30 2 350 311.6357 15.22012759
#> 31 2 351 287.6798 9.47954096
#> 32 2 352 337.8096 12.17935254
#> 33 2 369 298.4364 11.25740928
#> 34 2 370 282.3360 24.31494643
#> 35 2 371 308.4077 11.56134564
#> 36 2 372 318.3786 10.25454904
#> 37 3 308 355.4612 16.02163698
#> 38 3 309 227.1736 4.33310753
#> 39 3 310 241.4012 10.53718548
#> 40 3 330 303.0614 4.85415938
#> 41 3 331 286.1737 5.72360853
#> 42 3 332 313.9013 15.36090456
#> 43 3 333 306.0752 9.86804517
#> 44 3 334 295.8035 15.39322591
#> 45 3 335 250.4616 -0.43212266
#> 46 3 337 364.3758 17.16464726
#> 47 3 349 267.5221 12.34918763
#> 48 3 350 308.1251 17.16883479
#> 49 3 351 285.0292 11.45565020
#> 50 3 352 331.0759 12.72769679
#> 51 3 369 302.5843 12.91622415
#> 52 3 370 302.6116 15.49600950
#> 53 3 371 286.8786 10.37506480
#> 54 3 372 321.5099 10.03933892
#> 55 4 308 342.9693 19.93786371
#> 56 4 309 212.0790 5.48779241
#> 57 4 310 227.0907 10.41590543
#> 58 4 330 305.2046 8.76772308
#> 59 4 331 307.0714 2.19886267
#> 60 4 332 294.7439 13.06579305
#> 61 4 333 309.0650 6.80009240
#> 62 4 334 297.4442 11.16674191
#> 63 4 335 245.4160 -6.82376986
#> 64 4 337 364.3474 18.15414839
#> 65 4 349 274.5581 12.44377119
#> 66 4 350 307.9864 18.22179026
#> 67 4 351 295.9237 10.67980563
#> 68 4 352 320.6010 13.45480466
#> 69 4 369 302.7202 13.87352698
#> 70 4 370 309.5899 16.43483479
#> 71 4 371 298.1674 12.16441059
#> 72 4 372 327.7141 14.21748454
#> 73 5 308 323.1387 15.46521826
#> 74 5 309 215.8321 9.14482731
#> 75 5 310 234.6459 7.41613885
#> 76 5 330 303.6025 6.44360842
#> 77 5 331 319.2735 9.07547969
#> 78 5 332 287.2094 4.43771144
#> 79 5 333 314.5059 11.21421182
#> 80 5 334 282.1397 13.30730352
#> 81 5 335 243.1427 1.42542749
#> 82 5 337 372.3530 13.46887996
#> 83 5 349 290.9026 14.21782581
#> 84 5 350 336.3735 18.49525324
#> 85 5 351 287.4741 7.09153862
#> 86 5 352 331.0923 9.31987390
#> 87 5 369 287.5721 10.91419418
#> 88 5 370 304.4750 16.64945067
#> 89 5 371 292.2828 11.33262972
#> 90 5 372 318.4176 12.32330012
#> 91 6 308 344.0166 19.35790348
#> 92 6 309 230.1489 4.27489523
#> 93 6 310 223.7029 8.33715608
#> 94 6 330 319.6468 4.85510354
#> 95 6 331 325.5620 7.23012695
#> 96 6 332 311.8333 8.92968165
#> 97 6 333 298.9319 6.11163650
#> 98 6 334 294.7657 8.85510886
#> 99 6 335 241.3640 -1.04502586
#> 100 6 337 366.2395 24.01533546
#> 101 6 349 278.9958 14.06935127
#> 102 6 350 317.2571 18.04938821
#> 103 6 351 296.7064 8.35058486
#> 104 6 352 341.7083 15.59168523
#> 105 6 369 314.7045 11.20947645
#> 106 6 370 282.6083 9.55759146
#> 107 6 371 292.2798 5.15516686
#> 108 6 372 314.6697 10.35903083
#> 109 7 308 341.4652 22.58733647
#> 110 7 309 234.5074 6.61121240
#> 111 7 310 237.2918 8.15831353
#> 112 7 330 302.3899 0.88561733
#> 113 7 331 324.4390 9.45114953
#> 114 7 332 316.1415 8.33906640
#> 115 7 333 319.5547 10.52101106
#> 116 7 334 282.6366 11.12338188
#> 117 7 335 256.2644 1.60795713
#> 118 7 337 365.3119 21.13322234
#> 119 7 349 292.3287 14.17578209
#> 120 7 350 307.2386 14.44544926
#> 121 7 351 279.9613 5.70320943
#> 122 7 352 319.7275 13.73262503
#> 123 7 369 293.4945 7.90013722
#> 124 7 370 297.5748 16.75364203
#> 125 7 371 304.3577 10.57831793
#> 126 7 372 315.9154 11.94946394
#> 127 8 308 326.9347 20.18933459
#> 128 8 309 222.9010 6.50551570
#> 129 8 310 235.6415 7.76686947
#> 130 8 330 292.2893 6.49009609
#> 131 8 331 293.4353 7.37816680
#> 132 8 332 310.2404 13.97218199
#> 133 8 333 322.6715 12.58587552
#> 134 8 334 281.1703 16.55651938
#> 135 8 335 250.2035 7.30476665
#> 136 8 337 365.0367 14.65948581
#> 137 8 349 283.4785 13.32564936
#> 138 8 350 323.0141 16.86931849
#> 139 8 351 283.2815 9.75414429
#> 140 8 352 332.8866 16.23004763
#> 141 8 369 310.9696 11.06495901
#> 142 8 370 293.8579 17.93000570
#> 143 8 371 285.9850 8.43751714
#> 144 8 372 319.1772 16.45333161
#> 145 9 308 346.0509 19.48016687
#> 146 9 309 227.5443 6.19428446
#> 147 9 310 229.4871 10.63377863
#> 148 9 330 293.9346 3.23905669
#> 149 9 331 311.2277 11.20724656
#> 150 9 332 301.9799 12.57524048
#> 151 9 333 313.5163 10.92184705
#> 152 9 334 283.6213 9.29317524
#> 153 9 335 268.2440 1.66246722
#> 154 9 337 371.6352 10.23503306
#> 155 9 349 287.0067 12.37500746
#> 156 9 350 316.3085 17.64305155
#> 157 9 351 293.6602 6.94351654
#> 158 9 352 350.4292 13.07486360
#> 159 9 369 304.0832 15.15972162
#> 160 9 370 288.9211 12.53791190
#> 161 9 371 277.9855 9.53538073
#> 162 9 372 311.7656 10.68880341
#> 163 10 308 341.2708 23.59978843
#> 164 10 309 224.6703 3.57382722
#> 165 10 310 244.9173 6.08970132
#> 166 10 330 307.7568 7.50906522
#> 167 10 331 308.6323 1.49444369
#> 168 10 332 303.2943 14.00895238
#> 169 10 333 314.0696 11.88653338
#> 170 10 334 299.5341 11.47677708
#> 171 10 335 231.0234 2.13023953
#> 172 10 337 367.7746 15.28219417
#> 173 10 349 272.8836 16.72448596
#> 174 10 350 303.5257 9.83422526
#> 175 10 351 295.4144 9.00702577
#> 176 10 352 334.6084 11.68093589
#> 177 10 369 319.5737 5.24335781
#> 178 10 370 283.5683 17.34988490
#> 179 10 371 298.7678 12.31887479
#> 180 10 372 310.8551 6.61010872
#> 181 11 308 347.6577 22.65220462
#> 182 11 309 219.8993 9.20688174
#> 183 11 310 248.1434 8.50765997
#> 184 11 330 306.8474 1.99761676
#> 185 11 331 327.4854 7.49610686
#> 186 11 332 306.0385 7.88778735
#> 187 11 333 292.8048 8.60077969
#> 188 11 334 287.3902 11.70313394
#> 189 11 335 265.8404 3.34928919
#> 190 11 337 356.7030 16.63218825
#> 191 11 349 269.0511 8.19692602
#> 192 11 350 303.6300 13.22204736
#> 193 11 351 282.7582 10.80613110
#> 194 11 352 321.0505 13.39605411
#> 195 11 369 306.3262 10.61014367
#> 196 11 370 292.3537 14.72972976
#> 197 11 371 294.7023 13.49807552
#> 198 11 372 329.8807 11.98805299
#> 199 12 308 335.4558 22.18399572
#> 200 12 309 209.5003 4.01879909
#> 201 12 310 237.3709 9.02145313
#> 202 12 330 299.7801 10.84127724
#> 203 12 331 301.4123 3.80943484
#> 204 12 332 325.8746 11.19417328
#> 205 12 333 311.2512 7.03360979
#> 206 12 334 299.3816 12.41661494
#> 207 12 335 242.5615 -0.03764573
#> 208 12 337 374.9328 18.15157978
#> 209 12 349 267.7217 17.81835039
#> 210 12 350 303.7574 17.84995109
#> 211 12 351 299.4708 7.56128212
#> 212 12 352 346.5809 19.77469751
#> 213 12 369 314.5311 13.65282671
#> 214 12 370 283.8065 20.59139473
#> 215 12 371 294.8348 10.51448550
#> 216 12 372 317.4155 10.62616572
#> 217 13 308 358.1113 24.89141191
#> 218 13 309 211.4792 2.83498429
#> 219 13 310 233.6796 3.81059573
#> 220 13 330 308.4465 6.09280619
#> 221 13 331 320.5605 4.14560060
#> 222 13 332 315.5969 12.18165973
#> 223 13 333 314.6144 6.13156235
#> 224 13 334 293.6102 16.29276841
#> 225 13 335 258.0520 1.39206618
#> 226 13 337 349.4123 17.37467084
#> 227 13 349 282.0155 9.93812350
#> 228 13 350 309.4231 18.54185744
#> 229 13 351 292.0443 9.78520351
#> 230 13 352 335.2547 12.84937311
#> 231 13 369 315.0199 11.34712448
#> 232 13 370 281.6302 14.21220138
#> 233 13 371 301.8017 11.48677596
#> 234 13 372 323.2950 12.57277654
#> 235 14 308 347.8397 18.02221904
#> 236 14 309 233.4972 1.03775420
#> 237 14 310 242.8236 6.34137952
#> 238 14 330 308.6788 4.71207131
#> 239 14 331 304.7936 6.83315241
#> 240 14 332 309.2339 9.93757981
#> 241 14 333 313.2748 7.36284235
#> 242 14 334 305.3716 12.09122574
#> 243 14 335 256.2787 3.98192698
#> 244 14 337 383.2854 21.02458121
#> 245 14 349 275.2225 14.59319987
#> 246 14 350 325.7448 14.07899021
#> 247 14 351 279.9024 7.90240693
#> 248 14 352 329.5989 18.50671636
#> 249 14 369 315.6739 4.47617225
#> 250 14 370 278.9449 18.48933210
#> 251 14 371 298.5216 15.10827670
#> 252 14 372 318.1578 15.93902724
#> 253 15 308 323.0461 15.94360655
#> 254 15 309 221.3585 3.78231653
#> 255 15 310 228.2355 4.76409286
#> 256 15 330 298.1527 1.68720659
#> 257 15 331 302.9006 9.30479140
#> 258 15 332 300.1191 12.31456990
#> 259 15 333 317.3754 10.73994425
#> 260 15 334 277.5942 7.24639314
#> 261 15 335 253.8430 5.68305153
#> 262 15 337 378.2863 16.36791071
#> 263 15 349 287.2326 12.46506248
#> 264 15 350 312.2120 20.72260071
#> 265 15 351 280.5372 7.52521377
#> 266 15 352 327.5799 11.51567758
#> 267 15 369 291.9143 11.59866000
#> 268 15 370 299.7506 15.75820093
#> 269 15 371 292.3968 11.30075236
#> 270 15 372 327.3229 13.68476990
#> 271 16 308 346.1977 23.29096119
#> 272 16 309 210.7023 4.29628648
#> 273 16 310 241.5848 5.91971284
#> 274 16 330 308.4787 5.26072884
#> 275 16 331 317.0078 2.01816032
#> 276 16 332 309.1918 12.78882792
#> 277 16 333 319.4283 13.54103049
#> 278 16 334 294.4355 11.94082751
#> 279 16 335 259.0316 -0.02171594
#> 280 16 337 383.3306 22.88942810
#> 281 16 349 271.8507 10.38460305
#> 282 16 350 323.3234 19.73326458
#> 283 16 351 296.0060 5.20859001
#> 284 16 352 336.5053 15.99895657
#> 285 16 369 302.2511 12.40854530
#> 286 16 370 272.5784 19.82328151
#> 287 16 371 283.3146 5.74511671
#> 288 16 372 327.2521 13.39690110
#> 289 17 308 335.2901 17.98893805
#> 290 17 309 224.8838 6.91337212
#> 291 17 310 232.2632 5.34536501
#> 292 17 330 299.1314 5.79851434
#> 293 17 331 306.8887 3.86302511
#> 294 17 332 298.7432 9.65593158
#> 295 17 333 314.8442 8.10323210
#> 296 17 334 293.2873 12.26107293
#> 297 17 335 255.2739 4.46738013
#> 298 17 337 393.5447 18.30968698
#> 299 17 349 268.6022 14.59663718
#> 300 17 350 302.8697 21.17092409
#> 301 17 351 282.7424 8.10134394
#> 302 17 352 351.4722 14.58613303
#> 303 17 369 305.2914 12.72787911
#> 304 17 370 305.0190 14.01643243
#> 305 17 371 299.5960 12.85525575
#> 306 17 372 312.1804 12.39898959
#> 307 18 308 338.6575 19.60939002
#> 308 18 309 229.7206 3.16892825
#> 309 18 310 234.2153 5.66644867
#> 310 18 330 296.0288 -1.55095180
#> 311 18 331 302.8813 10.47957789
#> 312 18 332 306.1700 6.10299564
#> 313 18 333 310.3183 5.46629697
#> 314 18 334 283.2953 14.09553115
#> 315 18 335 250.9227 1.20610528
#> 316 18 337 365.5168 13.30758626
#> 317 18 349 272.5434 10.55470208
#> 318 18 350 305.0299 17.44890465
#> 319 18 351 298.7536 7.05689323
#> 320 18 352 318.3833 17.69577669
#> 321 18 369 315.1591 17.27332963
#> 322 18 370 290.2072 21.78585617
#> 323 18 371 297.1259 8.18212386
#> 324 18 372 309.8317 7.01270975
#> 325 19 308 348.4585 22.27535673
#> 326 19 309 229.6730 6.16888033
#> 327 19 310 237.5383 7.91970905
#> 328 19 330 319.6113 3.80601676
#> 329 19 331 303.6157 7.67456887
#> 330 19 332 302.3359 6.61635166
#> 331 19 333 308.1509 10.81022952
#> 332 19 334 311.3662 15.82327912
#> 333 19 335 244.1206 0.97446701
#> 334 19 337 388.2827 20.51813457
#> 335 19 349 285.2332 13.32742418
#> 336 19 350 323.0863 20.05872384
#> 337 19 351 289.3300 5.48438491
#> 338 19 352 333.1701 9.20300872
#> 339 19 369 313.1202 10.91983458
#> 340 19 370 273.3242 13.42797824
#> 341 19 371 298.6807 12.16093550
#> 342 19 372 320.1878 8.25259406
#> 343 20 308 333.5854 28.08409359
#> 344 20 309 231.9529 4.21897875
#> 345 20 310 234.1778 6.46506937
#> 346 20 330 309.5382 3.39235575
#> 347 20 331 304.4817 5.14676266
#> 348 20 332 307.1498 12.13370562
#> 349 20 333 297.2818 14.41258614
#> 350 20 334 297.0918 12.90761980
#> 351 20 335 258.3064 -2.75755110
#> 352 20 337 362.2099 14.80082768
#> 353 20 349 281.9191 14.54623853
#> 354 20 350 302.9644 17.66271185
#> 355 20 351 278.6308 7.91875176
#> 356 20 352 331.1181 10.73387698
#> 357 20 369 307.9894 14.08901835
#> 358 20 370 278.4048 15.33300874
#> 359 20 371 311.0351 10.72603547
#> 360 20 372 338.3062 7.21007404
#> 361 21 308 342.2465 17.11130231
#> 362 21 309 227.0016 7.02721309
#> 363 21 310 236.3128 11.11029994
#> 364 21 330 311.0279 0.23926328
#> 365 21 331 306.8469 8.66870715
#> 366 21 332 287.8043 10.34638187
#> 367 21 333 311.9622 8.01101285
#> 368 21 334 286.0188 12.68565755
#> 369 21 335 253.3595 0.79759160
#> 370 21 337 376.8475 18.47400470
#> 371 21 349 280.6304 10.99977972
#> 372 21 350 321.6189 20.59456581
#> 373 21 351 290.8117 6.17382843
#> 374 21 352 336.9908 15.37035485
#> 375 21 369 301.3858 13.17614843
#> 376 21 370 297.6046 11.01508227
#> 377 21 371 289.2423 5.81206077
#> 378 21 372 326.9194 8.12393754
#> 379 22 308 335.4401 17.06906435
#> 380 22 309 228.6792 5.76184862
#> 381 22 310 233.1068 5.34442146
#> 382 22 330 306.2613 4.70525344
#> 383 22 331 286.8394 8.78379339
#> 384 22 332 309.8780 8.32617880
#> 385 22 333 315.6480 7.30157261
#> 386 22 334 311.2674 11.48451071
#> 387 22 335 263.2429 4.25396931
#> 388 22 337 364.8506 14.12097595
#> 389 22 349 278.4798 11.12908082
#> 390 22 350 321.0889 16.52347082
#> 391 22 351 293.9199 3.01591558
#> 392 22 352 328.9795 12.76400485
#> 393 22 369 308.3297 8.28057004
#> 394 22 370 291.4620 14.23852732
#> 395 22 371 300.9596 9.85554263
#> 396 22 372 316.1865 13.22034658
#> 397 23 308 339.1286 19.48563033
#> 398 23 309 217.9621 3.54510354
#> 399 23 310 237.0612 8.68054660
#> 400 23 330 305.4788 4.36968380
#> 401 23 331 307.9708 8.57238279
#> 402 23 332 311.6822 8.33549364
#> 403 23 333 318.4000 6.79909859
#> 404 23 334 291.9974 5.13091475
#> 405 23 335 245.1525 0.12869603
#> 406 23 337 390.9388 11.56286528
#> 407 23 349 277.7066 10.40585306
#> 408 23 350 311.1509 19.40364667
#> 409 23 351 281.8974 4.07478896
#> 410 23 352 335.6503 8.44535995
#> 411 23 369 292.1709 9.61545287
#> 412 23 370 288.5924 11.93921185
#> 413 23 371 285.7155 9.01440017
#> 414 23 372 305.6664 14.26490706
#> 415 24 308 323.1753 16.61654824
#> 416 24 309 221.0979 1.07606373
#> 417 24 310 232.2864 4.69646244
#> 418 24 330 308.4107 5.60870160
#> 419 24 331 315.5031 7.15893941
#> 420 24 332 324.0506 9.40911649
#> 421 24 333 316.9993 5.27968343
#> 422 24 334 292.6491 14.19286087
#> 423 24 335 252.9418 4.24170293
#> 424 24 337 388.6829 14.99058566
#> 425 24 349 260.3096 15.64383211
#> 426 24 350 328.2453 20.08392092
#> 427 24 351 288.7174 7.67796166
#> 428 24 352 339.2961 14.09847349
#> 429 24 369 308.7808 9.37432269
#> 430 24 370 287.2833 17.09096612
#> 431 24 371 289.3059 10.65375384
#> 432 24 372 320.6715 9.64289999
#> 433 25 308 338.0348 19.13611613
#> 434 25 309 237.9016 2.12072611
#> 435 25 310 231.3158 11.33426735
#> 436 25 330 290.8301 7.61016046
#> 437 25 331 296.4592 7.61135718
#> 438 25 332 307.4752 8.79651914
#> 439 25 333 307.7714 10.31860976
#> 440 25 334 280.7865 12.04587102
#> 441 25 335 249.7664 3.22631187
#> 442 25 337 383.4059 15.09124944
#> 443 25 349 275.7748 13.05927419
#> 444 25 350 321.9059 20.65571833
#> 445 25 351 301.5994 6.42614275
#> 446 25 352 340.2362 15.13540819
#> 447 25 369 319.3925 8.79025093
#> 448 25 370 310.3025 18.05120279
#> 449 25 371 304.1258 6.64670421
#> 450 25 372 313.2080 10.84961504
#> 451 26 308 332.8784 16.22460179
#> 452 26 309 210.5569 4.60582431
#> 453 26 310 244.6384 5.51334190
#> 454 26 330 287.4886 6.02359011
#> 455 26 331 312.4973 7.62622536
#> 456 26 332 311.4584 2.06999940
#> 457 26 333 330.0975 10.11575520
#> 458 26 334 298.2656 8.24757064
#> 459 26 335 256.1615 -1.02599169
#> 460 26 337 372.4179 15.41236758
#> 461 26 349 276.8684 13.61776924
#> 462 26 350 311.8223 14.86063611
#> 463 26 351 283.4173 6.45196620
#> 464 26 352 330.6528 14.38197522
#> 465 26 369 305.9609 9.30203994
#> 466 26 370 295.1399 13.82424185
#> 467 26 371 310.4481 6.01782483
#> 468 26 372 314.9949 12.94718667
#> 469 27 308 334.3321 18.42169932
#> 470 27 309 222.9114 2.54494088
#> 471 27 310 242.6277 3.75722330
#> 472 27 330 300.8096 1.49123569
#> 473 27 331 301.7754 12.90266334
#> 474 27 332 304.4267 15.19007994
#> 475 27 333 312.9692 8.57559748
#> 476 27 334 290.9033 13.90369694
#> 477 27 335 244.7864 3.56227879
#> 478 27 337 365.1722 21.25432886
#> 479 27 349 270.1594 8.86177672
#> 480 27 350 324.7561 19.05771026
#> 481 27 351 282.4341 5.05924148
#> 482 27 352 335.0539 14.82855729
#> 483 27 369 303.3612 11.50123988
#> 484 27 370 293.0220 17.60454463
#> 485 27 371 285.3060 16.12457305
#> 486 27 372 305.6831 8.11509254
#> 487 28 308 338.7322 19.19570432
#> 488 28 309 221.0303 3.46231373
#> 489 28 310 233.5414 6.61443650
#> 490 28 330 309.7026 4.33282292
#> 491 28 331 307.3319 2.21757715
#> 492 28 332 292.8833 8.85363909
#> 493 28 333 308.1272 4.95794703
#> 494 28 334 281.4154 13.85324589
#> 495 28 335 253.5459 1.13558979
#> 496 28 337 359.2303 15.58295230
#> 497 28 349 276.5651 11.15979235
#> 498 28 350 314.5713 16.26407478
#> 499 28 351 302.6252 8.78385358
#> 500 28 352 333.1993 14.11121435
#> 501 28 369 302.7771 14.04890293
#> 502 28 370 301.0879 15.58714408
#> 503 28 371 293.5754 16.41206566
#> 504 28 372 318.1752 8.64965970
#>
#> $draws$popef.dispersion
#> (Intercept) Days_c
#> [1,] 1137.6249 23.38361
#> [2,] 1062.1152 27.25176
#> [3,] 947.2972 16.10941
#> [4,] 1064.6554 36.45147
#> [5,] 1126.5534 15.23363
#> [6,] 1318.9726 50.78292
#> [7,] 940.4095 25.73853
#> [8,] 1207.9720 46.20914
#> [9,] 1056.3420 21.48129
#> [10,] 1123.0794 33.24365
#> [11,] 1185.4269 32.60768
#> [12,] 1147.3777 38.07960
#> [13,] 1107.0121 31.77784
#> [14,] 895.2636 28.20770
#> [15,] 1176.4605 18.04722
#> [16,] 1419.9710 33.54756
#> [17,] 1587.0645 26.55614
#> [18,] 871.4135 42.14163
#> [19,] 1215.6899 29.67104
#> [20,] 1144.3437 34.79392
#> [21,] 2132.5297 22.78564
#> [22,] 882.5063 16.62091
#> [23,] 1555.0622 28.95669
#> [24,] 1109.5086 32.68807
#> [25,] 1375.5023 29.03299
#> [26,] 992.4001 19.25377
#> [27,] 1806.4240 29.36981
#> [28,] 1039.3785 22.31448
#>
#> $draws$popef.iters
#> (Intercept) Days_c
#> [1,] 13 15
#> [2,] 8 9
#> [3,] 11 19
#> [4,] 8 24
#> [5,] 23 25
#> [6,] 17 19
#> [7,] 37 13
#> [8,] 6 14
#> [9,] 6 11
#> [10,] 13 15
#> [11,] 19 8
#> [12,] 23 11
#> [13,] 4 23
#> [14,] 10 21
#> [15,] 16 11
#> [16,] 17 13
#> [17,] 8 25
#> [18,] 25 34
#> [19,] 12 8
#> [20,] 18 8
#> [21,] 39 19
#> [22,] 12 38
#> [23,] 30 15
#> [24,] 12 17
#> [25,] 13 4
#> [26,] 9 10
#> [27,] 27 10
#> [28,] 7 25
#>
#> $draws$groupef.iters
#> [1] 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
#>
#> $draws$popef.mode
#> $draws$popef.mode$`(Intercept)`
#> (Intercept)
#> 298.5079
#>
#> $draws$popef.mode$Days_c
#> (Intercept)
#> 10.37476
#>
#>
#> $draws$popef.init
#> $draws$popef.init$`(Intercept)`
#> (Intercept)
#> 298.5079
#>
#> $draws$popef.init$Days_c
#> (Intercept)
#> 10.37476If the pilot mean averages n_pilot independent chains,
its sampling distribution has covariance \(\Sigma/n_{\text{pilot}}\), so the squared
standardized start distance obeys \(D_0 =
\chi^2_{p,\alpha} / n_{\text{pilot}}\) at confidence level \(\alpha\). Bigger pilots give smaller \(D_0\) and therefore fewer main-stage sweeps
— but the pilot itself costs sweeps. The engine minimizes total cost
over n_pilot; the optimization is worked through in
vignette("Chapter-B03", package = "lmebayesCore").
Both stages leave a sweep history:
names(fit$sweep_history)
#> [1] "pilot" "main"
fit$sweep_history$pilot
#>
#> --- two-block [pilot stage summary: fixef by sweep (3 sweeps)] ---
#> Block 2 fixed effects (mode and chain stats after each sweep):
#> Random effect Covariate mode/sweep mean sd
#> --------------------------------------------------------------------------------------------
#> (Intercept) (Intercept) mode 298.5079
#> (Intercept) (Intercept) sweep 1 297.8788 10.0647
#> (Intercept) (Intercept) sweep 2 298.4804 7.7030
#> (Intercept) (Intercept) sweep 3 298.1352 7.5755
#> Days_c (Intercept) mode 10.3748
#> Days_c (Intercept) sweep 1 10.7468 1.4401
#> Days_c (Intercept) sweep 2 10.2859 1.2920
#> Days_c (Intercept) sweep 3 10.6915 1.6056
#>
#> Block 2 RE precision (1/tau^2; chain mean and SD after each sweep):
#> Random effect sweep mean sd
#> ------------------------------------------------------------
#> (Intercept) sweep 1 0.0008 0.0003
#> (Intercept) sweep 2 0.0008 0.0002
#> (Intercept) sweep 3 0.0009 0.0002
#> Days_c sweep 1 0.0392 0.0129
#> Days_c sweep 2 0.0372 0.0110
#> Days_c sweep 3 0.0378 0.0122fit$sweep_history$main
#>
#> --- two-block [main stage summary: fixef by sweep (3 sweeps)] ---
#> Block 2 fixed effects (mode and chain stats after each sweep):
#> Random effect Covariate mode/sweep mean sd
#> --------------------------------------------------------------------------------------------
#> (Intercept) (Intercept) mode 298.5079
#> (Intercept) (Intercept) sweep 1 298.6079 8.1706
#> (Intercept) (Intercept) sweep 2 299.0546 9.0954
#> (Intercept) (Intercept) sweep 3 298.3910 8.2584
#> Days_c (Intercept) mode 10.3748
#> Days_c (Intercept) sweep 1 10.5485 1.4322
#> Days_c (Intercept) sweep 2 10.3930 1.3979
#> Days_c (Intercept) sweep 3 10.3266 1.5139
#>
#> Block 2 RE precision (1/tau^2; chain mean and SD after each sweep):
#> Random effect sweep mean sd
#> ------------------------------------------------------------
#> (Intercept) sweep 1 0.0008 0.0002
#> (Intercept) sweep 2 0.0008 0.0002
#> (Intercept) sweep 3 0.0008 0.0002
#> Days_c sweep 1 0.0365 0.0125
#> Days_c sweep 2 0.0365 0.0121
#> Days_c sweep 3 0.0355 0.0120The pilot history starts at the ICM mode and moves; the main history starts at the pilot mean and barely moves. That contrast is the pilot doing its job.
summary(fit)
#> Call:
#> lmerb(formula = form, data = dat, pfamily_list = pf_ing, dispersion_ranef = ps$group.dispersion,
#> n = 300L, progbar = FALSE)
#>
#> Bayesian linear mixed model fit [300 draws, two-block Gibbs]
#> Formula: Reaction ~ Days_c + (1 + Days_c || Subject)
#>
#> Random effects (lmer reference; tau^2 sampled for non-dNormal components):
#> Groups Name Std.Dev.
#> Subject (Intercept) 37.533
#> Subject.1 Days_c 5.922
#> Residual 25.592
#> Number of obs: 180, groups: Subject, 18
#>
#> === Block 2 dispersion (RE variance tau^2_k) ===
#>
#> Prior and lmer reference:
#>
#> Prior 1/E[1/tau2] E[tau2] sqrt(E[tau2])
#> (Intercept) dIndependent_Normal_Gamma 931.9470 1346.1456 36.6899
#> Days_c dIndependent_Normal_Gamma 20.3346 29.3722 5.4196
#> Residual fixed 709.8260 709.8260 26.6426
#> disp_lower disp_upper lmer lmer SD
#> (Intercept) 632.3273 2998.2033 1408.7541 37.5334
#> Days_c 13.7971 65.4193 35.0716 5.9221
#> Residual NA NA 654.9399 25.5918
#>
#> Overview:
#> Post.Mode Post.Mean Post.Sd MC Error Mean SD Pr(Prior_tail)
#> (Intercept) 1346.1456 1391.4003 450.1460 25.9892 36.8477 0.4767
#> Days_c 29.3722 31.4504 10.4277 0.6020 5.5328 0.4867
#> Residual 709.8260 709.8260 NA NA 26.6426 NA
#> Cand/draw
#> (Intercept) 5.080
#> Days_c 4.783
#> Residual 1.000
#>
#> Distribution percentiles (tau^2):
#>
#> 1.0% 2.5% 5.0% Median 95.0% 97.5% 99.0%
#> (Intercept) 677.53 757.35 810.32 1328.48 2292.71 2506.01 2725.03
#> Days_c 14.84 15.49 16.84 29.57 49.38 56.05 58.96
#> Residual 709.83 709.83 709.83 709.83 709.83 709.83 709.83
#>
#> SD credible interval (sqrt(tau^2) draws; lmer SD for reference):
#>
#> 2.5% Median 97.5% lmer SD
#> (Intercept) 27.519 36.448 50.055 37.533
#> Days_c 3.935 5.437 7.486 5.922
#> Residual 26.643 26.643 26.643 25.592
#>
#> === Block 2: Level-2 fixed effects (hyperparameters) ===
#>
#> Prior and lmer reference:
#>
#> Prior Mean Prior.sd lmer lmer.se
#> (Intercept)::(Intercept) 298.51 90.05 298.51 9.050
#> Days_c::(Intercept) 0.00 15.38 10.47 1.546
#>
#> Overview:
#> fixef.mode fixef.means Post.Sd MC.Error Pr(Prior_tail)
#> (Intercept)::(Intercept) 298.508 298.391 8.258 0.477 0.457
#> Days_c::(Intercept) 10.375 10.327 1.514 0.087 <2e-16
#>
#> (Intercept)::(Intercept)
#> Days_c::(Intercept) ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> Distribution percentiles:
#>
#> 1.0% 2.5% 5.0% Median 95.0% 97.5% 99.0%
#> (Intercept)::(Intercept) 277.277 282.920 284.784 297.776 311.444 315.412 318.55
#> Days_c::(Intercept) 6.636 7.359 7.842 10.394 12.650 12.996 13.18
#>
#> === Block 1: Random effects (group-level) ===
#>
#> Summary of posterior mode (ranef.mode) across groups:
#>
#> Mean SD Min Q1 Median Q3 Max
#> (Intercept) 298.5079 36.4724 219.4041 290.8999 304.3215 314.6673 371.8535
#> Days_c 10.4463 5.0724 0.1222 7.1416 10.4194 12.8294 19.1842
#> MCMC.mean
#> (Intercept) 298.46
#> Days_c 10.43
#>
#> Per-group random effects: inspect fit$groupef.mode or fit$groupef,
#> or call summary(fit, groups = <level ids>) for selected groups.The Block 2 dispersion section is now informative.
Post.Mode, Post.Mean, Post.Sd,
and the percentile rows are real distributions rather than repeated
constants, disp_lower / disp_upper show the
truncation window, and Cand/draw reports how many
accept-reject candidates the envelope sampler needed per draw — a value
near 1 means the envelope is tight, larger values mean it is working
harder.
Sample them when the number of groups is small (a variance component from 18 subjects is barely identified), when the variance component is itself of interest, or when you are reporting intervals you want to be honest.
Fix them when there are many groups and \(\tau^2\) is well determined, when you want exact iid draws with no sweeps, or when speed matters and the sensitivity check shows it makes no difference.
The cheapest way to decide is to fit both, as above, and look at the ratio of posterior standard deviations.
vignette("Chapter-B05", package = "lmebayesCore") — the
ING Block 2 derivation and the truncation window.vignette("Chapter-B06", package = "lmebayesCore") — the
local rate diagnostic used when \(\tau^2\) is sampled.?pfamily_list,
?dIndependent_Normal_Gamma