Fitting¶

Functions for maximum-likelihood parameter estimation. The recommended entry point is fit_listener, which implements the full two-stage profile-likelihood procedure described in Barumerli et al.[1].

Two-stage fitting¶

bayesian_listener.fitting.fit_listener(sofa_path, obs_tbl, targets_coords, interpolation_method, subject_id=None, num_repetitions=200, num_repetitions_motor=200, num_grid_points=1, fix_sigma_ild=True, motor_estimation_seed=42, verbose=True)[source]¶

Run the full two-stage fit recommended in Barumerli et al.[1].

Estimate \(\hat{\kappa}_m\) (Eq. 12) via estimate_motor_noise.
Hold \(\hat{\kappa}_m\) fixed and jointly fit \(\sigma_{\mathrm{mon}}\) and \(\sigma_{\mathrm{prior}}\) (Eq. 13) via fit_listener_partial.

Parameters:

sofa_path (str) – Path to the participant’s SOFA file.
obs_tbl (pandas.DataFrame) – Behavioural observations. Must contain columns 'azi_response', 'ele_response', 'azi_target', 'ele_target'. Must additionally contain 'participant' when subject_id is given.
targets_coords (pyfar.Coordinates) – Target directions presented in the experiment.
interpolation_method ({'SH', 'SHMAX', 'barycentric', 'barumerli2023'}) – Forwarded to compute_template.
subject_id (str or None, default=None) – Participant identifier. None uses every row of obs_tbl.
num_repetitions (int, default=200) – Monte Carlo repetitions for the stage-2 likelihood (Eq. 13).
num_repetitions_motor (int, default=200) – Monte Carlo repetitions for the stage-1 motor-noise estimation (Eq. 12).
num_grid_points (int, default=1) – Initialisation grid size per parameter dimension before BADS.
fix_sigma_ild (bool, default=True) – If True, fix \(\sigma_{\mathrm{ild}} = 1.0\) dB; if False, fit it alongside the other free parameters.
motor_estimation_seed (int, default=42) – Seed forwarded to estimate_motor_noise.
verbose (bool, default=True) – Print stage-by-stage progress.

Returns:

Result mapping with keys sigma_itd, sigma_ild, sigma_spectral, kappa_motor, sigma_motor, sigma_prior, nll, n_trials, time_* (timing breakdown), success (bool) and on failure error (str).

Return type:

dict

Likelihood¶

bayesian_listener.fitting.negloglik(model, targets, responses, resp_targets_idx, sigmas_log, num_repetitions=200)[source]¶

Negative full-sphere log-likelihood for BADS optimisation (Eq. 13 of Barumerli et al.[1]).

For each observation, runs Monte Carlo inference with the supplied parameters, builds the von Mises–Fisher pdf around each MC prediction, and accumulates \(-\log p(\hat{\boldsymbol{\varphi}}^* \mid \boldsymbol{\varphi}, \boldsymbol{\theta})\).

Parameters:

model (BayesianListener) – Model instance whose parameters and target are overwritten in place.
targets (Barumerli2023) – Subset of the model template at the unique presented directions.
responses (pyfar.Coordinates) – Observed responses, shape (n_obs,).
resp_targets_idx (numpy.ndarray) – Integer mapping of each response to its target direction, shape (n_obs,).
sigmas_log (numpy.ndarray) –
Log-transformed parameter vector of shape (4,):
- sigmas_log[0] — \(\log \sigma_{\mathrm{ild}}\) (dB).
- sigmas_log[1] — \(\log \sigma_{\mathrm{mon}}\) (dB).
- sigmas_log[2] — \(\log \kappa_m\).
- sigmas_log[3] — \(\log \tau_{\mathrm{prior}}\) with \(\tau_{\mathrm{prior}} = 1 / \sigma_{\mathrm{prior}}^2\) (precision parametrisation for a better optimisation landscape).
num_repetitions (int, default=200) – Monte Carlo repetitions per trial in infer.

Returns:

Negative log-likelihood in nats.

Return type:

float

bayesian_listener.fitting.von_mises_loglik_mc(kappa, resp_lat, est_lat_mc)[source]¶

Negative log-likelihood of a von Mises with Monte Carlo predictions.

Implements the lateral-only likelihood of Eq. 12 of Barumerli et al.[1], approximated by averaging the von Mises pdf over n_mc Monte Carlo samples per observation.

Parameters:

kappa (float) – Von Mises concentration \(\kappa\) (positive).
resp_lat (numpy.ndarray) – Observed lateral angles, shape (n_obs,) in radians.
est_lat_mc (numpy.ndarray) – Monte Carlo model predictions, shape (n_obs, n_mc) in radians.

