Learning Representations for Independence Testing
Authors: Nathaniel Xu
Feng Liu
Danica J. Sutherland
Venue: Accepted by TMLR ()

Abstract: Many tools exist to detect dependence between random variables, a core question across a wide range of machine learning, statistical, and scientific endeavors. Although several statistical tests guarantee eventual detection of any dependence with enough samples, standard tests may require an exorbitant amount of samples for detecting subtle dependencies between high-dimensional random variables with complex distributions. In this work, we study two related ways to learn powerful independence tests. First, we show how to construct powerful statistical tests with finite-sample validity by using variational estimators of mutual information, such as the InfoNCE or NWJ estimators. Second, we establish a close connection between these variational mutual information-based tests and tests based on the Hilbert-Schmidt Independence Criterion (HSIC); in particular, learning a variational bound (typically parameterized by a deep network) for mutual information is closely related to learning a kernel for HSIC. Finally, we show how to, rather than selecting a representation to maximize the statistic itself, select a representation which can maximize the power of a test, in either setting; we term the former case a Neural Dependency Statistic (NDS). While HSIC power optimization has been recently considered in the literature, we correct some important misconceptions and expand to considering deep kernels. In our experiments, while all approaches can yield powerful tests with exact level control, optimized HSIC tests generally outperform the other approaches on difficult problems of detecting structured dependence.

Published in Transactions on Machine Learning Research (03/2026)
Learning Representations for Independence Testing
Nathanicl Xu
sunathan@cs.ubc...
University of British Columbia
Feng Liu University of Melbourne
feng, liul @unimell. edu.
Danica J. Sutherland
dsuth Beuthh@cs.ubc.oo
University of British Columbia Alberta Machine Intelligence Institute
Reviewed OpenReview: https: openreview. net/ per/formm? id= odDDuooXRssn
Abstract
Many tools exist to detect dependence between random variables, core question across wide range machine learning, statistical, and scientific endeavors. Although several statistical tests guarantee eventual detection of any dependence with enough samples, standard tests may require an exorbitant amount of samples for detecting subtle dependencies between high-dimensional random variables with complex distributions. In this work, we study two related ways to learn powerful independence tests. First, we show how to construet powerful statistical tests with finite-sample validity by using variational estimators of mutual information, such as the InfoNCE or NWJ estimators. Second, we establish a close connection between these variational mutual information-baseed tests and tests based on the Hilbert- Schmidt Independence Critcrion (HSIC); in particular, learning a variational bound (typically parameterized by a deep network) for mutual information is closely related learning a kernel for HSIC. Finally, we show how to, rather than selecting a representation to maximize the statistic itself, select representation which maximize the power of test, in either setting; term the former case Neural Dependency Statistic (NDS). While HSIC power optimization has been recently considered in the literature, we correct some important misconceptions and expand to considering deep kernels. In our experiments, while all approaches can yield powerful tests with exact level control, optimized HSIC tests generally outperform the other approaches on difficult problems of detecting structured dependence.
1 Introduction
Independence testing, the question of using paired samples to determine whether a random variable X and another Y are associated with one another or if they are statistically independent, is one of the most common tasks across scientific and data-based fields. Traditional methods make strong parametric assumptions, for instance assuming that X and Y are jointly normal so that dependence is characterized by covariance, and/or operate only in limited scttings, for instance the tabular sctting of the celebrated x² test or Fisher's exact test. Applying these approaches to high-dimensional continuous data difficult at best.
One characterization of independence is the Shannon mutual information (MI): this quantity is zero if two random variables are independent, and positive if they are dependent. Substantial effort has been made in est iml this quantity with a variely of estimators; see e.g. the broad list based on binning, nearest-neighhbrr, kernel density estimation, and so on implemented by Szabó (2014). work has focused on estimating variational bounds defined by deep networks (sce c.g. Poole al., 2019), which can ideally learn problem-specific structure. To our knowledge, this class of estimators has not yet been used to construct statisticnl tests of independence: in particular, if we the estimation algorithm and
digh-dimensional cases, however, recent

---

Published in Transactions on Machine Learning Research (03/2026)
estimate a lower bound on the MI of 0.12, does that mean that the variables are (perhaps weakly) dependent, or might that be due to random noise on the samples we saw? The first contribution of this paper is to construct and evaluate tests of this form.
Mutual information, though, is notoriously difficult to estimate (Paninski, 2003). If the question we want to ask "are X and Y dependent,"" we can also consider turning to a different characterization of dependence which may be statistically easier to estimate. The Hilbert-Schmidt Independence Criterion (HSIC), introduced by Gretton et al. (2005), measures the total cross-covariance between feature representations of X and Y, and can be estimated from samples efficiently. The construction supports even infinite- minniteeddmeesionnll features in reproducing kernel Hilbert space (RKHS), equivalent to choosing a kernel function. With appropriate choices of kernel (see Szabó & Sriperumbudur, 2018), the HSIC is also zero if and only if X TL Y. Székely et al. (2007) and Lyons (2013) separately proposed distance covariance tests, which measure the covariance of pairwise distances between X values with distances between Y values; this can be viewed as HSIC with a particular kernel (Sejdinovic et al., 2013). Our second contribution describes how HSIC, in fact, is also a lower bound on mutual information, and that. lests based on either are very closely related.
Alternatively, we can characterize independence problem are concerned with the question: given samples from P and Q, is P = Q? Under this framework, when we consider samples from the joint distribution Pazy and the product of the marginals Pa x Py (where X L Y by construction), the two-sample problem characterizes an equivalent independence problem between variables X and Y. One way we can measure the discrepancy between Pay and P. x Py is with the kernel-based Maximum- Mean Discrepancy (MMD) (Gretton et al., 2012a); doing so recovers HSIC.
two-sample one. Two-sample problems
Any reasonable choice of kernel, such the Gaussian kernel with unit bandwidth or the distance kernel that yields distance covariance, will eventually (with enough samples) be able to detect any fixed dependence. In practice, however, this scheme can perform extremely poorly; if, for instance, the data varies on a very different scale than the bandwidth of the Gaussian kernel, it may take exorbitant quantitics of data to achieve any reasonable test power. Thus, tests using Gaussian kernels often rely on the median hewristic to choose a kernel relevant to the data at hand: choosing a bandwidth based on the median pairwise distance among data points (Gretton et al., 2012a). While this is a reasonable first guess for many data types, there exist many two-sample testing problems where it can be dramatically better to instead select a bandwidth that optimizes measure of lest power (Sutherland et al., 2017). Beyond that, there are many problems where no Gaussian kernel performs well - for instance, many problems on natural image data or involving sparsity but Gaussian kernels applied to latent representations of such data do. As an extreme example in the dependence case, consider the construction X ~ (0,, 1), Y, ~ (N(, but then Y is obtained by replacing the fifteenth decimal place of Yo with that of X. Representations based on Euclidean distance would require an absurd number of samples to detect the dependence between X and Y, but a representation that extracts only the fifteenth decimal place will make the dependence obvious.
These considerations apply similarly to tests based on mutual information. In fact, we can view estimation of variational MI bound as learning a representation of the data that maximizes its mcasure of dependency. Thus, learning representations of the data that make any dependency more explicit is central to developing more powerful independence tests. Our main algorithmic contribution involves developing a scheme to learn these optimal representations, both for HSIC tests and for a class of tests closcly related to the variational MI-based tests. In both settings, this is based on maximizing an estimate of the asymptotic power of the test, the primary term of which is an estimate of the signal-to-noise ratio of the statistic: the estimator divided by its standard deviation. In both cases, we can efficiently, and differentiably, estimate this quantity; we also show via uniform convergence argument that optimizing the power estimate leads to representation which generalizes well. Using data splitting and permutation testing (as in e.g. Rindt et al., 2021), we obtain a test which is exactly valid and achieves high power. In experiments, we find that while both methods often work well, there are settings where the HSIC optimization scheme far outperforms the other categories of tests.
In both cases, the overall approach builds on prior work in kernel two-sample testing (Gretton et al., 2012b; Jitkrittum et al., 2016; Sutherland et al., 2017; Liu et al., 2020). We show, however, that performing the direct reduction to two-sample lesting and applying these techniques, in addition to breaking the assumptions of their theoretical results, yields notably less powerful test than our direct HSIC power optimization, due to
2

---

Published in Transactions on Machine Learning Research (03/2026)
HSCC
HSIC
-4-
NGS
HNONCE
MA]
C25TL
C2315
NMD-5
MSICO
-4-
HBKC
H3KCAgg
empiniall power
empircaal power
6.8
empinical power
0.0
9.6
0.8
0.6
04
0.
9.4
0.
0.2
0.0
500
1000
15:00
2000
500
1006
1500
2000
500
1000
1500
2000
test samples
test samples
test samples
(a) Samples X, ~ Pay
(b) Power vs. m; d
(c) Power vs. m; 10
(d) Power vs. m; d 15
Figure 1: Vanilla kernel-based HSIC tests struggle on high-dimensional data. (a) A bimodal Gaussian mixture, with visible dependence between X and Y. (b)-(d) Test power as we add an increasing number of independent noise dimensions. When the median heuristic (HSIC-M) confidently detects dependence between and with a reasonably large sample size, but struggles when detects dependence high dimensions with many fewer samples.
or 15. Our method (HSIC-D)
statistical properties of the estimators. We also note that Ren et al. (2024) recently proposed a closely related scheme for the HSIC power optimization setting; we argue, however, that several of their main algorithmic suggestions are not well-supported (particularly in Appendix F.5) and identify an apparent gap in the proof of their main theorem (see Footnote 10).
2 Tests Based On Variational Mutual Information Bounds
domain X x y, and Pz, Py be the corresponding marginal distributions of X and y. We observe m paired samples X = (I1, and Y Y= (y..... vm) jointly drawn from Pzy. We wish to conduct a null hypothesis significance test, with null hypothesis 5o Pay P, Py and alternative 51 Pyy + P. p.
Problem 2.1 (Independence testing). Let 7 (X, Y) ~ Pay
We wish to solve this problem without making strong (parametric) assumptions about the form of IP zy. P., Py- Most independence tests based estimating the "amount" of dependence between X and Y, equivalently the discrepancy between between I Pay and P. x P.y, Given nonnegative quantity which is zero if X 1 Y, we can reject Ho if our estimate is large enough that we are confident the true value is positive,
The most famous such quantity is the Shannon mutual information I(; can be defined as the Kullback-Liebler divergence of Pay from P,, x y; and zero if and only X T Y; equivalently, it is the amount by which knowledge of Y decreases the entropy of X. While many estimators exist based roughly on various forms of density estimation (see e.g. Szabó, 2014), as discussed above these can fail detect subtle or structured forms of dependence with reasonable numbers of samples. Variational estimators of mutual information give an opportunity for problem-specific representations. As an example, consider the following lower bound (van den Oord et al., 2018; Poole et al., 2019), called InfoNCE for its connection to noiso-contrastive estimation (Gutmann & Hyvärinen, 2012):
K.K NCE (X:Y)
E
f..K NCE
S I(X;Y)
where
f.K NCE
Σ log K
1
(as
Σ ellerss)
Here K is a batch size, which in our setting will typically be the total number available points. Each f X' x y R (and each K) leads to different lower bound; the largest NCE the tightest bound. In practice, users typically parameterize / as deep network and maximize INE on minibatches from a training set, providing an opportunity for ff to learn useful feature adapted to the problem at hand. There are a variety bounds of this general type; Poole et al. (2019) give a unified accounting of many. Different bounds yield different bias-variance tradeoffs, but run up against various limitations on the possibility of estimation (Song& & Ermon, 2020; McAllester & Stratos, 2020).
f..K
For independence testing, the key question is whether the true value I(X;Y) > 0. This is guaranteed if lower bound, such NCE for some particular and K, is positive at the population level, i.e. with the true
3

---

Published in Transactions on Machine Learning Research (03/2026)
expectation. To estimate this, by far the easiest approach is data splitting: choose f to maximize if.K (minibatches from) a training set, then evaluate f on the heldout test set.
INCE
How large should KCE be in order to be confident that NCE > 0? That is, what does the distribution of if.K NCE look like when Y( = and hence NCE S O? For a given S, we can answer this question with permutation testing, which estimates values under the mull hypothesis (Fay Ps x Fy) by randomly shuffling the test data, breaking dependence. To construct test with probability of false rejection at most a, can compute the empirical 1 - a quantile from this permuted set, as long as we include the original paired data in this shuffling (Hemerik & Goeman, 2018, Theorem 2). We reject the null hypothesis if this quantile is smaller than the test statistic
I
K
i NCE
Σ
log
( j=l
(1)
Written in this form, notice that the second term of ICE is permutationn-nnvarian changing each yj with Yoy for some permutation a changes the value of the first term, but does not change the value of the second. Put another way, the test statistic and each of its permuted versions are * ΣE f(z,, 3:) shifted by the same constant. Thus, although this second term plays vital role in selecting the critic function f, at test time the only thing that matters is whether the mean value of f(x, is higher for the true pairings than for random pairs. We call this term the neural dependency statistic (NDS) and denote it for batch of K paired samples X (f.. Tk) and Y = (yI,..., yx) by
K
1 (X, = K (a.
Σ
(2)
The same is true for the NWJ (Nguyen et al., 2010), DV (Donsker & Varadhan, 1983), Lıs, and Ia lower bounds discussed by Poole et al. (2019), as well as the MINE estimator (Belghazi et al., 2018); that is, at test time only the NDS statistic (2) matters.
2.1 Asymptotic test power
Typically, a mutual information lower bound is considered botter the population value of the bound is larger: that makes it tighter bound. This viewpoint, however, neglects the issue of statistical estimation of that i bound, e.g. the difference between INE and INCE The best statistical test, among tests with appropriate Type (false rejection) control, is the one with highest test power: the probability of correctly rejecting the mull hypothesis 5o when X L Y. Taking this into account, we examine the behavior of permutation test based on (1) or the many other mutual information bounds based only on the NDS (2). Each bound will choose a different f during training, but at test time only the NDS matters.
Let Ts, Ex_xP, f(x, y) and Ts, Ee,, f(z, y) be the population level statistics under the mull and alternative distributions. Assuming fixed critic chosen independently from test samples X and W (e.g. because of sample splitting), and that 0 < 7 Varr,y f(x co, the central limit theorem implies that (f To.) 4 NN(O 1). Then, using & for the standard normal edf, the rejection threshold rm satisfies TSi
/mInoo =
=1-6 - (rm Vml'sso + o(1).
Pr Do (vmf Tm)
Pr 50
Tno
Tso
Doo
4

---

Published in Transactions on Machine Learning Research (03/2026)
0.6
0.9
аsymp emarrical
0.9
empiical perm
empirical perm
0.8
0.5
power
power
0.7
0.7
0.4
0.6
0.3
1200 1400
test size m
test size m
test size m
(a) HDGM-4
(b) Simnsoid
(c) RatIn ABox
Figure 2: Power of NDS tests with a particular critic f (learned by the approach of Section 5) on three datasets described in Section For each problem, the blue line (asymp) is power described by the asymptotic formula (3) with the o(1) error terms set to zero, where we use 20 000 samples to estimate the population-level statistics T5 Tso and their standard deviations Th T50- The orange line (empirical) computes the power based on the simulated null distribution, i.c. the threshold estimated via the empirical CDF of To with test size m. The green line (perm) is the permutation test power. All three lines roughly agree with each other.
As is Lipschitz, this implies rm = /m + Tso T5o (1- a) + o1)). Using the asymptotic normality of T under 51, we then have that
Pr (vmit > rm) 51
= Pr Bi
VII
- Tss
I'
mmTs
Thi
T5i
v /ml Do
+ Tso
10-1(1 (1-a)) (1- + o(())
vnIt.
= =1-0
+
Tsy
To,
Too
Thi
( (1 a)) + +0(1)) o(1)
(3)
Thn
Tso
To confirm that these asymptotics give reasonable approximation in finite sample sizes, the asymptotic test power (3) for given critic function f roughly agrees with the empirical test power using con-asymptoiic rejection threshold estimated by obtaining fresh samples from the null distribution, i.e. the empirical quantile of / under the mull. (This would not be available in practical situations, but is syuthetic problemn where we can obtain arbitrary numbers of fresh samples.) Figure indeed confirms that (3) gives a good estimate of the empirical power even at moderate sample sizes. (Intriguingly, the power of permutation test consistently slightly higher than the test based the correct non-asymptotic threshold; will discuss this issue in Scction 7.)
check that
This argument shows that, as long as T5 > Ts. and T has finite positive variance under both Pay and P. Py, test based on I will almost surcly eventually reject any fixed altcrnative: the test consistent. How quickly in m it reaches high power, however, depends on the expression (3), which in particular for large m is dominated by the signal-to-noise ratio (SNR) given by (Ts, T/5)/TT.. When we choose a test. based on maximizing the value of mutual information lower bound, we maximize criterion, e.g., (1), which does not directly correspond to the test power. We will discuss instead maximizing the asymptotic power in Section 5.
3 Tests with the Hilbert-Schmide Independence Criterion (HSIC)
The Hilbert-Schmidt Independence Criterion (HSIC, Gretton et al., 2005) is also zero and only X and are independent when an appropriate kernel is chosen (Szab6 Sriperumbudur (2018) give precise conditions). Unlike the mutual information, HSIC is easy estimate from samples; for a fixed bounded kernel, the typical istimasors concentrate to the population value with deviation O,(1/m))
5

---

Published in Transactions on Machine Learning Research (03/2026)
To define HSIC, we first briefty review positive-defiite kernels. A (real-valued) kernel is function k: XxX I that can be expressed as the inner product between feature maps Φ X - F, k(x,x')= (0(x)) ((2')) Iº where F is any Hilbert space. An important special case is == and k(r,x') = x(z)- - where ((T)
extracts p-dimensional features of x. For every kernel function, there exists a unique reproducing kernel Hilbert space (RKHS), which consists of functions f X + R. The key reproducing properly of an RKHS F states that for any function f F and any point x E X, we have (f.((x))> = f(x).
Suppose have kernel k on X' with RKHS F and feature map Φ, as well another kernel I on y with RKHS G and feature map sb Let @ denote the outer product." The mrss-covarianes operator is
Czy
E (b(x) ( Ed))) & (M(V) (o(u)-E<(O) Evw))]: ay
V
for kernels with fiitee-limensional feature maps, this is exactly the standard jrooss)covvaiance matrix between the features of X and those of Y. Under mild integrability conditions on the kernel and the distributions, 2 the reproducing property shows that (f, Czyg) = Cov(f(X), independence is whether there exist any correlated "test functions" f and g. Thus, for rich enough choices of kernel - using universal & and suffices, but not necessary (Szabó Sriperumbudur, 2018) - we have that X - Y if and only if the operator Cuy 0. We can thus check whether the operator is zero, and hence whether X L Y, by checking the squared Hilbert-Schmidt norm of Cry, HSIC(X, Y) || lliss With fnite-dimensionall features, this is the squared Frobenius norm of the feature cross-covariance matrix.
for all f E F, 9 E G. One definition of
Another way to interpret HSIC is as distance between Fzy and P. x Fy, similarly to how the mutual information is the KL divergence between those same distributions.
Proposition 3.1 (Gretton el al., 2012a, Theorem 25). Let and be kernels on X and V, and defone a kernel on X x y by ((x,,, y), (x', x')) k(x, a x'), v/) with RKHS H. Then
(X, MMD, (P my: P. x xP))
sup E [S(X,))
E [((X, X~P.
SEN (X,Y)~Fav I1
Y'MP,
E (XY)) Y" Y""'
X'XXX") - 2X(XX + A(X, X')("" Y"")
Taking the last form and rearranging to save repeated computation yields two similar, popular estimators of HSIC. "The biased estimator" is HSIC(P P2gy , Pa Pa) for empirical distributions P:
HSIC, Y)
K HLH for (A. = BAy
(4)
m
where here K is the m m matrix with entries kı Iy), L similarly has entries 4g = [[ By). and is the "centering matrix" In - m 1n (so the estimator can be easily implemented without matrix multiplication). This estimator has ((1/m) bias, but is consistent.
The other common estimator, "the unbiased estimator," " is a U-statistic (Song et al., 2012):
1 (X, = m(m - 3) (K 1) + KT1,,,, I T L1m 21 KL
(5)
(m 1)(m - 2)
where K and L are m X m matrices whose diagonal entries are zero but whose off- off-diagonal entries agree with those of K or L. It is unbiased, E HSIC, Y) = HSIC(X, Y); it is also consistent, and can be computed in the same Omm) time as HSIČ,, without matrix multiplication. Either statistic can be used with a permutation test. to construct an independence test with finite-sample validily, in the same way described for mutual information bounds. The following subsection describes the behavior of such tests in the large-sample limit.
'Th Euclidean outer product ab is a matrix with JabT = a(bTb/); in Hilbert spaces, f 0g : G + F has [ " = f/g. d')9. 2It, suffices that / ()), v\l < oog this guaranteed regardless the distribution when & are bounded.
6

---

Published in Transactions on Machine Learning Research (03/2026)
1.0
100
0.9
0.55
0.8
power
0.5
power
0.93
pawer
10p111 emplics
0.03
Jomp amp real
037
0.7
0.5
test Size
size
test size
(a) HDGM-4
(b) Sinnsoid
(c) RatInABox
Figure 3: Power of HSIC tests based on kernels h and learned with the method of Section 5, evaluated on problems described in Section 6. For cach problem, the blue line (asymp) is power described by the asymptotic formula (6), where we use 10000 samples to estimate the population level statistics HSIC, its deviation )) and the threshold (1- a) via simulation. The orange line (empirical) computes the power based on the simulated null distribution, i.e. the threshold estimated via the empirical CDF of HSIC with test size m. The green line (perm) the permutation test power. All three lines ronghly agree.
3.1 Asymptotic test power
The unbiased HSIC estimator is asymptotically normal with /m - HSIC) & (O (3 (Proposi- tion A.I, in Appendix A). is possible, however, to have o3 = 0, which makes that result less useful; this occurs in particular whenever HSIC 0, and hence always true under Hp- In that case, it more informative to say instead that m HSIC,, converges in distribution (which is not the case when 2 > 0)- The distribution to which it converges is a mixture of shifted chi-squareds with complex dependence on Px, pg. k and I (Proposition A.2). Thus, if we consider test statistic m HSTC when HSIC = 0 this statistic has mean zero and standard deviation (()). When HSIC > O, though, the statistic has mean and standard deviation @(ym) + co. Thus, as m + oc, eventually the test will reject if HSIC > 0 (Rindt et al., 2021).
Recall that the optimal test one with highest. test power among those with appropriate level control. As for NDS-based tests, can describe this power asymptotically. Let v be the edf of the null distribution of m HSIC which depends on k, L, Pr, and Py. It follows that the rejection threshold I'm satisfies
1 a) m
Pr(SSIC rm) (mrm) + o(1) 50
So
Im
Using Proposition A.1, when as > 0 the asymptotic test power then satisfies
HSIC HSIC
Tm HSIC
Pr
HSIC Tm)
Pr
m
m
51
Th
TB.
(1 a) +(1) HSIC
vm
m
+ +o(1)
051
Vm HSIC
T v- - a) + o(1)) + o())
(6)
Vmon.
Figure 3 verifies that this expression lines up well with the empirical test power even at moderate test sizes. Additionally, (6) also tells that, as long HSIC and HSIC, has finite positive variance, an HSIC-based test will eventually reject with probability of one. The rate at which it increases in power is described by the gap between the signal-to-noise ratio HSIC /as, and a threshold-depeeddet term v- (1-(11))//5, a)/as, with the former dominating the latter. This gap increases at a faster rate than the one described for NDS-based tests (3); the ratio of the two terms is ((m) for HSIC, while it only @(ym) for NDS. (Appendix B gives another motivation for the power being dominated by this signal-to-noise ratio, without using the central limit theorem.)

