mixle.stats.processes.birth_death module¶

Evaluate, estimate, and sample from a general birth-death-sampling process.

Defines BirthDeathSamplingDistribution, BirthDeathSamplingSampler, BirthDeathSamplingAccumulatorFactory, BirthDeathSamplingAccumulator, BirthDeathSamplingEstimator, and BirthDeathSamplingDataEncoder.

A continuous-time linear birth-death process on a population n(t): each individual independently gives birth at rate birth_rate, dies at rate death_rate, and is sampled (observed through time, without removal) at rate sampling_rate. This is a general population/epidemic model; the fossilized birth-death model is the special case where sampling_rate is the fossilization rate. Pure birth-death is sampling_rate = 0.

Data type: one fully-observed trajectory (n0, T, events) – initial count n0, observation window length T, and a time-ordered list of (time, type) events with type in {0: birth, 1: death, 2: sampling}. The log-likelihood is

sum_events log n_i + n_b log(birth) + n_d log(death) + n_s log(sampling) - (birth+death+sampling) * I,

where n_i is the population just before event i and I = integral_0^T n(t) dt (n is piecewise constant between events). The MLE is closed-form: each rate is its event count divided by I (summed over trajectories).

Reference: Feller, An Introduction to Probability Theory and Its Applications, Vol. 1 (Wiley).

class BirthDeathSamplingDistribution(birth_rate, death_rate, sampling_rate=0.0, initial_population=1, horizon=10.0, name=None, keys=None)[source]

Bases: SequenceEncodableProbabilityDistribution

General linear birth-death-sampling process (fossilized birth-death is the sampling_rate>0 case).

Parameters:

birth_rate (float)
death_rate (float)
sampling_rate (float)
initial_population (int)
horizon (float)
name (str | None)
keys (str | None)

density(x)[source]

Probability density of one trajectory (see log_density).

Parameters:: x (Any)
Return type:: float

log_density(x)[source]

Log-likelihood of one fully-observed trajectory (n0, T, events).

Parameters:: x (Any)
Return type:: float

seq_log_density(x)[source]

Vectorized log-likelihood for an (N, 5) array of per-trajectory sufficient statistics.

Parameters:: x (ndarray)
Return type:: ndarray

sampler(seed=None)[source]

Return a BirthDeathSamplingSampler for this distribution.

Parameters:: seed (int | None)
Return type:: BirthDeathSamplingSampler

estimator(pseudo_count=None)[source]

Return a BirthDeathSamplingEstimator (closed-form rate MLE).

Parameters:: pseudo_count (float | None)
Return type:: BirthDeathSamplingEstimator

dist_to_encoder()[source]

Return the encoder for birth-death samples.

Return type:: BirthDeathSamplingDataEncoder

class BirthDeathSamplingSampler(dist, seed=None)[source]

Bases: DistributionSampler

Exact Gillespie simulation of birth-death-sampling trajectories on [0, horizon].

Parameters:

dist (BirthDeathSamplingDistribution)
seed (int | None)

sample(size=None)[source]

Draw one trajectory (n0, T, events) or a list of size trajectories.

Parameters:: size (int | None)
Return type:: Any

class BirthDeathSamplingAccumulator(name=None, keys=None)[source]

Bases: SequenceEncodableStatisticAccumulator

Accumulate weighted event counts, the time-integral of the population, and trajectory count.

Parameters:

name (str | None)
keys (str | None)

update(x, weight, estimate)[source]

Accumulate weighted event counts and exposure for one trajectory.

Parameters:

x (Any)
weight (float)
estimate (BirthDeathSamplingDistribution | None)

Return type:

None

initialize(x, weight, rng)[source]

Initialize the sufficient statistics with one weighted trajectory.

Parameters:

x (Any)
weight (float)
rng (RandomState | None)

Return type:

None

seq_update(x, weights, estimate)[source]

Accumulate weighted event counts and exposure from encoded trajectories.

Parameters:

x (ndarray)
weights (ndarray)
estimate (Any | None)

Return type:

None

seq_initialize(x, weights, rng)[source]

Initialize the sufficient statistics from encoded trajectories.

Parameters:

x (ndarray)
weights (ndarray)
rng (RandomState | None)

Return type:

None

combine(suff_stat)[source]

Merge serialized birth-death-sampling statistics into this accumulator.

Parameters:: suff_stat (tuple[float, float, float, float, float, float])
Return type:: BirthDeathSamplingAccumulator

value()[source]

Return event counts, exposure, trajectory count, and horizon total.

Return type:: tuple[float, float, float, float, float, float]

from_value(x)[source]

Restore the accumulator from serialized birth-death-sampling statistics.

Parameters:: x (tuple[float, float, float, float, float, float])
Return type:: BirthDeathSamplingAccumulator

scale(c)[source]

Scale accumulated sufficient statistics by a constant.

Parameters:: c (float)
Return type:: BirthDeathSamplingAccumulator

key_merge(stats_dict)[source]

Merge this accumulator into a keyed statistics dictionary.

Parameters:: stats_dict (dict[str, Any])
Return type:: None

key_replace(stats_dict)[source]

Replace this accumulator from a keyed statistics dictionary.

Parameters:: stats_dict (dict[str, Any])
Return type:: None

acc_to_encoder()[source]

Return an encoder for trajectory sufficient statistics.

Return type:: BirthDeathSamplingDataEncoder

class BirthDeathSamplingAccumulatorFactory(name=None, keys=None)[source]

Bases: StatisticAccumulatorFactory

Factory for BirthDeathSamplingAccumulator.

Parameters:

name (str | None)
keys (str | None)

make()[source]

Create an empty birth-death-sampling accumulator.

Return type:: BirthDeathSamplingAccumulator

class BirthDeathSamplingEstimator(name=None, keys=None)[source]

Bases: ParameterEstimator

Closed-form rate MLE: rate = (total events of that type) / integral_n.

Parameters:

name (str | None)
keys (str | None)

accumulator_factory()[source]

Return a factory for birth-death-sampling sufficient-statistic accumulators.

Return type:: BirthDeathSamplingAccumulatorFactory

estimate(nobs, suff_stat)[source]

Estimate per-capita rates from accumulated event counts and population time-integral.

Parameters:

nobs (float | None)
suff_stat (tuple[float, float, float, float, float, float])

Return type:

BirthDeathSamplingDistribution

class BirthDeathSamplingDataEncoder[source]

Bases: DataSequenceEncoder

Encode trajectories (n0, T, events) into an (N, 6) array of sufficient statistics + T.

seq_encode(x)[source]

Encode trajectories as sufficient-statistic rows.

Parameters:: x (Sequence[Any])
Return type:: ndarray