neighbayes.models.SARProbit

class neighbayes.models.SARProbit(formula=None, data=None, y=None, X=None, W=None, region_col=None, region_ids=None, mobs=None, priors=None, robust=False)[source]

Bayesian spatial probit with regional random effects.

A binary-response model in which the latent utility includes a spatially autoregressive regional random effect. Each observation \(i\) belongs to one of \(m\) regions; region-level effects \(a\) follow a SAR process on the region-level weights matrix \(W\), and observation-level disturbances are standard Normal (probit link).

\[y_i = \mathbb{1}[z_i > 0],\quad z_i = x_i'\beta + a_{r(i)} + \varepsilon_i,\quad \varepsilon_i \sim \mathcal{N}(0, 1),\]
\[a = \rho W a + u,\quad u \sim \mathcal{N}(0, \sigma_a^2 I_m),\]

so that \(a \sim \mathcal{N}(0, \sigma_a^2 (I_m - \rho W)^{-1} (I_m - \rho W)^{-T})\). The marginal choice probability is \(P(y_i = 1 \mid \beta, a) = \Phi(x_i'\beta + a_{r(i)})\).

Parameters:
formula : str, optional

Formula for the binary response model, e.g. "y ~ x1 + x2". Requires data and region_col.

data : pandas.DataFrame, optional

Data source used with formula mode.

y : array-like, optional

Binary dependent variable (0/1), required in matrix mode.

X : array-like or pandas.DataFrame, optional

Covariate matrix, required in matrix mode.

W : libpysal.graph.Graph or scipy.sparse matrix

Region-level m x m spatial weights matrix. Accepts a libpysal.graph.Graph (the modern libpysal graph API) or any scipy.sparse matrix. The legacy libpysal.weights.W object is not accepted directly; pass w.sparse or convert with libpysal.graph.Graph.from_W(w). W should be row-standardized; a UserWarning is raised if not.

region_col : str, optional

Region identifier column in data (formula mode).

region_ids : array-like, optional

Region identifier per observation (matrix mode).

mobs : array-like, optional

Region observation counts (m,) in sorted region order (matrix mode alternative to region_ids).

priors : dict, optional

Override default priors. Supported keys:

  • rho_lower (float, default -0.95): Lower bound of the Uniform prior on \(\rho\).

  • rho_upper (float, default 0.95): Upper bound of the Uniform prior on \(\rho\).

  • beta_mu (float, default 0.0): Normal prior mean for \(\beta\).

  • beta_sigma (float, default 1e6): Normal prior std for \(\beta\).

  • sigma_a_sigma (float, default 10.0): HalfNormal scale for the regional random-effect std \(\sigma_a\).

robust : bool, default False

Not supported. The probit link uses a Normal CDF; a Student-t analogue is not implemented. Setting robust=True raises.

Notes

This class follows the core semip_g structure (binary response with spatially dependent regional effects). It uses a standard probit link with unit observation-level variance and does not currently sample the v_i/ r heteroskedastic hierarchy from legacy semip_g.

Robust regression

robust=True is not supported for SARProbit. The probit link function uses a Normal CDF; a robust version would require a Student-t CDF link, which is not yet implemented. Use robust=True with Gaussian models (OLS, SAR, SEM, etc.) instead.

__init__(formula=None, data=None, y=None, X=None, W=None, region_col=None, region_ids=None, mobs=None, priors=None, robust=False)[source]

Methods

__init__([formula, data, y, X, W, ...])

fit([draws, tune, chains, random_seed, ...])

Draw samples from the posterior.

fitted_probabilities()

Return posterior mean fitted probabilities for observed data.

fitted_values()

Return fitted values at posterior mean parameters.

random_effects_mean()

Return posterior mean regional effects.

residuals()

Return residuals y - fitted_values.

spatial_diagnostics()

Run Bayesian LM specification tests and return a summary table.

spatial_diagnostics_decision([alpha, format])

Return a model-selection decision from Bayesian LM test results.

spatial_effects([return_posterior_samples])

Compute average marginal effects (AME) for the spatial probit model.

summary([var_names])

Return posterior summary table.

Attributes

inference_data

Return the ArviZ InferenceData from the most recent fit.

pymc_model

Return the PyMC model object built for the most recent fit.

fit(draws=2000, tune=1000, chains=4, random_seed=None, progressbar=True, **sample_kwargs)[source]

Draw samples from the posterior.

fitted_probabilities()[source]

Return posterior mean fitted probabilities for observed data.

fitted_values()[source]

Return fitted values at posterior mean parameters.

Returns:

Posterior-mean fitted values (on the model’s native scale; fixed-effects-transformed for panel models).

Return type:

np.ndarray

property inference_data : arviz.data.inference_data.InferenceData | None[source]

Return the ArviZ InferenceData from the most recent fit.

Returns:

The inference data object, or None if the model has not been fit yet.

Return type:

arviz.InferenceData or None

property pymc_model : pymc.model.core.Model | None[source]

Return the PyMC model object built for the most recent fit.

Returns:

The model object used by fit(), or None if the instance has not been fit yet.

Return type:

pymc.Model or None

random_effects_mean()[source]

Return posterior mean regional effects.

residuals()[source]

Return residuals y - fitted_values.

Returns:

Residual vector y - fitted_values on the same scale as fitted_values().

Return type:

np.ndarray

spatial_diagnostics()[source]

Run Bayesian LM specification tests and return a summary table.

Looks up the diagnostic suite registered for this model class and calls each test function on this fitted model, collecting the results into a tidy DataFrame. The set of tests depends on the model type — for example, an OLS model runs LM-Lag, LM-Error, LM-SDM-Joint, and LM-SLX-Error-Joint, while an SAR model runs LM-Error, LM-WX, and Robust-LM-WX. Panel models run the Panel--prefixed analogues (e.g. Panel-LM-Lag).

Requires the model to have been fit (.fit() called) and a spatial weights matrix W to have been supplied at construction time.

Returns:

DataFrame indexed by test name with columns:

Column

Description

statistic

Posterior mean of the LM statistic

median

Posterior median of the LM statistic

df

Degrees of freedom for the \(\chi^2\) reference

p_value

Bayesian p-value: 1 - chi2.cdf(mean, df)

ci_lower

Lower bound of 95% credible interval (2.5%)

ci_upper

Upper bound of 95% credible interval (97.5%)

The DataFrame has attrs["model_type"] (class name) and attrs["n_draws"] (total posterior draws) metadata.

Return type:

pandas.DataFrame

Raises:
  • RuntimeError – If the model has not been fit yet.

  • ValueError – If no spatial weights matrix W was supplied.

See also

spatial_diagnostics_decision

Model-selection decision based on the test results.

spatial_effects

Posterior inference for direct/indirect/total impacts.

Examples

>>> ols = OLS(formula="price ~ income + crime", data=df, W=w)
>>> ols.fit()
>>> ols.spatial_diagnostics()
                 statistic  median  df  p_value  ci_lower  ci_upper
LM-Lag                3.21    2.98   1    0.073      0.12      8.54
LM-Error              5.67    5.34   1    0.017      0.34     12.10
LM-SDM-Joint          7.89    7.12   4    0.096      1.23     18.32
LM-SLX-Error-Joint    6.45    5.98   4    0.168      0.89     15.67
spatial_diagnostics_decision(alpha=0.05, format='graphviz')[source]

Return a model-selection decision from Bayesian LM test results.

Implements the decision tree from Koley and Bera [2024] (the Bayesian analogue of the classical stge_kb procedure in Anselin et al. [1996]). Panel models use the Panel--prefixed test analogues and the panel decision specs, following Elhorst [2014]. The decision logic depends on the current model type and the pattern of significant tests:

From OLS (6-test decision tree):

  1. If only LM-Lag is significant → SAR.

  2. If only LM-Error is significant → SEM.

  3. If both are significant → use the Anselin–Florax / Koley–Bera robust pair: Robust-LM-Lag → SAR, Robust-LM-Error → SEM, both → SARAR. If neither robust test is significant, fall back to the lower raw p-value.

  4. If neither naive test is significant → OLS.

From SAR (3-test decision tree):

  • LM-Error significant → SARAR; LM-WX significant → SDM; Robust-LM-WX significant → SDM.

From SEM (2-test decision tree):

  • LM-Lag significant → SARAR; LM-WX significant → SDEM.

From SLX (4-test decision tree):

  • Robust-LM-Lag-SDM significant → SDM; Robust-LM-Error-SDEM significant → SDEM; both → MANSAR; neither → SLX.

From SDM: LM-Error-SDM significant → MANSAR; else SDM.

From SDEM: LM-Lag-SDEM significant → MANSAR; else SDEM.

Parameters:
alpha : float, default 0.05

Significance level for the Bayesian p-values.

format : {"graphviz", "ascii", "model"}, default "graphviz"

Output format. "model" returns the recommended-model name string. "ascii" returns an indented box-drawing rendering of the full decision tree with the chosen path highlighted. "graphviz" returns a graphviz.Digraph object that renders inline in Jupyter; if the optional graphviz package is not installed a UserWarning is issued and the ASCII rendering is returned instead.

Returns:

Recommended model name when format="model", an ASCII tree string when format="ascii", or a graphviz.Digraph when format="graphviz" (with ASCII fallback on missing dep).

Return type:

str or graphviz.Digraph

See also

spatial_diagnostics

Compute the Bayesian LM test statistics.

References

Koley and Bera [2024], Anselin et al. [1996], Elhorst [2014]

spatial_effects(return_posterior_samples=False)[source]

Compute average marginal effects (AME) for the spatial probit model.

For the spatial probit with SAR regional random effects

\[z_i = x_i'\beta + a_{r(i)} + \varepsilon_i, \quad a = \rho W a + u,\]

the marginal effect of covariate \(k\) on the choice probability is

\[\frac{\partial P(y_i=1)}{\partial x_{ik}} = \phi(x_i'\beta + a_{r(i)}) \, \beta_k,\]

where \(\phi(\cdot)\) is the standard-normal PDF.

Because the spatial autoregression enters only through the unobserved regional effects \(a\) (not through a spatial multiplier on \(x\)), there is no indirect effect of \(x_j\) on \(y_i\) for \(i \neq j\). The indirect column is therefore zero and total equals direct.

Parameters:
return_posterior_samples : bool, default False

If True, also return a dict of posterior draws for each effect type.

Returns:

  • pandas.DataFrame – One row per non-intercept covariate with columns direct, direct_ci_lower, direct_ci_upper, direct_pvalue, indirect_*, total_*.

  • dict, optional – Only returned when return_posterior_samples=True. Keys: direct_samples, indirect_samples, total_samples.

summary(var_names=None, **kwargs)[source]

Return posterior summary table.

Parameters:
var_names : list, optional

Variable names to include in the summary.

**kwargs

Additional arguments passed to arviz.summary().

Returns:

Posterior summary statistics.

Return type:

pandas.DataFrame