---

Published in Transactions on Machine Learning Research (03/2026)
3.2 MMD- independence test
Given that more attention has been recently paid to learning two-sample tests than dependence tests, we can consider reformulating the independence problem into a two-sample onc: are Pay and Ps. x Py the same distribution? A natural measure of difference between distributions Pay and Ps x Py MMD Pg x Py) with some kernel h, which becomes HSIC(X, Y) if h is product kernel. For sufficiently powerful h, characterized by Szabó & Sriperumbudur 2018, MMD ay, Ixx x P2y) 0 if and only if X IL Y. Thus it makes sense to use any consistent MMD estimator as the test statistic; we use the biased U-statistic estimator described in Gretton et al. (2012a) for simplicity.
Next we consider a testing procedure. The standard independence problem observes paired samples 2 = ((ar)3))..... .i..... (zm: Vm)) drawn i.id. from the joint distribution F Pay: and so to simulate samples from P, Fy we first shufile the Y samples by some permutation on {4,..,m), We define this action of permuting Z with a as oA = (ot)),,) (Iim Bo(m))) These samples are then used compute MMD, oZ). The following proposition shows that, with enough samples, this statistic is consistent estimator of the true MMD value.
Proposition 3.2. Suppose satisfies supzex wey li((,, y), (x, S v². The eslimator MMD (Z, where Z is independent of uniformly random o, satisfies with probability at least 1 that
17v2
- (Fzy P. X P
2
m
This result is proved in Appendix F.4.
Viewing this permuted MMD estimator as a function of the original paired samples, i.c. f(Z) == MMD oZ), we then perform an independence permutation test by shuffling only the samples in Z. Under the null hypothesis the distribution of Z is unchanged by shuffling the Y samples, and so an independence permutation test on f(Z) is guaranteed by Theorem of Hemerik & Goeman (2018) have correct level control, provided the original sample order is included in the permutations."
It's also worth noting that, since the MMD estimator is a two-sample statistic, the number of samples of I and Pa x Iy can be different. As such, we can also consider using multiple permutations a = {o),..., Øp) rather than just one, and then simulate the null samples as oZ Uø It turns out that we consider all 2 rotation permutations, MMD, (Z, Z)) is exactly the biased HSIC estimator HSIC, Y). We elaborate on this connection in Appendix F.4 and show that using more permutations improves estimation quality. The HSIC estimator, which effectively uses all permutations, is the minimum variance estimator for this permuted MMD statistic (Proposition F.2).
4 Connecting HSIC and Mutual Information Tests
4.1 HSIC as a lower bound on MI
Suppose the kernel on (x, y) pairs of Proposition 3.1 bounded: supeeN h((,, v), (x,,)) S v2 for v 2 0. Then
If v) Kf. on v) S lf (l w S &1lfIlH
"One might wonder whether can instcad perform two-sample permutation Consider the test statistic as function of the pooled samples f(Z, ož)) M2, σZ). The two-sample permutation test involves shuffling the order of 2 oZ, however, exchangeability of 2 72 under the null is guaranteed as variables in the first split deterministically tied to variables in the sccoud; thus two-sample permutation test may be valid. A workaround would be to independently split the data into Z.j and Z2 and then construct the pooled sample Z U oZ2, but this process more involved and potentially power.
8

---

Published in Transactions on Machine Learning Research (03/2026)
by the reproducing property and Cauchy-Schwarzz so ||fl 2 1/2 sup su أ)و Illf This implies that ( 11 1) S ( b0 I v), and so by Proposition 3.1,
HHSICXX,Y) & sup
E [XX (X
- E X-P y'mBx
= = my. Ps x P.),
(7)
where TV is the total variation distance between distributions (Sriperumbudur et al., 2012).* Applying standard bounds relating the total variation to the KL divergence, we obtain the following.
Proposition 4.1. In the setting of Proposition 3.1, suppose yyy (((x,)), (2) S sv² Then
1 HSSC(XX Y) S ST(XYV) and
1 HSIC(X,Y)) s <T(XYY). - -log (1 (1-11² 4v²
2v²
Proof. The first bound applics Pinsker's incquality, which relates total variation to KL, to (7). The second instead applies the bound of Bretagnolle & Huber (1978) (also see Canonne, 2022).
The second bound is tighter for large values of YY,; but both are monotonic in HSIC(X, Y)/r². We could thus consider, as in Section 2, choosing kernels k, to maximize a lower bound on I(X;Y), by maximizing HSIC(X,Y)) Y)/b². Indeed, maximizing HSIC has been used by previous applications in many areas (e.g. Blaschko & Gretton, 2009; Song al., 2012; Li et al., 2021; Dong et al., 2023).
4.2 Kernel-based tests and variational MI tests
Now consider test based on MMD, (PP my: xP,) using a kernel of the form h((2, v), -0)) for some real-valued function f. If fz fi (z) (2))),), this is an HSIC test with kernels (((2,22) = fi fz and [(,, y') = /2v')) Because Φn (x, y) = f(x, v) E R is a valid feat ure map, every function in H is of the form a = [(x, y) + a (x,))) with Ilaf = la]. By Proposition 3.1,
2
S(X.Y1))
(MD, my: P, x P,,)²=
sup a |ais1
E f(X,)) P..
E P.
( f(X,Y) PAx
E Y')
P. x
The plug-in estimator would yield the lest statistic
2
(M, Y)
1
fI.
m70
m
(as, v5)
(8)
m
m
which when -
corresponds exactly to HSIC.
Comparing (8) to (1) with K m, can that the main term Í- E, v.) is identical, which we called the NDS (2). The other term is permutation-invariane permutations, T. Thus, a permutation test basod on (8) asks how far the value of f for the true data from T, while test based on any variational MI estimator based on (2) is equivalent to asking how much the value of T exceeds T. The only difference is that (8) gives two-sided test, while (1) is one-sided test.
it the mean of i over all possible
Our usual test uses HSIC, of (5) instead HSIC,, but the difference in estimators is typically small. Thus, if we use deep kernels of the form k(r, x') = (x') and 1(y v') = g(v)g(s'), obtaining a test quite closely related to witness two-sample test (Kübler et al., 2022) used for independence, the HSIC test is nearly equivalent to the NCE test with separable critic function (z, v) + f(x)g(v). This relationship is roughly analogous to the relationship between classifier two-sample tests and MMD tests observed by Liu et al. (2020): while each test chooses a critic/kernel kernel in a different way, at test time they are essentially equivalent.
SSriperumbudrr et al. (2012) define the TV as twice the more common definition, which we here. 5i.i et al. (2021, Section 3.1) also found relationship between HSIC and INCE for categorical Y.
9

---

