Higher-order U-centering: ANOVA residualization and fast unbiased estimation

The unbiased sample versions of distance covariance and HSIC are fourth-order U-statistics, yet U-centering evaluates them from pairwise arrays in $O(n^2)$ operations. We show that U-centering is exactly the least-squares residual obtained after fitting additive endpoint effects to a symmetric hollow array. This interpretation explains the zero row sums and the denominator $n(n-3)$ through the residual degrees of freedom. We extend the construction to arrays indexed by $r$-subsets. Higher-order U-centering removes all effects involving fewer than $r$ sample labels, leaves zero $(r-1)$-way margins, and projects onto a residual space of dimension $\binom nr-\binom n{r-1}$. For two kernels with $r$ arguments, the normalized inner product of the centered arrays is unbiased for the pairing of their $r$th Hoeffding components. Although the corresponding direct estimator can involve up to $2r$ distinct observations, either subset-margin inversion or higher-order U-centering followed by a normalized inner product evaluates it in $O(n^r)$ operations for fixed $r$. The same calculation yields the classical unbiased Hoeffding variance-component estimators, with the highest component represented as a nonnegative residual mean square.