hessian_ci

Contents

hessian_ci#

class braintools.fit.hessian_ci(loss, best_params, *, confidence=0.95, kind='sse', n_data=None)#

Compute confidence intervals from the Hessian of a loss at its minimum.

The loss is expanded to second order around the optimum \(\hat\theta\) (in the constrained, physical parameter space), with Hessian \(H = \nabla^2 L(\hat\theta)\) taken over all p scalar entries of best_params:

  • kind='nll': loss is a negative log-likelihood and the Laplace approximation gives \(\mathrm{Cov} = H^{-1}\); intervals use the normal quantile.

  • kind='sse': loss is a sum of squared errors over n_data residuals. With \(s^2 = L(\hat\theta) / (n - p)\), \(\mathrm{Cov} = 2 s^2 H^{-1}\), and intervals use the Student-t quantile with \(n - p\) degrees of freedom. For a linear model this is exactly the ordinary-least-squares covariance \(s^2 (X^\top X)^{-1}\).

The interval of parameter \(i\) is \(\hat\theta_i \pm q\,\sqrt{\mathrm{Cov}_{ii}}\).

Parameters:

  • loss (Callable[[Any], Any]) – loss(params) -> scalar, twice differentiable, taking a pytree with the structure of best_params.

  • best_params (Any) – The optimum (e.g. FitResult.best_params). Leaves may be brainunit Quantities; the loss then receives Quantities and the results carry the same units.

  • confidence (float) – Two-sided confidence level in (0, 1).

  • kind (str) – Meaning of the loss, see above.

  • n_data (int | None) – Number of residuals in the sum of squares; required for kind='sse' and must exceed the number of parameters.

Returns:

(estimate, stderr, low, high), each with the structure of best_params.

Return type:

HessianCI

Raises:

  • ValueError – If kind is unknown, confidence is outside (0, 1), n_data is missing or not larger than the number of parameters for kind='sse', or best_params has no leaves.

  • TypeError – If confidence or n_data is traced, or the loss carries a physical unit.

See also

Fitter

Produces best_params.

Notes

The intervals are asymptotic and assume the optimum is an interior minimum where the loss is locally quadratic. A singular or indefinite Hessian (non-identifiable parameters, a saddle) gives non-finite or meaningless standard errors.

confidence and n_data fix the quantile, so they must be concrete; with them static, the function works under jax.jit and jax.vmap over best_params. It is not meant to be differentiated.

Examples

>>> import jax.numpy as jnp
>>> from braintools.fit import hessian_ci
>>> x = jnp.arange(6.)
>>> y = jnp.array([0.1, 1.9, 4.2, 5.8, 8.1, 9.9])
>>> sse = lambda p: jnp.sum((y - p['a'] - p['b'] * x) ** 2)
>>> best = {'a': jnp.asarray(0.05714286), 'b': jnp.asarray(1.9771429)}
>>> ci = hessian_ci(sse, best, n_data=6)
>>> print(round(float(ci.stderr['b']), 4), round(float(ci.low['b']), 3), round(float(ci.high['b']), 3))
0.0398 1.867 2.088
Contents