Published in Transactions on Machine Learning Research (03/2026)
5 Learning Representations for Independence Testing
Choosing a test based on maximizing the value of its statistic does not directly correspond to the test power. Instead, we would perhaps be better served by maximizing the power directly. Recall the asymptotic power expressions for NDS and HSTC-based tests given by (3) and (6), respectively:
Pr(f > rm) $ [m - - 5)
The (1 a) and Pr Pr(HSIC,, > rm) 5)
(1
Tp:
Tn
Vmos,
In both cases, as the sample size m grows, the power is dominated by the signal-to-noise ratio
SNR[T] = (T - TTo)//T
and SNR[HSIC] = HSIC fass-
Maximizing the SNR maximizes the limiting power of the test as m - oo. (Also see Appendix B.) Thus, maximizing estimates of the SNR (based on a finite number of samples n) will roughly maximize the power of a test an arbitrarily large sample size (m + oo); also see discussion by Sutherland & Deka (2019); Deka & Sutherland (2023).
Alternatively, can also choose to include the threshold-dependent term. This may be helpful, particularly in estimating the power at points where the number of samples is not overwhelming: for tests with a 50% probability of rejection, the two terms in the power expression are necessarily the same size. For NDS, the gap between the two terms only grows with vm, and additionally estimating Ts, is no problem when already estimating T thus we propose to include the threshold term choosing an NDS critic. For HSIC, however, the threshold-dependent term is both relatively less important as m grows and much more difficult to estimate. Therefore, we propose to optimize the SNR with the threshold term for NDS, and without it for HSIC. Our population level objectives, as functions of random variables X and Y, are then
JNDS(X, T: f) Tss (X,, Ts) (X
(X, f:
(1-
(10)
VTTS (X,Y;
and
IIICIX,Y;K.I.
JuSIC Jusio(X,Y; X, Th (X, Y;k,
(11)
Then, if we have a batch of paired observations X = (x,, Xm) and V = (Y,,.., Ym) sampled jointly from Pxy, we can approximate JNDS and JASIC by the estimators JRDS and Jusice which we define below.
NDS objective. We define J as the plug-in estimator of JNDS with variance regularization À >> 0. We estimate Ts, = Epyr Y)] by its sample mean f(X, Y) 1 Em f(, Yi), and Tho ExxxPy Y( by the V-statistic Ts. (x, Y) = m Em ΣØ= f(G,, Yj). The variance terms The and r are estimated with a regularized version of the sample variance,
2
m
m
m
+2 (X,
f(X,,)))
Σ
Σ
+1
m
m
i=1
j
2
1
f(X,Y)) - 1 E(())
2/ Y:
A ,/
+
m
where À > 0. Putting it all together, our complete NDS objective is
JiDS (x, Y; 7(x, Y: f) To f) tr (X
V Asy (X, Y; 1 (1-(1
(12)
xx,
'Re et al. (2024) proposed using moment-mmatched gamma approximation to the threshold, and claimed that not including this threshold can lead to catastrophically wrong kernel choice. We argue in Appendix F.5, however, that their argument is unjustified. We also find experimentally in Section 6 that including it can chouse contemporary approachus
kernels in practice. More common estimating the threshold use permutation testing cigendocomposition; both involve substantial computational overhead, and the latter particularly expensive to differentiate while the former can be quite difficult to uscfully estimate on the same data the HSIC estimate, observed for MMD by Deka & Sutherland (2023) and discussed further in Section
10

---

Published in Transactions on Machine Learning Research (03/2026)
HSIC objective. Jisic is defined in similar manner. We estimate HSIC by the unbiased estimator HSIC (X, Y) given in (5), and we approximate the variance term of according to
1) = 16 HSIC +
which is just a regularized version of the variance estimate stated in Proposition A.1. Following Song et al. (2012), this can be computed more efficiently with R (n((m=1)])= [|]] where the vector h is
((n-4)))²
(n 2)2 (K D) 1- (ī1) (1 Ł1) (1 K1)L1 (1 + (n ((1 (K 1)111 ( K
with o denoting elementwise multiplication on matrices, and 1 = (1,..., E RR.. Thus, the complete HSIC objective
Jusic (X, Y;k, HSIC (X, 1; T (X, 1;
(13)
Maximizing (12) or (13) then gives us a structured approach to select the critic or kernel yielding the strongest asymptotic test. Typically, we run some variant of a gradient-based optimization algorithm with respect to the parameters of the critic or kernel class, and do hyperparameter selection and early stopping based on the same objective on a validation set.
Critic & kernel architecture. We consider function families parameterized by deep networks, as this allows to learn representations of the data that more efficiently capture dependency. For critics, this parameterization is straightforward and is typically some problem-specific architecture like multilayer percep- trons or convolutional networks. For kernels, we first incorporate a featurizer that acts on individual inputs æ and and is then fed to some standard kernel function. When our feature mapping is parameterized by deep networks, we get the class of deep kernels (Wilson et al., 2016), which have been successfully used in two-sample testing (Sutherland et al., 2017; Liu et al., 2020; 2021) and many other settings (c.g. Li al., 2017; Arbel et al., 2018; Jean al., 2018; Li et al., 2021; Gao et al., 2021). We use the following deep kernels for X and Y:
(1-x)))=11(xxxx((( 4(,.0') (y- (,)) +ey (r.
(xxx(xx))))
Here fw and By deep networks with paramcters in w, 7, which extract relevant features from X' Y to feature space D. These features are then used inside a Gaussian kernel the space RD, to compute the baseline similarity between data points. We then take convex combination of that kernel with Gaussian kernel on the input space; the weight of this component is determined by a parameter c E (0, 1). The lengthscale of x and q as well as the mixture parameter e are included in the overall parameters, w or 7, and learned during the optimization process.
Overall representation learning algorithm. The overall procedure, writlen based on minibatch gradient ascent for simplicity, is shown in Algorithm 1. In practice we use AdamW (Loshchilov & Hutter, 2019), and draw minibatches in epochs; experimental details are given in Appendix E.1. The randomized p-value is potentially conservative in the case of ties; breaking ties randomly can give a slight improvement when ties are present. (Hemerik & Goeman, 2018).
Time complexity. Each training iteration is dominated by computing the objectives JADS and JASIc- Suppose K is the minibatch size. For NDS tests, Ez the cost. of evaluating critic fo y), then one training step costs O(K² E2). For HSIC tests, if Ex and Ey are the costs of computing embeddings f=( and (v), and L the cost of computing (_, x²) and (y v²) given the embeddings, then cach training
Usinn é >> provides "backup" the deep kernel, perhaps giving some signal early in optimization when the deep kernel features are not yet useful, and guaranteeing that the overall kernel characteristic.
11

---

Published in Transactions on Machine Learning Research (03/2026)
Algorithm 1 Independence testing with learned representations
Input: paired samples Sz (S); split the data into sz U 5z with (S). + : (5x -2: model parameters 0 <- 0g and tost statistic T(X, 003; various hyperparameters like À <- 10 etc.
# Phase I; train the paromaters on 5%
for T = 1,2....,Tan Tmax do
(x, W) t minibatch from 82 = sV);
(X, 1;0) compute SNR cstimate;
# as in Equation (12) or Equetion (13) # gradient ascent step
+ A V
# Phase 2: permatation test on S2 with learned representations permn < ();
# evabuate test statistic
for i= 2....,,,.. Alperm do
perm; <- T(Sx shufflo(S5):0)
# no need to shuffic both sets
Output: kw by perm,; p-value 1 Di1
yyerm (perm, 2 perm,)
uperm
iteration costs OKEx + KEy + K2T.). Typically, for practical values of K, Ex + Fy >> KL, so this cost. is "almost" linear in practice.
Theoretical analysis. While neither J'DS nor Jusic are biased estimates of the relevant population quantities, 9 the following Theorems 5.1 and 5.2 show that they are both uniformly bounded in probability.
Theorem 5.1 (Uniform convergence of JiDs). Let (So : A E Q) be a critic family with parameter space o satisfying Assumptions (A) (C) in Appendir C.2 which define a critic bound B, dimension D, and smoothness L. Suppose T 2 s² for some positive s under both Ho and h1, and X = O(n 11/3. Then
B + BBL+BVVD
(D)(()) - JNDD(() =Op 2 ee
B²
sup
Theorem 5.2 (Uniform convergence of JisIcc). Let {kw w O} and [4 : I3} be kernel families with parometer spaces Ω and I' satisfying Assumptions (A') ' to (C") in Appendia D.2 which define kernel dimensions
Da, Dr and smoothnesses Lk, Lt. Suppose oū, 2 s² > 0 under fh1 and X ((n-1/3) Then
J Jhsscow,)) Jusic(w.)) Öp 2
1
L + L1 + DD + VDD))
sup yER,yEr
A Appendices C and D state and prove non-asymptotic versions of the above results based on covering numbers; the assumptions and proof techniques are similar to those of Liu et al. (2020)10, but is slightly more involved for few reasons: 1) Uniform convergence of HSIC and its variance estimators now include two kernels, & and I, instead of one, and a somewhat different estimator form. 2) For NDS we also need to show convergence of a mull-to-alternative variance ratio, a term which is not present in the MMD or HSIC objectives.
Following Theorem C.13 (a standard result from van der Vaart, 1998), successfully maximizing these estaes will thus also maximize the population quantity, and consequently the asymptotic test power, for sufficiently large training set sizes n.
6 Experiments
The repository code for experiments.
contains implementations of our methods and
*Equation (13) could use block estimators (Zaremba et al., 2014) or incomplete U'statistics (Blom, 1976) to reduce O(K*L)
to O(KPL) for any < 2, at the cost of increased variance (see Ramdas et al., 2016).
Indeed, no unbiased estimator is likely to exist; see Appendix A of Deka & Sutherland 2023.
Ren al. (2024) stated similar result Theorem 5.2, but with seemingly incorrect proof: in their Appendix D.3, the term (r(n) (17)-r1 in their (63) should be bounded uniformly over kernel parameters, but they appeal to result about convergence of distribution functions for single distribution. For the uniform result over kernel parameters that their theorem statement claims, they would need additionally use cover or other approach (as done for the HSIC and the variance estimators), but none attempted.
12

---

Published in Transactions on Machine Learning Research (03/2026)
Baselines. We compare our HSIC-based (HSIC-D/Dx/O) and MI-based (NDS/IffNCC//WW) tests with various alternative methods. All tests are performed using permutation testing.
HSIC-D: HSIC using deep kernels on each space X' and V; simultaneously trained via Section 5.
HSIC-Dxx HSIC using tied deep kernel, ie. h l, and trained via Section 5.
HSIC-O: HSIC using Gaussian kernels, with each bandwidth parameter optimized via Section 5.
NDS: The neural dependency statistic (2) trained via Section 5.
InfoNCE (van den Oord et al., 2018; Poole et al., 2019): the statistic INCE as in (1).
NWJ (Nguyen et al., 2010; Poole et al., 2019): another mutual information bound statistic INWJ-
HSIC-M: HSIC using Gaussian kernel, with bandwidth selected via the median heuristic.
HSIC-Agg (Albort et al., 2021; Schrab et al., 2023): aggregating Gaussian kernels of various bandwidths. We use the complete U-statistic with default settings: B1 = 500, B2 = 500, B3 = 50, and uniform weights.
MMD-D; The method of Liu et al. (2020) applied to xx P. x Py, y, 11 with Gaussian kernel on X x y.
C2ST-S (Lopez-Paz & Oquab, 2017) / C2ST-L (Cheng Cloninger, 2022): Sign/logit: based classifier two-sample test for Pay vs P2 x py, with samples set up the same way as for MMD-D.
Datasets. We consider four informative datasets, where the true answers are known.
High-dimensional Gaussian mixture. The distribution HDGM-d has d total dimensions (divided between X and Y), but has dependence between only two of them:
0.5(-1)
0 0
[X Ma/2)./ (4/22) Yil
N Od 0.5(-1)*
02-2
Oa-2
0a-2
Ia-2
where the odd dimensions are taken to be from P. and even dimensions be from Py, Moreover, for d 2 4 the dependent variables X2 and Yla/aj are at different dimensions. We perform independence tests at dimensions 4, 8, 10, 20, 30, 40, and 50. This distribution is similar to one used by Liu et al. (2020).
Sinusoid (Sejdinovic et al., 2012). We sample from sinusoidally dependent data with distribution Pzy x +++ (Er( sin(Ey) on support X' x == [--,,,;; Higher frequencies E produce subtler departures from the uniform distribution, resulting in a harder independence problem; we use l 4. A visualization of this density is given in Figure 8.
RatInABox (George et al., 2024). RatInABox simulates hippocampal cells of rat in motion. In particular, we test for dependence between firing rates of grid cells and the rat's 's head direction. Grid cells respond near points in a grid covering the environment surface, and should be subtly connected to head direction because of the geometry of the "box" (Figure 7). We consider 8 grid cells, and simulate motion for 100 000 seconds, taking a measurement every seconds as our dataset.
Wine Quality (Paulo al., 2009). Details physicochemical propertics (e.g., sugar, pH, chlorides) of different lypes of red and white wines and their perceived quality an integer value from 1 to 10. We test. for dependency between residual sugar levels and quality.
Power vs. test size. We first compare how well methods identify dependency with a large-size training set, by comparing the rate at which the learned tests achieve perfect power (1.0) as the test size m increases. Training and validation sizes for each dataset are given in Appendix E.1, and results are visualized in Figure 4 (with comprehensive results in Figure 9). Overall, HSIC-D outperforms baselines, and is able to reach perfect power at smaller test sizes m. We note that HSIC-O performs reasonably well for simpler problems like Sinusoid, but struggles on harder problems like RatInABox. This suggests that Gaussian kernel the input space is well-suited for the task. On the other hand, kernels applied to optimal representations of the data (HSIC-D) are the most powerful. Surprisingly, optimizing the power of NDS tests is significantly weaker than directly maximizing InfoNCE or NWJ; we elaborate on this in the following section.
11 Speeifiially, we compare 2 ~ Pay against single shuffle of the samples oZ, detailed in Section 3.2. We this approach over the data-splitting approach since it both minimizes memory (compared to using large number of permutations) and performs slightly better in practice, as seen Figure 13.
13

---

Published in Transactions on Machine Learning Research (03/2026)
ISCC
HSCC
-4-
NZS
HIONCE
MA]
C25TL
C23T5
NMD-5
MSICO
HBKC
H3KCAgg
1.0
1.0
1.0
empinicl power
0.8
emprrical paweer
0.8
empirical power
0.8
empirical power
0.8
0.6
0.6
0.4
04
02
0.0
0.0
1000
1500
samples
samples
samples
tesples
(a) HDGM-10
(b) Sinusoid
(c) RatInABox
(d) Wine
Figure Empirical power vs sample size m for different datasets, when trained with a large training set. The average power computed over training runs, where the empirical power determined over 100 permutation tests. The shaded region covers one standard error from the mean.
Power vs. dataset size. A drawback to kernel selection via optimization is that we must hold out a split of the data for training. contrast, HSIC-M and HSIC-Agg are able to utilize the entirety the dataset for their test. To examine this tradc-off, we consider consistent data splitting at datascts of varying sizes. For the HDGM and Sinusoid problems we 7:2:1 train-val-test split, and for RatInABox and Wine 3:2 train-test split. Conversely, HSIC-M and HSIC-Agg use the entire dataset for testing. Results are shown in Figure 5. Non-splitting methods are able to take advantage of the additional test samples and outperform some data-splitting methods at smaller dataset sizes; for large dataset sizes the reverse is true. The validation set in this experiment used only for early stopping of the training process, and not for hyperparameter selection.
NZS
C25TL
HSK-AAg
1.0
1.0
Empinical power
empinical power
empinical power
0.8
empirical power
0.6
0.6
0.4
0.2
0.2
0.2
0.0
0.0
0.0
1000
dataset size
dataset size
dataset size
dataset size
(a) HDGM-10
(b) Sinusoid
(c) RatInABoo
(d) Wine
Figure Empirical power vs dataset size. HDGM and Sinusoid uses a consistent 7:2:1 train-val-tess split across all dataset sizes, while RatInABox and Wine maintain a 3:2 train-test split. HISC-Agg and HSIC-M do split the data. The shaded region covers one standard deviation over 5 training runs.
SNR with vs. without threshold. Our HSIC objective Jusic disregards the threshold-dependeet term, whereas the NDS objective JADS does not. The impact this term has during training, however, still unclear. We compare the strength of tests maximized with and without this term in Figure 15. It is important to note that the threshold estimate for HSIC is challenging (see also Footnote 6); the moment-mon matehed Gamma approximation, but with improved gradient estimation rather than the finite differences estimate they used (Appendix F.6). This estimate, in either form, introduces substantial additional code complexity and its own bias and variance that could each hurt overall performance. In general, we find that including the threshold for HSIC has no strong effect on simpler problems, but degrades the test power for harder problems. On the other hand, NDS with the threshold yiclds stronger tests. This is not too surprising as, unlike for HSIC, the null distribution here follows the CLT and is easily approximated given enough samples.
method based on
14

---

Published in Transactions on Machine Learning Research (03/2026)
1.0
2.6
power
power
0.6
power
sryp
0.0
200
300
500
200
300
400
500
6:0
800
1000
1200
1400
test size
test size
test size m
(a) HDGM-4
(b) Sinusoid
(c) RatInABox
Figure 6: Power of NDS tests derived from a learned InfoNCE critic f. For each problem, the blue line (asymp) is power described by the asymptotic formula (3); where we use very large sample size (n 20 000) to estimate the population level statistics T, To and their deviations , T The orange line (empirical) computes the power based on the simulated mull distribution, i.c. the threshold estimated via the empirical CDF of Th with test size The green line (perm) is the permutation test power.
7 Discussion
Independence tests with MMD. We find that MMD-D underperforms HSIC-D in almost all settings, despite the equivalence between MMD and HSIC specified in Proposition 3.1. This equivalence, while true at the population level, is not the case when considering finite sample estimates. These estimators exhibit different statistical properties, with some being more conducive to learning than others. In particular, we notice that MMD estimators based on single permutation of the samples (as described in Section 3.2) exhibit larger variance compared to HSIC ones. We prove this for the biased estimators in Appendix F.4, and observe the same phenomenon when training MMD-D and HSIC-D.
Permutation tests sometimes disagree with asymptotics. We hypothesized that maximizing an estimate of the asymptotic power would produce stronger tests than maximizing the test statistic directly. Although this seems to be true for HSIC (Appendix F.3), it is surprisingly not the case for MI-based statistios; The previous power versus test size experiments show that directly maximizing InfoNCE or NWJ outperforms maximizing the NDS asymptotic power, at least when using a permutation test. Why is this?
Figure 6 examines several power estimates for an NDS test based on an InfoNCE-learued critic. Notice that the asymptotic formula. (blue) agrees well with the simulated test power (orange), indicating that the asymptotics describe the test well. However, the permutation power (green) is significantly greater than what the asymptotics suggest. Recalling that the simulated test power is in fact using near-exact threshold based the empirical quantile of the test statistic under the null, this should be roughly the most powerful valid test with a innpee-indeeendent threshold. Permutation tests, however, use a cample-dependent threshold, which in this case apparently can yield a much more powerful test than the best sample-indeeenddent threshold. We explore this discrepancy further in Appendix F.7.
Thus, while maximizing the power of the NDS test maximizes the asymptotic power of a test with a data- independent threshold, it does not necessarily maximize the power of a permutation test. InfoNCE, and other variational MI estimators, seem to frequently choose critics whose sample-independent threshold tests are worse than those obtained by NDS maximization, but whose permutation tests are actually much more powerful. Exploring this connection more deeply is very intriguing area for future work.
We did not observe this discrepancy for HSIC; maximizing the sample-independentt asymptotics for HSIC seems yield permutation tests (roughly matching their asymptotics) which still outperform InfoNCE or similar NDS permutation tests. Overall, HSIC-D seems to be the most reliable method of those considered here.
15

---

Published in Transactions on Machine Learning Research (03/2026)
8 Conclusion and Future Work
Independence testing aims to see if two paired random variables are statistically independent. We explored two families of tests to address this problem. The first are lests based on variational mutual information estimators, which, to the best of our knowledge, we are the first to construct. The second are tests based on maximizing the asymptotic power, which was explored for the two-sample problem but not for independence testing. Our findings show that learning representations of the data via our proposed methods lead to powerful tests, with HSTC-based tests generally outperforming MI-based ones. Future work may look to extend this learning scheme to conditional independence testing and apply this to causal discovery. Meanwhile, it will also be valuable to investigate approaches for mitigating, or even removing, the data-splitting procodure as done by Biggs et al. (2023); Kübler et al. (2020) while maintaining the ability to learn rich data representations.
Acknowledgments
The authors would like thank Roman Pogodin for productive discussions. This work enabled in part by support provided by the Natural Sciences and Engineering Research Council of Canada (with grant numbers RGPIN-2021-0297 and ALLRP-57708-2022 the Canada CIFAR AI Chairs program, Calcul Québec, the BC DRI Group, and the Digital Research Alliance of Canada. This work was also supported by the Australian Research Council (ARC) with grant numbers DE240101089 and DP230101540.
References
Mélisande Albert, Béatrice Laurent, Amandine Marrel, and Anouar Mcynaoui. Adaptive test of independence based on HSIC measures, 2021.
Michael Arbel, Danica J. Sutherland, Mikolaj Binkowski, and Arthur Gretton. On gradient regularizers for MMD GANs. In NewrIPS, 2018.
Mohamed Ishmael Belghazi, Aristide Baratin, Sai Rajeshwar, Sherjil Ozair, Yoshua Bengio, Aaron Courville, and Devon Hjelm. Mutual information neural estimation. In ICML, 2018.
Felix Biggs, Antonin Schrab, and Arthur Gretton. MMD-FUSE: Learning and combining kernels for two-sample testing without data splitting, 2023.
Matthew B. Blaschko and Arthur Gretton. Learning taxonomies by dependence maximization. In NearIPS, 2009.
Gunnar Blom. Some properties of incomplete U-statistics. Biometrika, 63( 1976.
Jean Bretagnolle and Catherine Huber. Estimation des densités: risque minimax. Séminaire de probabilités de Strassoourg 12:342-363, 1978.
Jaroslaw Blasiok. An extended version of McDiarmid's inequality. MathOverflow, 2025. URL https: //mathoverflow.co net/ wwtet//////eeeeeeeeeeeeeeeeeeeee 9/491099.
Clément L. Canonne. short note on an inequality between kl and tv, 2022.
Xiuyuan Cheng and Alexander Cloninger. Classification logit two-sample testing by neural networks for differentiating near manifold densities. JEEE Transactions on Information Theory, 68(10):6631-660) October 2022. doi:
URL http::/dx.doo. wwt.:ooi..oooooooo 1tt//0..119/TIT..
Felipe Cucker and Steve Smalc. On the mathematical foundations of lcarning. Bulletin of the American Mathemalicul Suciely, 39(1):1 49, 2002.
Namrata Deka and Danica J. Sutherland. MMD-B-Fair: Learning fair representations with statistical testing. In AISTATS, 2023.
Ruijiang Dong, Feng Liu, Haoang Chi, Tongliang Liu, Mingming Gong, Gang Niu, Masashi Sugiyama, and Bo Han. Diversity-enhancing generative network for few-shot. hypothesis adaptation. In ICML, 2023.
16

---

Published in Transactions on Machine Learning Research (03/2026)
M. D. Donsker and S. R.S. Varadhan. Asymptotic evaluation of certain Markov process expectations for large time. IV. Communicattins Pure and Applied Mathemalics, 36(2):183-22 -222, 1983.
Ruize Gao, Feng Liu, Jingleng Zhang, Tongliang Liu Bo Han, Gang Niu, and Masashi Sugiyama. Maximum mean discrepancy test is aware of adversarial attacks. In ICML, 2021.
Tom M George, Mehul Rastogi, William de Cothi, Claudia Clopath, Kimberly Stachenfeld, and Caswell
Barry, Ratinabox, a toolkit for modelling locomotion and neuronal activity in continuous environments, eLife, 13:c85274, feb 2024. doi:
URL attps://dol. http::///ooooooo/////o wttp://...rr/.....////...// w
Arthur Gretton, Olivier Bousquet, Alexander J. Smola, and Bernhard Schölkopf. Mcasuring statistical dependence with hilbert-schmidt norms. In ALT, 2005.
Arthur Gretton, Kenji Fukumizu, Choon Hui Teo, Le Song, Bornhard Schölkopf, and Alexander J. Smola. A kernel statistical test of independence. In NewrIPS, 2007.
Arthur Gretton, Karsten M. Borgwardt, Malte J. Rasch, Bernhard Schölkopf, and Alexander Smola. A kernel (wo-sample test. Journal of Machine Learning Research, 13(2)):723 773, 2012a.
Arthur Gretton, Bharath K. Sriperumbudur, Dino Sejdinovic, Heiko Strathmann, Sivaraman Balakrishnan, Massimiliano Pontil, and Kenji Fukumizu. Optimal kernel choice for large-scale two-sample tests. In NewrlPS, 2012b.
Michael U. Gutmann and Aapo Hyvärinen. Noise-contrastive estimation of unnormalized statistical models, with applications to natural image statistics. Journal of Machine Learning Research, 13(11):307-30) 2012.
Jesse Hemerik and Jelle Goeman. Exact testing with random permutations. TEST, 27:811 825, 2018.
Neal Jean, Sang Michael Xie, and Stefano Ermon. Semi-supervised deep kernel learning: Regression with unlabeled data by minimizing predictive variance. In NeurIPS, 2018.
Willawal Jitkrittum, Zoltán Szabó, Kacper P Chwialkowski, and Arthur Gretton. Interpretable distribution features with maximum testing power. In NewrIPS, 2016.
Jonas M. Kübler, Wittawat Jitkrittum, Bernhard Schölkopf, and Krikamol Muandet. Learning kernel tests without data splitting. In NewrIPS, 2020.
Jonas M. Kübler, Wittawat Jitkrittum, Bernhard Schölkopf, and Krikamol Muandet. A witness two-sample test. In AISTATS, 2022.
Chun-Liang Li, WeinCheng Chang, Yu Cheng, Yiming Yang, and Barnabás Póczos. MMD GAN: Towards deeper understanding of moment matching network. In NewrIPS, 2017.
Yazhe Li, Roman Pogodin, Danica J. Sutherland, and Arthur Gretton. Self-supervised learning with kernel dependence maximization. In NewrIPS, 2021.
Feng Liu, Wenkai Xu, Jie Lu, Guangquan Zhang, Arthur Gretton, and Danica Sutherland. Learning deep kernels for non-parametric two-sample tests. In ICML, 2020.
Feng Liu, Wenkai Xu, Jie Lu, and Danica J. Sutherland. Meta two-sample testing: Learning kernels for testing with limited data. In NewrIPS, 2021.
David Lopez-Paz and Maxime Oquab. Revisiting classifier two-sample tests. In ICLR, 2017.
Ilya Loshchilov and Frank Hutter. Decoupled weight decay regularization. In ICLR, 2019.
Russell Lyons. Distance covariance in metric spaces. The Annals of Probability, 2013.
Andreas Maurer. A bernstein-lype inequality for functions of bounded interaction, 2017.
David McAllester and Karl Stratos. Formal limitations on the mcasurement of mutual information. In A AISTATS, 2020.
17

---

Published in Transactions on Machine Learning Research (03/2026)
R. J. Moore. Algorithm as 187: Derivatives of the incomplete gamma integral. Journal of the Royal Statistical Sociely. Series C (Applied Slatislics), 31 (3):330-333) 1982.
XuanLong Nguyen, Martin J. Wainwright, and Michael I. Jordan. Estimating divergence functionals and the likelihood ratio by convex risk minimization. IEEE Trans. Inf. Theor., 56(11):5847-5888 5861, 2010.
Liam Paninski. Estimation of entropy and mutual information. Neural Compulation, 15:1191 1253, 2003.
Cortez Paulo, Cerdeira A., Almeida F., Matos T., and Reis J. Wine quality. UCI Machine Learning Repository, 2009.
Ben Poolc, Sherjil Ozair, Aaron van den Oord, Alexander A. Alemi, and George Tucker. On variational bounds of mutual information, 2019.
Aaditya Ramdas, Sashank J. Reddi, Barnabas Poczos, Aarti Singh, and Larry Wasserman. Adaptivity and compatationn-sttiitii tradeoffs for kernel and distance based high dimensional two sample testing, 2015.
Yixin Ren, Yewei Xia, Hao Zhang, Jihong Guan, and Shuigeng Zhou. Learning adaptive kernels for statistical independence tests. In AISTATS, 2024.
David Rindt, Dino Sejdinovie, and David Steinsaltz. Consistency of permulation tests of independence using distance covariance, HSIC and dHSIC. Stat, 10(1), 2021.
Antonin Schrab, Ilmu Kim, Mélisande Albert, Béatrice Laurent, Benjamin Guedj, and Arthur Gretton. MMD aggregated two-sample test. Journal of Machine Learning Research, 24(194), 2023.
Dino Sejdinovic, Arthur Gretton, Bharath Sriperumbudur, and Kenji Fukumizu. Hypothesis testing using pairwise distances and associated kernels, 2012.
Dino Sejdinovic, Bharath Sriperumbudur, Arthur Gretton, and Kenji Fukumizu. Equivalence of distance- based and RKHS--baeed statistics in hypothesis testing. The Annals of Statistics, pp. 2263 2291, 2013.
Jiaming Song and Stelano Ermon. Understanding the limitations of variational mutual information estimators. In ICLR, 2020.
Le Song, Alex Smola, Arthur Gretton, Justin Bedo, and Karsten Borgwardt. Feature selection via dependence maximization. Journal of Machine Learning Research, 13(47):1393-14333 2012.
Bharath K. Sriperumbudur, Kenji Fukumizu, Arthur Gretton, Bernhard Schölkopf, and Gert. R. G. Lanckriet. On the empirical estimation integral probability metrics. Electronic Journal of Statistics, 6:1550 1599, 2012. doi: 00.1214/12-EJS222
Danica J. Sutherland and Namrata Deka. Unbiased estimators for the variance of MMD estimators, 2019.
Danica J. Sutherland, Hsiao-Yu Tung, Heiko Strathmann, Soumyajit De, Aaditya Ramdas, Alex Smola, and Arthur Gretton. Generative models and model criticism via optimized maximum mean discrepancy. In ICLR, 2017.
Zoltán Szabó. Information thcoretical estimators toolbox. Journal of Machine Learning Research, 15(9): 283 287, 2014.
Zoltán Szabó and Bharath K. Sriperunbudur. Characteristic and universal tensor product kernels. Journal of Machine Learning Research, 18(23)):1-29, 2018.
Gábor J Székely, Maria L Rizzo, and Nail K Bakirov. Mcasuring and testing dependence by correlation of distances. The Annols of Slatistics, 2007.
Aaron van den Oord, Yazhe Li, and Oriol Vinyals. Representation learning with contrastive predictive coding, 2018.
A. W. van der Vaart. Asymptotic Statistics. Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, 1998.
18

---

Published in Transactions on Machine Learning Research (03/2026)
Andrew Gordon Wilson, Zhiting Hu, Ruslan Salakhutdinov, and Eric P Xing. Deep kernel learning. In A ATSTA7S, 2016.
Wojeiech Zaremba, Arthur Gretton, and Matthew Blaschko. B-tests: Low variance kernel two-sample tests. In Newr/PS, 2014.
19

---

Published in Transactions on Machine Learning Research (03/2026)
A Asymptotics of HSIC Estimators
Proposition A.1 (Song et al., 2012, Theorem 5). Under the alternative hypothesis $1 : Pay F FP x Py, the unbiased estimator of HSIC is asymptotioally normal: with m samples,
VM(HSIC, HSIC) 4 NN(0,
(14)
The asymptotic variance of vm HSIC S can be consistently estimated from samples as 16 16116(R- HSIC where R = -
(n-1) Eijamety (0) h(..i.g.)) Here i. \ {i} denotes the set of all e-tuples drawn 1 E
(n=4))
without replacement from the set {I.....m) \ {i), and h(i.j. = 24 (s.t.or))) lea + luo 21 where the sum ranges over all 4! = 24 ways to assign the distinct indices j. to the four variables (s,t,u,) without replacement.
(1.3g.)))
The behavior of the estimator is different under the null hypothesis (X 1 Y); in this regime m HSIC (rather than vmHSIC,,) converges in distribution to something with complex dependence on Pz, Py, h and I. Proposition A.2 (Gretton et al., 2007, Theorem 2). Under the null hypothesis Ho Pay = P2, x Py, the U-statistic estimator of HSIC is degenerate. In this case mHSIC,, converges in distribution as
mHSIC, $ > (= 1),
1=1
iid where ze ~ N((,,) 1), and de are the solutions to the eigenvalue problem
Xebe(zz)= hejor dF
where hajor == lo(., is defoned in Proposition A.1, and Fig. denoles the probabilily measure with respect to variables 24, Zyz and Zr-
B Bernstein-based Justification for the SNR
For both the NDS and HSTC tests, our justifications of the SNR criterion for test power, (3) and (6), relied the asymptotic central limit theorem. A Berry-Esseen theorem would provide uniform bounds the o(l) terms that appear in these derivations, but one might expect for HSIC in particular that these might not be especially tight in some cases: since the null distribution is asymptotically a mixture of chi-squared variables (Proposition A.2), when the dependence under the current kernel is weak, the distribution may be further from normal.
We can also justify the same criterion, although less precisely, with a non-asymptotic argument. While we present it here for HSIC, it also applics in the same way to NDS. The Bernstein-type bound of Maurer (2017, Theorem 6) gives that, for kernels bounded as in Assumption (C"),
Pr HSIC HSIC| > 8) < 2exp
m
S 5.
20
2301
+
256
Solving for E, obtain
mme2
8
< log
2 51
2304
256
2
11 m-
2304
256
2
log
20
+
+
m
m
4
3
256
2
2
2
2304
2
3m
log
E
m
0
91
log
m
2 0.
12Typically n, but might want few (n) samples roughly estimate the power an m-sample with n, as done different context, by Sutherland Deka (2019); Deka & Sutherland (2023).
20

---

Published in Transactions on Machine Learning Research (03/2026)
This can be written 2² - B= C² for positive B, C; restricting to positive E, this occurs when 2 B + 2 2 +
and hence also occurs when e is larger than the (slightly larger but simpler) bound B + C. That is, with probability at least 1 ô, it holds that
256
2
2
2
2304
HSIC HSTC| v
log
log
+
3m
m
m
256
2
2
2304FR2
3m
log
555 V
m7 σ
log
+
m 4
Recall from the discussion preceding (6) that the rejection threshold I'm (V--(1 a1) + +(1) (where W
will vary with the kernel choice). Now, sufficient condition for the test rejecting is that both HSIC HSIC and HSIC - T'm Plugging in the deviation from above, we obtain power at least 8 as long
256
2304
HSTC > T' +
log
05
log 5
+
3m
m 4
V
2
log
2
8
+
1
v
'1 - a) + 0(1) + 48
1/2// m 4
m
256 absorbing
log 3 (L/m). Equivalently, the power condition holds if
3m
HSIC
48
2
2 1 1 LV--1(1-0 + ma
log m
OS
m 4
Shh
5 m
Jh
Now, we know that scaling k should not change the difficulty of the problem: indeed, it will scale HSIC, 0h1 11 - a); and Vuyin all equally. Thus, the only term in this inequality which is potentially sensitive to scaling is the final O term. From this point, then, while other changes in the kernel could potentially affect 41 a)/o5, in mcaningful ways, for reasonably large m and Øs, not too drastically small, the right-hand side is only weakly dependent on the kernel choice. It is thus reasonable to think that we can maximize the lent-hand side of the inequality, and ignore the right- hand side, to roughly maximize power. It seems highly unlikely, however, for the coefficient 48/Vm - to be exactly the right trade-off between HSIC and the kernel bound; this comes out of the details of the (conservative and somewhat loose) concentration inequalities to reach this point. To maximize our chance at having high power, then, it very reasonable to maximize the SNR.
C Uniform Convergence: NDS
C.1 Preliminaries
We start by defining the notation used in the following proofs. Let So be the critic function with parameters 0, and (X,, be samples drawn from the joint distribution, with X (Xi)==... and Y = (Yihe=
We define the population level statistics based on critic fo as
Ti, == E [f(()) and PxY
Too = E [o((X,)), PxxPy
with their corresponding finite sample estimates
1 n²
ni Jo(X,, Yi)
and To
Σ
So Yz).
n
We note that 7, just the NDS (2), and that To is the permutation invariant term in the numerator of Equation (12). For brevity, we often omit the subscripts 0, and/or the parameters 0 when it is uninformative clear from the context.
21

---

log 2 8 Published in Transactions on Machine Learning Research (03/2026)
We respectively denote the alternative and mull variances by
o Var (X, Y) Fxy
and S = Var [So(X, Y)],
Px Pxx
with corresponding finite sample estimates
2
00
fo(X,,
1
Σ
T5
fo(
71
n² j=1
and
2 :- 12 1=
fo(X,, -
Σ y=1y²=1
So So Xy)
We choose this notation, rather than 75, and 2h used in the main body, for distinguishability, As before, we often omit the subscripts when it is clear from the context.
C.2 Assumptions
Our main results assume the following
(A) The set of critic parameters @ lies in a Banach space of dimension D, where each parameter in the space is bounded by R:
C (0 : |00| S R}.
(B) The critic f is uniformly Lipschitz on X x y with respect to the parameter space O:
fo((,,)) (1r - L||0- 0'll
4 X, ),0 e.
(C) The critic is uniformly bounded:
sup sup ( y) B; # xe yEY
that such a exists is implied by the previous two assumptions.
C.3
Lemmas
Lemma C.1 (Concentration of Ri). For any critic f satisfying Assumption (C), with probability at least 1- 8 we have
h - < B
Proof. The result follows from McDiarmid's inequality. First, we show that 7; satisfies bounded differences.
Let ((()) denote the f, statistic but with sample Xg, Y, replaced by x', v. For any h E [n], we have
sup (7((k) .g
fi sup (f(a'.y') f(X,, Yc)) S n
1
2B
of
n
x
Xx.Y. Y
so that with probability 5,
1 n
2B
2
2
2
Σ
loog
B
log
2
n
n
1=1
22

---

Published in Transactions on Machine Learning Research (03/2026)
Lemma C.2 (Concentration of To). For any critic f satisfying Assumption (C), with probability at least - 5 we have
1 -Elool < 2B log
Proof. The result follows from McDiarmid's inequality. First, we show that To satisfics bounded differences.
Let 1O/) denote the To statistic but with sample Xs, Yr replaced by x', v/" For any k E [n], we have
sup (704) Pol sup x,y
Σ
((('))))-(Xxx88 Ki)
a'
x
J#k
Σ ( - ((x,;
Σ I(X,, Pt))
it
z=j=
2B - 1)
22B - 2B 1 +2B
n²
n²
that with probability 1 6,
2
2 5
2 j
1 - -Efo S 2B
log
< 2B
log
n2
2n
for all integers n >
Lemma C.3 (Bias of To). For any critic f satisfying Assumption (C), the bias of To is bounded by
2B
EEO - Tol -
n
Proof. Notice that
Eí =
1
Σ ((x,,) ((x,.)))
n
To +
-Ti 71
n1 iデj
72
i=j
so that
1 ETG Tol
- [Ti
Tol s
2B
n1
n
Lemma C.4 (Concentration of 6²). For any crilic f salisfying Assumplion (C), with probabilily at least 8 we have
2 8 V 2n log
8² EE 5B2
Proof. The result follows from McDiarmid's inequality. First, we show that 0² satisfies bounded differences. Let i i= I(x,, Y:)² be the 2nd moment estimate so that 6² = Ai2 - 17. Let o( (2, Fi,(s) denote the corresponding statistics 0², f0z fi but with sample Ng, Yn replaced by x', v/ For any l E [n], we have
sup 107 821 < sup [02,0K) A2 + sup x'
T. 1?1
a' .
.9' x
x,7
x
1f(2 11/2211112 - f(X,, Ya) + sup (i,a) f:1li + fil
sup
.v'
X
x
B2 4
5B3
23

---

Published in Transactions on Machine Learning Research (03/2026)
where we use the finite differences computed in the proof of Lemma C.1 to bound (Ti.(k) - fil. Thus, with probability at least. - 6,
1
2 8
6² - Eo² S 5B
log
2n
Lemma C.5 (Bias of o²). For any critic f satisfying Assumption (C), the bias of d² is bounded by
B²
EE - - 0²1
Proof. Let P = - f(X,, Y2)² and N2 = Rl = Exx f(X, v)2. It follows that
27
E0 = -Elū- - -P-1 -E E { H s(Xi, Yi)
112
E
ΣΣ 1(x, xXxCX)))
n1
122
1
E
Σ
f(x,, vi)? + Yi,,,,,,, Y3)
n
i=j
1
(mp2 n(n -1)7?)
n
M2
n²
n
which yields
B²
E
E
0²1
0²
n1
n
Lemma C.6 (Concentration of (2). For any crilic I salisfying Assumplion (C), with probabilily al least
5 we have
[ - < 10B2
log
2n 8
Proof. The result follows from McDiarmid's inequality. We first show that E satisfies bounded differences.
Let D2 = Eisemm f(X,, Y2)² be the 2nd moment estimate so that E2 = b2 T3. Let Ely (2, To, denote the corresponding statistics ? b2 To but with sample X2, Ye replaced by x',y'. For any & [m] we have
sup Eiey
El < sup s'.0' X
12,(0) - 022 + sup a'.y' X,7
To - -
sup g X
f(X,, 1 + sup
To,er - TollTo, + fol
n
e's3' X
i=
or
j=t
B (2n 1)
2B(2n - 1)
10B²
5B²
+
2B <
n²
n²
n²
It follows that with probability at least 1- 8,
[ E S
10B2
5B²
V
log
2
n
2n
Lemma C.7 (Bias of (2). For any critic f satisfying Assumption (C), the bias of Q² is bounded by
7
E - {2]
<
24

---

Published in Transactions on Machine Learning Research (03/2026)
Proof. Let V2 (f and #2 Exxy [f(XX ! be the second moment under the mill and
alternative distributions, and let D2 = nス Eis [(X;; 13)² and j //2 ((X,, Yi) be their estimates. First, we can write the bias as
5² = Elb - Po - (w2 13) - Ei 22 + |
To
Now since
Eig
nt
E
Σ (ix..))) + E(xx
n2
(unn + n(n - 1)12)
422
n2
1/2
171
we have that
|E 02
v2l = = -22l S B2/n
For the remaining term, note that
((xx,))) )
E
E
1
Σ
1 E n2
Σ idigi
f(X,, (3)((,)))))
n
in
1
E
Σ
Yx)
Σ
s(Xi, (x)s(Xx, Y.)
n²
(1.....))EEz
1
[( + (n² - (n²)) Ba (1
6
11
6
6
11
6 n3
$
T
B²,
n2
11
32
72
n2
so that
[[[ - 78 = 1 # + 6 (B² - T8n1 $ (
#1
+
n3
6
Therefore, putting it all together we get
7B2
11
6B
7B
B -² <
n22
71
C.4 Propositions
Proposition C.8 (Uniform convergence of fi). For any critic family {fo 0 E O} and parameter space θ satisfying Assumptions (A) to (C), we have with probability at least 1 - 5 that
2L
B
sup 6ee
n Tı,ol S nn
= log + 1-2 log =0 vn
Proof. We employ a covering-number argument on the parameter space. Let = (0: C θ be p-cover of O; i.e., any point in θ in the closed p-neighborhood of some point in U. Assumption (A) ensures the minimal cover exists with at most N (4R/p)" points (Cucker & Smale, 2002, Proposition 5). Now, for any 0 O let 0' be the center in U closest to 0, and then decompose the bound into
sup
Ta Tol < sup To To max To To + sup Ter Tol e
8ee
(a)
(b)
(c)
We first handle terms (a) and (c). From Assumption (B), we have
sup Top Tol S sup E. fxx, sup 1.100' - A|| S Lp.
#eB
This same bound applies for (a) by replacing E with its empirical expectation É..
Next, we find a high probability bound on (b). Recalling the result of Lemmn C.I, and then laking the union bound over all N points in U, we get that with probability at least 1 - 5
2
2 8 log + D Doo
4R p
max To Tal B
n
25

---

Published in Transactions on Machine Learning Research (03/2026)
Finally, putting it all together, we get with probability at least 1 - 6 that
2 [oo Tol 2Lpp + B n
2
4R
sup
log
D log
EE0
p
which, for p p= 1/V gives the desired result.
Proposition C.9 (Uniform convergence of To). For any critic family { : 0 } and paremeter space θ satisfying Assumptions (A) (C), we have with probability at least 1-8 8 that
2L Vn
2B
2B
2 + 2D 2Dlog(4R/m) 2 2og vn
sup Fo.0 To,o S 85
n
Proof. We employ a lovering-uumbe argument on the parameter space. Let 10.)i- S θ be p-cover of O; ie., any point in θ is in the closed p-neighborhood of some point in U. Assumption (A) ensures the minimal cover exists with at most N p) points (Cucker & Smale, 2002, Proposition 5). Now, for any
0 E O let 0' be the center in U closest to 0, and then decompose the bound into
sup To Tol I sup To 0
To max fTg Efol max ETw
To + sup To Tol
o
O0
(a)
(h)
We first handle terms (a) and (d). From Assumption (B), we have
sup Tor Tol I sup E (X, (x) Jo(X, P)1 S sup 1/10' - 011 s Lp, o0
oee
which is the same bound when considering to, instead of To.e-
Next, we find a high probability bound on (b). Recalling the result of Lemma C.2, and then taking the union bound over all N points in U, get that with probability at least 1 1-6
2 n
2 D log
4K
max
Tog Etol S 2B
log
BELL
Finally, combining everything and using Lemma C.3 to bound (c), get that with probability at least 1 5
2B
4R
sup fo Tol 5 2Lp +
+ 2B n V
log + D log
0
p
which, for p = 1/ VT, gives the desired result.
Proposition C.10 (Uniform convergence of g²). For any critic family { : 0 E Q} and parameter space θ salisfying Assumplions (A) to (C), we have with probability at least - 5 that
8BL B² n
5 2
2 + +2D log(o -(n 2 5
sup e=0
| o
n
Proof. We employ a covering-number argument on the parameter space. Let 4 = (0. C θ be p-cover of θ; i.e., any point in θ in the closed p-neighborhood of some point in U. Assumption (A) ensures the minimal cover exists with most N (4R/p)D points (Cucker Smale, 2002, Proposition 5). Now, for any O let 0' be the center in U closest 0, and then decompose the bound into
sup 0 6Ee
2
sup 102 S
max log Peu
E. a
max
E. E
sup l o
OE
e'ell
(6)
(c)
(d)
26

---

Published in Transactions on Machine Learning Research (03/2026)
We first handle terms
(a) and (d). From Assumption (B), we have
sup 0E8
o3 S sup AE
[ (X, r)² fo(X, 8) + sup BEB
(E ((E)))))-EEEEEET) V))²
$ sup E (X, - fo(X, x for (X, + fo(
sup E (X, - fo x) E (X, + fo(X,
8CB
sup LI|0' - 011 2B + sup 1110' - 011 4BLp, ac
#E0
which yields the same bound for (a) when replacing E with its empirical expectation E.
Next, we find a high probability bound on (b). Recalling the result of Lemma C.4, and then taking the union bound over all N points in U, we get that with probability at least 1 1-6
1 2n
2
4R ) + Dlog
max
|0
Eojl S 5B2
log J
DEl
Finally, combining (a), (b), (d), and using Lemma C.5 bound (c), we get that with probability at least
1-5 6
sup 0 - ou v 8BL B² + 5B
2
4R p
log 8
Dlog
EE0
2n
which, for T/, yields the desired result.
Proposition C.11 (Uniform convergence of (?). For any critic family {fo : 0 E ))} and parameter space θ salisfying Assumplions (A) to (C), we have with probabilily at least ô that
8BL Vn
7B2
10B²
2 2/n V 2 8 20222ooog44RB 0On
sup
+
n
Proof. We employ a eeverng-nuumbe argument on the parameter space. Let 4 {0i]i=I. S θ be p-cover of O; i.c., any point in θ is in the closed p-ncighborhood of some point in U. Assumption (A) ensures the minimal cover exists with at most N = (4 points (Cucker & Smale, 2002, Proposition 5). Now, for any 0 E O let 0' be the center in U closest to 0, and then decompose the bound into
sup K 2l S
sup 1 BE0
max B'ELL
Sø
E
max B'ELl
E Eg-
sup
S REE
cl
#Ee
First, we bound terms (a) and (d) using Assumption (B) and the covering-mumber number argument. Let To = [x,o(X, and Uy [X, Y)a] be the first and second moments under the null distribution, so that 5 = Ug TZ. It follows that
sup E.- 3 < sup law val sup T3 To AE
sup E S (X,Y) fo(X, Y
+ sup Tor o2e
| Tor Tol
4BL sup |10' - 0| 4BL.p
We obtain the same bound for (a) by replacing E with its V-statistic estimator.
Next, we bound (b) using Lemma C.6 and a uniform bound over all N points in the p-cover. We get that with probability at least 1 8
1
2 g log
4R p
max BELd
S0 E
E
10B2
+ D Doog
2n
27

---

Published in Transactions on Machine Learning Research (03/2026)
Finally, combining (a), (b), (d), and using Lemma C.7 to bound (c), get that with probability at least
1
7B²
4R D log
sup
E El
8BLp + 10B²
log
AEe
n
2n
p
which, for p = 1/ Vn, yields the desired result.
C.5
Main Results
Theorem C.12 (Uniform convergence of (Dg)). Let {fo θ} be critic family with parameter space θ satisfying Assumptions (A) to (C). Suppose of s² and E s2 for positive 8, and n-133 Then, with probability at least 1 - 6,
2B
16B
5
12
sup JUDDS JNDS! 65
+
1/3 V
2 log 8
+ 2D log log(4IlV)) + o(n -1/3),
1/3
1/3
and thus
jNnss JNDS! =Öp
1
B
BBB+BBVVD
sup 6c
1/2
Proof. Let Za denote the standard normal a-quantile. We start by making the following decomposition:
T1,0
To.0
T1,0 To,0
Eo
sup 8c0
JNDSI
sup #ce
2-a
21
Œ8.A
0g
Vnoa
Vn
fi.0
Ti.0
To
To.0
21-a
60.
E0
sup
sup AEe
sup
#
0g
ōg
Vn
BEN
0p
(a)
(b)
(c)
Define 0b oo + and notice that 0i 0 0², o + 2 52 + X and 0ū = d A A.
since it is the biased variance estimator; this implies that
Now, rewrite (a) (with subscripts A omitted) as
f
T
fi
i
T
sup 8EB
< sup 0:
+ sup 8EB
+ sup 8E0
σx
0 03
of d²
sup AEe
+ sup AE0
+ sup #E0 σ
+ Ti
θxO 86 σ
O2 +
B
BA
sup 102
0²
sup
T - Ti
AV + 1 Vxs2 A)
85B
s(82 X)
22 2 X
6
B
sup 102 - 02|
sup
Tı
Vxs²
BE0
8EB
where the first. term in the last line the slowest, growing upper bound X - 0.
28

---

Published in Transactions on Machine Learning Research (03/2026)
Employ Propositions C.8 and C.10 as uniform bounds on [fi - -Til and Jō² -o²l respectively, then take the union bound to give us, with probability at least 1 - 6,
ii
Ti
L
B
5B
4 2D +2D log(4R\\))
sup Ag0
2log 8
0
0
82 n
Xnv
28 Vnd
nd
BA
33
2L
B
+
22og
+ +2D 22Dlog(RR/n)
s/n
sVn
BX
L
2L
B
5B3
syn )
B
4 g 2D log m.)
+ 82
+
2 log
8
VnV
8 syn
sanv
2n2Vnd
The best-case (fastest to decay) asymptotic bound follows when X 1/3 which yields with probability at least 1 - 8,
K xx σ
Ti
B san 1/3 +
L
5B3
4
sup 6
1/n
2s 1/3
2 2
2D log(4R\/)
Term (b) can be handled in similar manner, except we use Proposition C.9 as a uniform bound on 17o -Tol- It turns out this gives us the same asymptotic bound as (a).
For term (c), we first define Ea, E + X. Since E is the biased estimator we have & 2 0. Also,
2 2 + x N x2 x and Êa.A = Eg + + 2 A so that
3
9
mm
sup
< sup
+ sup OE5
sup
+
sup
σ
OEE
ax
X 82
6 oR
o²
{2 g2
sup
+ sup AEG
|g|
sup
(Ê)
sup
O2 + 3
θ (bx + σx)
8E8
( A + σ)
o(Ê E)
B
BA
1
<
sup
182
0²|
sup
- 52
Avg2 +
1s82 12 + A)
s(32 A)
8
8
1
B
sup OEN
BA
+
1
sup
_2 - g2
2
Propositions C.10 and C.11 are uniform bounds on |@2 gives us, with probability at least 1 8,
0² and
[2 - -2 Then, taking the union bound
BX
B
1
8BL
B
7
B²
sup l
σx
σ
33
82
V
VnI
82
B
2
5B²
4
2logg 8
+ 2D log log(AR//))
82
V
$2
2
2
When À = 1/33 we have 1/VnA=n~'> and 1/(n//) n O(1). Then, multiplying by the factor z1 =1 vn in front of (c) gives us order 0(n-1/2) Finally, taking the union bound of (a), (b), and (c) gives us, with probability at least 1 - 8,
5/6 5/6
that the overall order of this bound is
2B
L
5B3
12 + 2D 1olog(4R/m + o(n - 1133.
sup JADs JNDS| < 6ee
3/3 +
+ 2/ V 2 log
1/3
5
29

---

Published in Transactions on Machine Learning Research (03/2026)
Theorem C.13 (Consistency [van der Vaart, 1998, Theorem 5.7]). Let In be the finite sample estimate of J with somples, as function of the parumeters 0, and suppose that Jn salisfies uniform convergence
sup (0) f O. ocê
If 0* a well-separated marimum of J, then any sequence of estimalors Ó, approximalely marimizing In, i.e.
Ja(0") In (0m) S or(1),
is consistent estimator, i.e. On & o"
D Uniform Convergence: HSIC
D.1 Preliminaries
We start by defining the notation used in the proofs. Let kernels has and ly be paramcterized by some w E Ω and y e r, with _ Ω I as the joint parameter space. The samples (X,, Yi) ~ Pay are drawn i.id. from the joint distribution. We denote the will often omit the kernel parameters w and when it is clear from the context.
gram matricos of hu and l, by K(w) and I respectively. We
Let be the HSIC test statistic and f to be its U-statistic estimator given by
Σ Hijor
i (1)4 ((.--.)))
where (n)) = ml/(n-
(n- k)! is the Pochhammer symbol and i7 is all possible 4-tuples drawn without replacement from 1 to H is the kernel gram matrix of the U-statistic defined by
1 (1i.)))
Hijer
Σ
Kab (Lib + Led - 2Lac),
4!
(a,b.c,d)
where sum represents all 4! combinations of tuples (a, b, c, d) that can be selected without replacement from r),
D.2 Assumptions
Our main uniform convergence results require the following assumptions.
(A') The set of kernel parameters Ω lies in Banach space of dimension Dn, and the set of kernel parameters I lies in Banach space of dimension Dr. Furthermore, each parameler space is bounded by Rg and Rr respectively, i.e.
(w |||| < Ra).
r T (r 1r Rr).
(B') Both kernels ls and are Lipschitz with respect to the parameter space: for all xr, xf & and w, w/ E Ω (2, x') - (x, s LKll - w'll,
and for all u. of e y and 7,7' Er
(v.v0') - by (v, 1/1 S Lall7 - 711-
(C") The kernels bw and Lo are uniformly bounded:
sup sup hew (z, x) $ vk.
we
sup sup L, (y y) S v. >et k
For the kernels in practice, ii 1 1.
30

---

Published in Transactions on Machine Learning Research (03/2026)
D.3 Lemmas
Lemma D.1. (Concentration of 0 ), For any kernels & and 1 salisfying Assureplion (A'), with probabilily at least - ò we have
2
- S 6144 V log 8"
Proof. We apply McDiarmid's inequality to o². First, we show that the variance estimator satisfies bounded differences. For convenience, denote (i,), 12 simply as (.r.. and (i, ls to be the set of 3-tuples drawn without replacement from that exclude the number k. Recall that
8² 16
Σ (k.i...)) (b,c,d)
Hajer Hibed fp
(n)a(nn 1)3
Let. F denote the kernel tensor H but with sample (Xe, Yr) replaced by (X), Y2) so that F agrees with H except at indices é, and let f/ and o² denote the HSIC and its variance estimators according to this updated sample set. The deviation is then
16
102 - (n)a(m- - 1)3
HeiggrHibce FeborF + 16] f/21.
(a.j.g.r) (h,r,d)
We bound the first term by noticing that A = Hijgr Hihed Fijor Picd is zero when none of the indices {i, j, 9, b, df is é. Let S = ((i.j...r r,,,,, d) (i.j.g,r) E i2 (b, d) E M7Mi, L E r.jj.grr,b,c d}} be the set of indices where A may be non-zero. By Assumption (A') we know that | S 32r Thus, we can bound the first term by
16
512 (n))(n 1)3
ΣTAI=
|S|
(n)a(n- 1)3- 8
512
9(n - 1))(n 2) 4(n - 1) + 3(n 1-)11 222(n - 1)3
(n)a(n - 1)s
rElidagrt
RE(he,,)
BE (j.g.)) and ee (b.r,d)
512
81920307
(Vn > 1).
We can bound the second term using M < 40 and the bounded difference result (15) from Proposition D.3:
16/0 - - 16 + // 1'' S 16 Sugn 32v 4096xZ0?
Combining these two terms, the maximal bounded difference for a² is
12288v
102 -
0 -
Finally, applying McDiarmid's inequality gives us, with probability at least 8,
2
66" s 61440Ru²
V
log 5
Lemma D.2 (Bias of ). For any kernels ls and satisfying Assumption (A'), the bias is bounded by
4608v2u?
E 2
0²1
n
31

---

Published in Transactions on Machine Learning Research (03/2026)
Proof. The expectation of the variance estimator is
1
E 16
E Hige
E ElHes Habed)
(n)a(n 1)a
(n)
(4.-9.)) (b,e,d)
(KJ.a.))
(a.b,e,d)
First, we can break down the left-hand sum into only terms of E1 H1577 by considering the cases where (,r,, d] are unique. Tet S {(i,j, possible indices of our left-hand sum. It follows that
r, g,
be the set of all c, d) (i.j,g,r) E (b, 13\{i]
EB. Hibed]
Σ (ii.a.bbbe
ElHtigr Hibea]
Σ 8/14
If all indices are unique, then the expectation HaHdd is equivalent to can bound the expectation by 16 via Assumption (A'). Thus, the bound on the left-hand sum is
otherwise, we
B[H Habed] & (n17 123 1167] + ((m)a(n- - 1)3 - (m)7 16UR
(B.-g.r) (b,c,a)
Similarly, we can break down the right-hand sum into only terms of H2 H5678] Let R (i,j,,rr d, (i.j.,.) 1n (a.bb,.c. e it} be the possible indices of our right-hand sum. We have that
Σ Billiggr Habed] =
Σ
Bllligr Habedl + Σ
Bl Habeel
((i.g.)) (a.b.c.d)
1915
S (n)s (25H122 + (nm)3 - (n)s) 16v2v?.
Now, using these two results and Assumption (A'), can compute bound the desired bias of 8²:
E 021- = (n)4(n- - 1)3
EjH Hstel] E H25
Σ
E Hijar Hadeal + 1211
(n))
(f...g.)) (b,n,d)
(L..g.))
(a.b.r,d)
(m)7
(n)s
< 16
1
16ugu1
B|H1234
E E|H H3878] - 16
(n) 1)3
(n)
o5 $100gu?
(n))
16 1
32
5512221
(Vn 2 4)
(n)4(n - 1)3
4608vRy?
n
D.4 Propositions
Proposition D.3 (Uniform convergence of HSIC,).. Under Assumplions (A') to (C'), we have that with probability at least 8,
SugM
LE
2 x + 2Da loggoooo///) + 2Dr log(ARr L + 2 2 log +
sup f10 nol S
Vn
ux
in
Proof. We use e-net arguments on both spaces Ω and r. Let (wi)Tn be arbitrarily placed centers with radius pa such that any point Ω satisfies min ||w will pa- Similarly, let (r be centers with
32

---

Published in Transactions on Machine Learning Research (03/2026)
radius pr satisfying min || --11 - pr for any ? Г. Assumption (A') ensures this is possible with at most To = (4Ra/pa) and Tr = (4Rr/pr) points respectively (Cucker & Smale, 2002, Proposition 5).
We can decompose the convergence bound into simpler components and tackle each component individually
sup he - mel S sup he forl
max
fiw 170 + sup me n7l.
E c'f T'E(71...)7p)
First, let us analyze Me m for any 0, O' E O. Recall that H123) where H1231
4 (a,bead)) Kab Lab +
1 E (1,2,344)
Led 2Loc). We have that
(1234) Σ
11(0)
11122
1234
k (L (L + +L 2L(22 k(m²)) (L + +L -2L(1t))
4!
(abed)
(1234)
[(0) (7) - (w)) (7)) 1
ke (
Lod K L + K (w
(T)
(w0)) (T')
4!
(abed)
From Assumption (A') we know that |Kabl S uy and |Labl S vi, and via Assumption (B') we notice that
k(w) (()) 1<(="
(((a) L(r)) -k(l) L(6))
L((()
ab
ab
ab k(w) k(w(() ab
k(w)
1,(7) -1((7') + [,(7²)
s ab
ab
s LLllll - r'll + in In |ww - w'll
S n Lipr + nLLpn
This expression true for all three components of |
H 1123al and so follows that
nno norl
E 1234
(0)
B(H
1234
E
H
1234
Lip + 4v1Lapa.
norl
(0) Gjrr
H
Σ
H (0)
H12dd I AugLipr + 4u Lnpa-
(n)a (6.igg-))ei
(n)4
(Ki.))))
Now, we study the random error function A = ター - Note that E = 0 since v is unbiased, and [Hsgrll I 42g2 via Assumption (A'). This f and hence A, satisfics bounded differences. Let F denote the kernel tensor but with sample (Xe, Y2) replaced by (?) so that agrees with except at indicies f, and let f/ T E(idaar)tty Fsjer be it's HSIC estimator.
For convenience, we denote (i.j.q, r) E i7 simply as (i, J, q,r), and (1, 4) k to be the set of 3-tuples drawn without replacement from i that exclude the number We can compute the maximal bounded difference A A'l = f//
In f/1l
Hijer Fejgr S (n)s (6.5.g.))
Σ
Haijr Feigr!
(15)
(n)) (k...g.))
Hejor Fejgr
+
Hisgr Fitgrl + Σ
Hijer Fijerl +
Σ
Hiigt
Fisgel
(n)s
(eg.r)e
32
(n)1 (nn 1)3 1.V8uv-
33

---

Published in Transactions on Machine Learning Research (03/2026)
Then, applying McDiarmid's inequality on A : 17 followed by union bound over the ToTr center pairs gives us, with probability at least 1 6, that
TT
max w'forr.. n' E ((2
1i0 170ll - 32
log
2n
8
16
2 log Th + 2 2og Tr
2 2log
Vii
16
2
4Rg
4Rr
2 2lo
+ +D log
+
2Dp log
vn
Pn
pr
Finally, we combine these results to get our uniform convergence bound:
16vgn
2
4Rg 2Dr
4Rr
sup lee
M/0 mo0l S Sun Lupr
Sun Lupa
2 log + +DD log
log
/ni
P
pr
Ik pat
La
2
2
4Ra
4Rr
pr +
2 2log
+ +2D log
+ 22r log
Vk
vi
vn
8
Pn
pr
Setting Pn Pr 1//n yields the desired result.
Proposition D.4 (Uniform convergence of 02,7) Under Assumptions (A') to (C"), we have that with
probability at least. 1 - 6,
2048vI Vn
Lk Li
9 2 + 2D log(4Rn + 2Dr log log(ARrVT) + +3/ log
sup jo lE0
UR
VI
4yn
Proof. We use e-net argument on both spaces Q and Т. Using the same construction as in Proposition D.3, we once again decompose our convergence bound:
sup dee
lo ogو s sup 100 0g. w WWW y'Elm....
max
160 og.
sup oo o²
0
First, let us analyze
10 - og for any 0, 0' E θ. Recall that
²= 16 [E
H133 H 1567| - 1)²). It follows that
10 o = 16 (()
(8) (()))
H(00) (1)))
E(H(
E
(1)
E
(8)))
1234
1567
1234
1234
5678
1231
< 16 E
(0)
1234
(0) H (0')
1567
11 (0')
1567
+ 16E
1234
(0)
5678
H
(0')
1234
5678
Under Assumptions (A') and (C') we know that |H1234| and
S 4unu and |11233
11 S Japr + Aur Lapas
(0)
(0')
H1 S (0) (r(0)
11 (8') + 1234
(0')
(0')
11
1234 1567
1234
1667
1234
1567
1567
1234
(1)4 (0)
(88)) + H (8') (
1234 H(
1234
1567
1567
< 32 LIPT - in L.P))
This expression true for both components of 1o - o and it follows that
1o ob S 10 + LLpp).
(16)
Similarly, replacing the expectations E H H1567] and 12 H5678] with the respective estimators (m)a(n-1)s E(sor)) Hijor ihed and (n)] Σ E(tigr)(abea)) Hijgr Habed give us the same bound
10 0b. 1024m vLLxpa).
(17)
34

---

Published in Transactions on Machine Learning Research (03/2026)
Next, using Lemma D.1 and Lemma D.2 followed by union bound over the T&Ir center combinations gives us, with probability at least. 1 - 8,
2
2TaIr
4608
max w Ee y (
10 0z. S 6144
log
(18)
8
2048vRu2
2
4R
4RD
9
31
2log
2
2D log
+ 2Dr log
+
Vn
6
P
pr
4/n
Finally, we combine Equations (16) to (18) to get
2048vgu7 3 22 2 + +2D log
4Ra
4Rr
9
Tok
Ia
pr) )
sup j0 8EB
+ 2Dr log
+
+
vn
patt
+
vn
Pn
pr
4/7i
ux
n
Setting Pn pr 1/ Vn gives us our desired uniform convergence bound.
D.5 Main Results
Theorem D.5 (Uniform convergence of f (isic). Under Assumptions (A') to (C"), let Q be the set
of kernel parameters E @ for which oî 2 s², and take X = n 1-133 Assume Vk. i1 2 1. Then, with probability at least 6,
9216
sup 00
00,
06
I
n/3
Vn
88 ) + 122888 + n/1
(
Tok (Th
La + v 2
+
2DDn lonlog((R + 2Dr log vn)
un is
and thuus, treating Vk, Ut constants,
170
Öp
La La Dr + vDrl)
sup A50
00.2
σg
1/3
Proof. Let 07 = 07 + x be our regularized variance estimator from which we can assume is positive. We start by decomposing
fi
sup #
< sup 8E9
+ sup 8E0
+ sup
θ0,A
09
T0,A
TA,
0H
0
09
[ipl
100
22
100
0
sup
+ sup
+ sup
170
0o. 00,
σ0.A + 00.
σ0
08.
0g
41gui\
sup |03
og
sup
SVi(S + VX
SV g + X(s
32
4Fr
2RRuA
sup |03
0z
sup
V0o
mol
2² 0
Proposition D.3 and Proposition D.4 show the uniform convergence of ne and ãg, from which we get that with probability at least 1 - 6, the error is at most
10
770
20r
184320Bu
81920 2/n
SUkUI syn
Tik
La
sup lc0
Ō0,A
Da
V/
V1
2457608u2
16
2200
2Da log(4Ravn)) 2Dr log(4Rr
vm
sVn
35

---

Published in Transactions on Machine Learning Research (03/2026)
Taking A n- 1-1/3 gives
100
7]0
2/
18432v8u2 s²n55/6
8192v 1/3
Sugu
(
La (L
+
sup 8
00,
00
sV
us
vi
24576ugu2
16 1/3
V2log
+ 2Dn log log(4Rav/) + 2Dr log(4R vn))
$ syn
Using Mk, Vt 2 1 we can slightly simplify our bound to
110
244776uu
16 syn
Lk
Li
sup 8c0
Ōg.A
1/3
Vk
vi
+ =
2log 2D2 log(4Rav/) 2Dr log (og(4Rr
18432v
+ 1/3
82n5/6
E Experimental Details
Figure 7: RatInABox ABox simulation environ- ment. The red dot the current position of the rat and the purple circles indicate the past trajectory over 5 seconds. The box is designed to have only single protruding wall.
Figure 8: Samples drawn from the sinusoidal problem with frequency f 1 (left) and f (right). We consider the latter frequency in our experiments.
E.1 Training & Test Details
We design the featurizers dw and ΦN of our deep kernels kn and l to be neural networks with ReLU activations. We avoid using normalization it may affect test power. Moreover, make the Gaussian bandwidth of both l and l learnable parameter, possible, we use similar neural network architectures for each deep learning based method. In general, we let the featurizer HSIC-D and MMD-D be identical up to a concatenation layer which concatenates and frame the problem as a two-sample test.. We construct the C2ST-S/L. classifier as the MMD-D fcaturizer plus a lincar layer classification head with scalar output, and we let C2ST-S/L, InfoNCE, NWJ, and NDS all use identical architectures for the classifier and critic. Detailed descriptions of each architecture are demonstrated in the following subsections.
well as the smoothing rate 6. To make comparisons as fair as
JDS//nfoNCE/NW. MMD-D, C2ST-S/I. are first trained on an identical split of the data, and then tested on the remaining split. In contrast, HSIC-M selects the median bandwidth based on the entire dataset, and evaluated on the test set. We train our models using the AdamW optimizer with a learning rate of 1c-4 over 1,000 epochs for HDGM and RatInABox, and 10,000 epochs for Sinusoid and Wine. We use a batch size of 512. All methods are implemented in PyTorch and trained NVIDIA A A100SXM4 GPU.
All optimization-basce methods DUSIC-D/Dx/2A
36

---

Published in Transactions on Machine Learning Research (03/2026)
Regarding the power vs. test size experiments, we use training size of 10,000 for HDGM & 30, 100,000 for HDGM > 30, and 2000 validation samples for all dimensions. For the Sinusoid problem, we train on 5,000 samples and use 1,000 for validation. RatInABox uses a training size of 4,000 samples without validation, and Wine uses 1,200 training samples without validation.
Once learned, each methods' empirical power evaluated on 100 test sets (S ,.... S 0t) Each test set contains m test samples S (4 y... I ) which are then used compute the average rejection rate under the mull via a permutation test. We use 500 permutations for cach test and with predetermined type-I error rate of 0.05.
E.2 Validation
Although validation sets were only used for early stopping in our experiments, it's certainly possible to use them for hyperparameter selection. The validity of permutation tests still hold as long as this validation set is separate from the test set. As for the validation criterion, a natural approach is to evaluate the same approximation of the asymptotic power used for our training objective. This avoids the need for relatively-expensive permutations, while providing a reasonable estimate of performance. Alternatively, we could also directly estimate the test statistic in cases where there isn't simple expression for asymptotic power, though this approach may be less reliable since the statistic is not necessarily correlated with test power. It is also possible to select hyperparameters through rross-validatin,, where the training data is split into folds: one fold for evaluating hyperparameters and the remaining folds for training the model. Again, as long as this data is separate from the test set, any permutation test would be well-controlled. The drawback to these validation-based approaches is that effective power is reduced since we need to construct this validation set from either the training or test data.
Ideally, we would like to perform hyperparameter selection while avoiding sample splitting, but this is a chal- lenging problem. One potential approach is to perform separate tests corresponding to each hyperparameter, and then combine these tests using a multiple test correction procedure like MMDAgg (Schrab et al., 2023) to ensure valid level control on the test data. However, this is potentially expensive if hyperparameters need to be optimized. It's also unclear if this would perform better than simply using a validation set due to the conservatism of the multiple test correction procedure. A related approach is MMDFuse (Biggs et al., 2023), which constructs an aggregate statistic by combining MMD estimates at varying bandwidths. This approach could be naturally extended to consider estimates based on different training sets and would naturally fit into optimization framework; though doing so in way that does not explode compute reduce critic quality is another challenge.
E.3 Architectures
In all experiments we consider deep kernels with Gaussian feature and smoothing kernels and q, where each bandwidth is a trainable parameter randomly initialized around 1.0. We let the smoothing weight also be a learnable parameter initialized to 0.01. No batch normalization is used and all hidden layers use ReLU activations. Dataset-specific designs are elaborated below.
High-dimensional Gaussian mixture. We use a feed-forward network for our deep kernel featurizer with latent dimensions 2d, 3d, and 2d. Details of each model is given in Table 1.
Sinusoid. The deep kernel fcaturizer is taken to be foed-forward network with widths 1x8x12x8. C2st, infoNCE, and NWJ use a similar architecture one with widths 2x8x12x8x1 which includes an additional scalar output layer.
We use a ford-forward fcaturizer with details given Table Unlike the previous two problems, the sample spaces X' and y are not equivalent, and so the deep featurizers for and I have different. architectures.
RatInABox.
37

---

Published in Transactions on Machine Learning Research (03/2026)
mndde
Impos
ISSC-D
Y
2- 1--44-6--
HDGM-1
MMDDDD
4448
15
CUST-8/L
X,
HDXM--3
MMD+D
X Yi
30 -+ 00
HSIC-D
Y
1448
CuST-S/L
X, Y1
30
BOGM-8
MMD+D
X
CIST-8/L
X
16
HIXM 40
MMD-D
X.
40 80
120
80
ISICCDD
10 15
C2ST-8/L
X. Yj
30
50
HIXI5C-10
MMD-D
X.
10
HSTC-D
X
75
CIST-8/1
X
20-20--30--0
IIDM so
MMD-D
X,
50 100
130
HSSC-D
Y
10-+20+30---
C2ST S/L
X. Yj
30 -> 100
130 - 200 1
IDKGM
MMD D
X] Y)
20
11
Table 1: Featurizer architectures used in deep kernels for HSIC-D, MMD-D, and classifier architecture used for C2ST-S/L on the HDGM problem. Brackets denote a sequence of linear layers with correspouding input. and output features.
method
input
network
X
[ [8 + 32 - 64 + 32
HSIC-D
2 22 4
MMD-D
[X, Y]
10 - 32 + 64
C2ST-S/L
X Y]
[ 10 + 32 64 32
Table 2: Featurizer architectures used in deep kernels for HSTC-D, MMD-D, and classifier architect used for C2ST-S/L on the RatInABox problem. Brackets denote sequence of linear layers with corresponding input and output features.
F Additional Experiments
F.1 High-Dimensional Gaussian Mixture
We provide comprehensive power versus test size results for the HDGM problem at all dimensions (2,4, 15,20 20} in Figure 9. Overall, our method HSIC-D/Dx achieves highest power the smallest test. sizes.
Additionally, we demonstrate the effectiveness of our method at various dimensions d/2 by examining the empirical test power at HDGM-d for E 84,8, 10, 20, 30, 40, 50} with fixed test sizes m. Results are shown in Figure 10. Again, HSIC-D exhibits the highest test power across all dimensions. When using small number of test samples (e.g. m = 100), the performance of HSIC-D slightly degrades with increasing dimension, whercas at larger test sample sizes consistently has ncar-perfect power.
F.2 Type-l Error
Table 3 shows that the type-T error rates for optimization-based tests well-controlled.
38

---

Published in Transactions on Machine Learning Research (03/2026)
ISCC HSCCOr
NZS
HnONCE MA]
C25TL
C25T5
NMD-5
MSICO
-6-
HBXCM
HSKC-Ag
1.0
1.0
empirical power
empirical power
0.8
empirical power
empirical power
04
0.4
24
0.2
02
02
22
0.0
0.0
0.0
0.0
500
1000
1500
2000
500
1000
1500
2000
500
1000
1.500
2000
500
1000
1500
2000
test samples
test samples
test samples
test samples
HDGM-4
(b) HDGM-8
(c) HDGM-10
(d) HDGM-20
1.0
empirical power
empirical power
empirical power
02
0.2
0.0
0.0
500
1000
2000
500
1.000
1500
2000
500
1800
1500
2000
test samples
test samples
test samples
HDGM-30
HDGM-40
(g) HDGM-50
Figure 9: Power vs test size m for the HDGM problem at dimensions d= (2, 4, 10, 15, 20). The average test power computed over 5 training runs, where the empirical power is determined over 100 permutation tests. The shaded region covers standard crror from the mcan
HSI.D
MBCOx
NUS
HONCC [M)
C25TL
C25T5
NMD-D
HBIOO
H6K.
HSK
1.0
1.0
1.0
1.0
power
power
empirical power
0.8
0.6
0.6
0.6
0.6
0.4
0.4
0.4
04
0.0
10
20
30
40
50
0.0
10
20
50
40
50
0.0
10
20
30
40
50
0.0
10
20
30
40
50
dimension
dimension
dimension
dimension
(a) m 100
(b) m 200
(c) 500
(d) m = 1000
Figure 10: Empirical power dimension across various test sample sizes HDGM. The shaded region covers one standard error over 5 training runs.
{100, 200, 500, 1000} for
F.3 Maximizing HSIC vs. its SNR
We examine the trade-off between directly optimizing HSIC versus its SNR objective J. The results on power versus test size are shown in Figure 11. Optimizing our proposed objective significantly outperforms optimizing HSIC for all problems.
F.4 Independence Testing with MMD
Let -= ((X, x) (Xmm Ym)) be a test set and Sm be the permutation group of [m] with elements E Sm- Suppose we take permutations o = foi.... Jp}, with each permutation sampled uniformly from 9m. We define the action of o on test samples as oZ = ((x,, Yoa() (o2())))(m))))))) and the action of a on the empirical distribution Pxxy the empirical distribution of the permuted samples, i.e. of xx ap Σn El 8x, XY
39

---

Published in Transactions on Machine Learning Research (03/2026)
Method
HDGM-4
HDGM-8
HDGM-100
HDGM-20
HDGM-30
HDGMM40
HDGM-50
Sinusoid
RatInABox
IISIC-D
0.043
0.043
0.050
0.050
0.062
0.057
0.052
0.050
0.048
MMD-D
0.048
0.055
0.040
0.053
0.048
0.048
0.055
0.054
0.050
C2ST-L
0.060
0.030
0.053
0.048
0.053
0.058
0.045
0.046
0.048
InfoNCE
0.046
0.046
0.046
0.054
0.044
0.050
0.048
0.048
0.045
NWJ
0.030
0.054
0.058
0.052
0.044
0.064
0.054
0.052
0.042
Table 3: Average type-I error rates under the mull distribution over 400 tests. We use 512 samples.
1.0
2.0
2.0
1.0
empirical power
0.8
empirical
0.6
empirical
0.6
0.6
empirical
0.4
8.
92
0.2
0.0
500
1000
500
1000
1000
1500
1000
1500
test samples
test samples
test samples
test Samples
(a) HDGM-10
(b) HDCM-30
(c)
HDGM-50
(d) RatInABox
Figure 11: Test power using deep kernels optimized for the approximate asymptotic test power J (red) versus optimizing just the test statistic HSIC (blue).
One way construct an independence criterion is by taking the MMD between Pxy and Ixx x Py. We do this for MMD-D by using the empirical distributions PxY and of respectively, which yields
MMD? (Pxy
Σ (.E(mm)
kaglis
Σ
keglia,))) + m²pp i.j
Σ
kigla,(0).0.0))
m²
m²p
J.C mm] (E
q.rE
We first show that even with one permutation, this yields consistent estimator.
Proposition 3.2. Suppose h satisfies suP li((,, y), (x. v)) < v². The eslimator MMD, A (Z, oZ)), Z independent of a uniformly random o, satisfies with probability at least 1 that
where
17v2
M - MMD2 (1Z,, Pr PF
2
Proof. First notice that
MMD mz: P.
PPay,
x x Py) (X Y), (X', (X,, (X,, (1(7)))<<((X.... +E<(E, (X", (X",Y")) -2 222((( -X) (), (X' Y"),
P2
44
2 while the estimator (2, oZ) is given by
1
h (ax, (1), (2.,
(((a,, Vo. (I -2
h((z,, wi), (I), B),))
m
m²
Ti
Ta
Ts
The first term Ti is a typical V-statistic:
m
1 R1,
h((x,, v:), (t), (z), + m² h (((,, (x),))))
Σ
U
m
m
40

---

Published in Transactions on Machine Learning Research (03/2026)
where have defined a U-statistic U1 m) Edes (((,,,,s) (zg, and remainder term R1 1 E; (((1,, vi), (zi wi)). For U1, we have immediately that EU, = E (E(((, (X',))) =: p1. Noting that both R1 and ir are necessarily in [0, r²], the overall error from this term is therefore
1 - M1l aal= (1- m (U - -1)) + 1 ( m1) s (1-
_
1 -v2 U1 - V1l +
m
m
Changing a single (ri, Ui) pair changes U1 by at most m(m-1] 2(m 1) . v2 = 20² inequality gives that with probability at least 1 81,
and so McDiarmid's
2v m
2
2
2
UU1 M1l
log
v²
log
m
2
81
m 81
and with probability at least 1-
1d1
2
v2
2
T mul s 1
2
log
+
v²
2
(19)
m
m
81
m
mm
81
mm
Turning to T3 next, we can use similar approach if we also take into account the random o. We can write T3
N(00
(1-
Ng m2
Σ
b((zi,))) vi), (1), + m²
Σ
h((2;, v.), (1,)
m² U3
Ra,
m
where we let N( be the (random) number of (i,j) pairs for which i, aj are all distinct, Us the mean of h((x,, vi), (Xj, )) for which these indices are distinct, and Rs the mean for which they are not. Note that Ipal, RR S v², regardless of the choice of a and Z, and so
No (U3 - -1a)+ + (1- 11_1" N2 m²
(1 U1 ual + (1 0 m²
(
N:O0)
[73 pal =
(Ra 13)
S
2
m²
m²
Fix the choice of a, but let E be random. Then EU/3 μ3, and changing a single (I, y4) pair changes the value of U3 by at most 1 3(mm 1 1). v?. Thus, applying McDiarmid's incquality conditionally on the choice
of a, with probability at least 8/3
3v2
/m(m
2 8/
UU3 M3l 1
log
V2
obtaining that
3v² m log
2
1
N(
ITg Msl <
+
2v
V2 mmm
8'3
m²
We will also need to show that N/ /m2 is nearly 1. We have that
EN = E
1HHi,j, = 3) Σ Pr(a, & (i.i)) = m(m 1) m - = (m - 1)(m - 2),
so that E (1-N /m² b/ 2. Moreover, if a and o' almost agree except that of ai and of = Tke then (Nv( -
v((') - < 2m: the only (i.j) which are potentially affected are those of the form (-, k) or (-,1). Using a version of McDiarmid's inequality for uniform distributions of permutations (Lemma F.1) then gives
us that with probability at least 1 - &2 over the choice of o,
(a) 2 (m 1) - 2) - my 8m log
8/3
41

---

Published in Transactions on Machine Learning Research (03/2026)
and
N(
3
8
<
+
log
m
m
m
m
0
thus with probability at least 8 + 83, we have that
u²
3 m 1
2
6
T3 x
M3
log
+
+
2 log
V2
8/
my m
8/3
y²
3
2
1
6
2 log
+ 4
2 log
Vm
2
8
8
m
To simplify this a little further, let 8/3 = 66 and
8/²-
1Ø3; then it holds with probability at least
63 that
6v²
I3's pal -
vm V
2
5g + vm
(20)
It remains to handle T2, which is similar to T3. Letting o be the number of (i, pairs for which i, o 0j are all distinct, we can similarly define U2 with mean #2 and R2 with [Rz] S u² so that
No (U2 - 122) + (1- N/ m²
No m
No
T2 - P2l =
(R2 U2)
S
U2
H2l +
2v2.
m²
m²
For fixed 0, changing single (I,, 34) pair changes the value of U2 by at most obtain like before that with probability least
4(m 1) v22 and we
N
8/2
4y² m
2
(
12/2 p2l <
log
+
2
2v²
mym
02
m²
Since when i i j any permutation satisfies that a; # oj, we have that
EN EPr(on,o $ (i.ji) = Σ
771 2 m
(m- - 2)(m - 3).
#$
i
A single transposition changes changes N( by no more than 4m, and so by Lemmna F.1 we obtain that with probability least 8/
2 2 8n
m²
m
m
Letting 8/-- = 66 and 8" = 160, it thus holds with probability at least 1 8n that
2v²
m
6
3 + v 2 log
[2 u2l <
2 log
Vm
m
8n
vm
m/m
h2
20v2
3
1
2 log
+
(21)
/ m
82
m
To combine (19), (20), and (21), will be convenient to use d1 = 1 and 6g = 53 = 8 since then 61/2 = 8/8, 82/3 = 6/8 and on + 8g + 83 = 8. Then we obtain that with probability at least - 8,
1702 Vm
T + T2 - 2T3 - (p) + 422 - 2|²i <
log
mm
42

---

Published in Transactions on Machine Learning Research (03/2026)
Lemma F.1 (McDiarmid's inequality for uniform permutations). Let f Sm - R, where Sm is the symmetric group of permalations on [m], satisfy that for every o' E Sm; ( S (i: [i:::: # d. Let S be random variable which is uniform on Sm. Then it holds with probability at least - ó that f(S) Es f(S) S /v 8m log
Proof (based on Blasiok, 2025). For i {O....,mi, m] and any permutation a, define the Doob martingale
X; = [If(S) S = o1,.. , = o.],
so that X. =
=== Es f(S) and Xm = f(o).
We show that |X; Xi+1l 2c. No matter the known values (o1,-..; oi), for any (0i+1,---,ou) on); identify uniquely paired (d'41...., differing in at most two positions, which implies that If(o)- - f(o')l & 2c. Because these pairs are equiprobable under the uniform distribution of S, this implies that any E [X++1] 1 Si+1 = j] differs from any E [X411 SS+1 = =j'l by at most 2c, and so the same holds on average over j.
The Azuma- Hoema--oo inequality then gives Pr(Xm Xo 5) S 2xxp (-=2/ (2m(22))) and the desired result follows by solving for
When considering more than permutation, this looks reminiscent of the biased HSIC estimator. Indeed, when we consider the set of n circular shifts d. (i.e. permulations where the order of the Ys are rolations of one another), then the two are equal with
MMD xx (Pxx, MMD
RxI (Px, Px Px & H(IY),
In practice, MMD-D underperforms against HSTC-D. We believe this might be because the variance of the permuted MMD estimator often higher than that of HSIC. The following proposition proves this for the biased estimators.
Proposition F.2. Suppose we have samples ((X,, (Xm, Ym)) drawn üd from Pxy, and let a (o1....,oo Jp} be a set of p permutations sampled uniformly from Sm, the permutation group of [m]. We define the action of a samples Z oZ {(X,, estimators we have
Then for the biased HSIC and MMD
Var HISCH < Var[MMD' (2,
Proof. First we note that
MMD (Z, oZ) m2 [JE[m]
Σ
ki
Σ kigleo,))
Σ
ka.lo,,
m²p
m²p²
(.E[(m] (E(p]
(JEmm) q.rE(p)
Taking the expectation of this estimator conditioned on the samples Z gives us Rao-Blackwellization. Now notice that for part of the second torm
E
= Kigh
1 1 E
Hoo = m
kasli
where we use the fact that Eo, = = Em 4, since the permutation of j is equally likely to be any number in [m]. Similarly part of the last term yields
1)
E
(a.lo,(0)..0)) = kiis
E
E Un Lo
Σ
kiglew
p
m2
g.relpl
(,uE[m]
43

---

Published in Transactions on Machine Learning Research (03/2026)
Therefore, the Rao- Bno-BBaakwellizatton is just the biased HSIC estimator
E MMD (2, 02 I HZI
(22)
Tt. follows from the law of total variance that
Var[MMD (Z,MMDD(Z,, (Z, = Var Var[E (2, 2]] + Var
MMD (Z) [ Var[HSIC(Z]].
The result follows.
24
MMD 02)
242-04
0.85
HSICiZ)
0.80
variance
1.30-04
0.75
1:60-4
0.70
0.65
(z,
1.0e01
HSIC(2)
0.0
50 35 permutations
15
45 50
15
25
25
# permutations
variance
(b) power
Figure 12: (a) compares the sample variance of the biased MMD estimator using p shuffles of the original sample. (b) examines how the number of shuffles p affects the empirical test power. All results based 50 samples from the Simusoid problem with frequency parameter l = 1. Both MMD and HSIC use Gaussian kernels k and I with bandwidth 1.
Figure 12 shows that the higher variance of the permuted MMD estimator negatively impacts its overall test. power compared to HSIC. Also, including more permutations seems to decrease the overall variance, although never lower than that of HSIC.
In practice, the kernels learned by MMD-D and HSTC-D may not correspond, and so Proposition F.2 is inapplicable. That said, we observe the same phenomenon in experiments. Figure 13 shows estimates of the asymptotic variance of MMD and HSIC along a training trajectory. For permuted MMD we consider both a single shuffling of the data (MlD-full), as well as a split shuffling (MMD-split) where we use half the data for our joint distribution sample, and the other half to permute for product-of-marinnal sample. We note that the initial variance of MMD-split is substantially higher than that of MMD-full, which is much higher than the variance of HSIC, MMD-split/fun /full also exhibit greater final variances, particularly at larger batch sizes.
mmd
mmsifal mms split
mmn--g
It
Variance
variaace
Seriance
epochs
epochs
epochs
(a) 128
(b) 512
1024
Figure 13: Estimates of the asymptotic variance of HSIC (blue), MMD with single permutation (orange), and MMD with split permutation (green) along training trajectory for HDGM-10 at sample sizes 128 (a), m 512 (b), and m = 1024 (c).
44

---

Published in Transactions on Machine Learning Research (03/2026)
F.5 SNR Pitfall
Recent work by Ren et al. (2024) argues so-called pitfall of the HSIC signal-to-noise ratio paradigm. They identify a corner case whereby, when the bandwidth of one of kernel & or approaches 0, ignoring the threshold term causes the SNR objective Jujo HSIC,/ / to differ from the true asymptotic power objective J HSIC, (rsic,rrmm) We argue that this is not a cause for concern, and in some cases, ignoring this is even preferred.
r/m)
/ by a factor of (m 1), where r denotes the asymptotic threshold.
First, the behavior of the SNR objective in the bandwidth limit tells us nothing about the actual global maximum. Their argument tells us that Jwjo = (mm - 1) when the bandwidth wx - and for a specific estimator = E[MHSTC.). This does not imply, however, that Juyo explodes as wx - since the SNR objective is optimized
fixed sample size m, and even in the bandwidth limit Ju/o can not be -cc. The latter observation is true because HSIC, is bounded and our variance estimates 0 are also bounded away from zero. Therefore even though Jwjo is high when wx Y 0, this value is not necessarily the global maximum. Second, their argument seems entirely dependent on the choice of estimators. For instance, if we use HSIC,, rather than HSIC,, and take the asymptotic threshold to be E|mHSIC., m HSIC(P, x PPy)) - O, get that Ju/o J in this regime.
As further evidence, we plot both Jw/ and Jwjo at varying bandwidths wx on the ISA dataset used in Ren el al. (2024). We use the same settings as they do: m = 250, d = 3, 0 = x 7/10, wy = 1.0. Results are shown in Figure 14. For very small bandwidths, Jujo does not explode and is not the global maximum. Moreover, the actual maximum agrees relatively well between Ju/ and Jw/o at larger sample sizes m. We were unable to reproduce their Figure 1.
Additionally, Figure 14 (a) shows that ignoring this threshold term may actually be preferred. Consider the limiting behavior the bandwidth wx goes to infinity. In this regime, the Gram and centered Gram matrices are respectively K wwo 11T and Kolw=xo O, which tells us that HSIC, Tr[K,L]/m 0. Since mHSIC& is degenerate its variance and 1 - quantile are zero, and so both Jus and Jw/o are also 0 in this bandwidth limit. 13 Now consider a difficult problem where the power at initialization is less than 0.5, meaning Juy is less than 0. Then the maximum of Jw/ erroneously this bandwidth limit.
10.0
SNR objective
SNR Bejeciive
20
5.0
lffrestredd
SNR
SNR
2.5
0.0
-13
-3
-15 -10
-10
-15 -10
Log-bandwidth (log mx)
Log bag-bandwieth (log mx))
Log bandwidth (lag 6xy)
Log-bandwidth (lag (x)
(a) 250
(b)
500
750
1000
Figure 14: Plots of Jwjo (blue) and Ju/ (orange) very small bandwidths wx with sample size m. For the threshold estimate, use the -5-quantile of the Gamma approximation. The objective maximum is indicated by stars.
F.6 SNR with Threshold
Although maximizing the HSIC signal-to noine objective JUSIC is sufficient for learning powerful tests, it may also be beneficial to estimate the threshold-dependentt term in the asymptotic power expression of Equation (6). Adopting the same notation as before, we define the SNR with threshold expression as
JuSIC ( X, (, HSIC(X,
V v a1)
σs (X, I; mos, (XYY
Recall that our variance estimate adds small positive constant for stability, preventing J from being undefined.
45

---

Published in Transactions on Machine Learning Research (03/2026)
where V is the problem-dependent edf described in Proposition A.2. Estimating JUSIC proves troublesome, largely due to Wp-1 depending on both the variables X, Y well as kernels k,1. Additionally, it needs to be differentiated with respect to those kernel parameters. We opt to approximate the null distribution & with a Gamma distribution FN. with shape and scale parameters
E Var [HSCC]
0 Var HSIC&]
and
[[SIC,
as suggested by Gretton et al. (2007). When performing gradient updates, we need to compute the gradients V.ryo (1 - a) and Vo a). This is typically not possible with auto-differentiatio,, as the inverse odf is not tractable. Instead, we propose a solution based on implicit differentiation: suppose Su, the Gamma pdf and FyJ a). Then
V.Fi (1- - aa) VuFuo(r) and
VOlu - aa) voFv,e())
Sua(r)
Su,e(r)
The gradient of the Gamma odf with respect to the scale 0 yields the simplified expression
1 T.T(v) (g) e-
VoFue(r)
The gradient (r) is more troublesome. We resort to approximating this gradient via series expansions, detailed in Moore (1982).
We compare the performance of HSIC and NDS tests maximized with and without the threshold-dependeess term in Figure 15. For HSIC, including this term scems to hurt the overall test power for harder problems. This is expected, since the null distribution is difficult to accurately estimate. We use two approxinations here: one being the Gamma approximation, and two being the series expansion used to estimate the derivative of the Gamma edf with respect to the shape parameter. On the other hand, NDS with the threshold term yields higher test power. Again, this is not too surprising as, unlike for HSIC, the mull distribution here follows the CLT and thus is easily approximated given enough samples.
power
nds w//
power
power
hsic w!
W-W
hsic-w/
nds
hsicww/o
test
best size
test size
test size
(a) NDS HDGM-4
(b) NDS & HDGM-10
(c) HSIC @ HDGM-4
(d) HSIC & HDGM-10
Figure 15: Power vs test size based on NDS and HSIC signal-to-noist objectives, both with (blue) and without (red) estimating the lhreshold-depeenden term. Our proposed methods NDS (12) and HSIC-D (13) correspond to nds-w/ and hsic-w/o, respectively. We consider both HDGM-8 and HDGM-10, and use the same training configuration as specified in Section
F.7 Why do Permutation Tests deviate from Asymptotics?
As Figure 6 suggests, the asymptotics (i.c. tests based on the asymptotic null and alternative distributions) may not explain permutation tests well. Figure 16 illustrates why this might be: we observe strong dependency between the sample statistic f(x, Y) and its rejection threshold, estimated via permutation.
For the InfoNCE critic this dependency is substantial, with a correlation between the two of 98%. While permutation tests still ensure the appropriate level, this coupling between the test statistic and permutation threshold is not accounted for in asymptotic analysis. Tf these tests do follow their asymptotics, then
46

---

Published in Transactions on Machine Learning Research (03/2026)
expect this correlation to be relatively small, since the threshold should always be reasonable estimate of the 1 - a quantile. This is certainly not true for InfoNCE, and explains why its permutation test power is so much higher than what its asymptotics suggest. As for why or how InfoNCE critics exhibit such strong coupling, we do not yet know; this very interesting area for future work.
power: 46%
power: 63%
correlation: 37%
correlation: 98%
Te% m)
i v
permutation threshold
permutation threshold
(a) NDS critic
(b) InfoNCE critic
Plot of the NDS statistic î evaluated on a sample (x, Y) from RatInABox, versus its rejection threshold estimated with 200 permutations. Each point corresponds a separate test set (X, Y) of size 1000, but is evaluated with the same critic function f. The results in (a) use critic that maximizes the NDS test power, in Section 5. The critic in (b) directly maximizes InfoNCE. The red line is y 2, which is the rejection boundary; points above this line (colored yellow) are tests that reject the mull hepothesis, while points below (purple) do not.
Figure 16:
47


Reviews:
[{"review_id": "F0bS8zdzqj", "reviewer": "Reviewer_7pdE", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths:\n- The method certainly appears to me correct, is well motivated, and well explained. \n- The switch to optimising the asymptotic power statistics is an interesting idea and it appears to lead to good results. \n- The experimental results are extensive and convincing\n- The paper includes a thorough theoretical analysis of their methods. \n\nWeaknesses:\n- There is little discussion of hyperparmeter optimisation for the learning algorithm which seems like it would require some thought for the best way to incorporate into the framework. For example, this could be done by using a validation set but then there is the question of if this should use a permutation test of the asymptotic statistic for validation. Alternatively some option of doing a hyperparameter sweep, permuting on the test set, and then doing multiple test corrections could be used. This is just to say the question of how to do hyperparameter optimisation seems an interesting and important one in this case and to the best of my knowledge it is not discussed.\n- Based on the above, the authors do discuss the use of the validation set but it is not totally clear how it is used. I would appreciate more discussion on this. \n- Not a weakness per say, but I think the method could potentially benefit from some form of sample splitting of the type more common in doubly robust estimation of causal effects, or k fold validation, and it would be nice to see some discussion of this. For example instead of having to split into one train/test set, the authors could use k different folds where the each fold consists of a different train/test split to train a model on and the final test statistic is calculated by evaluating each of the k functions on their test set. This would have the advantage of improving effective sample size.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper explores the problem of learning a witness function to test for dependence between two random variables. The paper aims to learn a function that provides the most powerful permutation test for dependence between two variables. In order to learn this function, they optimise asymptotic power expressions in the training set, which used as an optimisable proxy for the power of the final test. They use two forms of this test, one where the test statistic is the mean of the learned function (they term this the neural dependency statistic) and one using HSIC with deep kernels. They experimentally demonstrate strong performance of their method compared to alternatives.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "Based on the above, I would ask for the following:\n- Can the authors include a discussion or provide their thoughts on hyperparameter optimisation in their setting? \n- Can the authors add discussion on sample splitting or K-fold validation.", "broader_impact_concerns": "None", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "2imMWD4FyC", "reviewer": "Reviewer_Hm98", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths:\n1. Apart from the contributions mentioned above, the authors also provide a clear analysis of the results and comparison+background  with existing methods.\n2. The paper is well written and the contributions are clear.\n3. The authors highlight the difference between maximizing a statistic value vs maximizing its statistical power esp for MI based bounds.\n4. HSIC is shown to be a lower bound on mutual information.\n5. The authors justify the ``learning representations'' approach through Uniform Convergence arguments. They show that (under regularity assumptions)  the representations optimized by Algorithm 1 are maximizes a population SNR objective within a restricted class of kernels.\n\n\nWeaknesses:\n1. The main idea of optimizing deep kernels by maximizing the test power (SNR) rather than the statistic value has been studied by several papers including Liu et al. (2020) and others. Although this has been acknowledged in the paper, the main proof sections (Appendix B and C) are also acknowledged to be based on Liu et al. (2020). It would be more helpful for readers if the authors could provide a more detailed comparison of proof techniques with Liu et al. (2020) and may be with other related papers.  \n2. In section 3.1, the optimization objective is based on maximizing the SNR (HSIC/sigma), which is derived from a CLT-based Gaussian approximation. However, this approximation is least accurate when the signal is weak. I am curious if maximizing a Gaussian-based SNR is \"right thing to do\" when the true distribution of the statistic is closer to a weighted sum of Chi-squares (as mentioned in Prop. A.2). In this regime the learned representations are even more valuable.\n3. (Minor) The tests in this paper are based on data splitting between training and test sets (which might reduce the power for a fixed set of samples). Although I believe the analysis will become more complex if we try to remove the data splitting.\n4. (Minor) All learned methods are evaluated using sample splitting between representation learning and evaluating the test statistic. While this is acceptable from theoretical point of view, it might be informative to include a diagnostic experiment showing what happens if the same data are used for both training and testing, to better understand the practical impact of sample splitting and potential overfitting effects (such as seeing the percentage of times the outputs in the two cases agree).", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "The paper makes a number of contributions in the area of independence testing via variational mutual information (MI) bounds and Hilbert-Schmidt Independence Criterion (HSIC) bounds. The focus is mainly on statistical tests of independence for high-dimensional continuous variables where classical tests like chi-square or parametric assumptions are not applicable. The kernel based methods might not perform well due to poor representation of the data. The authors propose a method to optimize the representation of the data to make the dependency more explicit.  In short the key contributions are:\n\n1. The authors introduce a new testing technique based on variational MI estimators (e.g., InfoNCE, NWJ). Previously, only estimators based on variational bounds were of main focus. \n2. Establishes a connection between variational MI-based tests and the HSIC, showing that HSIC serves as a lower bound on MI.\n3. Proposes an algorithmic framework to learn optimal data representations by maximizing the Signal-to-Noise Ratio (SNR) of the test statistic, which directly increases the power of the test.\n4. Provides uniform convergence proofs (Theorems 5.1 and 5.2) for the SNR estimators, guaranteeing that representations optimized on finite samples reliably generalize to population-level test power.\n5. Corrections to Literature such as identifying and correcting mathematical errors in prior work (Ren et al., 2024) and expands the use of deep kernels in HSIC optimization.\n6. Theoretical results are supported by experiments on a range of synthetic and real world datasets. Empirical results suggest that HSIC-D dominates MI-based tests.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "See weaknesses 1 (critical) and 4 (minor).", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "9oRXqtM73X", "reviewer": "Reviewer_iko6", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "## Strengths\n- The focus on optimizing the signal-to-noise ratio rather than the raw statistic is well motivated and leads to more meaningful gains in test power.\n- The connection between MI-based statistics and HSIC is clearly articulated and helps unify two commonly used approaches in independence testing.\n- The empirical evaluation is thorough and shows that representation learning can significantly improve performance in high-dimensional settings.\n- The uniform convergence results provide useful theoretical support for the proposed objectives.\n\n## Weaknesses\nSeveral of the points discussed below are connected to specific suggestions in the Requested Changes section.\n- The practical impact of the MI–HSIC connection on test power remains unclear, as the bounds used offer limited interpretability.\n- Some assumptions in the theoretical analysis, including the structure of the parameter space and the choice of $\\lambda=\\Theta(n^{-1/3})$, are not fully justified.\n- Several arguments rely on CLTs without sufficiently explaining how they support the tail-probability approximations used in testing.\n- The assumption that certain population quantities are uniformly bounded away from zero may be restrictive in practice.\n- Experiments focus mostly on low-dimensional tabular data, leaving open how well the method scales to more complex modalities where deep architectures are most beneficial.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This work studies representation learning for independence testing using variational MI estimators and HSIC-based kernel methods. It analyzes their connections, proposes a framework that optimizes an estimate of asymptotic test power through the signal-to-noise ratio, and demonstrates its effectiveness through theoretical results and experiments.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "I do not separate major and minor concerns in this review, and all revision suggestions are grouped together in the Requested Changes section.\n\n---\n### References\n\n[1] Kim, Ilmun, and Aaditya Ramdas. “Dimension-agnostic inference using cross U-statistics.” Bernoulli, 2024.\n\n[2] Shekhar, Shubhanshu, et al. “A permutation-free kernel two-sample test.” Advances in Neural Information Processing Systems, 2022.\n\n[3] Gretton, Arthur, et al. “A kernel two-sample test.” Journal of Machine Learning Research, 2012.\n\n[4] Cai, Zhanrui, et al. “Asymptotic distribution-free independence test for high-dimension data.” Journal of the American Statistical Association, 2024.\n\n---\n## Requested Changes\n1. On page 4, the expression\n$$\n\\frac{1}{\\tau_{\\mathfrak{H}1}} \\sqrt{m}\\big(\\hat{T}-T_{\\mathfrak{H}_1}\\big) \\xrightarrow{d} \\mathcal{N}(0,1)\n$$\n\nrequires that the function $f$ used in $\\hat{T}$ does not depend on the evaluation samples $(Z_1, \\ldots, Z_m)$. A brief explanation regarding the need for sample splitting, or at least an explicit acknowledgment of this issue, would improve clarity.\n\nEmphasizing, as discussed in [1], that sample splitting is needed to remove the dependency between $f$ and the evaluation samples would help clarify why the stated asymptotic normality holds.\n\n2. $$\n\\frac{\\sqrt{m} \\hat{T}}{\\tau_{\\Im_1}}>r_m, \\quad m \\widehat{\\mathrm{HSIC}}_u>r_m\n$$\nwould align better with conventional forms of test statistics based on CLTs.\n\n3. The text should explicitly clarify that the quantities $a=a_m$ and $b=b_m$ depend on $m$, and that the relation $a \\sim b$ refers to the limit as $m \\to \\infty$. Since this passage is intended to define the notation $a \\sim b$, it would be helpful to state these points directly.\n\n4. Clarify why CLT-based approximations justify ratio convergence of tail probabilities. The current argument does not fully justify this; further explanation is needed, especially in the context of Equation (3).\n\n5. Clarify why the CLT is needed for discussing testing power. It would be helpful for the paper to explain why asymptotic normality is required to characterize power consistency, rather than directly showing that power converges to one as $(m \\to \\infty)$. The Proof of Theorem 8 in [2] may help illustrate how asymptotic normality contributes to such arguments.\n\n6. When constructing samples from $P_X P_Y$ via permutation, it is unclear whether this yields i.i.d. samples from $P_X P_Y$. Please clarify how this differs from standard permutation HSIC tests and under what assumptions the MMD-based test is valid.\n\n7. $\\widehat{\\operatorname{MMD}}^2_b$ first appears in Proposition 3.2, but its definition is provided only later in the text. Introducing the definition before this proposition would improve readability.\n\n8. Clarifying with respect to which distribution each expectation is taken in the definition of $\\mathrm{MMD}^2$  introduced above equation (8) would make the presentation clearer.\n\n9. Explain the role of the concentration bound in constructing valid tests. The paper provides a bound but does not tie it directly to hypothesis testing. Additional explanation, possibly similar to Corollary 9 in [3] would help the reader understand how to form a rejection region.\n\n10. Clarify the role of the parameter space $\\Theta$ in the theoretical analysis. It is unclear why this specific space must be considered and why uniform convergence over this space is required. Also, the relation between parameter-space constraints and conditions such as $\\lambda=\\Theta(n^{-1 / 3})$ should be discussed.\n\n11. Discuss the realism of the assumption, appearing in Theorem 5.1 and Theorem 5.2, that the population quantity governing the variance of the statistic is uniformly bounded away from zero. This condition may be too strong in practice. Commenting on possible relaxations or alternative assumptions would strengthen the theoretical contribution.\n\n12. Expand the empirical evaluation to more challenging representation-learning settings. Since the method is motivated by the flexibility of deep networks, evaluating it on higher-dimensional modalities (such as images or text) would strengthen the empirical contribution. In addition, incorporating comparisons in high-dimensional tabular settings, as explored in [4], would provide a useful benchmark for understanding the method’s performance in more classical data regimes.\n\n13. Clarify whether dimensions of $X$ and $Y$ must match. It would help to explicitly state whether the method requires equal-dimensional inputs or whether the framework supports mismatched dimensionalities.", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Decision: Accept (poster)