mixle.stats.univariate.discrete.logseries module

Logarithmic (log-series) distributions over positive integers.

Observations are positive integers k >= 1. A log-series distribution with shape p in (0, 1) has log-mass

log(P(k; p)) = k*log(p) - log(k) - log(-log(1 - p)), k = 1, 2, …,

a one-parameter exponential family used for over-dispersed positive counts (species abundance, word-frequency / type-token models). The mean is -p / ((1 - p) * log(1 - p)) (a value in (1, inf) that increases with p), which the estimator inverts.

The per-row score is linear in the encoded k and log(k) fields once the scalar normalizer log(-log(1 - p)) is precomputed, so the family gets generated NumPy, Torch, and Numba kernels.

Reference: Fisher, Corbet & Williams, ‘The relation between the number of species and the number of individuals…’, J. Animal Ecology (1943).

class LogSeriesDistribution(p, name=None, keys=None)

Bases: SequenceEncodableProbabilityDistribution

Logarithmic (log-series) distribution on k = 1, 2, … with shape parameter p in (0, 1).

Parameters:

• p (float)

• name (str | None)

• keys (str | None)

classmethod compute_capabilities()

Describe backend support for generated log-series kernels.

classmethod compute_declaration()

Return the structured compute declaration for log-series distributions.

static backend_legacy_sufficient_statistics(x, params, engine)

Return per-row (count, k) sufficient statistics in accumulator order.

Parameters:

• x (tuple[Any, Any])

• params (dict[str, Any])

• engine (Any)

Return type:

tuple[Any, …]

static exp_family_sufficient_statistics(x, engine)

Return the log-series sufficient statistic T(k) = (k,).

Parameters:

• x (tuple[Any, Any])

• engine (Any)

Return type:

tuple[Any, …]

static exp_family_natural_parameters(params, engine)

Return the log-series natural parameter eta = log(p).

Parameters:

• params (dict[str, Any])

• engine (Any)

Return type:

tuple[Any, …]

static exp_family_log_partition(params, engine)

Return the log-series log partition A = log(-log(1 - p)).

Parameters:

• params (dict[str, Any])

• engine (Any)

Return type:

Any

static exp_family_base_measure(x, engine)

Return the log-series base measure log h(k) = -log(k) for k >= 1, and -inf below.

-log(k) alone is +inf at an admitted k = 0 and NaN at a negative count, so the generated engine kernels returned +inf and NaN where every numpy route returns -inf (V01-F04). The encoder guarantees integers.

Parameters:

• x (tuple[Any, Any])

• engine (Any)

Return type:

Any

static exp_family_from_natural(eta)

Return the log-series with natural parameter eta = log(p) (so p = exp(eta)).

Parameters:

eta (Any)

Return type:

LogSeriesDistribution

static backend_log_density_from_params(k, log_k, log_p, log_norm, engine)

Engine-neutral log-series log-mass from explicit parameters (-inf for k < 1).

Linear in k and log k, which is +inf at an admitted k = 0 and NaN at a negative count, so the torch kernel returned both where the numpy scorer returns -inf (V01-F04).

Parameters:

• k (Any)

• log_k (Any)

• log_p (Any)

• log_norm (Any)

• engine (Any)

Return type:

Any

density(x)

Return the probability mass at a single observation.

Parameters:

x (int)

Return type:

float

log_density(x)

Return the log-mass at a single positive integer (or -inf off support).

Parameters:

x (int)

Return type:

float

seq_log_density(x)

Return vectorized log-mass values for sequence-encoded (k, log k) observations.

Parameters:

x (tuple[ndarray, ndarray])

Return type:

ndarray

backend_seq_log_density(x, engine)

Engine-neutral vectorized log-mass for encoded data.

Parameters:

• x (tuple[Any, Any])

• engine (Any)

Return type:

Any

classmethod backend_stacked_params(dists, engine)

Return stacked parameters for a homogeneous mixture kernel.

Parameters:

• dists (Sequence[LogSeriesDistribution])

• engine (Any)

Return type:

dict[str, Any]

classmethod backend_stacked_log_density(x, params, engine)

Return an (n, k) matrix of log-series log densities.

Parameters:

• x (tuple[Any, Any])

• params (dict[str, Any])

• engine (Any)

Return type:

Any

classmethod backend_stacked_sufficient_statistics(x, weights, params, engine)

Return stacked sufficient statistics using engine-resident arrays.

Parameters:

• x (tuple[Any, Any])

• weights (Any)

• params (dict[str, Any])

• engine (Any)

Return type:

tuple[Any, Any]

mean()

Mean E[X] = -p / ((1-p) log(1-p)).

Return type:

float

variance()

Variance Var[X] = -p (p + log(1-p)) / ((1-p)^2 log(1-p)^2).

Return type:

float

entropy()

Shannon entropy in nats, by exact summation of the standard series.

The log-series entropy -sum_k p_k log p_k has no closed form (Fisher, Corbet & Williams, 1943). The series is summed from k = 1 up to this distribution’s own quantile at 1 - 1e-16 (plus a safety margin), beyond which the remaining tail mass is far below double rounding.

Return type:

float

cdf(x)

Cumulative distribution function P(X <= x), support x >= 1.

Summed directly in bounded blocks (the series has no closed form); past _SERIES_CAP terms one minus the Euler-Maclaurin tail mass stands in, exact to double rounding by then and never below the partial sum at the cap, so the result is monotone in x and within [0, 1].

Parameters:

x (float)

Return type:

float

