lineage.states.StateDistributionCR

Competing-risks state distributions.

The Gamma/GaPhs emissions in this package treat a cell's fate (Bernoulli) and its phase duration (Gamma) as independent observations, and discard the duration of any cell that dies. That throws away every death time and, more subtly, mis-states the likelihood of a time-censored cell: "no event yet at time t" is written as P(division > t) when it should be P(division > t AND death > t).

Here each phase instead carries two latent clocks,

T_D ~ Gamma(a, s) division / transition T_X ~ Gamma(a_x, s_x) death

and we observe min(T_D, T_X) together with an indicator of which fired. The three likelihood cases are the standard competing-risks ones:

transition seen at t f_D(t) * S_X(t) death seen at t f_X(t) * S_D(t) censored at t S_D(t) * S_X(t)

The division probability is then derived, P(divide) = int f_D(t) S_X(t) dt, rather than fit as a free Bernoulli parameter, so the death fraction and the death timing are forced to agree. With the G1 death clock pinned to a constant hazard (shape 1) this costs no degrees of freedom relative to the Bernoulli/Gamma model.

Attributes

TIME_FLOOR

Classes

StateDistribution

One cell-cycle phase with competing division and death clocks.

StateDistributionPhase

G1 and G2 phases, each with its own pair of competing clocks.

Functions

_gamma_logpdf(t, a, scale)

log Gamma(a, scale) pdf, vectorized directly against scipy.special.

_gamma_logsf(t, a, scale)

log Gamma(a, scale) survival function, vectorized directly against scipy.special.

event_masks(x)

Split cells into the three competing-risks cases.

exponential_estimator(obs, events, weights, param_idx, K)

Weighted right-censored MLE for exponential scales, one per group.

fit_clocks(distributions, x_list, gammas_list, state_j)

Fit both clocks of one state, sharing each shape across conditions.

atonce_estimator(all_tHMMobj, x_list, gammas_list, phase)

M step across several conditions at once, matching the Gamma model's interface.

Module Contents

lineage.states.StateDistributionCR.TIME_FLOOR = 1e-10
lineage.states.StateDistributionCR._gamma_logpdf(t, a, scale)

log Gamma(a, scale) pdf, vectorized directly against scipy.special.

Parameters:

  • t (numpy.ndarray)

  • a (float)

  • scale (float)

Return type:

numpy.ndarray

lineage.states.StateDistributionCR._gamma_logsf(t, a, scale)

log Gamma(a, scale) survival function, vectorized directly against scipy.special.

Parameters:

  • t (numpy.ndarray)

  • a (float)

  • scale (float)

Return type:

numpy.ndarray

lineage.states.StateDistributionCR.event_masks(x)

Split cells into the three competing-risks cases.

x has the usual three per-phase columns [fate, duration, censoring], where fate is 1 for surviving the phase, 0 for dying in it, and NaN when unknown, and the censoring flag is 1 when the phase was seen through to its end.

Cells whose duration is negative have been masked for cross validation, and cells with a NaN duration never entered the phase; both are excluded everywhere.

Returns:

boolean masks for (division observed, death observed, censored)

Parameters:

x (numpy.ndarray)

Return type:

tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray]

lineage.states.StateDistributionCR.exponential_estimator(obs, events, weights, param_idx, K)

Weighted right-censored MLE for exponential scales, one per group.

With a constant hazard the MLE is total time at risk over events observed, which needs no iteration. A pseudocount keeps a group with no observed deaths finite.

Parameters:

  • obs (numpy.ndarray)

  • events (numpy.ndarray)

  • weights (numpy.ndarray)

  • param_idx (numpy.ndarray)

  • K (int)

lineage.states.StateDistributionCR.fit_clocks(distributions, x_list, gammas_list, state_j)

Fit both clocks of one state, sharing each shape across conditions.

distributions is either a single StateDistribution or a list of them, one per condition. The shape parameters are shared across conditions and the scales are free, mirroring how atonce_estimator treats the Gamma model.

Parameters:

  • x_list (list[numpy.ndarray])

  • gammas_list (list[numpy.ndarray])

  • state_j (int)

lineage.states.StateDistributionCR.atonce_estimator(all_tHMMobj, x_list, gammas_list, phase)

M step across several conditions at once, matching the Gamma model's interface.

