Source code for pybnf.priors.base
"""The ``Prior`` abstraction: a distribution family evaluated entirely in the
sampling space ``u`` (ADR-0010).
A ``Prior`` is **scale-agnostic** -- it knows nothing about ``theta`` or
``log10``. The owning ``FreeParameter`` holds the ``Scale`` and applies the
``theta <-> u`` transform, calling ``prior.logpdf(scale.forward(theta))`` and
``scale.inverse(prior.rvs(rng))``. This keeps each family's math pure and the
scale in one place (ADR-0003).
Concrete families (``Normal``, ``Uniform``, ...) live one-per-file and
self-register with ``@register_prior_family``. ``NoPrior`` is the first-class
null-object for ``var``/``logvar`` Simplex start points: a free parameter with a
scale but no distribution.
"""
from abc import ABC, abstractmethod
import numpy as np
from ..printing import PybnfError
[docs]
class Prior(ABC):
"""A distribution family operating in the sampling space ``u``.
Subclasses expose ``logpdf``/``rvs``/``ppf`` in ``u`` and report their
``support`` (in ``u``) and ``has_bounded_support``. ``frozen`` is the
underlying ``scipy.stats`` frozen distribution (or ``None`` for ``NoPrior``),
surfaced for ``FreeParameter._distribution`` back-compat.
"""
#: Whether the family has a proper distribution (``False`` only for ``NoPrior``).
has_prior = True
#: Whether the family's support is finite -- drives reflecting-bounds
#: eligibility and latin-hypercube participation. ``Uniform`` overrides.
has_bounded_support = False
#: The family's natural lower support endpoint in sampling space ``u`` (the lower
#: edge of ``support()``, a family constant independent of the distribution's
#: parameters). ``-inf`` for the doubly-unbounded families; the positive-support
#: families (gamma/exponential/chisquare/rayleigh, the half-* scale priors,
#: inv_gamma/weibull, and beta's ``[0,1]``) override to ``0.0``. The owning
#: ``FreeParameter``'s ``Scale.inverse`` maps it to the theta-space floor a one-sided
#: truncation measures bounds against (ADR-0047).
support_lo_u = -np.inf
#: How many config numbers the family's parameterization takes -- ``2`` for the
#: location/scale/bounds families; the one-parameter families (exponential/chisquare/
#: rayleigh, the half-* scale priors) override to ``1`` so the positional grammar admits a
#: single number (ADR-0010/#417); the three-parameter families (student_t) override to
#: ``3``. A ``n_params >= 3`` family is authored only through the new-era ``parameter:``
#: record -- the legacy positional ``*_var`` grammar carries at most two numbers, so
#: ``var_keyword_grammar`` omits it (ADR-0057).
n_params = 2
#: The config field names for the family's distribution parameters, in ``build()`` order
#: -- the new-era ``parameter:`` record names each one (ADR-0043), so a positional
#: ``p1 p2`` becomes ``mean: .. , sd: ..``. Concrete families override; the length must
#: match ``n_params`` (the record builds ``p1``/``p2``/``p3`` from the first three).
field_names = ('p1', 'p2')
#: The underlying scipy frozen distribution, or ``None``.
frozen = None
[docs]
@abstractmethod
def rvs(self, rng):
"""Draw one sample in sampling space ``u`` using ``rng``.
``rng`` is the caller's :class:`numpy.random.Generator`; it is passed to
scipy as ``random_state`` so prior sampling draws from the algorithm's
seeded Generator rather than NumPy's legacy global RNG.
"""
[docs]
@abstractmethod
def ppf(self, q):
"""Inverse CDF at quantile ``q`` (in ``[0, 1]``), in sampling space ``u``."""
[docs]
@abstractmethod
def support(self):
"""The ``(lo, hi)`` support in sampling space ``u`` (may be infinite)."""
[docs]
def logpdf_jax(self, u):
"""JAX-traceable log prior density at sampling-space ``u`` (ADR-0059).
The gradient-based reference sampler (``job_type = hmc``) composes its
target log-density entirely from these per-family JAX logpdfs plus the
model's JAX NLL, so ``jax.grad`` differentiates it. Every family in the
edition-2 catalog overrides this with a hand-written JAX density,
oracle-checked against its scipy ``logpdf`` (ADR-0059 item 4), so this base
implementation is the fallback for a family that has not supplied one -- it
raises a pointed error rather than silently producing a wrong target. ``u``
is a JAX scalar; the return is a JAX scalar."""
raise PybnfError(
"Prior family %r has no JAX log-density (logpdf_jax), so it cannot be "
"used with job_type = hmc (ADR-0059). The edition-2 prior catalog is "
"fully mapped to JAX; this error means a custom/unregistered family was "
"supplied -- implement its logpdf_jax (oracle-checked against its scipy "
"logpdf), or run a gradient-free sampler (am / dream / p_dream)."
% type(self).__name__)
[docs]
class FrozenPrior(Prior):
"""A :class:`Prior` backed entirely by a ``scipy.stats`` frozen distribution in the
sampling space ``u`` (ADR-0010).
Every location/scale/shape family (Normal, Laplace, Cauchy, Gamma, Exponential,
ChiSquare, Rayleigh) has the *same* density/sampling/inverse-CDF/support shape -- a thin
delegation to its frozen distribution -- so it lives here once. A subclass sets
``self.frozen`` in ``__init__`` (and ``has_bounded_support`` as a class attribute) and
provides the family-specific ``build`` classmethod. ``NoPrior`` (no distribution) and
``Uniform`` (a custom latin-hypercube ``ppf``) do not use this base.
"""
[docs]
class NoPrior(Prior):
"""The absence of a prior: a ``var``/``logvar`` Simplex start point.
Contributes nothing to the log prior, cannot be sampled, and has no support.
It still pairs with a ``Scale`` (``logvar`` is ``Log10``); the scale lives on
the ``FreeParameter``, not here.
"""
has_prior = False
has_bounded_support = False
field_names = ()
frozen = None
[docs]
@classmethod
def build(cls, p1, p2, scale, p3=None):
"""Factory matching the family ``build`` signature; ``p1``/``p2``/``p3``/``scale``
are ignored -- a no-prior parameter carries only a start value."""
return cls()
[docs]
def logpdf_jax(self, u):
"""A no-prior carrier contributes ``0`` to the HMC target (ADR-0059), the
JAX peer of :meth:`logpdf`. Returned as a JAX scalar so it sums cleanly
with the family logpdfs without forcing a host/device transfer."""
import jax.numpy as jnp
return jnp.asarray(0.0)