Nonparametric independence testing via mutual

arXiv:1711.06642v1 [stat.ME] 17 Nov 2017

information
Thomas B. Berrett and Richard J. Samworth
Statistical Laboratory, University of Cambridge
November 20, 2017

Abstract
We propose a test of independence of two multivariate random vectors, given a
sample from the underlying population. Our approach, which we call MINT, is based
on the estimation of mutual information, whose decomposition into joint and marginal
entropies facilitates the use of recently-developed efficient entropy estimators derived
from nearest neighbour distances. The proposed critical values, which may be obtained
from simulation (in the case where one marginal is known) or resampling, guarantee
that the test has nominal size, and we provide local power analyses, uniformly over
classes of densities whose mutual information satisfies a lower bound. Our ideas may
be extended to provide a new goodness-of-fit tests of normal linear models based on
assessing the independence of our vector of covariates and an appropriately-defined
notion of an error vector. The theory is supported by numerical studies on both
simulated and real data.

1

Introduction

Independence is a fundamental concept in statistics and many related fields, underpinning the
way practitioners frequently think about model building, as well as much of statistical theory.
Often we would like to assess whether or not the assumption of independence is reasonable,
for instance as a method of exploratory data analysis (Steuer et al., 2002; Albert et al., 2015;
Nguyen and Eisenstein, 2017), or as a way of evaluating the goodness-of-fit of a statistical
model (Einmahl and van Keilegom, 2008). Testing independence and estimating dependence
1

are well-established areas of statistics, with the related idea of the correlation between two
random variables dating back to Francis Galton’s work at the end of the 19th century (Stigler,
1989), which was subsequently expanded upon by Karl Pearson (e.g. Pearson, 1920). Since
then many new measures of dependence have been developed and studied, each with its own
advantages and disadvantages, and there is no universally accepted measure. For surveys of
several measures, see, for example, Schweizer (1981), Joe (1989), Mari and Kotz (2001) and
the references therein. We give an overview of more recently-introduced quantities below;
see also Josse and Holmes (2016).
In addition to the applications mentioned above, dependence measures play an important
role in independent component analysis (ICA), a special case of blind source separation, in
which a linear transformation of the data is sought so that the transformed data is maximally
independent; see e.g. Comon (1994), Bach and Jordan (2002), Miller and Fisher (2003) and
Samworth and Yuan (2012). Here, independence tests may be carried out to check the
convergence of an ICA algorithm and to validate the results (e.g. Wu, Yu and Li, 2009).
Further examples include feature selection, where one seeks a set of features which contains
the maximum possible information about a response (Torkkola, 2003; Song et al., 2012),
and the evaluation of the quality of a clustering in cluster analysis (Vinh, Epps and Bailey,
2010).
When dealing with discrete data, often presented in a contingency table, the independence
testing problem is typically reduced to testing the equality of two discrete distributions via
a chi-squared test. Here we will focus on the case of distributions that are absolutely continuous with respect to the relevant Lebesgue measure. Classical nonparametric approaches to
measuring dependence and independence testing in such cases include Pearson’s correlation
coefficient, Kendall’s tau and Spearman’s rank correlation coefficient. Though these approaches are widely used, they suffer from a lack of power against many alternatives; indeed
Pearson’s correlation only measures linear relationships between variables, while Kendall’s
tau and Spearman’s rank correlation coefficient measure monotonic relationships. Thus, for
example, if X has a symmetric distribution on the real line, and Y = X 2 , then the population quantities corresponding to these test statistics are zero in all three cases. Hoeffding’s
test of independence (Hoeffding, 1948) is able to detect a wider class of departures from
independence and is distribution-free under the null hypothesis but, as with these other
classical methods, was only designed for univariate variables. Recent work of Weihs, Drton
and Meinshausen (2017) has aimed to address some of computational challenges involved in
extending these ideas to multivariate settings.
2

Recent research has focused on constructing tests that can be used for more complex
data and that are consistent against wider classes of alternatives. The concept of distance
covariance was introduced in Székely, Rizzo and Bakirov (2007) and can be expressed as
a weighted L2 norm between the characteristic function of the joint distribution and the
product of the marginal characteristic functions. This concept has also been studied in
high dimensions (Székely and Rizzo, 2013; Yao, Zhang and Shao, 2017), and for testing
independence of several random vectors (Fan et al., 2017). In Sejdinovic et al. (2013) tests
based on distance covariance were shown to be equivalent to a reproducing kernel Hilbert
space (RKHS) test for a specific choice of kernel. RKHS tests have been widely studied
in the machine learning community, with early understanding of the subject given by Bach
and Jordan (2002) and Gretton et al. (2005), in which the Hilbert–Schmidt independence
criterion was proposed. These tests are based on embedding the joint distribution and
product of the marginal distributions into a Hilbert space and considering the norm of their
difference in this space. One drawback of the kernel paradigm here is the computational
complexity, though Jitkrittum, Szabó and Gretton (2016) and Zhang et al. (2017) have
recently attempted to address this issue. The performance of these methods may also be
strongly affected by the choice of kernel. In another line of work, there is a large literature
on testing independence based on an empirical copula process; see for example Kojadinovic
and Holmes (2009) and the references therein. Other test statistics include those based on
partitioning the sample space (e.g. Gretton and Györfi, 2010; Heller et al., 2016). These have
the advantage of being distribution-free under the null hypothesis, though their performance
depends on the particular partition chosen.
We also remark that the basic independence testing problem has spawned many variants.
For instance, Pfister et al. (2017) extend kernel tests to tests of mutual independence between a group of random vectors. Another important extension is to the problem of testing
conditional independence, which is central to graphical modelling (Lauritzen, 1996) and also
relevant to causal inference (Pearl, 2009). Existing tests of conditional independence include
the proposals of Su and White (2008), Zhang et al. (2011) and Fan, Feng and Xia (2017).
To formalise the problem we consider, let X and Y have (Lebesgue) densities fX on RdX
and fY on RdY respectively, and let Z = (X, Y ) have density f on Rd , where d := dX + dY .
Given independent and identically distributed copies Z1 , . . . , Zn of Z, we wish to test the null
hypothesis that X and Y are independent, denoted H0 : X ⊥⊥ Y , against the alternative
that X and Y are not independent, written H1 : X 6⊥⊥ Y . Our approach is based on
constructing an estimator Ibn = Ibn (Z1 , . . . , Zn ) of the mutual information I(X; Y ) between
3

X and Y . Mutual information turns out to be a very attractive measure of dependence in this
context; we review its definition and basic properties in Section 2 below. In particular, its
decomposition into joint and marginal entropies (see (1) below) facilitates the use of recentlydeveloped efficient entropy estimators derived from nearest neighbour distances (Berrett,
Samworth and Yuan, 2017).
The next main challenge is to identify an appropriate critical value for the test. In
the simpler setting where either of the marginals fX and fY is known, a simulation-based
approach can be employed to yield a test of any given size q ∈ (0, 1). We further provide
regularity conditions under which the power of our test converges to 1 as n → ∞, uniformly
over classes of alternatives with I(X; Y ) ≥ bn , say, where we may even take bn = o(n−1/2 ).
To the best of our knowledge this is the first time that such a local power analysis has
been carried out for an independence test for multivariate data. When neither marginal is
known, we obtain our critical value via a permutation approach, again yielding a test of the
nominal size. Here, under our regularity conditions, the test is uniformly consistent (has
power converging to 1) against alternatives whose mutual information is bounded away from
zero. We call our test MINT, short for Mutual Information Test; it is implemented in the R
package IndepTest (Berrett, Grose and Samworth, 2017).
As an application of these ideas, we are able to introduce new goodness-of-fit tests of
normal linear models based on assessing the independence of our vector of covariates and
an appropriately-defined notion of an error vector. Such tests do not follow immediately
from our earlier work, because we do not observe realisations of the error vector directly;
instead, we only have access to residuals from a fitted model. Nevertheless, we are able to
provide rigorous justification, again in the form of a local analysis, for our approach. It seems
that, when fitting normal linear models, current standard practice in the applied statistics
community for assessing goodness-of-fit is based on visual inspection of diagnostic plots such
as those provided by the plot command in R when applied to an object of class lm. Our aim,
then, is to augment the standard toolkit by providing a formal basis for inference regarding
the validity of the model. Related work here includes Neumeyer (2009), Neumeyer and Van
Keilegom (2010), Müller, Schick and Wefelmeyer (2012), Sen and Sen (2014) and Shah and
Bühlmann (2017).
The remainder of the paper is organised as follows: after reviewing the concept of mutual
information in Section 2.1, we explain in Section 2.2 how it can be estimated effectively
from data using efficient entropy estimators. Our new tests, for the cases where one of the
marginals is known and where they are both unknown, are introduced in Sections 3 and 4
4

respectively. The regression setting is considered in Section 5, and numerical studies on both
simulated and real data are presented in Section 6. Proofs are given in Section 7.
The following notation is used throughout. For a generic dimension D ∈ N, let λD and
k · k denote Lebesgue measure and the Euclidean norm on RD respectively. If f = dµ/dλD
and g = dν/dλD are densities on RD with respect to λD , we write f  g if µ  ν. For
z ∈ RD and r ∈ [0, ∞), we write Bz (r) := {w ∈ RD : kw − zk ≤ r} and Bz◦ (r) := Bz (r) \ {z}.
We write λmin (A) for the smallest eigenvalue of a positive definite matrix A, and kBkF for
the Frobenius norm of a matrix B.

2
2.1

Mutual information
Definition and basic properties

Retaining our notation from the introduction, let Z := {(x, y) : f (x, y) > 0}. A very natural
measure of dependence is the mutual information between X and Y , defined to be
Z
f (x, y) log

I(X; Y ) = I(f ) :=
Z

f (x, y)
dλd (x, y),
fX (x)fY (y)

when f  fX fY , and defined to be ∞ otherwise. This is the Kullback–Leibler divergence
between the joint distribution of (X, Y ) and the product of the marginal distributions, so is
non-negative, and equal to zero if and only if X and Y are independent. Another attractive
feature of mutual information as a measure of dependence is that it is invariant to smooth,
invertible transformations of X and Y . Indeed, suppose that X := {x : fX (x) > 0} is an
open subset of RdX and φ : X → X 0 is a continuously differentiable bijection whose inverse
φ−1 has Jacobian determinant J(x0 ) that never vanishes on X 0 . Let Φ : X × Y → X 0 × Y be

given by Φ(x, y) := φ(x), y , so that Φ−1 (x0 , y) = φ−1 (x0 ), y) also has Jacobian determinant

J(x0 ). Then Φ(X, Y ) has density g(x0 , y) = f φ−1 (x0 ), y |J(x0 )| on Φ(Z) and X 0 = φ(X)

has density gX 0 (x0 ) = fX φ−1 (x0 ) |J(x0 )| on X 0 . It follows that
g(x0 , y)
dλd (x0 , y)
I φ(X); Y =
g(x , y) log
0
gX 0 (x )fY (y)
Φ(Z)

Z
f φ−1 (x0 ), y
0

=
g(x , y) log
dλd (x0 , y) = I(X; Y ).
−1
0
fX φ (x ) fY (y)
Φ(Z)


Z

0

5

This means that mutual information is nonparametric in the sense of Weihs, Drton and
Meinshausen (2017), whereas several other measures of dependence, including distance covariance, RKHS measures and correlation-based notions are not in general. Under a mild
assumption, the mutual information between X and Y can be expressed in terms of their
joint and marginal entropies; more precisely, writing Y := {y : fY (y) > 0}, and provided
that each of H(X, Y ), H(X) and H(Y ) are finite,
Z

Z

Z

f log f dλd −

I(X; Y ) =
Z

fX log fX dµdX −
X

fY log fY dµdY
Y

=: −H(X, Y ) + H(X) + H(Y ).

(1)

Thus, mutual information estimators can be constructed from entropy estimators.
Moreover, the concept of mutual information is easily generalised to more complex situations. For instance, suppose now that (X, Y, W ) has joint density f on Rd+dW , and let
f(X,Y )|W (·|w), fX|W (·|w) and fY |W (·|w) denote the joint conditional density of (X, Y ) given
W = w and the conditional densities of X given W = w and Y given W = w respectively.
When f(X,Y )|W (·, ·|w)  fX|W (·|w)fY |W (·|w) for each w in the support of W , the conditional
mutual information between X and Y given W is defined as
Z
I(X; Y |W ) :=

f (x, y, w) log
W

f(X,Y )|W (x, y|w)
dλd+dW (x, y, w),
fX|W (x|w)fY |W (y|w)

where W := {(x, y, w) : f (x, y, w) > 0}. This can similarly be written as
I(X; Y |W ) = H(X, W ) + H(Y, W ) − H(X, Y, W ) − H(W ),
provided each of the summands is finite.
Finally, we mention that the concept of mutual information easily generalises to situations
with p random vectors. In particular, suppose that X1 , . . . , Xp have joint density f on Rd ,
where d = d1 +. . .+dp and that Xj has marginal density fj on Rdj . Then, when f  f1 . . . fp ,
and writing Xp := {(x1 , . . . , xp ) ∈ Rd : f (x1 , . . . , xp ) > 0}, we can define
Z
I(X1 ; . . . ; Xp ) :=

f (x1 , . . . , xp ) log
Xp

=

p
X

