BetaLS

Beta family for GAMLSS with mean-precision parameterisation.

Usage

Source

BetaLS()

BetaLS extends the plain Beta family by letting both the mean mu and the precision phi vary smoothly with covariates, rather than treating precision as a single estimated constant. Use it for responses strictly between 0 and 1 (rates, proportions, fractions) where not only the typical level but also how tightly the response clusters around that level changes across the range of the predictors — for example, a proportion that becomes more variable in some regions of the covariate space and more tightly concentrated in others. The mean uses the logit link (as in Beta) and the precision uses the log link, keeping mu in (0, 1) and phi > 0.

Notes

The two parameters use distinct link functions:

g_{\mu}(\mu) = \log\!\left(\frac{\mu}{1-\mu}\right), \qquad g_{\phi}(\phi) = \log(\phi).

If a = mu * phi and b = (1 - mu) * phi, the response follows y ~ Beta(a, b) with density

f(y \mid \mu, \phi) = \frac{y^{a-1}(1-y)^{b-1}}{B(a, b)}, \qquad a = \mu\phi,\ \ b = (1-\mu)\phi.

As in Beta, larger phi concentrates the distribution more tightly around mu (Var(Y) = mu(1-mu) / (1+phi)), but here phi is itself modeled as a smooth function of covariates rather than a single scalar.

Examples

Fit a GAMLSS where both the mean and precision of a proportion response vary smoothly:

import numpy as np
import whittaker as wk
from scipy.special import expit

rng = np.random.default_rng(0)
n = 300
x = np.linspace(0, 2 * np.pi, n)
mu = expit(np.sin(x))
phi = 10.0 + 15.0 * np.abs(np.cos(x))
y = rng.beta(mu * phi, (1 - mu) * phi)

data = {"x": x, "y": y}

model = wk.GAMLSS(
    formulas={"mu": "y ~ s(x)", "phi": "y ~ s(x)"},
    family=wk.BetaLS(),
)
model.fit(data)
print(model.summary())
GAMLSS fit summary
========================================
Family: BetaLS(mu=logit, phi=log)
N obs: 300
Global deviance: -548.1336
AIC: -521.9058
BIC: -473.3348
Log-likelihood: 274.0668
Converged: True (6 iterations)

--- mu ---
  EDF total: 6.37
  Smooth 1: edf = 5.37

--- phi ---
  EDF total: 6.74
  Smooth 1: edf = 5.74

Attributes

Name Description
parameter_names Names of the distributional parameters modeled by BetaLS.

parameter_names

Names of the distributional parameters modeled by BetaLS.

parameter_names: tuple[str, …]

Methods

Name Description
d2l_dtheta2() Expected Fisher information for mu or phi.
dl_dtheta() First derivative of the Beta log-likelihood with respect to mu or phi.
initialize() Starting values for mu and phi from the raw response.
link() Apply the link function for mu or phi.
link_derivative() Derivative of the link function for mu or phi.
link_inverse() Apply the inverse link for mu or phi.
log_likelihood() Total Beta log-likelihood at the current mu and phi.
simulate() Simulate response values from a Beta distribution with the given mu, phi.

d2l_dtheta2()

Expected Fisher information for mu or phi.

Usage

Source

d2l_dtheta2(
    param,
    y,
    params,
)

With a = mu * phi, b = (1 - mu) * phi, and trigamma function \psi_1:

-E\!\left[\frac{\partial^2 \ell}{\partial \mu^2}\right] = \phi^2\big(\psi_1(a) + \psi_1(b)\big), \qquad -E\!\left[\frac{\partial^2 \ell}{\partial \phi^2}\right] = \mu^2 \psi_1(a) + (1-\mu)^2 \psi_1(b) - \psi_1(\phi).

Parameters

param: str

Either "mu" or "phi".

y: NDArray

Observed response values, shape (n,) (unused, kept for interface consistency).

params: dict[str, NDArray]
Current estimates with keys "mu" and "phi", each of shape (n,).

Returns

NDArray
Elementwise working weights, shape (n,).

dl_dtheta()

First derivative of the Beta log-likelihood with respect to mu or phi.

Usage

Source

dl_dtheta(
    param,
    y,
    params,
)

With a = mu * phi, b = (1 - mu) * phi, y* = logit(y), and mu* = psi(a) - psi(b) (where \psi is the digamma function), the score functions are

\frac{\partial \ell}{\partial \mu} = \phi\,(y^{*} - \mu^{*}), \qquad \frac{\partial \ell}{\partial \phi} = \mu\,(y^{*} - \mu^{*}) + \psi(\phi) - \psi(b) + \log(1 - y).

Parameters

param: str

Either "mu" or "phi".

y: NDArray

Observed response values in (0, 1), shape (n,).

params: dict[str, NDArray]
Current estimates with keys "mu" and "phi", each of shape (n,).

Returns

NDArray
Elementwise first derivatives, shape (n,).

initialize()

Starting values for mu and phi from the raw response.

Usage

Source

initialize(y)

mu is initialized at the observed values themselves (clipped to [0.01, 0.99] to stay strictly inside the unit interval), and phi at a method-of-moments estimate of the precision derived from the sample mean and variance of y (clamped to at least 1.0), constant across observations.

Parameters

y: NDArray
Observed response values in (0, 1), shape (n,).

Returns

dict[str, NDArray]
{"mu": y_safe.copy(), "phi": <constant array of the moment-based precision>}.

log_likelihood()

Total Beta log-likelihood at the current mu and phi.

Usage

Source

log_likelihood(
    y,
    params,
)

Evaluated using the shape parameterization a = mu * phi, b = (1 - mu) * phi.

Parameters

y: NDArray

Observed response values in (0, 1), shape (n,).

params: dict[str, NDArray]
Current estimates with keys "mu" and "phi", each of shape (n,).

Returns

float
The log-likelihood summed over all observations.

simulate()

Simulate response values from a Beta distribution with the given mu, phi.

Usage

Source

simulate(
    params,
    rng,
)

Converts to the shape parameterization (a = mu * phi, b = (1 - mu) * phi) expected by rng.beta.

Parameters

params: dict[str, NDArray]

Current estimates with keys "mu" and "phi", each of shape (n,).

rng: np.random.Generator
A NumPy random generator used to draw the simulated values.

Returns

NDArray
Simulated response values in (0, 1), shape (n,).