Returns:

Negative log-likelihood in nats.

Return type:

float

Examples

>>> rng = np.random.default_rng(0)
>>> resp = rng.normal(scale=0.1, size=20)
>>> mc   = resp[:, None] + rng.normal(scale=0.1, size=(20, 50))
>>> nll  = von_mises_loglik_mc(50.0, resp, mc)
>>> bool(np.isfinite(nll))
True

bayesian_listener.fitting.fit_kappa_ml(resp_lat, est_lat_mc)[source]¶

Fit the von Mises concentration \(\kappa\) by 1-D bounded ML search.

Wraps scipy.optimize.minimize_scalar (method='bounded', Brent’s method) over the bracket \(\kappa \in [0.1, 1000]\).

Parameters:

resp_lat (numpy.ndarray) – Observed lateral angles, shape (n_obs,) in radians.
est_lat_mc (numpy.ndarray) – Monte Carlo model predictions, shape (n_obs, n_mc) in radians.

Returns:

Maximum-likelihood concentration \(\hat{\kappa}\).

Return type:

float

Examples

>>> rng = np.random.default_rng(0)
>>> resp = rng.normal(scale=0.05, size=200)
>>> mc   = resp[:, None] + rng.normal(scale=0.05, size=(200, 100))
>>> kappa = fit_kappa_ml(resp, mc)
>>> bool(0.1 <= kappa <= 1000.0)
True

Parameter conversion¶

bayesian_listener.fitting.sigma_to_kappa(sigma)[source]¶

Convert a circular standard deviation in degrees to a von Mises concentration.

Solves \(I_1(\kappa)/I_0(\kappa) = \exp(-\sigma^2/2)\) for \(\kappa\) via scipy.optimize.brentq on the Bessel-ratio identity. Falls back to the asymptotic approximation \(\kappa \approx 1/(2(1-R))\) for very small \(\sigma\).

Parameters:: sigma (float) – Circular standard deviation in degrees. Saturates to a near-uniform \(\kappa = 10^{-6}\) for sigma >= 180.
Returns:: Von Mises concentration \(\kappa\).
Return type:: float

Examples

>>> bool(abs(kappa_to_sigma(sigma_to_kappa(15.0)) - 15.0) < 1e-3)
True

bayesian_listener.fitting.kappa_to_sigma(kappa)[source]¶

Convert a von Mises concentration to a circular standard deviation in degrees.

Computes \(\sigma = \sqrt{-2 \log(I_1(\kappa)/I_0(\kappa))}\) and converts to degrees.

Parameters:: kappa (float) – Von Mises concentration \(\kappa\) (positive). Saturates to 180° for \(\kappa < 10^{-6}\).
Returns:: Circular standard deviation in degrees.
Return type:: float

Examples

>>> round(float(kappa_to_sigma(50.0)), 1)
8.1

Default parameters and bounds¶

bayesian_listener.fitting.DEFAULT_PARAMS = mappingproxy({'sigma_itd': 0.569, 'sigma_ild': 1.0, 'sigma_spectral': 10.0, 'kappa_motor': 15.105536409469574, 'sigma_prior': 40.0})¶: Read-only mapping of default noise-parameter values. Used as the fallback in fit_listener_partial for any parameter that is neither in params_to_fit nor fixed_params. Pass fixed_params to override individual values; do not mutate this object directly.

bayesian_listener.fitting.DEFAULT_PARAM_BOUNDS = mappingproxy({'sigma_ild': mappingproxy({'lb': 0.1, 'plb': 0.5, 'pub': 3.0, 'ub': 50.0}), 'sigma_spectral': mappingproxy({'lb': 0.1, 'plb': 1.0, 'pub': 10.0, 'ub': 50.0}), 'kappa_motor': mappingproxy({'lb': 0.8156, 'plb': 2.6852, 'pub': 131.8139, 'ub': 821.2015}), 'tau_prior': mappingproxy({'lb': 3.1e-05, 'plb': 0.0004, 'pub': 0.04, 'ub': 1.0})})¶: Read-only mapping of BADS optimiser search bounds for each free parameter. Keys are 'sigma_ild', 'sigma_spectral', 'kappa_motor', and 'tau_prior' (precision \(\tau = 1/\sigma_{\mathrm{prior}}^2\), used internally for a better optimisation landscape). Each value is a read-only mapping with keys 'lb', 'plb', 'pub', 'ub' (lower / plausible-lower / plausible-upper / upper bounds). 'kappa_motor' bounds are von Mises concentrations (higher = less noise).