Statistics

The statistical primitives shared by Combinations, Selectors and Stability monitoring. Every carver, selector and monitoring metric routes its arithmetic through this module, so a formula stated here is the one that runs.

Pearson’s \(\chi^2\)

AutoCarver.stats.pearson_chi2(observed: ndarray, *, guard_zero_expected: bool = False) → float

Pearson \(\chi^2\) of a (B, C) observed contingency table.

Replicates scipy.stats.chi2_contingency() defaults: expected frequencies via the outer product of marginals divided by N, with Yates correction iff the table is exactly 2x2 (matches scipy’s own threshold).

\[\chi^2 = \sum_{i, j} \frac{(O_{ij} - E_{ij})^2}{E_{ij}}, \qquad E_{ij} = \frac{n_{i.}\, n_{.j}}{n}\]

where \(n_{i.}\) and \(n_{.j}\) are the row and column marginals and \(n\) the grand total. When the table is exactly \(2 \times 2\), Yates’ continuity correction shrinks \(|O_{ij} - E_{ij}|\) by \(0.5\) before squaring (matches scipy’s own threshold for applying it).

guard_zero_expected replaces the 0/0 of an all-zero row or column with 0 instead of nan. The selector kernels need it (they build tables by bincount, which can produce empty rows); the combination evaluators must not use it — they shift every cell by +tol beforehand, and changing that would break bit-exactness.

Cramér’s \(V\) and Tschuprow’s \(T\)

AutoCarver.stats.cramerv_tschuprowt(chi2: float, n_obs: float, n_rows: int, n_cols: int, tol: float) → tuple[float, float]

Cramér’s V and Tschuprow’s T from a chi² computed on an (n_rows, n_cols) table.

\[V = \sqrt{\frac{\chi^2}{N (\min(B, K) - 1)}}, \qquad T = \sqrt{\frac{\chi^2}{N \sqrt{(B-1)(K-1)}}}\]

with \(B\) = n_rows, \(K\) = n_cols, \(N\) = n_obs. Both are NaN when their denominator vanishes (mirrors the binary/ordinal None-on-degenerate convention).

For n_cols == 2, T is instead derived from the (already rounded) V via \(T = V / \sqrt[4]{B - 1}\) — the exact expression the binary combination evaluator’s closed form uses. Both formulas are mathematically identical at K=2, but only computing it this way guarantees the binary and multiclass evaluators agree bit-for-bit (independent sqrt/pow sequences are not guaranteed to round identically) — pinned by the K=2 parity test.

AutoCarver.stats.cramerv_tschuprowt_unrounded(chi2: float, n_obs: float, n_mod_x: float, n_mod_y: float) → tuple[float, float]

Selector-side V / T: no tol quantisation, T from the raw chi² (not from V).

\[V = \sqrt{\frac{\chi^2}{N (\min(n_x, n_y) - 1)}}, \qquad T = \sqrt{\frac{\chi^2}{N \sqrt{(n_x-1)(n_y-1)}}}\]

with \(n_x\), \(n_y\) the two features’ modality counts. Different normalisation from cramerv_tschuprowt() — n_obs here is the non-missing pair count and there is no rounding — so this is kept as a separate function rather than merged into it; do not “simplify” the two into one.

Kruskal-Wallis’ \(H\)

AutoCarver.stats.tie_correction(values: ndarray) → float

Kruskal-Wallis tie factor; 1.0 when N < 2.

\[C_{tie} = 1 - \frac{\sum_t (t^3 - t)}{N^3 - N}\]

where the sum runs over each group of \(t\) tied values and \(N\) is the total sample size. Matches scipy.stats.tiecorrect (which takes ranks; ties are the same either way).

AutoCarver.stats.h_from_rank_sums(rank_sums: ndarray, counts: ndarray, n_obs: float, tie_corr: float) → float

Tie-corrected H from per-group rank sums and counts. nan when tie_corr == 0.

\[H = \frac{1}{C_{tie}} \left[\frac{12}{N(N+1)} \sum_g \frac{R_g^2}{n_g} - 3(N+1)\right]\]

where \(R_g\) and \(n_g\) are group g’s rank sum and count, \(N\) the total sample size, and \(C_{tie}\) the tie_correction() factor. nan also whenever any group is empty (0/0 propagates through the sum) — matches scipy.stats.kruskal, which rejects an empty sample.

Wilson frequency confidence bound

AutoCarver.stats.wilson_upper_bound(count: ndarray | int | float, nobs: int, alpha: float) → ndarray | float

