The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.
A typical Bayesian SAE workflow with hbsaems v1.1.0
follows seven steps, each backed by a single primary function. This
vignette walks through the canonical pipeline.
For deeper coverage of topics beyond this introductory workflow, see the articles on the package website:
https://madsyair.github.io/hbsaems/articles/
The articles cover:
data("data_fhnorm")
str(data_fhnorm, max.level = 1)
#> 'data.frame': 100 obs. of 9 variables:
#> $ y : num 6.98 9.01 6.22 10.47 10.62 ...
#> $ D : num 0.921 2.643 1.704 0.659 1.017 ...
#> $ x1 : num -0.9026 -0.7527 -0.2606 0.0759 0.2118 ...
#> $ x2 : num 0.712 -0.878 -0.607 -1.077 -0.502 ...
#> $ x3 : num 1.431 -0.318 -1.029 -1.071 -0.496 ...
#> $ theta_true: num 8.5 9.73 8.54 9.68 11.86 ...
#> $ u : num -0.8531 -0.0145 -1.2505 -0.5989 1.5845 ...
#> $ regency : chr "regency_001" "regency_002" "regency_003" "regency_004" ...
#> $ province : chr "province_01" "province_01" "province_01" "province_01" ...
#> - attr(*, "datalabel")= chr "data_fhnorm"
#> - attr(*, "var.labels")= chr [1:9] "" "" "" "" ...
head(data_fhnorm[, c("regency", "province", "y", "D", "x1", "x2", "x3")], 4)
#> regency province y D x1 x2 x3
#> 1 regency_001 province_01 6.982249 0.9213655 -0.9026345 0.7116302 1.4305563
#> 2 regency_002 province_01 9.005801 2.6432183 -0.7526745 -0.8779797 -0.3182451
#> 3 regency_003 province_01 6.216509 1.7035092 -0.2605613 -0.6066989 -1.0293040
#> 4 regency_004 province_01 10.469687 0.6586857 0.0759456 -1.0766939 -1.0711613data_fhnorm is a simulated Fay-Herriot dataset: 100
regencies nested within 5 provinces, with a known sampling variance
D per area and three covariates.
# This is the production call. Replace `chains`, `iter` with your
# usual values; the lighter settings below are a quick demonstration.
model_prior <- hbm(
formula = brms::bf(y ~ x1 + x2 + x3),
data = data_fhnorm,
re = ~ (1 | regency), # area-level random effect (Fay-Herriot)
sampling_variance = "D", # KNOWN sampling variance D_i
sample_prior = "only",
prior = c(
brms::prior(normal(0, 1), class = "b"),
brms::prior(normal(10, 5), class = "Intercept"),
brms::prior(normal(0, 2), class = "sd")
),
chains = 2, iter = 1000
)
pc <- prior_check(model_prior,
data = data_fhnorm,
response_var = "y")
plot(pc)The sampling_variance = "D" argument is the
defining feature of a Fay-Herriot model: the sampling
variance is treated as known from the design (it is not
estimated), so brms pins the observation-level to . Omitting this
argument makes the model unidentified because the residual and the area
random-effect would compete to explain the same variance, typically
producing divergent transitions and very wide credible intervals.
model <- hbm(
formula = brms::bf(y ~ x1 + x2 + x3),
data = data_fhnorm,
re = ~ (1 | regency),
sampling_variance = "D",
control = list(adapt_delta = 0.99),
chains = 4, iter = 4000, warmup = 2000, cores = 4
)
summary(model)Two settings deserve attention:
sampling_variance = "D" – the column
D in data_fhnorm holds the known sampling
variance for each area; hbsaems translates this into an
offset on the observation-level standard deviation so brms / Stan never
tries to estimate it.adapt_delta = 0.99 – Fay-Herriot
likelihoods can produce mild funnel geometry between and the area random
effects when many areas have small . As of v1.1.0 hbsaems
already raises the default from brms’s 0.8 to
0.95; pushing it further to 0.99 buys an even
more conservative leapfrog step and eliminates the occasional divergent
transition without re-parameterising the model.For the remainder of the vignette we use a tiny model (very few iterations) so that the diagnostic and prediction chunks below have a real object to operate on. In your own analysis, use the full production settings shown above, not the toy settings here.
# Mini fit -- iter = 200, chains = 1 -- for vignette demonstration only.
# Do NOT use these settings for inference: the chains have not
# converged at this length and the posterior will be biased.
fit_demo <- suppressWarnings(
hbm(
formula = brms::bf(y ~ x1 + x2 + x3),
data = data_fhnorm,
re = ~ (1 | regency),
sampling_variance = "D",
chains = 1,
iter = 200,
warmup = 100,
refresh = 0,
seed = 1
)
)# Operate on the mini fit_demo above (NOT a substitute for production
# diagnostics on full chains).
diag <- convergence_check(fit_demo)
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> No divergences to plot.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
is_converged(fit_demo)
#> [1] FALSE
summary(diag)
#>
#> ===== Convergence Diagnostics Summary =====
#>
#> Divergent transitions: 0 (0.00%)
#> E-BFMI (min over chains): 0.930
#> Max-treedepth hit rate: 0.00%
#>
#> R-hat:
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 0.9904 0.9961 1.0030 1.0126 1.0175 1.1435
#>
#> Bulk ESS:
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 43.85 131.49 171.06 160.25 200.00 200.00
#>
#> Tail ESS:
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 21.15 64.12 76.15 75.61 90.20 116.78
#>
#> Geweke test (Z-scores):
#> [[1]]
#>
#> Fraction in 1st window = 0.1
#> Fraction in 2nd window = 0.5
#>
#> b_Intercept b_x1
#> 0.439675 0.675512
#> b_x2 b_x3
#> -0.110359 -1.329367
#> sd_regency__Intercept Intercept
#> -0.351228 0.518596
#> r_regency[regency_001,Intercept] r_regency[regency_002,Intercept]
#> 1.482982 0.578184
#> r_regency[regency_003,Intercept] r_regency[regency_004,Intercept]
#> 0.064829 -1.694880
#> r_regency[regency_005,Intercept] r_regency[regency_006,Intercept]
#> -0.086819 -1.550906
#> r_regency[regency_007,Intercept] r_regency[regency_008,Intercept]
#> -0.225457 -1.634221
#> r_regency[regency_009,Intercept] r_regency[regency_010,Intercept]
#> -0.573915 0.352415
#> r_regency[regency_011,Intercept] r_regency[regency_012,Intercept]
#> 0.329977 -1.152600
#> r_regency[regency_013,Intercept] r_regency[regency_014,Intercept]
#> -0.184377 0.994029
#> r_regency[regency_015,Intercept] r_regency[regency_016,Intercept]
#> -0.734035 -0.459930
#> r_regency[regency_017,Intercept] r_regency[regency_018,Intercept]
#> 1.866025 -0.857966
#> r_regency[regency_019,Intercept] r_regency[regency_020,Intercept]
#> -0.910251 -0.124925
#> r_regency[regency_021,Intercept] r_regency[regency_022,Intercept]
#> -0.638350 -1.228459
#> r_regency[regency_023,Intercept] r_regency[regency_024,Intercept]
#> -2.227878 0.009446
#> r_regency[regency_025,Intercept] r_regency[regency_026,Intercept]
#> -0.570614 0.706367
#> r_regency[regency_027,Intercept] r_regency[regency_028,Intercept]
#> -0.153144 0.666069
#> r_regency[regency_029,Intercept] r_regency[regency_030,Intercept]
#> -0.675936 0.019550
#> r_regency[regency_031,Intercept] r_regency[regency_032,Intercept]
#> 0.157159 0.068315
#> r_regency[regency_033,Intercept] r_regency[regency_034,Intercept]
#> -1.454501 -1.038916
#> r_regency[regency_035,Intercept] r_regency[regency_036,Intercept]
#> -1.488091 -1.626934
#> r_regency[regency_037,Intercept] r_regency[regency_038,Intercept]
#> 0.353194 0.786427
#> r_regency[regency_039,Intercept] r_regency[regency_040,Intercept]
#> -1.751016 0.171068
#> r_regency[regency_041,Intercept] r_regency[regency_042,Intercept]
#> 0.223230 -1.462420
#> r_regency[regency_043,Intercept] r_regency[regency_044,Intercept]
#> -1.376935 -0.370675
#> r_regency[regency_045,Intercept] r_regency[regency_046,Intercept]
#> -1.063761 -0.341718
#> r_regency[regency_047,Intercept] r_regency[regency_048,Intercept]
#> 0.665556 -1.282668
#> r_regency[regency_049,Intercept] r_regency[regency_050,Intercept]
#> -1.338728 -0.110090
#> r_regency[regency_051,Intercept] r_regency[regency_052,Intercept]
#> -0.292501 -0.653364
#> r_regency[regency_053,Intercept] r_regency[regency_054,Intercept]
#> 0.761260 -0.265161
#> r_regency[regency_055,Intercept] r_regency[regency_056,Intercept]
#> -0.781514 -1.031518
#> r_regency[regency_057,Intercept] r_regency[regency_058,Intercept]
#> -1.326952 -0.026405
#> r_regency[regency_059,Intercept] r_regency[regency_060,Intercept]
#> 0.556657 1.998405
#> r_regency[regency_061,Intercept] r_regency[regency_062,Intercept]
#> -0.076548 -0.893495
#> r_regency[regency_063,Intercept] r_regency[regency_064,Intercept]
#> 1.830265 1.063698
#> r_regency[regency_065,Intercept] r_regency[regency_066,Intercept]
#> -1.320349 -2.125862
#> r_regency[regency_067,Intercept] r_regency[regency_068,Intercept]
#> 0.086755 0.826772
#> r_regency[regency_069,Intercept] r_regency[regency_070,Intercept]
#> -2.425150 -1.021848
#> r_regency[regency_071,Intercept] r_regency[regency_072,Intercept]
#> -0.921946 -0.061025
#> r_regency[regency_073,Intercept] r_regency[regency_074,Intercept]
#> 1.882800 -2.702286
#> r_regency[regency_075,Intercept] r_regency[regency_076,Intercept]
#> -1.796570 0.426268
#> r_regency[regency_077,Intercept] r_regency[regency_078,Intercept]
#> 1.079833 0.048222
#> r_regency[regency_079,Intercept] r_regency[regency_080,Intercept]
#> -0.276625 -0.434901
#> r_regency[regency_081,Intercept] r_regency[regency_082,Intercept]
#> 0.373386 -1.262929
#> r_regency[regency_083,Intercept] r_regency[regency_084,Intercept]
#> -1.747763 0.778250
#> r_regency[regency_085,Intercept] r_regency[regency_086,Intercept]
#> 0.027099 0.458924
#> r_regency[regency_087,Intercept] r_regency[regency_088,Intercept]
#> 1.064617 -1.144627
#> r_regency[regency_089,Intercept] r_regency[regency_090,Intercept]
#> -0.268693 0.005188
#> r_regency[regency_091,Intercept] r_regency[regency_092,Intercept]
#> 0.625051 0.007719
#> r_regency[regency_093,Intercept] r_regency[regency_094,Intercept]
#> -1.305631 -0.166948
#> r_regency[regency_095,Intercept] r_regency[regency_096,Intercept]
#> -0.715510 0.749165
#> r_regency[regency_097,Intercept] r_regency[regency_098,Intercept]
#> -1.672694 -2.087388
#> r_regency[regency_099,Intercept] r_regency[regency_100,Intercept]
#> -0.694299 -1.305847
#> lprior lp__
#> -0.052890 -1.287801
#> z_1[1,1] z_1[1,2]
#> 1.025456 0.657983
#> z_1[1,3] z_1[1,4]
#> -0.007537 -1.689116
#> z_1[1,5] z_1[1,6]
#> -0.185234 -1.505061
#> z_1[1,7] z_1[1,8]
#> -0.189963 -1.297634
#> z_1[1,9] z_1[1,10]
#> -0.611099 0.524402
#> z_1[1,11] z_1[1,12]
#> 0.250522 -1.047651
#> z_1[1,13] z_1[1,14]
#> -0.179149 1.020036
#> z_1[1,15] z_1[1,16]
#> -0.693134 -0.585218
#> z_1[1,17] z_1[1,18]
#> 2.100592 -0.860425
#> z_1[1,19] z_1[1,20]
#> -1.204257 -0.143581
#> z_1[1,21] z_1[1,22]
#> -0.639256 -1.453809
#> z_1[1,23] z_1[1,24]
#> -2.248067 -0.017670
#> z_1[1,25] z_1[1,26]
#> -0.369461 0.547181
#> z_1[1,27] z_1[1,28]
#> -0.231902 0.524318
#> z_1[1,29] z_1[1,30]
#> -0.648019 -0.013206
#> z_1[1,31] z_1[1,32]
#> 0.128067 0.077570
#> z_1[1,33] z_1[1,34]
#> -1.271948 -0.932256
#> z_1[1,35] z_1[1,36]
#> -1.468568 -1.599269
#> z_1[1,37] z_1[1,38]
#> 0.301756 0.989599
#> z_1[1,39] z_1[1,40]
#> -1.740444 0.216227
#> z_1[1,41] z_1[1,42]
#> 0.416618 -1.274363
#> z_1[1,43] z_1[1,44]
#> -1.133581 -0.041793
#> z_1[1,45] z_1[1,46]
#> -1.152071 -0.337986
#> z_1[1,47] z_1[1,48]
#> 0.762589 -1.349635
#> z_1[1,49] z_1[1,50]
#> -1.280454 -0.225432
#> z_1[1,51] z_1[1,52]
#> -0.303364 -0.562894
#> z_1[1,53] z_1[1,54]
#> 0.639453 -0.043853
#> z_1[1,55] z_1[1,56]
#> -0.564568 -1.157184
#> z_1[1,57] z_1[1,58]
#> -1.152060 -0.099215
#> z_1[1,59] z_1[1,60]
#> 0.935305 1.811214
#> z_1[1,61] z_1[1,62]
#> -0.092633 -0.761527
#> z_1[1,63] z_1[1,64]
#> 1.720061 1.020080
#> z_1[1,65] z_1[1,66]
#> -1.753982 -2.212438
#> z_1[1,67] z_1[1,68]
#> 0.202992 0.872749
#> z_1[1,69] z_1[1,70]
#> -2.673858 -1.243029
#> z_1[1,71] z_1[1,72]
#> -0.910726 -0.115829
#> z_1[1,73] z_1[1,74]
#> 1.713765 -2.525547
#> z_1[1,75] z_1[1,76]
#> -1.877583 0.737048
#> z_1[1,77] z_1[1,78]
#> 1.104834 -0.037961
#> z_1[1,79] z_1[1,80]
#> -0.308535 -0.377443
#> z_1[1,81] z_1[1,82]
#> 0.281008 -1.239269
#> z_1[1,83] z_1[1,84]
#> -1.444329 0.529094
#> z_1[1,85] z_1[1,86]
#> 0.062075 0.368325
#> z_1[1,87] z_1[1,88]
#> 1.414949 -1.230895
#> z_1[1,89] z_1[1,90]
#> -0.334978 -0.094349
#> z_1[1,91] z_1[1,92]
#> 0.474061 -0.059493
#> z_1[1,93] z_1[1,94]
#> -1.482723 0.052240
#> z_1[1,95] z_1[1,96]
#> -0.933766 0.583424
#> z_1[1,97] z_1[1,98]
#> -1.645358 -2.388447
#> z_1[1,99] z_1[1,100]
#> -0.720579 -1.331146
hbm_warnings(fit_demo)
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> Warning: The ESS has been capped to avoid unstable estimates.
#> [1] "R-hat > 1.1: model may not have converged."The full set of diagnostics on your production fit would also include trace plots, R-hat / ESS tables, and pp-checks. These all follow the same calling convention:
model2 <- hbm(brms::bf(y ~ x1 + x2), data = data_fhnorm,
re = ~ (1 | regency),
sampling_variance = "D",
control = list(adapt_delta = 0.99),
chains = 4, iter = 4000)
model3 <- hbm(brms::bf(y ~ x1), data = data_fhnorm,
re = ~ (1 | regency),
sampling_variance = "D",
control = list(adapt_delta = 0.99),
chains = 4, iter = 4000)
model_compare(model, model2)
model_compare_all(full = model, medium = model2, simple = model3)The simplest call uses the data the model was fit on – no
newdata argument needed:
# In-sample prediction (default): operates on the data stored
# inside the fitted object.
est <- sae_predict(fit_demo)
summary(est)
#>
#> ===== Small Area Estimation Summary =====
#>
#> Areas : 100
#> Overall RSE : 7.2 %
#>
#> Predictions:
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 6.336 8.938 9.769 9.855 10.637 13.537
#>
#> RSE by area:
#> Min. 1st Qu. Median Mean 3rd Qu. Max.
#> 4.416 5.981 7.062 7.196 8.171 11.258
head(est$result_table, 5)
#> Prediction SD RSE_percent
#> 1 8.131821 0.6843505 8.415711
#> 2 9.340465 0.7375581 7.896375
#> 3 8.286397 0.7710245 9.304701
#> 4 10.378301 0.7482991 7.210228
#> 5 10.442262 0.7160301 6.857041To predict at fresh data points, pass newdata. When you
do, any internal offset columns the original fit needed
(e.g. .hbsaems_sigma_fixed = sqrt(D) for the Fay-Herriot
sugar) are repopulated automatically when the new data frame has the
same number of rows as the training data – the standard “predict at the
same areas with updated covariates” use case.
For model averaging, fit several candidate models, then either let
model_average() weight them by stacking weights (LOO-based)
or supply your own weights:
# Stacking weights (default, derived from LOO)
avg_stacked <- model_average(model, model2, model3)
# User-specified weights -- must sum to 1
avg_manual <- model_average(model, model2, weights = c(0.7, 0.3))
# The returned object is an `hbsae_results`, identical in shape to
# what sae_predict() returns, so all post-processing helpers
# (sae_transform, sae_scale, sae_filter, sae_aggregate) work on it.
plot(avg_stacked, type = "predictions")All sae_* post-processing helpers operate on the
hbsae_results object returned by either
sae_predict() or model_average():
sessionInfo()
#> R version 4.6.0 (2026-04-24)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.3 LTS
#>
#> Matrix products: default
#> BLAS: /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
#> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.26.so; LAPACK version 3.12.0
#>
#> locale:
#> [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
#> [3] LC_TIME=en_US.UTF-8 LC_COLLATE=C
#> [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
#> [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
#> [9] LC_ADDRESS=C LC_TELEPHONE=C
#> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
#>
#> time zone: Asia/Jakarta
#> tzcode source: system (glibc)
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] hbsaems_1.1.0
#>
#> loaded via a namespace (and not attached):
#> [1] tidyselect_1.2.1 dplyr_1.2.1 farver_2.1.2
#> [4] loo_2.9.0 S7_0.2.2 fastmap_1.2.0
#> [7] TH.data_1.1-3 tensorA_0.36.2.1 rpart_4.1.27
#> [10] digest_0.6.39 estimability_1.5.1 lifecycle_1.0.5
#> [13] StanHeaders_2.32.10 survival_3.8-6 processx_3.9.0
#> [16] magrittr_2.0.5 posterior_1.7.0 compiler_4.6.0
#> [19] rlang_1.2.0 sass_0.4.10 tools_4.6.0
#> [22] yaml_2.3.12 knitr_1.51 bridgesampling_1.2-1
#> [25] pkgbuild_1.4.8 curl_7.1.0 plyr_1.8.9
#> [28] RColorBrewer_1.1-3 multcomp_1.4-28 abind_1.4-8
#> [31] withr_3.0.2 purrr_1.2.2 nnet_7.3-20
#> [34] grid_4.6.0 stats4_4.6.0 jomo_2.7-6
#> [37] xtable_1.8-8 mice_3.19.0 inline_0.3.21
#> [40] ggplot2_4.0.3 emmeans_1.11.2 scales_1.4.0
#> [43] iterators_1.0.14 MASS_7.3-65 dichromat_2.0-0.1
#> [46] cli_3.6.6 mvtnorm_1.3-7 rmarkdown_2.31
#> [49] reformulas_0.4.4 generics_0.1.4 otel_0.2.0
#> [52] RcppParallel_5.1.11-2 rstudioapi_0.17.1 reshape2_1.4.5
#> [55] minqa_1.2.8 cachem_1.1.0 rstan_2.32.7
#> [58] stringr_1.6.0 splines_4.6.0 bayesplot_1.15.0
#> [61] parallel_4.6.0 matrixStats_1.5.0 brms_2.23.0
#> [64] vctrs_0.7.3 boot_1.3-32 V8_8.2.0
#> [67] glmnet_5.0 Matrix_1.7-5 sandwich_3.1-1
#> [70] jsonlite_2.0.0 callr_3.7.6 mitml_0.4-5
#> [73] foreach_1.5.2 jquerylib_0.1.4 tidyr_1.3.2
#> [76] glue_1.8.1 pan_1.9 nloptr_2.2.1
#> [79] codetools_0.2-20 ps_1.9.3 distributional_0.7.0
#> [82] stringi_1.8.7 gtable_0.3.6 shape_1.4.6.1
#> [85] QuickJSR_1.9.2 lme4_2.0-1 tibble_3.3.1
#> [88] pillar_1.11.1 htmltools_0.5.9 Brobdingnag_1.2-9
#> [91] R6_2.6.1 Rdpack_2.6.6 evaluate_1.0.5
#> [94] lattice_0.22-9 rbibutils_2.4.1 backports_1.5.1
#> [97] broom_1.0.13 bslib_0.10.0 rstantools_2.6.0
#> [100] Rcpp_1.1.1-1.1 coda_0.19-4.1 gridExtra_2.3
#> [103] nlme_3.1-169 checkmate_2.3.4 xfun_0.57
#> [106] zoo_1.8-14 pkgconfig_2.0.3These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.