mixle.analysis.coverage module

Sampling completeness, species richness, and diversity from frequency counts.

“How much of the population have I actually seen?” The same estimand recurs across fields: unseen probability mass in language models, undiscovered species in ecology, unobserved rare classes in machine learning. The answer comes from the frequencies of frequencies – especially how many items were seen exactly once (singletons) or twice (doubletons), which carry the signal about what is still missing.

• turing_coverage() / good_turing() – sample coverage and the Good–Turing discounting that reallocates probability mass to unseen items (Good 1953; Gale & Sampson’s Simple Good–Turing).

• chao1() / chao2() – nonparametric lower-bound richness estimators from abundance (Chao1) or replicated incidence (Chao2) data, with standard errors and log-normal CIs.

• ace() / ice() – abundance/incidence coverage-based richness estimators (rare-species corrected).

• hill_numbers() – the unified diversity profile (q=0 richness, q=1 exp-Shannon, q=2 inverse Simpson).

• rarefaction_curve() – expected richness as a function of sample size (Hurlbert interpolation), the basis for coverage-standardised comparison.

Counts must be a one-dimensional vector of finite non-negative integer abundances per species (validated by _abund(), which raises rather than silently coercing values or flattening axes); incidence inputs must be a finite binary (species, sites) 0/1 matrix with at least one site (validated by _incidence()). Estimators raise CoverageInsufficientDataError rather than issuing contradictory numerical answers when no observations are present.

exception CoverageInsufficientDataError

Bases: ValueError

The requested coverage/diversity estimand is not identified by an empty sample.

turing_coverage(counts)

Turing’s sample-coverage estimate and the complementary unseen probability mass.

C = 1 - f1 / n where f1 is the number of singletons and n the total count: the estimated probability that the next observation is a previously seen species. 1 - C = f1/n is the Good–Turing estimate of the total probability of all unseen species.

Returns:

{'coverage', 'unseen_mass', 'n', 'f1'}.

Parameters:

counts (ndarray)

Return type:

dict[str, float | int]

good_turing(counts)

Simple Good–Turing smoothed probabilities (Gale & Sampson 1995).

Reallocates probability from seen to unseen items using the frequencies of frequencies. Empirical Turing estimates r* = (r+1) N_{r+1}/N_r are used for small r and a smoothed log-linear fit S(r) takes over once the two diverge (the Gale switch), giving stable discounts in the sparse tail.

The log-linear fit needs at least two distinct abundance classes to be identifiable. With exactly one – e.g. an all-singleton sample like [1, 1, 1], a canonical Good–Turing input, not a malformed one – there is nothing to fit a slope against, so smoothing is skipped and the seen mass is reallocated by raw frequency instead (the Simple Good–Turing small-support fallback; previously this raised a LinAlgError out of np.polyfit). With no observed counts at all, coverage is not identifiable and the common CoverageInsufficientDataError is raised.

Parameters:

counts (ndarray) – per-species abundances (zeros ignored).

Returns:

{'p0', 'proba', 'r_star', 'r', 'insufficient_evidence', 'reason'} – p0 is the total probability assigned to unseen species; proba are the smoothed probabilities of the input species (aligned to the positive entries of counts, summing to 1 - p0); r_star / r are the discounted and raw frequencies for the distinct abundance classes. insufficient_evidence is always False for a returned estimate; empty samples raise CoverageInsufficientDataError.

Return type:

dict[str, ndarray | float | bool | str]

chao1(counts, *, ci_level=0.95)

Chao1 nonparametric richness estimator from abundance data (bias-corrected).

S_chao1 = S_obs + f1 (f1 - 1) / (2 (f2 + 1)) (Chao 1984, bias-corrected form), a lower bound on total richness driven by the singleton (f1) and doubleton (f2) counts. Returns a standard error and a log-normal confidence interval for the number of undetected species (Chao 1987), so the interval respects S_chao1 >= S_obs.

Returns:

{'estimate', 'observed', 'f1', 'f2', 'se', 'ci_low', 'ci_high'}.

Parameters:

• counts (ndarray)

• ci_level (float)

Return type:

dict[str, float]

chao2(incidence, *, ci_level=0.95)

Chao2 richness estimator from replicated incidence (presence/absence) data.

Parameters:

• incidence (ndarray) – finite binary (n_species, n_sites) 0/1 matrix (validated by _incidence(); every entry must be exactly 0 or 1). At least _MIN_REPLICATED_SITES sampling units are required – Chao2 estimates unseen richness from replication across units, so one unit does not identify it (CoverageInsufficientDataError).

• ci_level (float) – confidence level for the log-normal interval.

Returns:

{'estimate', 'observed', 'q1', 'q2', 'se', 'ci_low', 'ci_high', 'sites'} where q1/q2 are the numbers of species found in exactly one / two sites.

Return type:

dict[str, float]

ace(counts, *, rare_threshold=10)

ACE: Abundance-based Coverage Estimator of richness (Chao & Lee 1992).

Splits species into abundant (> rare_threshold) and rare, estimates sample coverage from the rare group’s singletons, and corrects for the coefficient of variation of the rare abundances.

Returns:

{'estimate', 'observed', 's_rare', 's_abund', 'c_ace'}.

Raises:

CoverageInsufficientDataError – when the estimated sample coverage c_ace is zero (every rare individual is a singleton) – the 1 / c_ace correction is undefined there and unseen richness is unbounded, so no finite estimate is reported.

Parameters:

• counts (ndarray)

• rare_threshold (int)

Return type:

dict[str, float]

ice(incidence, *, rare_threshold=10)

ICE: Incidence-based Coverage Estimator of richness (the Chao–Lee estimator for incidence data).

Parameters:

• incidence (ndarray) – finite binary (n_species, n_sites) 0/1 matrix (validated by _incidence(); every entry must be exactly 0 or 1). At least _MIN_REPLICATED_SITES sampling units are required, for the same identifiability reason as chao2().

• rare_threshold (int) – species found in <= rare_threshold sites are treated as infrequent.

Returns:

{'estimate', 'observed', 's_infreq', 's_freq', 'c_ice'}.

Raises:

CoverageInsufficientDataError – for a single-unit design, or when the estimated incidence coverage c_ice is zero (every infrequent species seen in exactly one unit) – the correction term s_infreq / c_ice is undefined there.

Return type:

dict[str, float]

hill_numbers(counts, q=(0.0, 1.0, 2.0))

Hill numbers (effective number of species) of order q.

The unified diversity profile: q=0 is observed richness, q=1 is the exponential of Shannon entropy, q=2 is the inverse Simpson concentration. Larger q weights common species more, so the profile D(q) summarises evenness as well as richness.

Parameters:

• counts (ndarray) – per-species abundances.

• q (float | ndarray) – a scalar order or an array of orders.

Returns:

Array of Hill numbers, one per requested order (scalar input still returns a length-1 array).

Return type:

ndarray

rarefaction_curve(counts, sizes=None)

Individual-based rarefaction: expected richness when subsampling m individuals (Hurlbert).

E[S(m)] = sum_i (1 - C(n - x_i, m) / C(n, m)) – the expected number of species seen in a random subsample of m of the n individuals. Used to compare richness between samples at a common sample size (or coverage).

Parameters:

• counts (ndarray) – per-species abundances.

• sizes (ndarray | None) – subsample sizes m to evaluate; defaults to 1 .. n. Must be exact integers in [0, n] (validated by _rarefaction_sizes()); the result is returned in the same order as sizes (not sorted or deduplicated).

Returns:

{'sizes', 'expected_richness'}.

Return type:

dict[str, ndarray]