How to set priors

Every model in neighbayes takes a priors= argument. This guide covers what the defaults are, how to override them, how the package stops you from overriding them wrongly, and the one case where leaving them alone actually changes an answer you might report.

The short version: the defaults are weakly informative and scaled to your data, and you should usually leave them there. Override when you have real prior information, when you need to bound the spatial parameter, or when a referee asks for a sensitivity analysis.

import contextlib
import io

import pandas as pd

from neighbayes.dgp import simulate_sar
from neighbayes.models import OLS, SAR, SLX
from neighbayes.models.priors import SARPriors, SEMPriors

data = simulate_sar(n=20, rho=0.5, seed=3)
y, X, W = data["y"], data["X"], data["W_graph"]
print(f"n = {len(y)}, k = {X.shape[1]}")
n = 400, k = 2

Three ways to pass priors

priors= accepts None, a plain dict, or a typed dataclass. All three end up in the same place — the dict and None forms are coerced to the dataclass at construction.

m_default = SAR(y=y, X=X, W=W)  # None: defaults
m_dict = SAR(y=y, X=X, W=W, priors={"rho_lower": 0.0})  # dict of overrides
m_typed = SAR(y=y, X=X, W=W, priors=SARPriors(rho_lower=0.0))  # typed dataclass

print(m_dict.priors_obj.rho_lower, m_typed.priors_obj.rho_lower)
0.0 0.0

Prefer the typed form in code you intend to keep. It gives you autocomplete on the supported keys, and a type checker catches a bad field before you run anything. The dict form is convenient in a notebook.

Inspect what a model resolved

Two public attributes. priors_obj is the dataclass; priors is a dict of the keys that carry a concrete value.

print(f"class : {type(m_default.priors_obj).__name__}")
print(f"object: {m_default.priors_obj}")
print()
print(f"dict  : {m_default.priors}")
class : SARPriors
object: SARPriors(beta_mu=None, beta_sigma=None, sigma2_alpha=2.0, sigma2_beta=None, sigma_sigma=10.0, nu=4.0, rho_lower=-1.0, rho_upper=1.0)

dict  : {'sigma2_alpha': 2.0, 'sigma_sigma': 10.0, 'nu': 4.0, 'rho_lower': -1.0, 'rho_upper': 1.0}

type(model.priors_obj) is also how you discover which dataclass a given model wants — there are 39 of them and the mapping follows the model hierarchy, so looking it up beats guessing.

What the defaults actually are

Three fields come back as None above: beta_mu, beta_sigma, and sigma2_beta. None here does not mean “flat” — it means resolved from the data at construction.

Coefficients: a data-scaled weakly-informative normal

Following Gelman et al. (2008), each coefficient gets a normal prior scaled to the response and to its own covariate:

Column

\(\mu\)

\(\sigma\)

intercept

\(\overline{y}\)

\(2.5\,\mathrm{sd}(y)\)

slope \(j\)

\(0\)

\(2.5\,\mathrm{sd}(y) / \mathrm{sd}(x_j)\)

Dividing by \(\mathrm{sd}(x_j)\) is what makes this behave sensibly regardless of the units your covariates happen to be in — a prior that is weak for income in dollars would be crushing for income in millions.

Noise variance: a scale-aware inverse gamma

\(\sigma^2 \sim \mathrm{InverseGamma}(\alpha, \beta)\) with \(\alpha = 2\) (finite mean, weakly informative) and \(\beta\) resolved to \(\mathrm{Var}(y)\), so the prior mean sits at roughly the observed variance. Following LeSage (2009), the conjugate inverse gamma is what keeps the \(\sigma^2\) Gibbs block in closed form — the Gaussian models use it for both NUTS and Gibbs so that the two paths target exactly the same posterior.

Spatial parameter: bounded uniform

\(\rho \sim \mathrm{Uniform}(\texttt{rho\_lower}, \texttt{rho\_upper})\), default \((-1, 1)\); SEMPriors spells the same fields lam_lower / lam_upper.

