fit_bayesian

fit_bayesian(
    model,
    distribution=None,
    sampler='DEMCz',
    chains=4,
    iterations=3000,
    warmup=None,
    output_length=10000,
    thinning_interval=-1,
    seed=12345,
    point_estimator=None,
    credible_level=0.9,
    **knobs,
)

Bayesian MCMC fit.

Fit a model with a Bayesian MCMC analysis and return a fit carrying the raw chains, the posterior summary (mean/median/sd/credible interval, R-hat, effective sample size), and the usual Bayesian goodness-of-fit scalars. Wraps the shared C++ BayesianAnalysis ported from USACE-RMC RMC.BestFit – the same estimator :func:~corehydropy.univariate_analysis and :func:~corehydropy.estimation_diagnostics build on.

Parameters

Name Type Description Default
model Model or array_like See :func:fit_mle. required
distribution str See :func:fit_mle. None
sampler ('DEMCz', 'DEMCzs', 'ARWMH', 'NUTS') MCMC sampler. BayesianAnalysis can only construct these four; RWMH, HMC, and SNIS are real MCMC samplers but need :func:~corehydropy.mcmc_sample instead. "DEMCz"
chains int Number of parallel Markov chains. 4
iterations int Number of post-warmup MCMC iterations, per chain. 3000
warmup int Number of warmup (burn-in) iterations. Defaults to max(50, iterations // 2) when omitted: BayesianAnalysis’s own class default (1500) would otherwise silently trip its sampler’s warmup <= iterations / 2 guard for any iterations below about 3000. None
output_length int Number of posterior draws retained (thinned down from the raw chains) for the summary and any downstream uncertainty quantification. 10000
thinning_interval int MCMC thinning interval; -1 keeps the sampler’s own default. -1
seed int PRNG seed for the sampler (fixed for reproducibility – a seeded call returns identical draws in R and Python). 12345
point_estimator ('PosteriorMean', 'PosteriorMode') Which posterior summary .parameters reports. None (default) leaves BayesianAnalysis’s own class default, which is "PosteriorMean" – so by default .parameters is bit-identical to .posterior_summary["mean"]. Pass "PosteriorMode" to report the MAP point instead. "PosteriorMean"
credible_level float Width of the posterior credible interval .posterior_summary’s lower/upper report, between 0 and 1. The default is BayesianAnalysis’s own class default. :meth:Fit.confint reuses these bounds when its level matches and otherwise re-runs the identical seeded chain, so setting it here is how a bare confint() avoids that rebuild. 0.9
**knobs Sampler-specific tuning knobs. Passing a knob the chosen sampler does not use is an error: "DEMCz" accepts jump, jump_threshold, noise; "DEMCzs" additionally accepts snooker_threshold; "ARWMH" accepts scale, beta; "NUTS" accepts max_tree_depth. {}

Returns

Name Type Description
Fit A fit with .method == "BayesianAnalysis". .draws is a 3-D array [iteration, chain, parameter]; .posterior is the thinned draw matrix the analyses consume, (.posterior_rows, len(.parameters)), its columns aligned with .parameter_names. .posterior_summary has one entry per metric (mean, median, sd, lower, upper, rhat, ess), each an array aligned with .parameter_names. .map and .posterior_mean are the two posterior point estimates. .mean_log_likelihood has one entry per iteration. .acceptance_rates has one entry per chain. .warmup records the warmup actually used and .credible_level the width .posterior_summary’s bounds were computed at. .dic, .waic, .looic are the usual Bayesian goodness-of-fit scalars, .looic_se the standard error of .looic, and .waic_pd/.loo_pd their effective parameter counts. See :func:fit_diagnostics for leverage/influence diagnostics off a Bayesian fit.

Raises

Name Type Description
ValueError When a setting falls outside the range BayesianAnalysis accepts: chains between 4 and 20, iterations between 100 and 1,000,000, warmup between 50 and 100,000, output_length between 100 and 1,000,000, seed not negative, and credible_level strictly between 0 and 1.

See Also

fit_mle, fit_map, fit_gmm, fit_diagnostics, univariate_analysis

Examples

>>> from corehydropy import fit_bayesian, model_univariate
>>> peaks = [12500, 15300, 8900, 22100, 18700, 14200, 9800, 28500, 17400, 11600]
>>> f = fit_bayesian(model_univariate("Normal", peaks), sampler="DEMCz", iterations=200,
...                  output_length=500, seed=12345)
>>> f.posterior_summary["mean"]