f (x1 , . . . , xp )
dλd (x1 , . . . , xp )
f1 (x1 ) . . . fp (xp )

H(Xj ) − H(X1 , . . . , Xp ),

j=1

6

with the second equality holding provided that each of the entropies is finite. The tests
we introduce in Sections 3 and 4 therefore extend in a straightforward manner to tests of
independence of several random vectors.

2.2

Estimation of mutual information

For i = 1, . . . , n, let Z(1),i , . . . , Z(n−1),i denote a permutation of {Z1 , . . . , Zn } \ {Zi } such that
kZ(1),i − Zi k ≤ . . . ≤ kZ(n−1),i − Zi k. For conciseness, we let
ρ(k),i := kZ(k),i − Zi k
denote the distance between Zi and the kth nearest neighbour of Zi . To estimate the
unknown entropies, we will use a weighted version of the Kozachenko–Leonenko estimator
(Kozachenko and Leonenko, 1987). For k = k Z ∈ {1, . . . , n − 1} and weights w1Z , . . . , wkZ
P
satisfying kj=1 wjZ = 1, this is defined as

 d
n
k
ρ(j),i Vd (n − 1)
1 XX Z
d,wZ
Z
b
b
Hn = Hn,k (Z1 , . . . , Zn ) :=
w log
,
n i=1 j=1 j
eΨ(j)
where Vd := π d/2 /Γ(1 + d/2) denotes the volume of the unit d-dimensional Euclidean ball
and where Ψ denotes the digamma function. Berrett, Samworth and Yuan (2017) provided
b nZ
conditions on k, w1Z , . . . , wZ and the underlying data generating mechanism under which H
k

is an efficient estimator of H(Z) (in the sense that its asymptotic normalised squared error
risk achieves the local asymptotic minimax lower bound) in arbitrary dimensions. With
b nX and H
b nY of H(X) and H(Y ) defined analogously as functions of (X1 , . . . , Xn )
estimators H
and (Y1 , . . . , Yn ) respectively, we can use (1) to define an estimator of mutual information by
bX + H
bY − H
bZ.
Ibn = Ibn (Z1 , . . . , Zn ) = H
n
n
n

(2)

Having identified an appropriate mutual information estimator, we turn our attention in the
next two sections to obtaining appropriate critical values for our independence tests.

7

3

The case of one known marginal distribution

In this section, we consider the case where at least one of fX and fY in known (in our experience, little is gained by knowledge of the second marginal density), and without loss of
generality, we take the known marginal to be fY . We further assume that we can generate in(b)

dependent and identically distributed copies of Y , denoted {Yi

: i = 1, . . . , n, b = 1, . . . , B},

independently of Z1 , . . . , Zn . Our test in this setting, which we refer to as MINTknown (or
MINTknown(q) when the nominal size q ∈ (0, 1) needs to be made explicit), will reject H0
for large values of Ibn . The ideal critical value, if both marginal densities were known, would
therefore be
Cq(n) := inf{r ∈ R : PfX fY (Ibn > r) ≤ q}.
(b)

Using our pseudo-data {Yi

: i = 1, . . . , n, b = 1, . . . , B}, generated as described above, we

define the statistics

(b)
Ibn(b) = Ibn (X1 , Y1 ), . . . , (Xn , Yn(b) )
for b = 1, . . . , B. Motivated by the fact that these statistics have the same distribution as
(n)
Ibn under H0 , we can estimate the critical value Cq by

(n),B
b
Cq
:= inf r ∈ R :


B
1 X
1 (b) ≤ q ,
B + 1 b=0 {Ibn ≥r}

(0)
(B)
(0)
the (1 − q)th quantile of {Ibn , . . . , Ibn }, where Ibn := Ibn . The following lemma justifies this

critical value estimate.
Lemma 1. For any q ∈ (0, 1) and B ∈ N, the MINTknown(q) test that rejects H0 if and only
bq(n),B has size at most q, in the sense that
if Ibn > C
sup


bq(n),B ≤ q,
P Ibn > C

sup

k∈{1,...,n−1} (X,Y ):I(X;Y )=0

where the inner supremum is over all joint distributions of pairs (X, Y ) with I(X; Y ) = 0.
An interesting feature of MINTknown, which is apparent from the proof of Lemma 1, is
b X in (1), either on the original data, or on the pseudothat there is no need to calculate H
n
(b)
(b)
data sets {(X1 , Y1 ), . . . , (Xn , Yn ) : b = 1, . . . , B}. This is because in the
(b)
b X appears on both sides of
of the event {Ibn ≥ Ibn } into entropy estimates, H
n

decomposition
the inequality,

so it cancels. This observation somewhat simplifies our assumptions and analysis, as well as
8

reducing the number of tuning parameters that need to be chosen.
The remainder of this section is devoted to a rigorous study of the power of MINTknown
that is compatible with a sequence of local alternatives (fn ) having mutual information
In → 0. To this end, we first define the classes of alternatives that we consider; these
parallel the classes introduced by Berrett, Samworth and Yuan (2017) in the context of
entropy estimation. Let Fd denote the class of all density functions with respect to Lebesgue
measure on Rd . For f ∈ Fd and α > 0, let
Z

kzkα f (z) dz.

µα (f ) :=
Rd

Now let A denote the class of decreasing functions a : (0, ∞) → [1, ∞) satisfying a(δ) =
o(δ − ) as δ & 0, for every  > 0. If a ∈ A, β > 0, f ∈ Fd is m := (dβe − 1)-times
differentiable and z ∈ Z, we define ra (z) := {8d1/2 a(f (z))}−1/(β∧1) and
kf (t) (z)k
,
Mf,a,β (z) := max max
t=1,...,m
f (z)



kf (m) (w) − f (m) (z)k
sup
.
f (z)kw − zkβ−m
w∈Bz◦ (ra (z))

The quantity Mf,a,β (z) measures the smoothness of derivatives of f in neighbourhoods of z,
relative to f (z) itself. Note that these neighbourhoods of z are allowed to become smaller
when f (z) is small. Finally, for Θ := (0, ∞)4 × A, and θ = (α, β, ν, γ, a) ∈ Θ, let
Fd,θ



:= f ∈ Fd : µα (f ) ≤ ν, kf k∞ ≤ γ, sup Mf,a,β (z) ≤ a(δ) ∀δ > 0 .
z:f (z)≥δ

Berrett, Samworth and Yuan (2017) show that all Gaussian and multivariate-t densities
(amongst others) belong to Fd,θ for appropriate θ ∈ Θ.
Now, for dX , dY ∈ N and ϑ = (θ, θY ) ∈ Θ2 , define
n
o
FdX ,dY ,ϑ := f ∈ FdX +dY ,θ : fY ∈ FdY ,θY , fX fY ∈ FdX +dY ,θ
and, for b ≥ 0, let
n
o
FdX ,dY ,ϑ (b) := f ∈ FdX ,dY ,ϑ : I(f ) > b .
Thus, FdX ,dY ,ϑ (b) consists of joint densities whose mutual information is greater than b. In
Theorem 2 below, we will show that for a suitable choice of b = bn and for certain ϑ ∈ Θ2 ,
the power of the test defined in Lemma 1 converges to 1, uniformly over FdX ,dY ,ϑ (b).
Before we can state this result, however, we must define the allowable choices of k and
9

the weight vectors. Given d ∈ N and θ = (α, β, γ, ν, a) ∈ Θ let


α−d
4(β ∧ 1)
2α
,
,
τ1 (d, θ) := min
5α + 3d
2α
4(β ∧ 1) + 3d
and





d
d/4
τ2 (d, θ) := min 1 −
, 1−
2β
bd/4c + 1

Note that mini=1,2 τi (d, θ) > 0 if and only if both α > d and β > d/2. Finally, for k ∈ N, let
W

(k)


:=

T

k

w = (w1 , . . . , wk ) ∈ R :

k
X

wj

j=1
k
X