Overriding

Set a scalar and it broadcasts to every coefficient; set a vector of length \(p\) and it applies per coefficient, in design-matrix column order.

m_scalar = SAR(y=y, X=X, W=W, priors={"beta_sigma": 5.0})
m_vector = SAR(y=y, X=X, W=W, priors={"beta_mu": [0.0, 1.0], "beta_sigma": [10.0, 0.5]})

print("scalar :", m_scalar.priors_obj.beta_sigma)
print("vector :", m_vector.priors_obj.beta_mu, m_vector.priors_obj.beta_sigma)
scalar : 5.0
vector : [0.0, 1.0] [10.0, 0.5]

Bounding the spatial parameter

The common substantive override. If theory rules out negative spatial dependence — diffusion processes, contagion — say so, rather than letting the sampler spend draws in a region you would not believe anyway.

m_positive = SAR(y=y, X=X, W=W, priors={"rho_lower": 0.0})
idata_positive = m_positive.fit(
    draws=1000, tune=500, chains=4, random_seed=1, progressbar=False
)

print(f"rho posterior minimum: {float(idata_positive.posterior['rho'].min()):.4f}")
print("(the prior bound is a hard constraint, not a nudge)")
rho posterior minimum: 0.3603
(the prior bound is a hard constraint, not a nudge)
/home/runner/micromamba/envs/test/lib/python3.14/site-packages/neighbayes/_logdet/_jax.py:188: ComplexWarning: Casting complex values to real discards the imaginary part
  W_arr = np.asarray(W, dtype=np.float64)

Mistakes fail loudly

A misspelled key is silently ignored in most libraries. Here it raises at construction, with the allowed keys in the message.

try:
    SAR(y=y, X=X, W=W, priors={"rho_lo": 0.0})
except TypeError as err:
    print(f"TypeError: {err}")
TypeError: Unknown prior key(s) for SARPriors: ['rho_lo']. Allowed keys: ['beta_mu', 'beta_sigma', 'nu', 'rho_lower', 'rho_upper', 'sigma2_alpha', 'sigma2_beta', 'sigma_sigma'].
try:
    SAR(y=y, X=X, W=W, priors=SEMPriors())
except TypeError as err:
    print(f"TypeError: {err}")
TypeError: priors must be None, a dict, or a SARPriors instance; got SEMPriors.

The second case is why the typed form is worth using: SEMPriors and SARPriors differ only in whether the spatial bounds are named lam_* or rho_*, which is exactly the kind of thing that is easy to get wrong and hard to notice.

The case where the default matters

Warning

Wide priors are not “uninformative” for model comparison Setting beta_sigma to a large value to be “objective” reproduces the legacy near-improper prior, and it will change your Bayes factors by many nats. This is the Bartlett–Lindley paradox: the marginal likelihood integrates the likelihood against the prior, so spreading the prior over a wider region dilutes the evidence for the model carrying more coefficients — regardless of what the data say.

Compare a non-spatial OLS against SLX, which adds spatially lagged covariates, under the default prior and under a deliberately diffuse one:

FIT_KW = dict(draws=1000, tune=500, chains=4, random_seed=1, progressbar=False)


def log_bf_slx_vs_ols(priors):
    from neighbayes.diagnostics.bayesfactor import bayes_factor_compare_models

    with (
        contextlib.redirect_stdout(io.StringIO()),
        contextlib.redirect_stderr(io.StringIO()),
    ):
        ols = OLS(y=y, X=X, priors=priors)
        ols.fit(**FIT_KW)
        slx = SLX(y=y, X=X, W=W, priors=priors)
        slx.fit(**FIT_KW)
        table = bayes_factor_compare_models(
            [ols, slx], model_labels=["OLS", "SLX"], log=True
        )
    return float(table.loc["SLX", "OLS"])