quantile(q)

Inverse CDF F^{-1}(q): the smallest k with P(X <= k) >= q.

Summed directly in bounded blocks, never through scipy’s generic discrete ppf. That search brackets the answer by doubling an upper bound until cdf(bound) >= q, and the summed CDF can saturate one ULP short of a q this close to one (numpy’s SIMD exp and log round the last bit differently per CPU: the sum reached q = 1 - 1e-16 at k = 5210 on arm64 and emulated x86-64 and never did on AVX-512 hosted CI runners), so the bound doubled without limit and the vectorized CDF sums grew until the 16 GB runner was OOM-killed – through entropy(), which calls this at exactly that q. Here a block whose terms no longer move the running sum, or the term cap, ends the scan, and the Euler-Maclaurin tail mass answers for whatever lies beyond: the smallest k past the scan whose tail mass is at most 1 - q.

Parameters:

q (float)

Return type:

float

sampler(seed=None)

Return a sampler for drawing observations from this distribution.

Parameters:

seed (int | None)

Return type:

LogSeriesSampler

estimator(pseudo_count=None)

Return an estimator for fitting this distribution from data.

Parameters:

pseudo_count (float | None)

Return type:

LogSeriesEstimator

dist_to_encoder()

Return the data encoder used by this distribution for vectorized methods.

Return type:

LogSeriesDataEncoder

enumerator()

Return an enumerator over k = 1, 2, … in descending probability order.

Return type:

LogSeriesEnumerator

class LogSeriesEnumerator(dist)

Bases: DistributionEnumerator

Enumerate log-series support values in descending probability order.

The log-series pmf p^k / (k * -log(1 - p)) is strictly decreasing in k (each step multiplies by p * k / (k + 1) < 1), so value order 1, 2, 3, … IS descending-probability order. The iterator is infinite.

Parameters:

dist (LogSeriesDistribution)

class LogSeriesSampler(dist, seed=None)

Bases: DistributionSampler

Draw iid log-series observations.

Parameters:

• dist (LogSeriesDistribution)

• seed (int | None)

sample(size=None, *, batched=True)

Draw size iid positive integers (an int when size is None).

Parameters:

• size (int | None)

• batched (bool)

Return type:

int | ndarray

class LogSeriesAccumulator(keys=None)

Bases: SequenceEncodableStatisticAccumulator

Accumulate weighted count and sum for log-series estimation.

Parameters:

keys (str | None)

update(x, weight, estimate)

Accumulate weighted count and total for one positive integer.

Parameters:

• x (int)

• weight (float)

• estimate (LogSeriesDistribution | None)

Return type:

None

initialize(x, weight, rng)

Initialize statistics from one observation.

Parameters:

• x (int)

• weight (float)

• rng (RandomState | None)

Return type:

None

supported_rows(x)

Encoded rows this law can be fitted on: finite and at least one.

The encoder admits out-of-support observations so that a mixture can encode a whole batch against every component (P02-F03); this is the predicate that keeps them out of the sufficient statistics, and the one a latent model’s initialization consults before it hands this component any responsibility. Without it the initializer could not tell that this component does not own a row, seeded it from one anyway, and the refusal below then blamed the DATA for a mixture the model can express – a Poisson-plus-Geometric fit on counts with zeros was refused outright (Q01-F02/Q02-F02).

Return type:

ndarray

seq_update(x, weights, estimate)

Accumulate weighted count and total from encoded observations.

Parameters:

• x (tuple[ndarray, ndarray])

• weights (ndarray)

• estimate (LogSeriesDistribution | None)

Return type:

None

seq_initialize(x, weights, rng)

Initialize statistics from encoded observations.

Parameters:

• x (tuple[ndarray, ndarray])

• weights (ndarray)

• rng (RandomState | None)

Return type:

None

combine(suff_stat)

Merge another log-series sufficient-statistic tuple.

Parameters:

suff_stat (tuple[float, float])

Return type:

LogSeriesAccumulator

value()

Return accumulated count and integer total.

Return type:

tuple[float, float]

from_value(x)

Replace accumulator contents from a sufficient-statistic tuple.

Parameters:

x (tuple[float, float])

Return type:

LogSeriesAccumulator

acc_to_encoder()

Return the encoder used by this accumulator.

Return type:

LogSeriesDataEncoder

class LogSeriesAccumulatorFactory(keys=None)

Bases: StatisticAccumulatorFactory

Factory for LogSeriesAccumulator.

Parameters:

keys (str | None)

make()

Create a fresh log-series accumulator.

Return type:

LogSeriesAccumulator

class LogSeriesEstimator(pseudo_count=None, suff_stat=None, name=None, keys=None)

Bases: ParameterEstimator

Maximum-likelihood estimator for the log-series shape p (inverts the mean -> p relation).

Parameters:

• pseudo_count (float | None)

• suff_stat (float | None)

• name (str | None)

• keys (str | None)

accumulator_factory()

Return an accumulator factory for log-series count statistics.

Return type:

LogSeriesAccumulatorFactory

estimate(nobs, suff_stat)

Estimate the shape parameter by inverting the weighted mean.

Parameters:

• nobs (float | None)

• suff_stat (tuple[float, float])

Return type:

LogSeriesDistribution

class LogSeriesDataEncoder

Bases: DataSequenceEncoder

Encode log-series observations as (k, log k) pairs.

seq_encode(x)

Encode observations as integer values and log-values.

Parameters:

x (Sequence[int])

Return type:

tuple[ndarray, ndarray]