Statistics and comparison brackets

Two-group comparisons#

Each two-group test takes two samples (a, b), orients signs as a − b, and returns one reporting row: a dict ready to save as a CSV in the plot’s _stats dissection.

import polars as pl
import fluxplot as fp

# sd_values and sleep_values are the two samples from your analysis.
row = fp.stats.welch_hedges(sd_values, sleep_values, names=("SD", "sleep"))
pl.DataFrame([row]).write_csv("plots/_dissections/app_gapdh/_stats/welch_ttest.csv")

The tests are two-sided unless you pass alternative=. A row’s keys, in the order of fp.stats.REPORT_COLUMNS, are:

• Test and result: sig_test_used, test_statistic_value, p-value, p_corrected_holm, dof.
• Effect size: effect_size_method, effect_size_value, effect_size_95_CI, and the interval as numbers in effect_size_ci_low and effect_size_ci_high.
• Samples: n_a, n_b, n_total, and groups (the names compared).
• Also recorded: alternative, p_corrected_bh, dof_error.

Each reports this effect size with a 95% confidence interval:

The CIs target the population effect size, so they are never multiplied by Hedges’ J. Both rank-based CIs stay non-degenerate at complete separation (delta or r = ±1), which is common at n = 6. Rank tests report dof = None. fluxplot’s test suite (tests/test_stats.py) checks each test against SciPy and against the equation that defines its interval.

Multiple comparisons#

A row on its own has p_corrected_holm == p-value (a family of one). Holm’s step-down correction needs every p-value in the family, so it is a separate pass over the rows of the comparisons that belong together (e.g. every measure tested on the same animals in one figure); extra keys such as a measure name ride along:

# sd and sleep map each measure name to its samples.
rows = [dict(measure=m, **fp.stats.welch_hedges(sd[m], sleep[m])) for m in measures]
# Fill p_corrected_holm across the family; the row order is kept.
rows = fp.stats.holm(rows)

fp.stats.holm_adjusted(p) does the same for a bare array of p-values; fp.stats.bh(rows) / bh_adjusted(p) fill p_corrected_bh with Benjamini–Hochberg (false-discovery-rate) values for a screen of many measures. NaN p-values pass through both untouched.

Three or more groups#

The omnibus tests return one row, in which groups lists every group and F tests carry dof and dof_error. The post-hoc tests return one row per pair. pairwise runs any two-group test above over the pairs of a {name: sample} family and corrects across them.

Notes on the tests:

• rm_anova reports Greenhouse–Geisser-corrected degrees of freedom and p.
• tukey_hsd reports Tukey’s family-wise p as it is, without a further correction.
• games_howell uses the studentized range on Welch’s degrees of freedom.
• dunn is Dunn’s rank-sum test after Kruskal–Wallis, corrected with adjust="holm" by default.
• pairwise corrects across the pairs with adjust="holm" (the default) or "bh".

The group-based functions take the samples as positional arguments with optional names=, for example fp.stats.welch_anova(ctl, drug, sham, names=("ctl", "drug", "sham")). The repeated-measures functions take a long table and its column names: fp.stats.rm_anova(table, subject, within, dv). Exact signatures are in the API reference. fluxplot’s test suite (tests/test_stats_multi.py) checks them against pingouin and scikit-posthocs reference values.

Draw the comparisons#

fp.brackets(ax, rows, positions=…) draws one significance bracket per post-hoc row and stacks them automatically: shortest span first, each bracket one step above the data it spans and above every bracket it overlaps in x (multiplicative steps on a log axis), so nothing crosses. positions maps group name → x; a glowbar / fluxbox result provides it, and gb.brackets(rows) is the one-liner:

# ctl, drug and sham are the samples; df["group"] holds the same three names.
rows = fp.stats.pairwise(fp.stats.welch_hedges,
                         {"ctl": ctl, "drug": drug, "sham": sham})
gb = fp.glowbar(data=df, x="group", y="value", ax=ax)

# Stars (***, **, *, ns) from p_corrected_holm:
gb.brackets(rows)
# Or "p = 0.003" labels from the BH-corrected p-values, without the ns pairs:
gb.brackets(rows, label="p", ns=False, p_column="p_corrected_bh")

label is "stars", "p", "both" or a callable on the row. thresholds, top, step and tip shape the stack, and extra keywords are passed to fp.significance_bracket. Each bracket’s manifest overlay carries between: [a, b], p and a stats block — test, statistic, raw and corrected p, correction, effect size and its CI, sizes — so the figure states exactly which test each star came from.