Upper bound of the two-sided Wilson score interval for count / nobs.

\[\text{upper} = \frac{\hat{p} + \frac{z^2}{2n} + z\sqrt{\frac{\hat{p}(1-\hat{p})}{n} + \frac{z^2}{4n^2}}} {1 + \frac{z^2}{n}}\]

where \(\hat{p} = \text{count}/n\), \(n\) = nobs, and \(z\) is the two-sided normal quantile for alpha.

Parameters:

• count (array-like or scalar) – Observed successes. Accepts integer counts or float counts (e.g. weighted/aggregated frequencies).

• nobs (int) – Number of trials. Must be >= 0; returns 1.0 when nobs == 0 so callers treat empty samples as non-significant.

• alpha (float) – Two-sided significance level (e.g. 0.05 for a 95% interval).

Returns:

Wilson upper bound, same shape as count.

Return type:

array-like or scalar

AutoCarver.stats.is_significantly_below(count: ndarray | int | float, nobs: int, min_freq: float, alpha: float) → ndarray | bool

Whether the observed proportion count / nobs is significantly below min_freq.

A modality is significantly below min_freq when the Wilson upper bound of its observed proportion is strictly below min_freq.

Ridit scores

AutoCarver.stats.ridits_from_counts(counts: Series) → dict

Scores a count-marginal’s own levels, as a {level: ridit} dict.

Convenience wrapper over ridit_scores_for_levels() keeping the original (non-float-cast) level keys, so the result is directly y.map-able (the carver’s pre-sort scale).

AutoCarver.stats.ridit_scores_for_levels(levels, reference_counts: Series) → ndarray

Ridits of arbitrary numeric levels against a fixed train count-marginal.

\[\text{ridit}(j) = F(j^-) + \frac{f_j}{2}\]

where \(f_j\) is level \(j\)’s train frequency and \(F(j^-)\) the cumulative train frequency of all strictly lower levels. A level unseen in the reference gets \(F(j^-)\) alone (zero mass at that level).

Parameters:

• levels (iterable of numbers) – Levels to score (e.g. a crosstab’s columns) — need not all appear in the reference.

• reference_counts (pd.Series) – Train count-marginal, indexed by level (value_counts of the train target, or a train crosstab’s column totals). Order does not matter.

Returns:

One ridit per queried level: F(j-1) + f_j/2 for reference levels; a level unseen in the reference gets P_train(y < level) (the natural CDF extension: zero mass at that level), so tables carrying extra levels stay well-defined.

Return type:

np.ndarray

Rank association of an ordered table

For an ordered contingency table \((r \times c)\) — \(r\) feature groups (rows) × \(c\) ordinal target levels (cols), both ascending — three rank-association statistics are built from the same pair counts:

• \(C\) — concordant pairs (both members order the same way on the feature and on the target);

• \(D\) — discordant pairs (members order oppositely);

• \(P_0 = n(n-1)/2\) — all pairs, with \(n\) the number of observations;

• \(T_X\), \(T_Y\) — pairs tied on the feature / on the target (equal row / equal column); \(P_0 - T_X\) and \(P_0 - T_Y\) are the pairs untied on each margin;

• \(m = \min(r', c')\) — the smaller of the number of non-empty grouped rows \(r'\) and target levels \(c'\).

The three measures are monotone-comparable transforms of \(C - D\). Each is None for a degenerate table (its denominator vanishes), mirroring the continuous evaluator’s None convention. Parity against scipy.stats.kendalltau() (tau-b) and scipy.stats.somersd() is pinned by tests/combinations/ordinal/test_ordinal_associations.py and the property suite tests/properties/combinations/test_ordinal_combinations_properties.py.

AutoCarver.stats.concordant_minus_discordant(values: ndarray) → float

Concordant minus discordant pairs \(C - D\) of an ordered table.

\[C - D = \sum_{i,j} n_{ij} \left(\sum_{k>i,\, l>j} n_{kl} - \sum_{k>i,\, l<j} n_{kl}\right)\]

values is the (r, c) cell-count array with rows / columns already ascending. Computed in closed form from the table’s cumulative cell sums.

AutoCarver.stats.rank_associations(values: ndarray) → dict[str, float | None]

Kendall’s tau-b, tau-c and Somers’ D D(Y|X) for an ordered table.

values is the (r, c) cell-count array with rows = X (feature groups) and columns = Y (target levels), both already in ascending order. Each measure is None when its denominator vanishes (degenerate table), mirroring the continuous evaluator’s None convention.

