# GAMLSSFamily


Abstract base class for GAMLSS distributional families.


Usage

``` python
GAMLSSFamily()
```


Ordinary [Family](Family.md#whittaker.Family) subclasses (used by [GAM](GAM.md#whittaker.GAM)) model only the mean `mu` of the response as a function of covariates, treating any other distributional parameters (e.g. the variance or dispersion) as constant across observations. A [GAMLSSFamily](GAMLSSFamily.md#whittaker.GAMLSSFamily), used by [GAMLSS](GAMLSS.md#whittaker.GAMLSS) instead of [GAM](GAM.md#whittaker.GAM), generalizes the single-parameter [Family](Family.md#whittaker.Family) abstraction to the location-scale-shape setting: instead of one response parameter with a mean link, it defines a full response distribution with `K` named parameters \theta_1, \ldots, \theta_K (e.g. location `mu` and scale `sigma`), each with its own link function g_k and its own additive predictor \eta_k = g_k(\theta_k). This is useful whenever more than the mean of the response changes systematically with covariates -- for example, when the spread (heteroscedasticity), skew, or zero-inflation probability also varies across the range of the predictors.

[GAMLSS](GAMLSS.md#whittaker.GAMLSS) fits every parameter's additive predictor jointly by alternating penalized IRLS updates across parameters (the RS algorithm), which relies on each family supplying the per-parameter score `dl_dtheta`, (expected) Fisher information `d2l_dtheta2`, link/inverse-link/link-derivative, log-likelihood, initial values, and a `simulate` method. Subclasses must implement all of the abstract methods below.

Whittaker ships with the following concrete GAMLSS families:

- [GaussianLS](GaussianLS.md#whittaker.GaussianLS) -- location-scale Gaussian: identity link for the mean, log link for the standard deviation.
- [GammaLS](GammaLS.md#whittaker.GammaLS) -- location-scale Gamma: log link for both the mean and the coefficient of variation.
- [BetaLS](BetaLS.md#whittaker.BetaLS) -- mean-precision Beta: logit link for the mean, log link for the precision.
- [ZeroInflatedPoisson](ZeroInflatedPoisson.md#whittaker.ZeroInflatedPoisson) -- Poisson mean plus a zero-inflation probability, for count data with excess zeros.
- [ZeroInflatedNegativeBinomial](ZeroInflatedNegativeBinomial.md#whittaker.ZeroInflatedNegativeBinomial) -- overdispersed counts with excess zeros, combining [NegativeBinomial](NegativeBinomial.md#whittaker.NegativeBinomial)-style overdispersion with zero-inflation.


## Examples

GAMLSS families are passed to [GAMLSS](GAMLSS.md#whittaker.GAMLSS), and each distributional parameter gets its own formula:


``` python
import numpy as np
import whittaker as wk

rng = np.random.default_rng(0)
n = 300
x = np.linspace(0, 2 * np.pi, n)
mu = np.sin(x)
sigma = 0.2 + 0.3 * np.abs(np.cos(x))
y = rng.normal(mu, sigma)

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

# Model both the mean and the standard deviation as smooth functions of x
model = wk.GAMLSS(
    formulas={"mu": "y ~ s(x)", "sigma": "y ~ s(x)"},
    family=wk.GaussianLS(),
)
model.fit(data)
print(model.summary())
```


    GAMLSS fit summary
    ========================================
    Family: GaussianLS(mu=identity, sigma=log)
    N obs: 300
    Global deviance: 271.2831
    AIC: 297.3004
    BIC: 345.4816
    Log-likelihood: -135.6416
    Converged: True (4 iterations)

    --- mu ---
      EDF total: 6.63
      Smooth 1: edf = 5.63

    --- sigma ---
      EDF total: 6.38
      Smooth 1: edf = 5.38


## Attributes

| Name | Description |
|----|----|
| [parameter_names](#parameter_names) | Names of the distributional parameters modeled by this family. |

------------------------------------------------------------------------


### parameter_names


Names of the distributional parameters modeled by this family.


`parameter_names: tuple[str, …]`


Every concrete [GAMLSSFamily](GAMLSSFamily.md#whittaker.GAMLSSFamily) names its parameters explicitly, e.g. `('mu', 'sigma')` for [GaussianLS](GaussianLS.md#whittaker.GaussianLS) and [GammaLS](GammaLS.md#whittaker.GammaLS), or `('mu', 'phi')` for [BetaLS](BetaLS.md#whittaker.BetaLS). [GAMLSS](GAMLSS.md#whittaker.GAMLSS) uses this tuple to determine which formulas it expects (one per name) and the order in which it cycles through parameters during the RS fitting algorithm.


## Methods

| Name | Description |
|----|----|
| [d2l_dtheta2()](#d2l_dtheta2) | Negative expected second derivative of the log-likelihood with respect to `param`. |
| [dl_dtheta()](#dl_dtheta) | First derivative of the log-likelihood with respect to `param`. |
| [initialize()](#initialize) | Starting values for all distributional parameters given `y`. |
| [link()](#link) | Apply the link function for `param`: `eta = g(theta)`. |
| [link_derivative()](#link_derivative) | Derivative of the link for `param`: `d(eta)/d(theta)`. |
| [link_inverse()](#link_inverse) | Apply the inverse link for `param`: `theta = g^-1(eta)`. |
| [log_likelihood()](#log_likelihood) | Full log-likelihood evaluated at the given parameter values. |
| [simulate()](#simulate) | Simulate response values from the distribution. |

------------------------------------------------------------------------


### d2l_dtheta2()


Negative expected second derivative of the log-likelihood with respect to `param`.


Usage

``` python
d2l_dtheta2(
    param,
    y,
    params,
)
```


Returns the (expected) Fisher information -E\[\partial^2 \ell / \partial \theta^2\] for `param`, evaluated elementwise at the current parameter estimates. Must return positive values, since it is used directly as a working weight during fitting: it both scales the working pseudo-response and forms the diagonal weight matrix for the penalized IRLS update of `param`'s additive predictor.


#### Parameters


`param: str`  
Name of the distributional parameter to differentiate with respect to, one of `parameter_names`.

`y: NDArray`  
Observed response values, shape `(n,)`.

`params: dict[str, NDArray]`  
Current estimates of all distributional parameters, keyed by name, each of shape `(n,)`.


#### Returns


`NDArray`  
Elementwise working weights, shape `(n,)`, guaranteed positive.


------------------------------------------------------------------------


### dl_dtheta()


First derivative of the log-likelihood with respect to `param`.


Usage

``` python
dl_dtheta(
    param,
    y,
    params,
)
```


Returns \partial \ell / \partial \theta evaluated elementwise at the current parameter estimates. This score forms the basis of the working pseudo-response used by [GAMLSS](GAMLSS.md#whittaker.GAMLSS)'s RS algorithm when updating the additive predictor for `param`, holding the other distributional parameters fixed.


#### Parameters


`param: str`  
Name of the distributional parameter to differentiate with respect to, one of `parameter_names`.

`y: NDArray`  
Observed response values, shape `(n,)`.

`params: dict[str, NDArray]`  
Current estimates of all distributional parameters, keyed by name, each of shape `(n,)`.


#### Returns


`NDArray`  
Elementwise first derivatives, shape `(n,)`.


------------------------------------------------------------------------


### initialize()


Starting values for all distributional parameters given `y`.


Usage

``` python
initialize(y)
```


Produces a reasonable initial guess for every parameter in `parameter_names` from the raw response alone, before any covariate information is used. [GAMLSS](GAMLSS.md#whittaker.GAMLSS) uses these as the starting point for the first RS iteration.


#### Parameters


`y: NDArray`  
Observed response values, shape `(n,)`.


#### Returns


`dict[str, NDArray]`  
Initial values for each parameter, keyed by name, each of shape `(n,)`.


------------------------------------------------------------------------


### link()


Apply the link function for `param`: `eta = g(theta)`.


Usage

``` python
link(
    param,
    values,
)
```


Maps a distributional parameter from its natural scale (e.g. `mu` on the real line, `sigma > 0`) onto the unconstrained scale of the additive predictor `eta`, so that [GAMLSS](GAMLSS.md#whittaker.GAMLSS) can model it with an ordinary smooth linear predictor.


#### Parameters


`param: str`  
Name of the distributional parameter being linked, one of `parameter_names`.

`values: NDArray`  
Parameter values on the natural scale, shape `(n,)`.


#### Returns


`NDArray`  
Linked values `eta = g(theta)`, shape `(n,)`.


------------------------------------------------------------------------


### link_derivative()


Derivative of the link for `param`: `d(eta)/d(theta)`.


Usage

``` python
link_derivative(
    param,
    values,
)
```


Used by the RS algorithm to convert the score `dl_dtheta` into a working pseudo-response on the scale of the additive predictor, via the chain rule \partial \ell / \partial \eta = (\partial \ell / \partial \theta)\\ (d\eta/d\theta)^{-1}.


#### Parameters


`param: str`  
Name of the distributional parameter being linked, one of `parameter_names`.

`values: NDArray`  
Parameter values on the natural scale at which to evaluate the derivative, shape `(n,)`.


#### Returns


`NDArray`  
Elementwise derivatives `d(eta)/d(theta)`, shape `(n,)`.


------------------------------------------------------------------------


### link_inverse()


Apply the inverse link for `param`: `theta = g^-1(eta)`.


Usage

``` python
link_inverse(
    param,
    eta,
)
```


Maps the additive predictor back onto the natural scale of the parameter, enforcing any constraints (e.g. positivity for a scale parameter, or the unit interval for a probability) implied by the chosen link.


#### Parameters


`param: str`  
Name of the distributional parameter being linked, one of `parameter_names`.

`eta: NDArray`  
Additive predictor values, shape `(n,)`.


#### Returns


`NDArray`  
Parameter values on the natural scale, shape `(n,)`.


------------------------------------------------------------------------


### log_likelihood()


Full log-likelihood evaluated at the given parameter values.


Usage

``` python
log_likelihood(
    y,
    params,
)
```


Sums the per-observation log-density of the response distribution across all observations, using the current estimate of every distributional parameter. Used for convergence checks during fitting and for computing information criteria such as AIC/BIC in model summaries.


#### Parameters


`y: NDArray`  
Observed response values, shape `(n,)`.

`params: dict[str, NDArray]`  
Current estimates of all distributional parameters, keyed by name, each of shape `(n,)`.


#### Returns


`float`  
The total log-likelihood summed over all observations.


------------------------------------------------------------------------


### simulate()


Simulate response values from the distribution.


Usage

``` python
simulate(
    params,
    rng,
)
```


Draws one random sample per observation from the response distribution at the given per-observation parameter values. Used for posterior/parametric-bootstrap style simulation from a fitted [GAMLSS](GAMLSS.md#whittaker.GAMLSS) model.


#### Parameters


`params: dict[str, NDArray]`  
Distributional parameter values to simulate from, keyed by name, each of shape `(n,)`.

`rng: np.random.Generator`  
A NumPy random generator (e.g. `numpy.random.default_rng()`) used to draw the simulated values.


#### Returns


`NDArray`  
Simulated response values, shape `(n,)`.