Γ(j + 2`/d)
= 0 for ` = 1, . . . , bd/4c
Γ(j)


wj = 1 and wj = 0 if j ∈
/ {bk/dc, b2k/dc, . . . , k} .

j=1

Thus, our weights sum to 1; the other constraints ensure that the dominant contributions
to the bias of the unweighted Kozachenko–Leonenko estimator cancel out to sufficiently
high order, and that the corresponding jth nearest neighbour distances are not too highly
correlated.
Theorem 2. Fix dX , dY ∈ N, set d = dX + dY and fix ϑ = (θ, θY ) ∈ Θ2 with
n
o
min τ1 (d, θ) , τ1 (dY , θY ) , τ2 (d, θ) , τ2 (dY , θY ) > 0.
∗
∗
and k ∗ = kn∗ denote any deterministic sequences of positive integers
, kY∗ = kY,n
Let k0∗ = k0,n

with k0∗ ≤ min{kY∗ , k ∗ }, with k0∗ / log5 n → ∞ and with

max

k∗
nτ1 (d,θ)−

,

kY∗
nτ1 (dY ,θY )−

,

k∗
nτ2 (d,θ)

,

kY∗
nτ2 (dY ,θY )


→0

(k )

for some  > 0. Also suppose that wY = wY Y ∈ W (kY ) and w = w(k) ∈ W (k) , and that
lim supn max(kwk, kwY k) < ∞. Then there exists a sequence (bn ) such that bn = o(n−1/2 )
and with the property that for each q ∈ (0, 1) and any sequence (Bn∗ ) with Bn∗ → ∞,
inf

inf

inf

∗ k ∈{k ∗ ,...,k ∗ } f ∈F
Bn ≥Bn
Y
dX ,dY ,ϑ (bn )
0
Y
k∈{k0∗ ,...,k∗ }

b(n),Bn ) → 1.
Pf (Ibn > C
q

Theorem 2 provides a strong guarantee on the ability of MINTknown to detect alternatives,
10

uniformly over classes whose mutual information is at least bn , where we may even have
bn = o(n−1/2 ).

4

The case of unknown marginal distributions

We now consider the setting in which the marginal distributions of both X and Y are
unknown. Our test statistic remains the same, but now we estimate the critical value by
permuting our sample in an attempt to mimic the behaviour of the test statistic under
H0 . More explicity, for some B ∈ N, we propose independently of (X1 , Y1 ), . . . , (Xn , Yn ) to
simulate independent random variables τ1 , . . . , τB uniformly from Sn , the permutation group
(b)
(b)
(b)
(b)
of {1, . . . , n}, and for b = 1, . . . , B, set Z := (Xi , Yτ (i) ) and I˜n := Ibn (Z1 , . . . , Zn ). For
i

q ∈ (0, 1), we can now estimate
C̃q(n),B

(n)
Cq

b

by


:= inf r ∈ R :


B
1 X
1 (b) ≤ q ,
B + 1 b=0 {I˜n ≥r}

(0)
(n),B
where I˜n := Ibn , and refer to the test that rejects H0 if and only if Ibn > C̃q
as

MINTunknown(q).
Lemma 3. For any q ∈ (0, 1) and B ∈ N, the MINTunknown(q) test has size at most q:
sup


P Ibn > C̃q(n),B ≤ q.

sup

k∈{1,...,n−1} (X,Y ):I(X;Y )=0
(n),B
Note that Ibn > C̃q
if and only if
B

1 X
1 (b)
≤ q.
B + 1 b=0 {I˜n ≥Ibn }

(3)

This shows that estimating either of the marginal entropies is unnecessary to carry out
(b)
b Z − H̃n(b) , where H̃n(b) := H
b d,wZ (Z1(b) , . . . , Zn(b) ) is the weighted
the test, since I˜n − Ibn = H
n
n,k
Kozachenko–Leonenko joint entropy estimator based on the permuted data.
We now study the power of MINTunknown, and begin by introducing the classes of
marginal densities that we consider. To define an appropriate notion of smoothness, for

11

z ∈ {w : g(w) > 0} =: W, g ∈ Fd and δ > 0, let

rz,g,δ :=

δeΨ(k)
Vd (n − 1)g(z)

1/d
.

(4)

Now, for A belonging to the class of Borel subsets of W, denoted B(W), define
Mg (A) := sup sup
δ∈(0,2] z∈A

Z

1
d
g(z)
Vd rz,g,δ

g dλd − 1 .
Bz (rz,g,δ )

Both rz,g,δ and Mg (·) depend on n and k, but for simplicity we suppress this in our notation.
Let φ = (α, µ, ν, (cn ), (pn )) ∈ (0, ∞)3 × [0, ∞)N × [0, ∞)N =: Φ and define
GdX ,dY ,φ


:= f ∈ FdX +dY : max{µα (fX ), µα (fY )} ≤ µ, max{kfX k∞ , kfY k∞ } ≤ ν,

Z
fX fY dλd ≤ pn ∀n ∈ N .
∃Wn ∈ B(X × Y) s.t. MfX fY (Wn ) ≤ cn ,
Wnc

In addition to controlling the αth moment and uniform norms of the marginals fX and
fY , the class Gd,φ asks for there to be a (large) set Wn on which this product of marginal
densities is uniformly well approximated by a constant over small balls. This latter condition
is satisfied by products of many standard parametric families of marginal densities, including
normal, Weibull, Gumbel, logistic, gamma, beta, and t densities, and is what ensures that
nearest neighbour methods are effective in this context.
The corresponding class of joint densities we consider, for φ = (α, µ, ν, (cn ), (pn )) ∈ Φ, is
Hd,φ


:= f ∈ Fd : µα (f ) ≤ µ, kf k∞ ≤ ν,
Z
∃Zn ∈ B(Z) s.t. Mf (Zn ) ≤ cn ,


f dλd ≤ pn ∀n ∈ N .

Znc

In many cases, we may take Zn = {z : f (z) ≥ δn }, for some appropriately chosen sequence
(δn ) with δn → 0 as n → ∞. For instance, suppose we fix d ∈ N and θ = (α, β, ν, γ, a) ∈ Θ.
Then, by Berrett, Samworth and Yuan (2017, Lemma 12), there exists such a sequence (δn ),
as well as sequences (cn ) and (pn ), where
d

ka(k/(n − 1)) β∧1
δn =
log(n − 1),
n−1

β∧1

β∧1
15 2 d d3/2
cn =
log− d (n − 1),
7 d + (β ∧ 1)

12

α

for large n and pn = o((k/n) α+d − ) for every  > 0, such that Fd,θ ⊆ Hd,φ with φ =
(α, µ, ν, (cn ), (pn )) ∈ Φ. We are now in a position to state our main result on the power of
MINTunknown.
Theorem 4. Let dX , dY ∈ N, let d = dX +dY and fix φ = (α, µ, ν, (cn ), (pn )) ∈ Φ with cn → 0
∗
∗
and pn = o(1/ log n) as n → ∞. Let k0∗ = k0,n
and k1∗ = k1,n
denote two deterministic

sequences of positive integers satisfying k0∗ ≤ k1∗ , k0∗ / log2 n → ∞ and (k1∗ log2 n)/n → 0.
Then for any b > 0, q ∈ (0, 1) and any sequence (Bn∗ ) with Bn∗ → ∞ as n → ∞,
inf

inf

inf

∗ k∈{k ∗ ,...,k ∗ } f ∈G
Bn ≥Bn
dX ,dY ,φ ∩Hd,φ :
0
1
I(f )≥b

Pf (Ibn > C̃q(n),Bn ) → 1

as n → ∞.
Theorem 4 shows that MINTunknown is uniformly consistent against a wide class of alternatives.

5

Regression setting

In this section we aim to extend the ideas developed above to the problem of goodness-offit testing in linear models. Suppose we have independent and identically distributed pairs
(X1 , Y1 ), . . . , (Xn , Yn ) taking values in Rp × R, with E(Y12 ) < ∞ and E(X1 X1T ) finite and
positive definite. Then
β0 := argmin E{(Y1 − X1T β)2 }
β∈Rp

is well-defined, and we can further define i := Yi − XiT β0 for i = 1, . . . , n. We show in the
proof of Theorem 6 below that E(1 X1 ) = 0, but for the purposes of interpretability and
inference, it is often convenient if the random design linear model
Yi = XiT β0 + i ,

i = 1, . . . , n,

holds with Xi and i independent. A goodness-of-fit test of this property amounts to a test
of H0 : X1 ⊥⊥ 1 against H1 : X1 6⊥⊥ 1 . The main difficulty here is that 1 , . . . , n are not
observed directly. Given an estimator βb of β0 , the standard approach for dealing with this
problem is to compute residuals b
i := Yi − XiT βb for i = 1, . . . , n, and use these as a proxy
for 1 , . . . , n . Many introductory statistics textbooks, e.g. Dobson (2002, Section 2.3.4),
Dalgaard (2002, Section 5.2) suggest examining for patterns plots of residuals against fitted
13

values, as well as plots of residuals against each covariate in turn, as a diagnostic, though it
is difficult to formalise this procedure. (It is also interesting to note that when applying the
plot function in R to an object of type lm, these latter plots of residuals against each covariate
in turn are not produced, presumably because it would be prohibitively time-consuming to
check them all in the case of many covariates.)
The naive approach based on our work so far is simply to use the permutation test of
Section 4 on the data (X1 , b
1 ), . . . , (Xn , b
n ). Unfortunately, calculating the test statistic Ibn
on permuted data sets does not result in an exchangeable sequence, which makes it difficult
to ensure that this test has the nominal size q. To circumvent this issue, we assume that
the marginal distribution of 1 under H0 has E(1 ) = 0 and is known up to a scale factor
σ 2 := E(21 ); in practice, it will often be the case that this marginal distribution is taken to
be N (0, σ 2 ), for some unknown σ > 0. We also assume that we can sample from the density
fη of η1 := 1 /σ; of course, this is straightforward in the normal distribution case above. Let
X = (X1 · · · Xn )T , Y = (Y1 , . . . , Yn )T , and suppose the vector of residuals b
 := (b
1 , . . . , b
n )T is
computed from the least squares estimator βb := (X T X)−1 X T Y . We then define standardised
residuals by ηbi := b
i /b
σ , for i = 1, . . . , n, where σ
b2 := n−1 kb
k2 ; these standardised residuals
are invariant under changes of scale of  := (1 , . . . , n )T . Suppressing the dependence of
our entropy estimators on k and the weights w1 , . . . , wk for notational simplicity, our test
statistic is now given by

b np (X1 , . . . , Xn ) + H
b n1 (b
b np+1 (X1 , ηb1 ), . . . , (Xn , ηbn ) .
Iˇn(0) := H
η1 , . . . , ηbn ) − H
Writing η := /σ, we note that

1/2

n
ηi − XiT (X T X)−1 X T η
1
1
Tb
T b
i − Xi (β − β0 ) =
,
ηbi = (Yi − Xi β) =
σ
b
σ
b
kη − X T (X T X)−1 X T ηk
(b)

(b)

whose distribution does not depend on the unknown β0 or σ 2 . Let {η (b) = (η1 , . . . , ηn )T :
b = 1, . . . , B} denote independent random vectors, whose components are generated independently from fη . For b = 1, . . . , B we then set sb(b) := n−1/2 kb
η (b) k and, for i = 1, . . . , n,
let
(b)

ηbi :=

1  (b)
ηi − XiT (X T X)−1 X T η (b) .
(b)
sb

We finally compute

(b)
b np (X1 , . . . , Xn ) + H
b n1 (b
b np+1 (X1 , ηb1(b) ), . . . , (Xn , ηbn(b) ) .
Iˇn(b) := H
η1 , . . . , ηbn(b) ) − H
14

Analogously to our development in Sections 3 and 4, we can then define a critical value by
Čq(n),B


:= inf r ∈ R :


B
1 X
1 (b) ≤ q .
B + 1 b=0 {Iˇn ≥r}
(0)

(B)

Conditional on X and under H0 , the sequence (Iˇn , . . . , Iˇn ) is independent and identically
distributed and unconditionally it is an exchangeable sequence under H0 . This is the crucial
observation from which we can immediately derive that the resulting test has the nominal
size.
Lemma 5. For each q ∈ (0, 1) and B ∈ N, the MINTregression(q) test that rejects H0 if
(0)
(n),B
and only if Iˇn > Čq
has size at most q:
sup


P Iˇn > Čq(n),B ≤ q.

sup

k∈{1,...,n−1} (X,Y ):I(X;Y )=0
(0)
(b)
As in previous sections, we are only interested in the differences Iˇn − Iˇn for b = 1, . . . , B,
b np (X1 , . . . , Xn ) terms cancel out, so these marginal entropy
and in these differences, the H

estimators need not be computed.
In fact, to simplify our power analysis, it is more convenient to define a slightly modified
test, which also has the nominal size. Specifically, we assume for simplicity that m := n/2
is an integer, and consider a test in which the sample is split in half, with the second half of
the sample used to calculate the estimators βb(2) and σ
b2 of β0 and σ 2 respectively. On the
(2)

first half of the sample, we calculate
ηbi,(1) :=

Yi − XiT βb(2)
σ
b(2)

for i = 1, . . . , m and the test statistic

p+1
p
1
bm
bm
bm
(b
η1,(1) , . . . , ηbm,(1) ) − H
(X1 , ηb1,(1) ), . . . , (Xm , ηbm,(1) ) .
I˘n(0) := H
(X1 , . . . , Xm ) + H
(b)

Corresponding estimators {I˘n : b = 1, . . . , B} based on the simulated data may also be
computed using the same sample-splitting procedure, and we then obtain the critical value
(n),B

C̆q

in the same way as above. The advantage from a theoretical perspective of this
2
approach is that, conditional on βb(2) and σ
b(2)
, the random variables ηb1,(1) , . . . , ηbm,(1) are
independent and identically distributed.
15

To describe the power properties of MINTregression, we first define several densities:
for γ ∈ Rp and s > 0, let fηbγ,s and fηbγ,s
b1γ,s := (η1 − X1T γ)/s and
(1) denote the densities of η
(1),γ,s

ηb1

(1)

γ,s
γ,s
be the densities of (X1 , ηb1γ,s )
:= (η1 −X1T γ)/s respectively; further, let fX,b
η and fX,b
η (1)
(1),γ,s

and (X1 , ηb1

) respectively. Note that imposing assumptions on these densities amounts

to imposing assumptions on the joint density f of (X, ). For Ω := Θ × Θ × (0, ∞) × (0, 1) ×
∗
denote the class of joint
(0, ∞) × (0, ∞), and ω = (θ1 , θ2 , r0 , s0 , Λ, λ0 ), we therefore let Fp+1,ω

densities f of (X, ) satisfying the following three requirements:
(i)
o n
o
n
γ,s
γ,s
fηb : γ ∈ B0 (r0 ), s ∈ [s0 , 1/s0 ] ∪ fηb(1) : γ ∈ B0 (r0 ), s ∈ [s0 , 1/s0 ] ⊆ F1,θ1
and
n
o n
o
γ,s
γ,s
fX,b
:
γ
∈
B
(r
),
s
∈
[s
,
1/s
]
∪
f
:
γ
∈
B
(r
),
s
∈
[s
,
1/s
]
⊆ Fp+1,θ2 .
0 0
0
0
0 0
0
0
η
X,b
η (1)
(ii)
n
o
sup max E log2 fηbγ,1 (η1 ) , E log2 fηbγ,1
(η
)
≤ Λ,
1
(1)

(5)

o
n

(1),γ,1 
γ,1
2
2
≤ Λ.
sup max E log fη ηb1 , E log fη ηb1

(6)

γ∈B0 (r0 )

and
γ∈B0 (r0 )

(iii) Writing Σ := E(X1 X1T ), we have λmin (Σ) ≥ λ0 .
The first of these requirements ensures that we can estimate efficiently the marginal entropy
of our scaled residuals, as well as the joint entropy of these scaled residuals and our covariates.
The second condition is a moment condition that allows us to control |H(η1 − X1T γ) − H(η1 )|
(and similar quantities) in terms of kγk, when γ belongs to a small ball around the origin.
To illustrate the second part of this condition, it is satisfied, for instance, if fη is a standard
normal density and E(kX1 k4 ) < ∞, or if fη is a t density and E(kX1 kα ) < ∞ for some α > 0;
the first part of the condition is a little more complicated but similar. The final condition is
very natural for random design regression problems.
(0)
(B)
By the same observation on the sequence (I˘n , . . . , I˘n ) as was made regarding the
(0)
(B)
sequence (Iˇn , . . . , Iˇn ) just before Lemma 5, we see that the sample-splitting version of the
MINTregression(q) test has size at most q.
Theorem 6. Fix p ∈ N and ω = (θ1 , θ2 , r0 , s0 , Λ, λ0 ) ∈ Ω, where the first component of θ2 is
16

α2 ≥ 4 and the second component of θ1 is β1 ≥ 1. Assume that
n
o
min τ1 (1, θ1 ) , τ1 (p + 1, θ2 ) , τ2 (1, θ1 ) , τ2 (p + 1, θ2 ) > 0.
∗
∗
and k ∗ = kn∗ denote any deterministic sequences of positive integers
, kη∗ = kη,n
Let k0∗ = k0,n

with k0∗ ≤ min{kη∗ , k ∗ }, with k0∗ / log5 n → ∞ and with

max

k∗
nτ1 (p+1,θ2 )−

,

kη∗
nτ1 (1,θ1 )−

,

(kη )

for some  > 0. Also suppose that wη = wη

k∗
nτ2 (p+1,θ2 )

,

kη∗
nτ2 (1,θ1 )


→0

∈ W (kη ) and w = w(k) ∈ W (k) , and that

lim supn max(kwk, kwη k) < ∞. Then for any sequence (bn ) such that n1/2 bn → ∞, any
q ∈ (0, 1) and any sequence (Bn∗ ) with Bn∗ → ∞,
inf

inf

inf

∗ k ∈{k ∗ ,...,k ∗ } f ∈F ∗
Bn ≥Bn
η
η
0
1,p+1 :I(f )≥bn
k∈{k0∗ ,...,k∗ }

Pf (I˘n > C̆q(n),Bn ) → 1.

Finally in this section, we consider partitioning our design matrix as X = (X ∗ X ∗∗ ) ∈
Rn×(p0 +p1 ) , with p0 + p1 = p, and describe an extension of MINTregression to cases where
we are interested in testing the independence between  and X ∗ . For instance, X ∗∗ may
consist of an intercept term, or transformations of variables in X ∗ , as in the real data
example presented in Section 6.3 below. Our method for simulating standardised residual
vectors {b
η (b) : b = 1, . . . , B} remains unchanged, but our test statistic and corresponding
null statistics become
b p (X ∗ , . . . , X ∗ ) + H
b 1 (b
b p+1 (X ∗ , ηb1 ), . . . , (X ∗ , ηbn )
bn ) − H
I¯n(0) := H
n
1
n
n η1 , . . . , η
n
1
n




(b)
b p (X ∗ , . . . , X ∗ ) + H
b 1 (b
b p+1 (X ∗ , ηb1(b) ), . . . , (X ∗ , ηb(b) ) , b = 1, . . . , B.
I¯n(b) := H
bn(b) ) − H
n
1
n
n η1 , . . . , η
n
1
n n
(0) (1)
(B)
The sequence (I¯n , I¯n , . . . , I¯n ) is again exchangeable under H0 , so a p-value for this mod-

ified test is given by
B

1 X
1 (b) (0) .
B + 1 b=0 {I¯n ≥I¯n }

17

6

Numerical studies

6.1

Practical considerations

For practical implementation of the MINTunknown test, we consider both a direct, data-driven
approach to choosing k, and a multiscale approach that averages over a range of values of
k. To describe the first method, let K ⊆ {1, . . . , n − 1} denote a plausible set of values of
k. For a given N ∈ N and independently of the data, generate τ1 , . . . , τ2N independently
b (j) be the (unweighted)
and uniformly from Sn , and for each k ∈ K and j = 1, . . . , 2N let H
n,k

Kozachenko–Leonenko joint entropy estimate with tuning parameter k based on the sample
{(Xi , Yτj (i) ) : i = 1, . . . , n}. We may then choose
b
k := sargmin
k∈K

N
X

b (2j) − H
b (2j−1) )2 ,
(H
n,k
n,k

j=1

where sargmin denotes the smallest index of the argmin in the case of a tie. In our simulations
below, which use N = 100, we refer to the resulting test as MINTauto.
For our multiscale approach, we again let K ⊆ {1, . . . , n−1} and, for k ∈ K, let b
h(0) (k) :=
b n,k denote the (unweighted) Kozachenko–Leonenko entropy estimate with tuning parameter
H
k based on the original data (X1 , Y1 ), . . . , (Xn , Yn ). Now, for b = 1, . . . , B and k ∈ K,
(b)
we let b
h(b) (k) := H̃n,(k) denote the Kozachenko–Leonenko entropy estimate with tuning
P
(b)
(b)
b
h(b) (k)
parameter k based on the permuted data Z , . . . , Zn . Writing h̄(b) := |K|−1
1

k∈K

for b = 0, 1, . . . , B, we then define the p-value for our test to be
B

1 X
1 (0) (b) .
B + 1 b=0 {h̄ ≥h̄ }
By the exchangeability of (h̄(0) , h̄(1) , . . . , h̄(B) ) under H0 , the corresponding test has the nominal size. We refer to this test as MINTav, and note that if K is taken to be a singleton
set then we recover MINTunknown. In our simulations below, we took K = {1, . . . , 20} and
B = 100.

6.2

Simulated data

To study the empirical performance of our methods, we first compare our tests to existing
approaches through their performance on simulated data. For comparison, we present cor18

responding results for a test based on the empirical copula process decribed by Kojadinovic
and Holmes (2009) and implemented in the R package copula (Hofert et al., 2017), a test
based on the HSIC implemented in the R package dHSIC (Pfister and Peters, 2017), a test
based on the distance covariance implemented in the R package energy (Rizzo and Szekely,
2017) and the improvement of Hoeffding’s test described in Bergsma and Dassios (2014) and
implemented in the R package TauStar (Weihs, Drton and Leung, 2016). We also present
results for an oracle version of our tests, denoted simply as MINT, which for each parameter value in each setting, uses the best (most powerful) choice of k. Throughout, we took
q = 0.05 and n = 200, ran 5000 repetitions for each parameter setting, and for our comparison methods, used the default tuning parameter values recommended by the corresponding
authors. We consider three classes of data generating mechanisms, designed to illustrate
different possible types of dependence:
(i) For l ∈ N and (x, y) ∈ [−π, π]2 , define the density function
fl (x, y) =

1
{1 + sin(lx) sin(ly)}.
π2

This class of densities, which we refer to as sinusoidal, are identified by Sejdinovic et al.
(2013) as challenging ones to detect dependence, because as l increases, the dependence
becomes increasingly localised, while the marginal densities are uniform on [−π, π] for
each l. Despite this, by the periodicity of the sine function, we have that the mutual
information does not depend on l: indeed,
Z πZ π

1
{1 + sin(lx) sin(ly)} log 1 + sin(lx) sin(ly) dx dy
I(fl ) = 2
π −π −π
Z πZ π
1
= 2
(1 + sin u sin v) log(1 + sin u sin v) du dv ≈ 0.0145.
π −π −π
(ii) Let L, Θ, 1 , 2 be independent with L ∼ U ({1, . . . , l}) for some parameter l ∈ N,
Θ ∼ U [0, 2π], and 1 , 2 ∼ N (0, 1). Set X = L cos Θ + 1 /4 and Y = L sin Θ + 2 /4.
For large values of l, the distribution of (X/l, Y /l) approaches the uniform distribution
on the unit disc.
(iii) Let X,  be independent with X ∼ U [−1, 1],  ∼ N (0, 1), and for a parameter ρ ∈
[0, ∞), let Y = |X|ρ .
For each of these three classes of data generating mechanisms, we also consider a corre19

sponding multivariate setting in which we wish to test the independence of X and Y when
X = (X1 , X2 )T , Y = (Y1 , Y2 )T . Here, (X1 , Y1 ), X2 , Y2 are independent, with X1 and Y1 having the dependence structures described above, and X2 , Y2 ∼ U (0, 1). In these multivariate
settings, the generalisations of the TauStar test were too computationally intensive, so were
omitted from the comparison.
The results are presented in Figure 1(a). Unsurprisingly, there is no uniformly most
powerful test among those considered, and if the form of the dependence were known in
advance, it may well be possible to design a tailor-made test with good power. Nevertheless,
Figure 1(a) shows that, especially in the first and second of the three settings described
above, the MINT and MINTav approaches have very strong performance. In these examples,
the dependence becomes increasingly localised as l increases, and the flexibility to choose
a smaller value of k in such settings means that MINT approaches are particularly effective.
Where the dependence is more global in nature, such as in setting (iii), other approaches
may be better suited, though even here, MINT is competitive. Interestingly, the multiscale
MINTav appears to be much more effective than the MINTauto approach of choosing a single
value of k; indeed, in setting (iii), MINTav even outperforms the oracle choice of k.

6.3

Real data

In this section we illustrate the use of MINTregression on a data set comprising the average minimum daily January temperature, between 1931 and 1960 and measured in degrees Fahrenheit, in 56 US cities, indexed by their latitude and longitude (Peixoto, 1990).
Fitting a normal linear model to the data with temperature as the response and latitude
and longitude as covariates leads to the diagnostic plots shown in Figures 2(a) and 2(b).
Based on these plots, it is relatively difficult to assess the strength of evidence against the
model. Centering all of the variables and running MINTregression with kη = 6, k = 3 and
B = 1785 = b105 /56c yielded a p-value of 0.00224, providing strong evidence that the normal
linear model is not a good fit to the data. Further investigation is possible via partial regression: in Figure 2(c) we plot residuals based on a simple linear fit of Temperature on Latitude
against the residuals of another linear fit of Longitude on Latitude. This partial regression
plot is indicative of a cubic relationship between Temperature and Longitude after removing
the effects of Latitude. Based on this evidence we fitted a new model with a quadratic and
cubic term in Longitude added to the first model. The p-value of the resulting F -test was
< 2.2 × 10−16 , giving overwhelming evidence in favour of the inclusion of the quadratic and
cubic terms. We then ran our extension of MINTregression to test the independence of 
20

Figure 1: Power curves for the different tests considered, for settings (i) (left), (ii) (middle) and (iii)
(right). The marginals are univariate (top) and bivariate (bottom).

Figure 2: Plots of residuals against fitted values (left) and square roots of the absolute values of the
standardised residuals against fitted values (middle) for the US cities temperature data. The right
panel gives a partial regression plot after removing the effect of latitude from both the response
and longitude.

21

and (Longitude, Latitude), as described at the end of Section 5 with B = 1000. This resulted
in a p-value of 0.0679, so we do not reject this model at the 5% significance level.

7

Proofs

7.1

Proofs from Section 3
(1)

(B)

Proof of Lemma 1. Under H0 , the finite sequence (Ibn , Ibn , . . . , Ibn ) is exchangeable. If fol(1)
(B)
lows that under H0 , if we split ties at random, then every possible ordering of Ibn , Ibn , . . . , Ibn
(1)
(B)
is equally likely. In particular, the (descending order) rank of Ibn among {Ibn , Ibn , . . . , Ibn },
denoted rank(Ibn ), is uniformly distributed on {1, 2, . . . , B + 1}. We deduce that for any

k ∈ {1, . . . , n − 1} and writing PH0 for the probability under the product of any marginal
distributions for X and Y ,

PH0


bq(n),B = PH0
Ibn > C



B

1 X
1 (b)
≤q
B + 1 b=0 {Ibn ≥Ibn }



 bq(B + 1)c
≤ PH0 rank(Ibn ) ≤ q(B + 1) =
≤ q,
B+1
as required.
(X,Y )
Let V (f ) := Var log fX f(X)f
. The following result will be useful in the proof of TheoY (Y )

rem 2.
Lemma 7. Fix dX , dY ∈ N and ϑ ∈ Θ2 . Then
V (f )
< ∞.
1/4
f ∈FdX ,dY ,ϑ (0) I(f )
sup

Proof of Lemma 7. For x ∈ R we write x− := max(0, −x). We have

E log

f (X, Y )
fX (X)fY (Y )

 

Z

fX (x)fY (y)
f (x, y) log
dλd (x, y)
f (x, y)
{(x,y):f (x,y)≤fX (x)fY (y)}


Z
fX (x)fY (y)
≤
f (x, y)
− 1 dλd (x, y)
f (x, y)
{(x,y):f (x,y)≤fX (x)fY (y)}
Z
Z 
= sup
fX fY −
f ≤ {I(f )/2}1/2 ,
=

−

A∈B(Rd )

A

A

22

where the final inequality is Pinsker’s inequality. We therefore have that


f (X, Y )
f (X, Y ) 1/2
f (X, Y ) 3/2
= E log
V (f ) ≤ E log
log
fX (X)fY (Y )
fX (X)fY (Y )
fX (X)fY (Y )

1/2
1/2 
f (X, Y )
f (X, Y ) 3
≤ E log
E log
fX (X)fY (Y )
fX (X)fY (Y )
 1/2 
1/2


f (X, Y ) 3
f (X, Y )
E log
= I(f ) + 2E log
fX (X)fY (Y ) −
fX (X)fY (Y )



1/2
1/2
3
.
E log fX (X)fY (Y ) + E| log f (X, Y )|3
≤ 2 I(f ) + 21/2 I(f )1/2
2

By Berrett, Samworth and Yuan (2017, Lemma 11(i)), we conclude that
V (f )
< ∞.
1/4
f ∈FdX ,dY ,ϑ :I(f )≤1 I(f )
sup

The result follows from this fact, together with the observation that
sup

V (f ) ≤ 2

f ∈FdX ,dY ,ϑ

sup
f ∈FdX ,dY ,ϑ



E log2 f (X, Y ) + E log2 fX (X)fY (Y )



< ∞,

(7)

where the final bound follows from another application of Berrett, Samworth and Yuan
(2017, Lemma 11(i)).
Proof of Theorem 2. Fix f ∈ FdX ,dY ,ϑ (bn ). We have by two applications of Markov’s inequality that
1
{1 + Bn Pf (Ibn(1) ≥ Ibn )}
q(Bn + 1)

1
1
≤ Pf |Ibn − Ibn(1) − I(f )| ≥ I(f ) +
q
q(Bn + 1)
1
1
≤
Ef [{Ibn − Ibn(1) − I(f )}2 ] +
.
2
qI(f )
q(Bn + 1)

b(n),Bn ) ≤
Pf (Ibn ≤ C
q

(8)

It is convenient to define the following entropy estimators based on pseudo-data, as well as
(1)

oracle versions of our entropy estimators: let Zi

23

(1)T T

:= (XiT , Yi

(1)

(1)

) , Z (1) := (Z1 , . . . , Zn )T

(1)

(1)

and Y (1) := (Y1 , . . . , Yn )T and set
b Y,(1) := H
b Y,wY (Y (1) )
H
n
n,k
n
1X
∗,Y
Hn := −
log fY (Yi ),
n i=1

b Z,(1) := H
b Z,wZ (Z (1) ),
H
n
n,k
n
X
1
∗
Hn := −
log f (Xi , Yi ),
n i=1
n

Hn∗,(1)
Writing In∗ :=

1
n

n

1X
1X
(1)
(1)
log fX (Xi )fY (Yi ), Hn∗,Y,(1) := −
log fY (Yi ).
:= −
n i=1
n i=1

Pn

i=1

i ,Yi )
log fXf(X(Xi )f
we have by Berrett, Samworth and Yuan (2017, TheoY (Yi )

rem 1) that
sup

sup

∗}
kY ∈{k0∗ ,...,kY
k∈{k0∗ ,...,k∗ }

f ∈FdX ,dY ,ϑ

=

nEf {(Ibn − Ibn(1) − In∗ )2 }



 2
Y
Y,(1)
(1)
∗,Y
∗
∗,Y,(1)
∗,(1)
b
b
b
b
+ Hn
sup nEf Hn − Hn − Hn + Hn − Hn − Hn − Hn

sup

∗ } f ∈F
kY ∈{k0∗ ,...,kY
dX ,dY ,ϑ
k∈{k0∗ ,...,k∗ }

≤ 4n

sup

sup

∗ } f ∈F
kY ∈{k0∗ ,...,kY
dX ,dY ,ϑ
k∈{k0∗ ,...,k∗ }

n
b Y − H ∗,Y )2 + (H
b n − H ∗ )2
Ef (H
n
n
n
o
Y,(1)
∗,Y,(1) 2
(1)
∗,(1) 2
b
b
+ (Hn − Hn
) + (Hn − Hn ) → 0.


2
It follows by Cauchy–Schwarz and the fact that E In∗ − I(f )
= V (f )/n that
2n :=

sup

sup

∗ } f ∈F
kY ∈{k0∗ ,...,kY
dX ,dY ,ϑ
k∈{k0∗ ,...,k∗ }

≤

sup

sup

∗ } f ∈F
kY ∈{k0∗ ,...,kY
dX ,dY ,ϑ
k∈{k0∗ ,...,k∗ }

nEf


2
Ibn − Ibn(1) − I(f )
− V (f )

n

1/2 o
nEf {(Ibn − Ibn(1) − In∗ )2 } + 2 nEf {(Ibn − Ibn(1) − In∗ )2 }V (f )

→ 0,

(9)
1/2

where we use (7) to bound V (f ) above. Now consider bn := max(n n−1/2 , n−4/7 log n).

24

By (8), (9) and Lemma 7 we have that
inf

inf

inf

∗ k ∈{k ∗ ,...,k ∗ } f ∈F
Bn ≥Bn
Y
dX ,dY ,ϑ (bn )
0
Y
k∈{k0∗ ,...,k∗ }

≥1−

Pf (Ibn > Cq(n),Bn )

(1)
E[{Ibn − Ibn − I(f )}2 ]
1
−
2
∗
∗ } f ∈F
qI(f )
q(Bn + 1)
kY ∈{k0∗ ,...,kY
dX ,dY ,ϑ (bn )

sup

sup

k∈{k0∗ ,...,k∗ }

≥1−

V (f ) + 2n
1
−
2
q(Bn∗ + 1)
f ∈FdX ,dY ,ϑ (bn ) nqI(f )

≥1−

1
V (f )
1
n
−
→ 1,
sup
−
q
q(Bn∗ + 1)
q log7/4 n f ∈FdX ,dY ,ϑ (0) I(f )1/4

sup

as required.

7.2

Proofs from Section 4
(0)

(1)

(B)

Proof of Lemma 3. We first claim that (Ibn , Ibn , . . . , Ibn ) is an exchangeable sequence un(b)
der H0 . Indeed, let σ be an arbitrary permutation of {0, 1, . . . , B}, and recall that Ibn is
computed from (Xi , Yτb (i) )ni=1 , where τb is uniformly distributed on Sn . Note that, under H0 ,
d

we have (Xi , Yi )ni=1 = (Xi , Yτ (i) )ni=1 for any τ ∈ Sn . Hence, for any A ∈ B(RB+1 ),

PH0 (Ibn(σ(0)) , . . . , Ibn(σ(B)) ) ∈ A
X

1
bn(σ(0)) , . . . , Ibn(σ(B)) ) ∈ A|τ1 = π1 , . . . , τB = πB
=
P
(
I
H
0
(n!)B π ,...,π ∈S
n
1
B
X

1
bn(0) , . . . , Ibn(B) ) ∈ A|τ1 = πσ(1) π −1 , . . . , τB = πσ(B) π −1
P
(
I
=
H
0
σ(0)
σ(0)
(n!)B π ,...,π ∈S
n
1
B

= PH0 (Ibn(0) , . . . , Ibn(B) ) ∈ A .
By the same argument as in the proof of Lemma 1 we now have that
bq(B + 1)c
≤ q,
PH0 (Ibn > C̃q(n),B ) ≤
B+1
as required.
b n(1) given just after (3). We give the proof of Theorem 4 below,
Recall the definition of H
followed by Lemma 8, which is the key ingredient on which it is based.

25

Proof of Theorem 4. We have by Markov’s inequality that
Pf (Ibn ≤ C̃q(n),Bn ) ≤

1
1
1
b Z ≥ H̃ (1) ) +
{1 + Bn Pf (Ibn ≤ I˜n(1) )} ≤ Pf (H
.
n
n
q(Bn + 1)
q
q(Bn + 1)
(10)

But

bZ ≥ H
b (1) ) = Pf H
b Z − H(X, Y ) ≥ H
b (1) − H(X) − H(Y ) + I(X; Y )
Pf (H
n
n
n
n




1
b (1) − H(X) − H(Y )| ≥ 1 I(X; Y )
b Z − H(X, Y )| ≥ I(X; Y ) + P |H
≤ P |H
n
n
2
2

2
b Z − H(X, Y )| + E|H
b (1) − H(X) − H(Y )| .
≤
E|H
(11)
n
n
I(X; Y )
The theorem therefore follows immediately from (10), (11) and Lemma 8 below.
Lemma 8. Let dX , dY ∈ N, let d = dX + dY and fix φ = (α, µ, ν, (cn ), (pn )) ∈ Φ with cn → 0
∗
∗
and pn = o(1/ log n) as n → ∞. Let k0∗ = k0,n
and k1∗ = k1,n
denote two deterministic

sequences of positive integers satisfying k0∗ ≤ k1∗ , k0∗ / log2 n → ∞ and (k1∗ log2 n)/n → 0.
Then
sup

sup

k∈{k0∗ ,...,k1∗ }

f ∈GdX ,dY ,φ

sup
k∈{k0∗ ,...,k1∗ }

b (1) − H(X) − H(Y )| → 0,
Ef |H
n

b Z − H(X, Y )| → 0
sup Ef |H
n

f ∈Hd,φ

as n → ∞.
Proof. We prove only the first claim in Lemma 8, since the second claim involves similar
arguments but is more straightforward because the estimator is based on an independent
and identically distributed sample and no permutations are involved. Throughout this proof,
we write a . b to mean that there exists C > 0, depending only on dX , dY ∈ N and
φ ∈ Φ, such that a ≤ Cb. Fix f ∈ GdX ,dY ,φ , where φ = (α, µ, ν, (cn ), (pn )). Write ρ(k),i,(1)
(1)

for the distance from Zi
(1)
ξi

:= e

−Ψ(k)

Vd (n −

(1)

(1)

to its kth nearest neighbour in the sample Z1 , . . . , Zn

1)ρd(k),i,(1)

and

so that
n

X
(1)
b (1) = 1
H
log ξi .
n
n i=1
Recall that Sn denotes the set of all permutations of {1, . . . , n}, and define Snl ⊆ Sn to be the
26

set of permutations that fix {1, . . . , l} and have no fixed points among {l + 1, . . . , n}; thus
P
r
|Snl | = (n−l)! n−l
r=0 (−1) /r!. For τ ∈ Sn , write Pτ (·) := P(·|τ1 = τ ), and Eτ (·) := E(·|τ1 = τ ).
Let ln := blog log nc. For i = 1, . . . , n, let
(1)

A0i := {The k nearest neighbours of Zi
(1)

(1)

are among Zl+1 , . . . , Zn(1) },
(1)

and let Ai := A0i ∩ {Zi

(1)

∈ Wn }. Using the exchangeability of ξ1 , . . . , ξn , our basic

decomposition is as follows:
b n(1) − H(X) − H(Y )|
EH

ln   X  X
l
n
1 X
1 X
n
1
(1)
(1)
≤
Eτ | log ξi | + Eτ
log ξi − H(X) − H(Y )
n! l=0 l
n i=1
n i=l+1
l
τ ∈Sn
 
n
1 X n X
b n(1) − H(X) − H(Y )
+
Eτ H
n! l=l +1 l
l
τ ∈Sn

n

ln 
X



l
n X 1X

n
X

1
(1)
(log ξi )1Ai − H(X) − H(Y )
l
n i=1
n i=l+1
l
l=0
τ ∈Sn

 
n
n


1 X
1 X n X
(1)
b (1) | + |H(X) + H(Y )| .
+
Eτ | log ξi |1{(Aci \(A0i )c )∪(A0i )c } +
Eτ |H
n
n i=l+1
n! l=l +1 l
l
≤

1
n!

(1)

Eτ | log ξi | + Eτ

n

τ ∈Sn

(12)
The rest of the proof consists of handling each of these four terms.
Step 1: We first show that
sup

sup

k∈{1,...,n−1} f ∈GdX ,dY ,φ

 
n

1 X n X
b n(1) | + |H(X) + H(Y )| → 0
Eτ |H
n! l=l +1 l
l
n

(13)

τ ∈Sn

as n → ∞. Writing ρ(k),i,X for the distance from Xi to its kth nearest neighbour in the
sample X1 , . . . , Xn and defining ρ(k),i,Y similarly, we have
max{ρ2(k),i,X , ρ2(k),τ1 (i),Y } ≤ ρ2(k),i,(1) ≤ ρ2(n−1),i,X + ρ2(n−1),τ1 (i),Y .

27

Using the fact that log(a + b) ≤ log 2 + | log a| + | log b| for a, b > 0, we therefore have that

 2

ρ(k),i,(1)
Vd (n − 1)
d
≤ log
+ d| log ρ(k),i,X | + log 2
eΨ(k)
2
ρ(k),i,X


Vd (n − 1)
d
≤ log
+ 2d| log ρ(k),i,X | + log 2 + d| log ρ(n−1),i,X | + d| log ρ(n−1),τ1 (i),Y |
Ψ(k)
e
2
(14)


Vd (n − 1)
d
≤ log
+ log 2
Ψ(k)
e
2


(1)
| log ξi |

+ 4d max max{− log ρ(1),j,X , − log ρ(1),j,Y , log ρ(n−1),j,X , log ρ(n−1),j,Y }.
j=1,...,n

(15)

Now, by the triangle inequality, a union bound and Markov’s inequality,
n
o


hn 
o i
E max log ρ(n−1),j,X − log 2 ≤ E log max kXj k ≤ E log max kXj k
j=1,...,n
j=1,...,n
j=1,...,n
+
Z ∞ 

=
P max kXj k ≥ eM dM
0

j=1,...,n



1
1
max 0, log n + log E(kX1 kα ) + nE(kX1 kα ) exp − log n − log E(kX1 kα )
α
α
1
≤
1 + log n + max(0, log µ) ,
(16)
α
≤

and the same bound holds for E maxj=1,...,n log ρ(n−1),j,Y . Similarly,
hn
E

min log ρ(1),j,X

j=1,...,n

o i
−

Z

∞


min ρ(1),j,X ≤ e−M dM
j=1,...,n
0
Z
−1
≤ 2dX log n + n(n − 1)VdX kfX k∞

=



P

∞

2d−1
X

≤

2d−1
X

and the same bound holds for E

log n +

e−M dX dM
log n

d−1
X VdX ν,


minj=1,...,n log ρ(1),j,Y

(17)

−

once we replace dX with dY .

By (15), (16), (17) and Stirling’s approximation,
 
n  
n
X
X
n
n
(1)
b
max
sup E(|Hn |1{τ1 ∈Snl } ) . log n
P(τ1 ∈ Snl )
k=1,...,n−1 f ∈Gd ,d ,φ
l
l
X Y
l=ln +1
l=ln +1



ln +1
n
log n
e
log n
log n
(n − ln − 1)! =
≤p
→0
≤
n! ln + 1
(ln + 1)!
2π(ln + 1) ln + 1

28

(18)

as n → ∞. Moreover, for any f ∈ GdX ,dY ,φ and (X, Y ) ∼ f , writing  = α/(α + 2dX ),
Z
|H(X)| ≤
X




kfX k∞
fX | log fX | ≤
fX (x) log
dx + | log kfX k∞ |
fX (x)
X
Z
kfX k∞
≤
fX (x)1− dx + | log kfX k∞ |

X
Z

ν
1−
α − 1−

(1 + kxk )
dx
≤ (1 + µ)

RdX

 α dX

VdX µ (α + dX )α+dX
1
+ max log ν , log
,
α
αα ddXX
Z

(19)

where the lower bound on kfX k∞ is due to the fact that VdαX µα (fX )dX kfX kα∞ ≥ αα ddXX /(α +
dX )α+dX by Berrett, Samworth and Yuan (2017, Lemma 10(i)). Combining (19) with the
corresponding bound on H(Y ), as well as (18), we deduce that (13) holds.
Step 2: We observe that by (15), (16) and (17),

sup

sup

k∈{1,...,n−1} f ∈GdX ,dY ,φ

ln   X
l
ln log n
n
1X
1 X
(1)
→ 0.
Eτ | log ξi | .
n! l=0 l
n i=1
n
l

(20)

τ ∈Sn

Step 3: We now show that
sup

sup

k∈{k0∗ ,...,k1∗ } f ∈GdX ,dY ,φ

ln   X
n

1 X
n
1 X
(1)
Eτ | log ξi |1{(Aci \(A0i )c )∪(A0i )c } → 0.
n! l=0 l
n i=l+1
l

(21)

τ ∈Sn

We use an argument that involves covering Rd by cones; cf. Biau and Devroye (2015,
Section 20.7). For w0 ∈ Rd , z ∈ Rd \ {0} and θ ∈ [0, π], define the cone of angle θ centred at
w0 in the direction z by

C(w0 , z, θ) := w0 + w ∈ Rd \ {0} : cos−1 (z T w/(kzkkwk)) ≤ θ ∪ {0}.
By Biau and Devroye (2015, Theorem 20.15), there exists a constant Cπ/6 ∈ N depending
(1)

only on d such that we may cover Rd by Cπ/6 cones of angle π/6 centred at Z1 . In each
(1)

cone, mark the k nearest points to Z1

(1)

(1)

(1)

among Z2 , . . . , Zn . Now fix a point Zi

not marked and a cone that contains it, and let

29

(1)
(1)
Zi1 , . . . , Zik

that is

be the k marked points in this

cone. By Biau and Devroye (2015, Lemma 20.5) we have, for each j = 1, . . . , k, that
(1)

(1)

(1)

(1)

kZi − Zij k < kZi − Z1 k.
(1)

Thus, Z1

(1)

is not one of the k nearest neighbours of the unmarked point Zi , and only the
(1)

marked points, of which there are at most kCπ/6 , may have Z1

as one of their k near-

est neighbours. This immediately generalises to show that at most klCπ/6 of the points
(1)

(1)

Zl+1 , . . . , Zn

(1)

(1)

may have any of Z1 , . . . , Zl

among their k nearest neighbours. Then,

by (15), (16) and (17), we have that, uniformly over k ∈ {1, . . . , n − 1},
1
Eτ
max
sup
ln
l n
f ∈GdX ,dY ,φ τ ∈∪l=0 Sn

X
n

1

(1)
| log ξi | (A0i )c


.

i=l+1

Now, for x ∈ X and r ∈ [0, ∞), let hx (r) :=

R
Bx (r)

kln
log n.
n

fX , and, for s ∈ [0, 1), let h−1
x (s) :=

Γ(n)
sk−1 (1
Γ(k)Γ(n−k)

inf{r ≥ 0 : hx (r) ≥ s}. We also write Bk,n−k (s) :=

(22)

− s)n−k−1 for s ∈ (0, 1).

Then, by (4) and (13) in Berrett, Samworth and Yuan (2017), there exists C > 0, depending
only on dX and φ, such that, uniformly for k ∈ {1, . . . , n − 1},

max Eτ | log ρ(k),i,X |1{Z (1) ∈W c }
n
i
i=l+1,...,n


Z 1
Z
VdX (n − 1)h−1
x (s)
log
fX (x)fY (y)
Bk,n−k (s) ds dλd (x, y)
=
eΨ(k)
0
Wnc
Z
Z 1

.
fX (x)fY (y)
log n − log s − log(1 − s) + log(1 + Ckxk) Bk,n−k (s) ds dλd (x, y).
Wnc

0

(23)
Now, given  ∈ (0, 1), by Hölder’s inequality, and with K := max{1, (α−1) log 2, C α }/(α),
Z

Z
fX (x)fY (y) log(1 + Ckxk) dλd (x, y) ≤ K

Wnc

Z
≤ K

α

fX (x)fY (y)(1
Wnc
 Z

fX (x)fY (y)(1 + kxk ) dλd (x, y)

+ kxkα ) dλd (x, y)
1−

fX (x)fY (y) dλd (x, y)

Wnc

Wnc

≤ K (1 + µ) p1−
n .

(24)

We deduce from (23) and (24) that for each  ∈ (0, 1) and uniformly for k ∈ {1, . . . , n − 1},
max

sup

l+1=1,...,n f ∈Gd


Eτ | log ρ(k),i,X |1{Z (1) ∈W c } . pn log n + K pn1− .
i

X ,dY ,φ

30

n

From similar bounds on the corresponding terms with ρ(k),i,X replaced with ρ(n−1),i,X and
ρ(n−1),τ1 (i),Y , we conclude from (14) that for every  ∈ (0, 1) and uniformly for k ∈ {1, . . . , n−
1},
max


(1)
Eτ | log ξi |1{Z (1) ∈W c } . pn log n + K pn1− .

sup

i=l+1,...,n f ∈Gd

(25)

n

i

X ,dY ,φ

The claim (21) follows from (22) and (25).
Step 4: Finally, we show that
sup

sup

k∈{k0∗ ,...,k1∗ } f ∈GdX ,dY ,φ

ln   X
n
1 X
1 X
n
(1)
Eτ
(log ξi )1Ai − H(X) − H(Y ) → 0.
n! l=0 l
n
l
i=l+1

(26)

τ ∈Sn

To this end, we consider the following further decomposition:
n
1 X
(1)
(log ξi )1Ai − H(X) − H(Y )
Eτ
n i=l+1
n
n
X
1 X
1
(1)
(1) 
(1) 
≤
Eτ log ξi fX fY (Zi ) 1Ai + Eτ
log fX fY (Zi ) 1{(Aci \(A0i )c )∪(A0i )c }
n i=l+1
n i=l+1
n
n − l
1 X
(1)
log fX fY (Zi ) −
H(X) + H(Y )
+ Eτ
n i=l+1
n

+

l
H(X) + H(Y ) .
n

(27)

We deal first with the second term in (27). Now, by the arguments in Step 3,


n
X
Cπ/6 kl
1
(1) 
(1)
Eτ
log fX fY (Zi ) 1(A0i )c ≤
Eτ
max | log fX fY (Zi )|
i=l+1,...,n
n i=l+1
n
 



Cπ/6 k1∗ l
E max | log fX (Xi )| + E max | log fY (Yi )| .
≤
i=l+1,...,n
i=l+1,...,n
n

(28)

Now, for any s > 0, by Hölder’s inequality,
α

P fX (X1 ) ≤ s ≤ s 2(α+dX )

Z

α
1− 2(α+d

fX (x)

X)

dx

x:f (x)≤s

≤s

α
2(α+dX )

α
1− 2(α+d

(1 + µ)

Z
X)

RdX

1
dx
α
(1 + kxk )(α+2dX )/α

31

α
 2(α+d

X)

α

=: K̃dX ,φ s 2(α+dX ) .

X)
Writing t∗ := − 2(α+d
log(n − l), we deduce that for n ≥ 3,
α




E max | log fX (Xi )| ≤ log ν − E
min


log fX (Xi )
i=l+1,...,n
i=l+1,...,n
Z t∗
αt
2(α + dX )
e 2(α+dX ) dt +
≤ log ν + (n − l)K̃dX ,φ
log n
α
−∞
2(α + dX )
2(α + dX )
= log ν +
log n.
(29)
K̃dX ,φ +
α
α

Combining (29) with the corresponding bound on E maxi=l+1,...,n | log fY (Yi )| , which is obtained in the same way, we conclude from (28) that
sup

sup

k∈{k0∗ ,...,k1∗ } f ∈GdX ,dY ,φ

n
X
k ∗ ln log n
1
(1) 
log fX fY (Zi ) 1(A0i )c . 1
Eτ
→ 0.
n i=l+1
n

(30)


Moreover, given any  ∈ 0, (α + 2d)/α , set 0 := α/(α + 2d). Then by two applications of
Hölder’s inequality,
n
X
1
(1) 
Eτ
log fX fY (Zi ) 1{Z (1) ∈W c }
n
i
n i=l+1
Z
0 Z
ν
0
fX fY (z)| log fX fY (z)| dz ≤ νpn + 0
≤
fX fY (z)1− dz
 Wnc
Wnc


Z
0

ν
2d/(α+2d)
≤ νpn + 0 p1−
f
f
(z)
dz
X Y
 n
Wnc

Z
α/(α+2d) 
0
ν  1− 
1
2d/(α+2d)
α−1
. p1−
≤ νpn + 0 pn
1 + max(1, 2 )µ
dz
n .
α )2d/α