AutoCarver.stats.rank_associations_from_counts(cd: float, n: float, untied_on_feature: float, untied_on_target: float, m: int) → dict[str, float | None]

Assembles tau-b, tau-c and Somers’ D from pre-computed pair counts.

\[\tau_b = \frac{C - D}{\sqrt{(P_0 - T_X)(P_0 - T_Y)}}, \qquad \tau_c = \frac{2 \, m \, (C - D)}{n^2 \, (m - 1)}, \qquad D(Y \mid X) = \frac{C - D}{P_0 - T_X}\]

with cd the \(C - D\) count, untied_on_feature \(P_0 - T_X\), untied_on_target \(P_0 - T_Y\), and m the smaller of the number of non-empty rows and columns. tau_b matches scipy.stats.kendalltau(), tau_c applies Stuart’s rectangular-table correction and somersd is the original asymmetric Somers’ D D(Y|X).

Shared by the closed form (rank_associations()) and the ordinal combination DP so both produce bit-identical values. Each measure is None when its denominator vanishes.

Kendall/Stuart’s \(\tau_c\) (ordinal default)

Stuart’s tau-c applies a \(\min(r, c)\) correction tailored to rectangular tables — exactly our shape (few feature groups × many target levels):

\[\tau_c = \frac{2 \, m \, (C - D)}{n^2 \, (m - 1)}.\]

Because the denominator depends only on \((n, m)\) and not on how observations distribute across groups, its magnitude stays comparable across combinations with different group counts. It self-balances toward fewer, robust modalities, only adding one when a split is genuinely discriminative — like Tschuprow’s T and the Kruskal effect sizes. This is the default for OrdinalCarver.

Kendall’s \(\tau_b\)

Kendall’s tau-b normalises \(C - D\) by the geometric mean of the two margins’ untied pairs:

\[\tau_b = \frac{C - D}{\sqrt{(P_0 - T_X)(P_0 - T_Y)}}.\]

It is bit-exact with the tau-b variant of scipy.stats.kendalltau() on the grouped table and tends to retain more modalities on smoothly monotone signals than \(\tau_c\).

Somers’ D

The original asymmetric Somers’ D D(Y|X) — concordant minus discordant pairs over pairs untied on the feature \(X\):

\[D(Y \mid X) = \frac{C - D}{P_0 - T_X}.\]

It matches scipy.stats.somersd(table).statistic. Being asymmetric it leans strongly toward the coarsest split (its maximum over groupings is typically two modalities); offered for users who specifically want raw Somers’ D rather than the self-balancing Kendall taus.

Correspondence analysis

AutoCarver.stats.fit_ca_axis(xtab: DataFrame, tol: float = 1e-10) → CAAxis

Fits the correspondence-analysis first axis of a crosstab.

xtab is an (n_rows, K) count crosstab (rows = modalities/groups, columns = target classes). Standardizes the table (row/column mass normalization), takes its SVD, and returns the (sign-fixed) first right singular vector plus the column mass vector needed to project any row’s own profile onto it (see ca_row_scores()).

AutoCarver.stats.ca_row_scores(xtab: DataFrame, axis: CAAxis) → Series

Projects each row of xtab onto a fixed CAAxis.

\[\text{score}_i = \sum_k \frac{p_{ik} - c_k}{\sqrt{c_k}}\, v_{1k}\]

where \(p_{ik}\) is row \(i\)’s own profile (proportions across columns), \(c_k\) the fixed (training) column mass, and \(v_1\) the fitted first right singular vector. Only the row’s own profile and the fixed (training) column masses / axis are needed, so this is well-defined for any row set sharing xtab’s columns — including a dev-sample grouping the axis was never fit on, or a carver’s grouped candidate table.

Falls back to (deterministic) descending-frequency scoring when axis.degenerate (encoded as -row_total so ascending sort still yields frequency-descending order).

Raises:

ValueError – When xtab doesn’t carry exactly the classes the axis was fit on — typically a target class present in a later sample but unseen at fit time. The axis is fixed by construction, so such a table cannot be projected onto it (see AutoCarver.stability).

class AutoCarver.stats.CAAxis(col_mass: ndarray, v1: ndarray, degenerate: bool = False)

A fixed correspondence-analysis first axis, reusable to score new rows.

degenerate is True when the training table carried too little structure to define a meaningful axis (fewer than 3 rows, fewer than 2 columns, or an ~zero first singular value); callers then fall back to a frequency-based order (see ca_row_scores()).