jwspecabund.forward

Bayesian forward model for emission-line abundance determination.

Implements the forward modelling approach of Cullen et al. (2025): free parameters (ionic abundances, T_e, n_e) are sampled via MCMC or nested sampling. For each parameter vector, CEL emissivities from PyNEB and the Aller (1984) Hbeta recombination formula predict line flux ratios that are compared to the observed ratios via a Gaussian likelihood.

References

  • Aller L. H., 1984, Physics of Thermal Gaseous Nebulae (Reidel)

  • Cullen et al. (2025) — forward modelling for EXCELS galaxies

Functions

forward_model(line_fluxes, line_errors, *[, ...])

Run a Bayesian forward model on emission-line fluxes.

hbeta_emissivity_aller84(Te)

Hbeta volume emissivity coefficient using the Aller (1984) Case B formula.

jwspecabund.forward.forward_model(line_fluxes, line_errors, *, sampler='emcee', n_walkers=32, n_steps=5000, n_burn=1000, n_live=500, seed=42, progress=True)[source]

Run a Bayesian forward model on emission-line fluxes.

Free parameters are log(T_e), log(n_e), and one log(ionic abundance) per detected ion. The model predicts line/Hbeta flux ratios using PyNEB CEL emissivities and the Aller (1984) Hbeta formula, and the Gaussian likelihood compares predictions to observations.

Parameters:
  • line_fluxes (dict) – {line_name: flux}. Must include "HBETA".

  • line_errors (dict) – {line_name: flux_err}.

  • sampler (str) – "emcee" (default) or "dynesty".

  • n_walkers (int) – Number of walkers for emcee (default 32).

  • n_steps (int) – Total MCMC steps for emcee (default 5000).

  • n_burn (int) – Burn-in steps to discard for emcee (default 1000).

  • n_live (int) – Number of live points for dynesty (default 500).

  • seed (int) – Random seed (default 42).

  • progress (bool) – Show progress bar (default True).

Returns:

Result dictionary with keys:

  • param_names, param_labels — parameter metadata

  • samples — raw posterior samples, shape (n, n_dim)

  • Te, Te_err, Te_posterior

  • ne, ne_err, ne_posterior

  • OH, OH_err, OH_posterior

  • NO, NO_err, NO_posterior (if nitrogen detected)

  • ionic — dict of median log-ionic abundances

  • ionic_posteriors — dict of log-ionic posterior arrays

  • ratio_names — observable labels

Return type:

dict

Examples

>>> import jwspecabund
>>> fluxes = {
...     "HBETA": 1.0, "OIII_5007": 5.0,
...     "OIII_4363": 0.05, "NeIII_3869": 0.3,
... }
>>> errors = {k: 0.05 * v for k, v in fluxes.items()}
>>> res = jwspecabund.forward_model(fluxes, errors, n_steps=2000)
>>> print(f"12+log(O/H) = {res['OH']:.2f}")
jwspecabund.forward.hbeta_emissivity_aller84(Te)[source]

Hbeta volume emissivity coefficient using the Aller (1984) Case B formula.

\[\varepsilon_{H\beta} = \alpha_{H\beta}^{\rm eff}(T)\;h\nu_{H\beta}\]

where \(\alpha_{H\beta}^{\rm eff} = 3.03\times10^{-14}\,(T/10^4)^{-0.874}\) cm3 s-1.

Unlike PyNEB’s recombination tables, this formula is valid at all temperatures — no upper-bound limitation at 30,000 K.

Parameters:

Te (float) – Electron temperature in K.

Returns:

\(\varepsilon_{H\beta}\) in erg cm3 s-1.

Return type:

float