Parameters:

  • all_tHMMobj (list)

  • x_list (list)

  • gammas_list (list[numpy.ndarray])

  • phase (Literal['all', 'G1', 'G2'])

class lineage.states.StateDistributionCR.StateDistribution(gamma_a=7.0, gamma_scale=4.5, death_a=1.0, death_scale=40.0, fixed_death_shape=True)

One cell-cycle phase with competing division and death clocks.

params is [bern_p, gamma_a, gamma_scale, death_a, death_scale]. The first three entries keep the meaning they have in StateDistribution, so downstream figure code that indexes them positionally continues to work; bern_p is now derived from the two clocks rather than fit.

Parameters:

  • gamma_a (float)

  • gamma_scale (float)

  • death_a (float)

  • death_scale (float)

  • fixed_death_shape (bool)

atonce_estimator
fixed_death_shape = True
params
property div_clock
property death_clock
division_probability()

P(the division clock fires first) = int f_D(t) S_X(t) dt.

quad evaluates the integrand at ~100-200 points per call, and this runs once per state per M step. Using the frozen rv_continuous machinery for each point (as div.pdf / death.sf do) costs orders of magnitude more than the arithmetic itself, so the integrand is written directly against scipy.special instead.

Return type:

float

rvs(size, rng=None)

Draw min(T_D, T_X) and record which clock fired.

Parameters:

size (int)

dist(other)

Wasserstein distance between the division clocks of two states.

Kept on the division clock alone so that the number is comparable with the Bernoulli/Gamma model's.

Return type:

float

dof()

Two division-clock parameters plus the death clock's scale, and its shape when that is not pinned. The Bernoulli is derived, so it is not counted.

Return type:

int

logpdf(x)

Competing-risks log likelihood of each cell's phase observation.

Called once per state on every E step, over the whole lineage array, so it is vectorized directly against scipy.special rather than going through the frozen rv_continuous machinery (self.div_clock.logpdf etc.), which pays generic validation/broadcasting overhead on every call regardless of array size.

Parameters:

x (numpy.ndarray)

Return type:

numpy.ndarray

estimator(x, gammas)

Weighted M step for a single condition.

The two clocks separate in the complete-data likelihood, so each is just a weighted right-censored Gamma fit over every timed cell.

Parameters:

  • x (numpy.ndarray)

  • gammas (numpy.ndarray)

censor_lineage_array(censor_condition, tree, obs, states, desired_experiment_time=2000000000000.0)

Applies censoring to array representation directly.

Parameters:

  • censor_condition (int)

  • tree (scipy.sparse.csr_array)

  • obs (numpy.ndarray)

  • states (numpy.ndarray)

Return type:

tuple[scipy.sparse.csr_array, numpy.ndarray, numpy.ndarray]

class lineage.states.StateDistributionCR.StateDistributionPhase(gamma_a1=7.0, gamma_scale1=3.0, gamma_a2=14.0, gamma_scale2=6.0, death_scale1=40.0, death_a2=3.0, death_scale2=20.0)

G1 and G2 phases, each with its own pair of competing clocks.

params is [bern_p1, bern_p2, a1, s1, a2, s2, death_a1, death_s1, death_a2, death_s2]. The leading six entries match StateDistribution exactly.

Parameters:

  • gamma_a1 (float)

  • gamma_scale1 (float)

  • gamma_a2 (float)

  • gamma_scale2 (float)

  • death_scale1 (float)

  • death_a2 (float)

  • death_scale2 (float)

atonce_estimator
G1
G2
params
_sync()

Mirror the sub-distributions' parameters into the flat params array.

rvs(size, rng=None)
Parameters:

size (int)

dist(other)
Return type:

float

dof()
Return type:

int

logpdf(x)
Parameters:

x (numpy.ndarray)

Return type:

numpy.ndarray

estimator(x, gammas)
Parameters:

  • x (numpy.ndarray)

  • gammas (numpy.ndarray)

censor_lineage_array(censor_condition, tree, obs, states, desired_experiment_time=2000000000000.0)
Parameters:

  • censor_condition (int)

  • tree (scipy.sparse.csr_array)

  • obs (numpy.ndarray)

  • states (numpy.ndarray)

Return type:

tuple[scipy.sparse.csr_array, numpy.ndarray, numpy.ndarray]