rows = {
    "default (data-scaled)": log_bf_slx_vs_ols(None),
    "beta_sigma = 1e6": log_bf_slx_vs_ols({"beta_sigma": 1e6}),
}
pd.DataFrame({"log BF (SLX vs OLS)": rows}).round(3)
log BF (SLX vs OLS)
default (data-scaled) 29.196
beta_sigma = 1e6 17.803

Same data, same models, same sampler. The only thing that changed is how much prior mass sits in regions the data rule out — and the evidence swings by more than ten log units.

The direction is the tell. SLX carries the extra spatially-lagged covariates, so it is the model that pays when the prior is spread thin: widening beta_sigma dilutes the likelihood mass that its additional coefficients can claim, and its Bayes factor falls. Push the prior wide enough and the richer model loses on prior width alone, whatever the data say. That is the paradox, and it does not go away by making the prior wider still — it gets worse.

Neither number is “wrong”; they answer different questions. But the diffuse one is not the neutral choice it looks like. If you compare models with Bayes factors, keep the default and say so, or report the comparison across a range of prior widths.

Sensitivity as routine practice

Refitting across a grid of prior widths costs little and answers the question a reader will ask. Do it on the smallest dataset you plan to report, because that is where a prior can still move the answer — at \(n = 400\) above, nothing reasonable would.

small = simulate_sar(n=7, rho=0.5, seed=3)  # 7 x 7 grid = 49 observations
ys, Xs, Ws = small["y"], small["X"], small["W_graph"]

rows = []
for width in (0.25, 1.0, 5.0, None):
    label = "default" if width is None else f"beta_sigma={width:g}"
    priors = None if width is None else {"beta_sigma": width}
    idata = SAR(y=ys, X=Xs, W=Ws, priors=priors).fit(**FIT_KW)
    rows.append(
        pd.Series(
            {
                "rho mean": float(idata.posterior["rho"].mean()),
                "beta[0] mean": float(idata.posterior["beta"][..., 0].mean()),
                "beta[1] mean": float(idata.posterior["beta"][..., 1].mean()),
            },
            name=label,
        )
    )
pd.DataFrame(rows).round(3)
/home/runner/micromamba/envs/test/lib/python3.14/site-packages/neighbayes/_logdet/_jax.py:188: ComplexWarning: Casting complex values to real discards the imaginary part
  W_arr = np.asarray(W, dtype=np.float64)
/home/runner/micromamba/envs/test/lib/python3.14/site-packages/neighbayes/_logdet/_jax.py:188: ComplexWarning: Casting complex values to real discards the imaginary part
  W_arr = np.asarray(W, dtype=np.float64)
/home/runner/micromamba/envs/test/lib/python3.14/site-packages/neighbayes/_logdet/_jax.py:188: ComplexWarning: Casting complex values to real discards the imaginary part
  W_arr = np.asarray(W, dtype=np.float64)
/home/runner/micromamba/envs/test/lib/python3.14/site-packages/neighbayes/_logdet/_jax.py:188: ComplexWarning: Casting complex values to real discards the imaginary part
  W_arr = np.asarray(W, dtype=np.float64)
rho mean beta[0] mean beta[1] mean
beta_sigma=0.25 0.459 1.229 1.315
beta_sigma=1 0.459 0.937 2.058
beta_sigma=5 0.459 0.917 2.106
default 0.459 0.917 2.107

Two things to read off this table.

The coefficients shrink when the prior is tight and then stop. At beta_sigma=0.25 the slope is pulled hard toward zero — the prior is narrow relative to the scale of the data. By beta_sigma=5 the posterior has stopped caring, and the default sits on that plateau. Reaching the plateau is what tells you an estimate is data-driven rather than an artefact of the prior; it is the thing worth reporting, and it is cheap to check.

rho does not move at all. beta_sigma is the prior on the coefficients; the spatial parameter has its own, and shrinking \(\beta\) does not touch it. If you want to probe sensitivity in \(\rho\), vary rho_lower / rho_upper — varying the wrong hyperparameter and finding nothing is not a sensitivity analysis.

See also