Skip to content

Edge Statistics

The edge statistic itself: one vectorised OLS GLM covering Pearson and Spearman and their covariate-controlled forms, plus ranks, residualisation and the exact t-distribution p-values of the edge coefficient.

p-values are exact. The pipeline does not compute one per edge — that is a CPU call per edge and permutation — but thresholds the t statistic against the exact critical value of each p-threshold (critical_t), which selects the same edges.

This is what UnivariateEdgeSelection(selection_statistic=..., selection_input=...) dispatches to — see Confound control for which combination means what.

The edge-level statistics.

One vectorised OLS GLM underlies every edge statistic CCCPM offers -- Pearson, Spearman, point-biserial and their partial variants are all the same fit with different preprocessing of the inputs (see edge_statistics). This module is that fit, plus the two transforms it composes with: rank conversion and confound residualisation, and the exact t-distribution p-values of its coefficient.

The fit stays on the device and yields a t statistic per edge. p-values are exact (scipy.special.stdtr) and computed on the CPU, which is too slow to do for every edge x permutation in a pipeline run -- so the pipeline never does. It thresholds on the t statistic instead, against the exact critical value of each p-threshold (see critical_t and PThreshold.select), which selects exactly the edges an exact p-value would.

Split out of edge_selection.py so that the maths sits apart from the policy that consumes it -- which edges to keep, and how the selectors are configured.

correlations_and_pvalues(X, Y_perms, correlation_type='pearson', confounds=None)

edge_statistics with exact p-values in place of the t statistic.

For direct use and for tests. The pipeline itself never materialises a p-value per edge -- see critical_t.

Returns: r_matrix: (N_features, N_perms) — semi-partial correlation (effect size) p_matrix: (N_features, N_perms) — exact two-sided p-value of the edge coefficient

critical_t(alpha, df)

The |t| above which a two-sided p-value falls below alpha.

p < alpha iff |t| > critical_t(alpha, df), because the two-sided p-value is strictly decreasing in |t| for fixed df. One scalar scipy call per threshold; the comparison itself stays on the device.

edge_statistics(X, Y_perms, correlation_type='pearson', confounds=None)

Univariate edge selection as a vectorised OLS GLM, batched over permutations.

For each edge (column of X) and each (permuted) target (column of Y_perms) this fits, in one batched linear-algebra pass on CPU/GPU, the linear model

Text Only
target ~ intercept [+ confounds] + edge

and returns the edge effect and the t statistic of its coefficient. By the Frisch–Waugh–Lovell theorem the coefficient only requires residualising the edge on the confounds — the target is never residualised to obtain the coefficient (its raw values drive it; this is what we want when the target is the thing we ultimately predict). This single path is mathematically identical to:

  • Pearson correlation (continuous target, no confounds)
  • point-biserial correlation (binary 0/1 target, no confounds)
  • the partial-correlation / coefficient F-test (confounds present)
  • Spearman, when X and Y are rank-transformed first.

A binary (0/1) target needs no special handling: regressing it on an edge is the linear-probability model, whose coefficient test equals the point-biserial correlation. (Its homoskedastic p-values are the conventional point-biserial ones; OLS standard errors are not heteroskedasticity-robust, which is standard for this kind of screening filter.)

The reported r is the semi-partial correlation (confounds removed from the connectome edge only, not the target). Its sign and the p-value are those of the regression coefficient, so the choice of semi-partial vs partial affects only the reported effect-size magnitude, never which edges are selected.

Args: X: (N_samples, N_features) — continuous edge values (fixed across perms) Y_perms: (N_samples, N_perms) — target(s), one column per permutation correlation_type: 'pearson' (linear / point-biserial) or 'spearman' (ranks) confounds: optional (N_samples, N_confounds); when given, the edge effect controls for these covariates.

Returns: r_matrix: (N_features, N_perms) — semi-partial correlation (effect size) t_matrix: (N_features, N_perms) — t statistic of the edge coefficient df: int — its degrees of freedom, N - 2 - N_confounds

get_residuals(data, confounds)

Regresses out 'confounds' from 'data' using OLS and returns the residuals.

Args: data: (..., N_samples) or (N_samples, ...) The target data (can be X or Y). Can be numpy or torch. confounds: (N_samples, N_confounds). Can be numpy or torch.

Returns: residuals: Same shape and type as data

t_pvalues(t, df)

Exact two-sided p-values of t statistics with df degrees of freedom.

Computed with scipy on the CPU and returned on t's device and dtype. Exact, but ~4 s for a 35,778-edge x 1000-permutation batch, which is why the pipeline thresholds with critical_t instead of calling this.

torch_rankdata(data, dim=-1)

Computes ranks of the data along a given dimension. Equivalent to scipy.stats.rankdata (method='ordinal') but fully vectorized on GPU.