(1
+
kzk
d
R
(31)

From (30) and (31), we find that
n
X
1
(1) 
sup
sup
max
Eτ
log fX fY (Zi ) 1{(Aci \(A0i )c )∪(A0i )c } → 0,
ln
l n
Sn
k∈{k0∗ ,...,k1∗ } f ∈GdX ,dY ,φ τ ∈∪l=0
i=l+1

(32)

which takes care of the second term in (27).
(1)

For the third term in (27), since Zi

 (j)
and Zi : j ∈
/ {i, τ (i), τ −1 (i)} are independent,

32

we have
 X

n
n
1 X
n − l
(1)
(1)
1/2 1
H(X) + H(Y ) ≤ Varτ
log fX fY (Zi ) −
log fX fY (Zi )
Eτ
n i=l+1
n
n i=l+1



31/2 
31/2
1/2
1/2
2
2
1/2
(1)
+ E log fY (Y1 )
.
≤ 1/2 Varτ log fX fY (Zn ) ≤ 1/2 E log fX (X1 )
n
n
By a very similar argument to that in (19), we deduce that
n
1 X
n − l
(1)
max Eτ
sup
sup
H(X) + H(Y )
log fX fY (Zi ) −
ln
l
n i=l+1
n
Sn
k∈{k0∗ ,...,k1∗ } f ∈GdX ,dY ,φ τ ∈∪l=0

→ 0, (33)

which handles the third term in (27).
For the fourth term in (27), we simply observe that, by (19),
sup

sup

max

ln
l
Sn
k∈{k0∗ ,...,k1∗ } f ∈GdX ,dY ,φ τ ∈∪l=0

ln
|H(X) + H(Y )| → 0.
n

(34)

Finally, we turn to the first term in (27). Recalling the definition of rz,fX fY ,δ in (4), we
define the random variables
Ni,δ :=

X

1kZ (1) −Z (1) k≤r

l+1≤j≤n

j

i

.
(1)
Z
,fX fY ,δ
i

j6=i

We study the bias and the variance of Ni,δ . For δ ∈ (0, 2] and z ∈ Wn , we have that
(1)
Eτ (Ni,δ |Zi

Z
= z) ≥ (n − ln − 3)

fX fY (w) dw
Bz (rz,fX fY ,δ )

d
≥ (n − ln − 3)Vd rz,f
f ,δ fX fY (z)(1 − cn )
 X Y

Ψ(k) n − ln − 3
= δ(1 − cn )e
.
n−1

(35)

Similarly, also, for δ ∈ (0, 2] and z ∈ Wn ,
(1)
Eτ (Ni,δ |Zi

Z
= z) ≤ 2 + (n − ln − 3)

fX fY (w) dw
Bz (rz,fX fY ,δ )

≤ 2 + δ(1 + cn )e

33

Ψ(k)




n − ln − 3
.
n−1

(36)

Note that if j2 ∈
/ {j1 , τ (j1 ), τ −1 (j1 )} then for δ ∈ (0, 2] and z ∈ Wn ,
Covτ



1

(1)
1

kZj −zk≤rz,fX fY ,δ

, 1

(1)
2

kZj −zk≤rz,fX fY ,δ

(1)
|Zi



= z = 0.

Also, for j ∈
/ {i, τ (i), τ −1 (i)} we have
Varτ



1

(1)
|Zi

(1)
kZj −zk≤rz,fX fY ,δ



Z

=z ≤

fX fY (w) dw ≤
Bz (rz,fX fY ,δ )

δ(1 + cn )eΨ(k)
.
n−1

When j ∈ {i, τ (i), τ −1 (i)} we simply bound the variance above by 1/4 so that, by Cauchy–
Schwarz, when n − 1 ≥ 2(1 + cn eΨ(k) ),
(1)

Varτ (Ni,δ |Zi

= z) ≤ 3(n − ln − 3)

δ(1 + cn )eΨ(k)
+ 1.
n−1

(37)

1/2

Letting n,k := max(k −1/2 log n, cn ), there exists n0 ∈ N, depending only on dX , dY , φ and
k0∗ , such that for n ≥ n0 and all k ∈ {k0∗ , . . . , k1∗ }, we have n,k ≤ 1/2 and



 k

n,k
Ψ(k) n−ln −3
Ψ(k) n−ln −3
min (1+n,k )(1−cn )e
−k , k−2−(1−n,k )(1+cn )e
≥
.
n−1
n−1
2
We deduce from (37) and (35) that for n ≥ n0 ,


(1)
(1)
Pτ {ξi fX fY (Zi ) − 1}1Ai ≥ n,k = Pτ Ai ∩ {ρ(k),i,(1) ≥ rZ (1) ,fX fY ,1+ }
n,k
i
Z
(1)
(1)
fX fY (z)Pτ (Ni,1+n,k ≤ k|Zi = z) dz
≤ Pτ (Ni,1+n,k ≤ k, Zi ∈ Wn ) =
Wn
(1)
Varτ (Ni,1+n,k |Zi = z)
fX fY (z)
(1)
{Eτ (Ni,1+n,k |Zi = z) − k}2
Wn

Z
≤

≤

dz

(1+n,k )(1+cn )eΨ(k)
+1
n−1

2.
n −3
cn )eΨ(k) n−l
−
k
n−1

3(n − ln − 3)
(1 + n,k )(1 −

(38)

Using very similar arguments, but using (36) in place of (35), we also have that for n ≥ n0 ,
(1)
(1)
Pτ {ξi fX fY (Zi )

− 1}1Ai ≤ −n,k ≤ 

(1−n,k )(1+cn )eΨ(k)
+1
n−1
 2.
n −3
n,k )(1 + cn )eΨ(k) n−l
n−1

3(n − ln − 3)



k − 2 − (1 −

34

(39)

Now, by Markov’s inequality, for k ≥ 3, i ∈ {l + 1, . . . , n}, τ ∈ Snl and δ ∈ (0, 2],
Pτ




(1)
(1)
(1)
ξi fX fY (Zi ) ≤ δ ∩ Ai = Pτ (Ni,δ ≥ k, Zi ∈ Wn )
Z
Z
n−l−3
≤
fX fY (z)
k−2
Wn
Bz (rz,f
≤

fX fY (w) dw dz
)
X fY ,δ

n−l−3
δeΨ(k)
(1 + cn )
.
k−2
n−1

(40)

Moreover, for any i ∈ {l + 1, . . . , n}, τ ∈ Snl and t > 0, by two applications of Markov’s
inequality,
Pτ




(1)
(1)
(1)
ξi fX fY (Zi ) ≥ t ∩ Ai ≤ Pτ (Ni,t ≤ k, Zi ∈ Wn )
Z
Z
n−l−3
≤
fX fY (z)
fX fY (w) dw dz
n − l − k − 3 Wn
Bz (rz,fX fY ,t )c
Z
µ + kzkα
n − ln − 3
α−1
max{1, 2 }
fX fY (z) α
≤
dz
n − ln − k − 3
rz,fX fY ,t
Wn

α/d
n − ln − 3
α−1
2α/d α/d n − 1
max{1, 2 }2µν
Vd
≤
t−α/d .
Ψ(k)
n − ln − k − 3
e
(41)


Writing sn := log2 (n/k) log2d/α n , we deduce from (38), (39), (40) and (41) that for n ≥ n0 ,
Z ∞
 2 (1)
 (1)

1/2
(1) 
(1)
Eτ log ξi fX fY (Zi ) 1Ai =
Pτ ξi fX fY (Zi ) ≤ e−s
∩ Ai ds
0
Z sn
Z ∞
 (1)

1/2
(1)
2
+
Pτ ξi fX fY (Zi ) > es
∩ Ai ds
+ log (1 + n,k ) +
log2 (1+n,k )
Ψ(k) Z ∞

≤

e
n−l−3
(1 + cn )
k−2
n−1
+ sn

sn

e−s

1/2

ds + log2 (1 + n,k )

0

12(n − ln − 3)

(1+n,k )(1+cn )eΨ(k)
n−1
k 2 2n,k

+1


α/d
n − ln − 3
α−1
2α/d α/d n − 1
+
max{1, 2 }2µν
Vd
n − ln − k − 3
eΨ(k)

Z

∞

e−s

1/2 α/d

ds.

sn

We conclude that
sup

sup

sup

max


(1)
(1) 
max Eτ log2 ξi fX fY (Zi ) 1Ai < ∞.

n S l i=l+1,...,n
n∈N k∈{k0∗ ,...,k1∗ } f ∈GdX ,dY ,φ τ ∈∪ll=1
n

35

(42)

Finally, then, from (38), (39) and (42),
n
1 X
(1)
(1) 
(1)
(1) 
Eτ log ξi fX fY (Zi ) 1Ai ≤ max Eτ log ξi fX fY (Zi ) 1Ai
i=l+1,...,n
n i=l+1
h 
i1/2

(1)
(1) 
(1)
(1) 
≤ 2n,k + max Pτ log ξi fX fY (Zi ) 1Ai ≥ 2n,k Eτ log2 ξi fX fY (Zi ) 1Ai
i=l+1,...,n
h 
i1/2

(1)
(1) 
(1)
(1)
≤ 2n,k + max Pτ ξi fX fY (Zi ) − 1 1Ai ≥ n,k Eτ log2 ξi fX fY (Zi ) 1Ai
,
i=l+1,...,n

so
sup

sup

max

ln
l
Sn
k∈{k0∗ ,...,k1∗ } f ∈GdX ,dY ,φ τ ∈∪l=1

n
1 X
(1)
(1) 
Eτ log ξi fX fY (Zi ) 1Ai → 0.
n i=l+1

(43)

as n → ∞. The proof of claim (26) follows from (27), together with (32), (33), (34) and (43).
The final result therefore follows from (12), (13), (20), (21) and (26).

7.3

Proof from Section 5

T
T T
Proof of Theorem 6. We partition X = (X(1)
X(2)
) ∈ R(m+m)×p and, writing η (0) := /σ,

T
(b)
(b)
T (0)
T
η(2) ,
X(2) )−1 X(2)
let η (b) = (η(1) )T , (η(2) )T ∈ Rm+m , for b = 0, . . . , B. Further, let γ
b := (X(2)
(1)

T
T
b(2) /σ. Now define the events
η(2) and sb := σ
X(2) )−1 X(2)
γ
b(1) := (X(2)

n
o
Ar0 ,s0 := max(kb
γ k, kb
γ (1) k) ≤ r0 , sb ∈ [s0 , 1/s0 ] , sb(1) ∈ [s0 , 1/s0 ]
and A0 := {λmin (n−1 X T X) > 21 λmin (Σ)}. Now

P(I˘n(0) ≤ I˘n(1) ) ≤ P {I˘n(0) ≤ I˘n(1) } ∩ Ar0 ,s0 ∩ A0 + P(Ac0 ) + P(A0 ∩ Acr0 ,s0 ).
For the first term in (44), define functions h, h(1) : Rp → R by
(1)

h(b) := H(η1 − X1T b) and h(1) (b) := H(η1 − X1T b).

36

(44)

Writing R1 and R2 for remainder terms to be bounded below, we have
b m (b
b m (X1 , ηb1,(1) ), . . . , (Xm , ηbm,(1) )
I˘n(0) − I˘n(1) = H
η1,(1) , . . . , ηbm,(1) ) − H




(1)
(1)
b m (b
b m (X1 , ηb(1) ), . . . , (Xm , ηb(1) )
−H
η1,(1) , . . . , ηbm,(1) ) + H
1,(1)
m,(1)
= I(X1 ; 1 ) + h(b
γ ) − h(0) − h(1) (b
γ (1) ) + h(1) (0) + R1
= I(X1 ; 1 ) + R1 + R2 .

(45)

To bound R1 on the event Ar0 ,s0 , we first observe that for fixed γ, γ (1) ∈ Rp and s, s(1) > 0,
η1 − X1T γ
H
s





 (1)



(1)
η1 − X1T γ
η1 − X1T γ (1)
η1 − X1T γ (1)
− H X1 ,
−H
+ H X1 ,
s
s(1)
s(1)

 (1)


η

 η (1)

η1 − X1T γ (1)
η1 − X1T γ
1
1
−H
X1 − H
+
H
X
=H
1
s
s
s(1)
s(1)


= H(η1 ) − H(η1 |X1 ) + h(γ) − h(0) − h(1) (γ (1) ) + h(1) (0)
= I(X1 ; 1 ) + h(γ) − h(0) − h(1) (γ (1) ) + h(1) (0).

(46)

Now, for a density g on Rd , define the variance functional
Z
v(g) :=

g(x){log g(x) + H(g)}2 dx.

Rd

If (b
γ, γ
b(1) , sb, sb(1) )T are such that Ar0 ,s0 holds, then conditional on (b
γ, γ
b(1) , sb, sb(1) ), we have
(1)

(1)

γ
b,b
s
γ
b ,b
s
fηbγb,bs ∈ F1,θ1 , fX,b
η ∈ Fp+1,θ2 , fηb(1)

(1)

(1)

γ
b ,b
s
∈ F1,θ1 and fX,b
η (1)

∈ Fp+1,θ2 . It follows by (46),

Berrett, Samworth and Yuan (2017, Theorem 1 and Lemma 11(i)) that
lim sup n1/2
n→∞

sup

sup Ef (|R1 |1Ar0 ,s0 )

∗
kη ∈{k0∗ ,...,kη∗ } f ∈Fp+1,ω
∗
∗
k∈{k0 ,...,k }

≤ lim sup n1/2
n→∞

sup

∗
kη ∈{k0∗ ,...,kη∗ } f ∈Fp+1,ω
k∈{k0∗ ,...,k∗ }

≤ 4 sup v(g) + 4
g∈F1,θ1

sup

sup



E(R12 1Ar0 ,s0 )

v(g) < ∞.

1/2

(47)

g∈Fp+1,θ2

Now we turn to R2 , and study the continuity of the functions h and h(1) , following the
approach taken in Proposition 1 of Polyanskiy and Wu (2016). Write a1 ∈ A for the fifth

37

component of θ1 , and for x ∈ R with fη (x) > 0, let
ra1 (x) :=

1
.
8a1 fη (x)

Then, for |y − x| ≤ ra1 (x) we have by Berrett, Samworth and Yuan (2017, Lemma 12) that
|fη (y) − fη (x)| ≤

15
15
fη (x)a1 (fη (x))|y − x| ≤ fη (x).
7
56

Hence
| log fη (y) − log fη (x)| ≤

120
a1 (fη (x))|y − x|.
41

When |y − x| > ra1 (x) we may simply write
| log fη (y) − log fη (x)| ≤ 8{| log fη (y)| + | log fη (x)|}a1 (fη (x))|y − x|.
Combining these two equations we now have that, for any x, y such that fη (x) > 0,
| log fη (y) − log fη (x)| ≤ 8{1 + | log fη (y)| + | log fη (x)|}a1 (fη (x))|y − x|.
By an application of the generalised Hölder inequality and Cauchy–Schwarz, we conclude
that
fη (η1 − X1T γ)
fη (η1 )
 2
2 1/2 
 1/4
1/2  
≤ 8 Ea1 (fη (η1 ))
E 1 + | log fη (η1 − X1T γ)| + | log fη (η1 )|
E(|X1T γ|4

1/2 
1/2  
1/4
≤ 16kγk Ea21 (fη (η1 ))
E 1 + log2 fη (η1 − X1T γ) + log2 fη (η1 )
E(kX1 k4 )
.

E log

(48)
We also obtain a similar bound on the quantity E log

fηbγ,1 (η1 −X1T γ)
fηbγ,1 (η1 )
d

when kγk ≤ r0 . Moreover,

for any random vectors U, V with densities fU , fV on R satisfying H(U ), H(V ) ∈ R, and
writing U := {fU > 0}, V := {fV > 0}, we have by the non-negativity of Kullback–Leibler
divergence that
Z
H(U ) − H(V ) =

Z
fV log fV −

V

Z
fU log fV −

U

fU log
U

38

fU
fV (V )
≤ E log
.
fV
fV (U )

Combining this with the corresponding bound for H(V ) − H(U ), we obtain that


fU (U )
fV (V )
.
|H(U ) − H(V )| ≤ max E log
, E log
fU (V )
fV (U )
Since fη , fηbγ,1 ∈ F1,θ1 when kγk ≤ r0 , we may apply (48), the second part of Berrett, Samworth and Yuan (2017, Proposition 9), (5), (6) and the fact that α2 ≥ 4 to deduce that
sup

sup

∗
f ∈Fp+1,ω

γ∈B0◦ (r0 )

|h(γ) − h(0)|
kγk



fηbγ,1 (η1 − X1T γ)
1
fη (η1 − X1T γ)
≤ sup
sup
max E log
, E log
< ∞.
∗
fη (η1 )
fηbγ,1 (η1 )
f ∈Fp+1,ω
γ∈B0◦ (r0 ) kγk
(49)
Similarly,
|h(1) (γ) − h(1) (0)|
< ∞.
∗
kγk
f ∈Fp+1,ω
γ∈B0◦ (r0 )
sup

sup

(50)

We now study the convergence of γ
b and γ
b(1) to 0. By definition, the unique minimiser of the
function
R(β) := E{(Y1 − X1T β)2 } = E(21 ) − 2(β − β0 )T E(1 X1 ) + (β − β0 )T Σ(β − β0 )
is given by β0 , and R(β0 ) = E(21 ). Now, for λ > 0,

R β0 + λE(1 X1 ) = E(21 ) − 2λkE(1 X1 )k2 + λ2 E(1 X1 )T ΣE(1 X1 ).
If E(1 X1 ) 6= 0 then taking λ =

kE(1 X1 )k2
E(1 X1 )T ΣE(1 X1 )


we have R β0 + λE(1 X1 ) < E(21 ), a

contradition. Hence E(1 X1 ) = 0. Moreover,
p
p
n
n
n
 X

X
1 X
1 X
1 X
1/2
E
Xi ηi ≤
E
Xij ηi ≤
Var
Xij ηi
m i=m+1
m j=1 i=m+1
m j=1
i=m+1
p
1 X
21/2 p
4 1/4
= 1/2
Var1/2 (X1j η1 ) ≤ 1/2 E(η14 )1/4 max E(X1j
) .
j=1,...,p
m j=1
n

39

(51)

It follows that
n
o

21/2 p
T
T (0)
4 1/4
E kb
γ k1A0 = E (X(2)
X(2) )−1 X(2)
η(2) 1A0 ≤
E(η14 )1/4 max E(X1j
) .
j=1,...,p
λ0 n1/2

(52)


Similar arguments yield the same bound for E kb
γ (1) k1A0 . Hence, from (49), (50) and (52),
we deduce that
lim sup n1/2
n→∞

sup

sup Ef (|R2 |1Ar0 ,s0 ∩A0 ) < ∞.

(53)

∗
kη ∈{k0∗ ,...,kη∗ } f ∈Fp+1,ω
k∈{k0∗ ,...,k∗ }

We now bound P(Ac0 ). By the Hoffman–Wielandt inequality, we have that
T
T
{λmin (m−1 X(2)
X(2) ) − λmin (Σ)}2 ≤ km−1 X(2)
X(2) − Σk2F .

Thus
P(Ac0 )



−1

≤ P km

T
X(2)
X(2)

p

X
1
4
− ΣkF ≥ λmin (Σ) ≤
Var(X1j X1l )
2
mλ2min (Σ) j,l=1

≤

8p2
4
max E(X1j
).
nλ2min (Σ) j=1,...,p

(54)

Finally, we bound P(A0 ∩ Acr0 ,s0 ). By Markov’s inequality and (52),


1
21/2 p
4 1/4
P {kb
γ k ≥ r0 } ∩ A0 ≤ E kb
γ k1A0 ≤
E(η14 )1/4 max E(X1j
) .
j=1,...,p
r0
r0 λ0 n1/2

(55)


The same bound also holds for P {kb
γ (1) k ≥ r0 } ∩ A0 . Furthermore, writing P(2) :=
T
T
X(2) (X(2)
X(2) )−1 X(2)
, note that
n
1
1 X 2
1 (0)T
1 (0) 2
(0) 2
(0)
kη(2) k − 1 − kP(2) η(2) k ≤
(ηi − 1) + η(2) P(2) η(2) ,
|b
s − 1| =
m
m
m i=m+1
m
2

40

so by Chebychev’s inequality, Markov’s inequality and (51), for any δ > 0

P {|b
s2 − 1| > δ} ∩ A0





n
 λ δ 1/2 
δ
1 T (0)
1 X 2
0
∩ A0 + P
X η
>
∩ A0
≤P
(η − 1) >
m i=m+1 i
2
m (2) (2)
2
≤

8
2p
4 1/4
) .
E(η14 ) +
E(η14 )1/4 max E(X1j
1/2 1/2
2
j=1,...,p
1/2
nδ
δ λ0 n

(56)


The same bound holds for P {|(b
s(1) )2 − 1| > δ} ∩ A0 . We conclude from Markov’s inequality, (44), (45), (47), (53), (54), (55) and (56) that
sup

sup

sup

∗ k ∈{k ∗ ,...,k ∗ } f ∈F ∗
Bn ≥Bn
η
η
0
1,p+1 :I(f )≥bn
k∈{k0∗ ,...,k∗ }

≤

Pf (I˘n ≤ C̆q(n),Bn )

1
+
sup
sup
Pf (I˘n(1) ≥ I˘n(0) )
∗
+ 1) kη ∈{k0∗ ,...,kη∗ } f ∈F1,p+1
:I(f )≥bn

q(Bn∗

k∈{k0∗ ,...,k∗ }

n

1
˘(1) ≥ I˘(0) } ∩ Ar0 ,s0 ∩ A0
+
sup
sup
P
{
I
≤
f
n
n
∗
q(Bn∗ + 1) kη ∈{k0∗ ,...,kη∗ } f ∈F1,p+1
:I(f )≥bn
k∈{k0∗ ,...,k∗ }

o
+ Pf (A0 ∩ Acr0 ,s0 ) + Pf (Ac0 )
n1
1
≤
+
sup
Ef (|R1 + R2 |1Ar0 ,s0 )
sup
∗
q(Bn∗ + 1) kη ∈{k0∗ ,...,kη∗ } f ∈F1,p+1
:I(f )≥bn bn
k∈{k0∗ ,...,k∗ }

o
+ Pf (A0 ∩ Acr0 ,s0 ) + Pf (Ac0 ) → 0,
as required.
Acknowledgements: Both authors are supported by an EPSRC Programme grant. The
first author was supported by a PhD scholarship from the SIMS fund; the second author is
supported by an EPSRC Fellowship and a grant from the Leverhulme Trust.

References
Albert, M., Bouret, Y., Fromont, M. and Reynaud-Bouret, P. (2015) Bootstrap and permutation tests of independence for point processes. Ann. Statist., 43, 2537–2564.

41

Bach, F. R. and Jordan, M. I. (2002) Kernel independent component analysis. J. Mach.
Learn. Res., 3, 1–48.
Bergsma, W. and Dassios, A. (2014) A consistent test of independence based on a sign
covariance related to Kendall’s tau. Bernoulli, 20, 1006–1028.
Berrett, T. B., Grose, D. J. and Samworth, R. J. (2017) IndepTest:
parametric

independence

tests

based

on

entropy

estimation.

non-

Available

at

https://cran.r-project.org/web/packages/IndepTest/index.html.
Berrett, T. B., Samworth, R. J. and Yuan, M. (2017) Efficient multivariate entropy estimation via k-nearest neighbour distances. https://arxiv.org/abs/1606.00304v3.
Biau, G. and Devroye, L. (2015) Lectures on the Nearest Neighbor Method. Springer, New
York.
Comon, P. (1994) Independent component analysis, a new concept?. Signal Process., 36,
287–314.
Dalgaard, P. (2002) Introductory Statistics with R. Springer-Verlag, New York.
Dobson, A. J. (2002) An Introduction to Generalized Linear Models. Chapman & Hall,
London.
Einmahl, J. H. J. and van Keilegom, I. (2008) Tests for independence in nonparametric
regression. Statistica Sinica, 18, 601–615.
Fan, J., Feng, Y. and Xia, L. (2017) A projection based conditional dependence measure with applications to high-dimensional undirected graphical models. Available at
arXiv:1501.01617.
Fan, Y., Lafaye de Micheaux, P., Penev, S. and Salopek, D. (2017) Multivariate nonparametric test of independence. J. Multivariate Anal., 153, 189–210.
Gibbs, A. L. and Su, F. E. (2002) On choosing and bounding probability metrics. Int. Statist.
Review, 70, 419–435.
Gretton A., Bousquet O., Smola A. and Schölkopf B. (2005) Measuring Statistical Dependence with Hilbert-Schmidt Norms. Algorithmic Learning Theory, 63–77.

42

Gretton, A. and Györfi, L. (2010) Consistent nonparametric tests of independence. J. Mach.
Learn. Res., 11, 1391–1423.
Heller, R., Heller, Y., Kaufman, S., Brill, B. and Gorfine, M. (2016) Consistent distributionfree K-sample and independence tests for univariate random variables. J. Mach. Learn.
Res., 17, 1–54.
Hoeffding, W. (1948) A non-parametric test of independence. Ann. Math. Statist., 19, 546–
557.
Hofert, M., Kojadinovic, I., Mächler, M. and Yan, J. (2017) copula:
variate Dependence with Copulas.

Multi-

R Package version 0.999-18. Available from

https://cran.r-project.org/web/packages/copula/index.html.
Jitkrittum, W., Szabó, Z. and Gretton, A. (2016) An adaptive test of independence with
analytic kernel embeddings. Available at arXiv:1610.04782.
Joe, H. (1989) Relative entropy measures of multivariate dependence. J. Amer. Statist.
Assoc., 84, 157–164.
Josse, J. and Holmes, S. (2016) Measuring multivariate association and beyond. Statist.
Surveys, 10, 132–167.
Kojadinovic, I. and Holmes, M. (2009) Tests of independence among continuous random
vectors based on Cramér–von Mises functionals of the empirical copula process. J. Multivariate Anal., 100, 1137–1154.
Kozachenko, L. F. and Leonenko, N. N. (1987) Sample estimate of the entropy of a random
vector. Probl. Inform. Transm., 23, 95–101.
Lauritzen, S. L. (1996) Graphical Models. Oxford University Press, Oxford.
Mari, D. D. and Kotz, S. (2001) Correlation and Dependence. Imperial College Press, London.
Miller, E. G. and Fisher, J. W. (2003) ICA using spacings estimates of entropy. J. Mach.
Learn. Res., 4, 1271–1295.
Müller, U. U., Schick, A. and Wefelmeyer, W. (2012) Estimating the error distribution
function in semiparametric additive regression models. J. Stat. Plan. Inference, 142,
552–566.
43

Neumeyer, N. (2009) Testing independence in nonparametric regression. J. Multivariate
Anal., 100, 1551–1566.
Neumeyer, N. and Van Keilegom, I. (2010) Estimating the error distribution in nonparametric multiple regression with applications to model testing. J. Multivariate Anal., 101,
1067–1078.
Nguyen, D. and Eisenstein, J. (2017) A kernel independence test for geographical language
variation. Comput. Ling., to appear.
Pearl, J. (2009) Causality. Cambridge University Press, Cambridge.
Pearson, K. (1920) Notes on the history of correlation. Biometrika, 13, 25–45.
Peixoto, J. L. (1990) A property of well-formulated polynomial regression models. The
American Statistician, 44, 26–30.
Pfister, N., Bühlmann, P., Schölkopf, B. and Peters, J. (2017) Kernel-based tests for joint
independence. J. Roy. Statist. Soc., Ser. B, to appear.
Pfister,

N.

and

Peters,

J.

(2017)

Schmidt Independence Criterion.

dHSIC: Independence Testing via Hilbert
R

Package

version

2.0.

Available

at

https://cran.r-project.org/web/packages/dHSIC/index.html.
Polyanskiy, Y. and Wu, Y. (2016) Wasserstein continuity of entropy and outer bounds for
interference channels. IEEE Trans. Inf. Theory, 62, 3992–4002.
Rizzo, M. L. and Szekely, G. J. (2017) energy:
Inference via the Energy of Data.

E-Statistics:

Multivariate

R Package version 1.7-2. Available from

https://cran.r-project.org/web/packages/energy/index.html.
Samworth, R. J. and Yuan, M. (2012) Independent component analysis via nonparametric
maximum likelihood estimation. Ann. Statist., 40, 2973–3002.
Schweizer, B. and Wolff, E. F. (1981) On nonparametric measures of dependence for random
variables. Ann. Statist., 9, 879–885.
Sejdinovic, D., Sriperumbudur, B., Gretton, A. and Fukumizu, K. (2013) Equivalence of
distance-based and RKHS-based statistics in hypothesis testing. Ann. Statist., 41, 2263–
2291.
44

Sen, A. and Sen, B. (2014) Testing independence and goodness-of-fit in linear models.
Biometrika, 101, 927–942.
Shah, R. D. and Bühlmann, P. (2017) Goodness of fit tests for high-dimensional linear
models. J. Roy. Statist. Soc., Ser. B, to appear.
Song, L., Smola, A., Gretton, A., Bedo, J. and Borgwardt, K. (2012) Feature selection via
dependence maximization. J. Mach. Learn. Res., 13, 1393–1434.
Steuer, R., Kurths, J., Daub, C. O., Weise, J. and Selbig, J. (2002) The mutual information:
detecting and evaluating dependencies between variables. Bioinformatics, 18, 231–240.
Stigler, S. M. (1989) Francis Galton’s account of the invention of correlation. Stat. Sci., 4,
73–86.
Su, L. and White, H. (2008) A nonparametric Hellinger metric test for conditional independence. Econometric Theory, 24, 829–864.
Székely, G. J., Rizzo, M. L. and Bakirov, N. K. (2007) Measuring and testing dependence
by correlation of distances. Ann. Statist., 35, 2769–2794.
Székely, G. J. and Rizzo, M. L. (2013) The distance correlation t-test of independence in
high dimension. J. Multivariate Anal., 117, 193–213.
Torkkola, K. (2003) Feature extraction by non-parametric mutual information maximization.
J. Mach. Learn. Res., 3, 1415–1438.
Vinh, N. X., Epps, J. and Bailey, J. (2010) Information theoretic measures for clusterings
comparison: variants, properties, normalisation and correction for chance. J. Mach. Learn.
Res., 11, 2837–2854.
Weihs, L., Drton, M. and Leung, D. (2016) Efficient computation of the Bergsma–Dassios
sign covariance. Comput. Stat., 31, 315–328.
Weihs,
ances:

L.,

Drton,

M.

and

Meinshausen,

N.

(2017)

Symmetric

rank

covari-

a generalised framework for nonparametric measures of dependence.

https://arxiv.org/abs/1708.05653.
Wu, E. H. C., Yu, P. L. H. and Li, W. K. (2009) A smoothed bootstrap test for independence
based on mutual information. Comput. Stat. Data Anal., 53, 2524–2536.
45

Yao, S., Zhang, X. and Shao, X. (2017) Testing mutual independence in high dimension via
distance covariance. J. Roy. Statist. Soc., Set. B, to appear.
Zhang, K., Peters, J., Janzing, D. and Schölkopf, B. (2011) Kernel-based conditional independence test and application in causal discovery. https://arxiv.org/abs/1202.3775.
Zhang, Q., Filippi, S., Gretton, A. and Sejdinovic, D. (2017) Large-scale kernel methods for
independence testing. Stat. Comput., 27, 1–18.

46

