Measuring and testing the dependency between two random variables is a fundamental problem in statistics. The classical Pearson’s correlation coefficient can be inaccurate and even misleading when nonlinear dependency exists [6]. propose the novel measure–distance correlation–which is exactly zero if and only if two random variables are independent. A limitation is that if the distance correlation is implemented based on its original definition, the corresponding computational complexity can be as high as O(n ²), which is not desirable when n is large.

We review the definition of the distance correlation in [6]. Let us consider two random variables X ∈ R p , Y ∈ R q , p ≥ 1 , q ≥ 1 . Let the complex-valued functions ϕ _X,Y (⋅), ϕ _X (⋅), and ϕ _Y (⋅) be the characteristic functions of the joint density of X and Y, the density of X, and the density of Y, respectively. For any function ϕ, we denote | ϕ | 2 = ϕ ϕ ̄ , where ϕ ̄ is the conjugate of ϕ; in words, |ϕ| is the magnitude of ϕ at a particular point. For vectors, let us use |⋅| to denote the Euclidean norm. In [6], the definition of distance covariance between random variables X and Y is V 2 ( X , Y ) = ∫ R p + q | ϕ X , Y ( t , s ) − ϕ X ( t ) ϕ Y ( s ) | 2 c p c q | t | p + 1 | s | q + 1 d t d s , (2.1) where two constants c _p and c _q have been defined at the end of Section 1. The distance correlation is defined as R 2 ( X , Y ) = V 2 ( X , Y ) V 2 ( X , X ) V 2 ( Y , Y ) . The following property has been established in the aforementioned paper.

Theorem 2.1

Suppose X ∈ R p , p ≥ 1 and Y ∈ R q , q ≥ 1 are two random variables, the following statements are equivalent:

1) X is independent of Y;

Given sample (X ₁, Y ₁), …, (X _n , Y _n ), we can estimate the distance covariance by replacing the population characteristic function with the sample characteristic function: for i = − 1 , t ∈ R p , s ∈ R q , we define ϕ ̂ X ( t ) = 1 n ∑ j = 1 n e i X j t t , ϕ ̂ Y ( s ) = 1 n ∑ j = 1 n e i Y j t s ,  and  ϕ ̂ X , Y ( t , s ) = 1 n ∑ j = 1 n e i X j t t + i Y j t s . Consequently one can have the following estimator for V 2 ( X , Y ) : V n 2 ( X , Y ) = ∫ R p + q | ϕ ̂ X , Y ( t , s ) − ϕ ̂ X ( t ) ϕ ̂ Y ( s ) | 2 c p c q | t | p + 1 | s | q + 1 d t ⋅ d s . (2.2) Note that the above formula is convenient to define a quantity, however, is not convenient for computation, due to the integration on the right-hand side. In the literature, other estimates have been introduced and will be presented in the following.

The paper [10] gives an equivalent definition for the distance covariance between random variables X and Y: V 2 ( X , Y ) = E [ d ( X , X ′ ) d ( Y , Y ′ ) ] = E [ | X − X ′ ‖ Y − Y ′ | ] − 2 E [ | X − X ′ ‖ Y − Y ″ | ] + E [ | X − X ′ | ] E [ | Y − Y ′ | ] , (2.3) where the double centered distance d(⋅, ⋅) is defined as d ( X , X ′ ) = | X − X ′ | − E X [ | X − X ′ | ] − E X ′ [ | X − X ′ | ] + E [ | X − X ′ | ] , where E X , E X ′ and E are expectations over X, X′ and (X, X′), respectively.

Motivated by the above definition, one can give an unbiased estimator for V 2 ( X , Y ) . The following notations will be utilized: for 1 ≤ i, j ≤ n, a i j = | X i − X j | , b i j = | Y i − Y j | , a i ⋅ = ∑ l = 1 n a i l , b i ⋅ = ∑ l = 1 n b i l , a ⋅ ⋅ = ∑ k , l = 1 n a k l , and b ⋅ ⋅ = ∑ k , l = 1 n b k l . (2.4) It has been proven [8, 28] that Ω n ( X , Y ) = 1 n ( n − 3 ) ∑ i ≠ j a i j b i j − 2 n ( n − 2 ) ( n − 3 ) ∑ i = 1 n a i ⋅ b i ⋅ + a ⋅ ⋅ b ⋅ ⋅ n ( n − 1 ) ( n − 2 ) ( n − 3 ) (2.5) is an unbiased estimator of V 2 ( X , Y ) . In addition, a fast algorithm has been proposed [8] for the aforementioned sample distance covariance in the univariate cases with complexity order O(n log  n) and storage O(n). We list the result below for reference purpose.

Theorem 2.2

(Theorem 3.2 & Corollary 4.1 in [8]). Suppose X ₁, …, X _n and Y 1 , … , Y n ∈ R . The unbiased estimator Ω _n defined in (2.5) can be computed by an O(n log  n) algorithm.

In addition, as a byproduct, the following result is established in the same paper.

Corollary 2.3

The quantity a ⋅ ⋅ b ⋅ ⋅ n ( n − 1 ) ( n − 2 ) ( n − 3 ) = ∑ k , l = 1 n a k l ∑ k , l = 1 n b k l n ( n − 1 ) ( n − 2 ) ( n − 3 ) can be computed by an O(n log  n) algorithm.

We will use the above result in our test of independence. However, as far as we know, in the multivariate cases, there does not exist any work on the fast algorithm of the order of complexity O(n log  n). This paper will fill in this gap by introducing an order O(nK log  n) complexity algorithm in multivariate cases.

Ref. [6] proposed an independence test using the distance covariance. We summarize it below as a theorem, which serves as a benchmark. Our test will be aligned with the following one, except that we introduced a new test statistic, which can be more efficiently computed, and it has comparable asymptotic properties with the test statistic that is used below.

Theorem 2.4

([6], Theorem 6). For potentially multivariate random variables X and Y , a prescribed level α _s , and sample size n, one rejects the independence if and only if n V n 2 ( X , Y ) S 2 > ( Φ − 1 ( 1 − α s / 2 ) ) 2 , where V n 2 ( X , Y ) has been defined in (2.2), Φ(⋅) denote the cumulative distribution function of the standard normal distribution and S 2 = 1 n 4 ∑ i , j = 1 n | X i − X j | ∑ i , j = 1 n | Y i − Y j | . Moreover, let α(X, Y, n) denote the achieved significance level of the above test. If E [ | X | + | Y | ] < ∞ , then for all 0 < α _s < 0.215, one can show the following: lim n → ∞ α ( X , Y , n ) ≤ α s ,  and  sup X , Y lim n → ∞ α ( X , Y , n ) : V ( X , Y ) = 0 = α s . Note that the quantity V n 2 ( X , Y ) that is used above as in [6] differs from the one that will be used in our proposed method. As mentioned, we use the above as an illustration for distance-based tests of independence, as well as the theoretical/asymptotic properties that such a test can achieve.

We consider how to use a fast algorithm for univariate random variables to compute or approximate the sample distance covariance of random vectors. The main idea works as follows: first, projecting the multivariate observations on some random directions; then, using the fast algorithm to compute the distance covariance of the projections; at last, averaging distance covariances from different projecting directions.

More specifically, our estimator can be computed as follows. For potentially multivariate X 1 , … , X n ∈ R p and Y 1 , … , Y n ∈ R q , let K be a predetermined number of iterations, we do:

1) For each k (1 ≤ k ≤ K), randomly generate u _k and v _k from Uniform ( S p − 1 ) and Uniform ( S q − 1 ) , respectively. Here S p − 1 and S q − 1 are the unit spheres in R p and R q , respectively.

Note that samples u k t X and v k t Y are now univariate.

3) Utilize the fast (i.e., order O(n log  n)) algorithm that was mentioned in Theorem 2.2 to compute for the unbiased estimator in Eq. 2.5 with respect to u k t X and v k t Y . Formally, we denote

See Algorithm 1 in the Supplementary Appendix for a stand-alone presentation of the above method. In the light of Theorem 2.2, we can handily declare the following.

Theorem 3.1

For potentially multivariate X 1 , … , X n ∈ R p and Y 1 , … , Y n ∈ R q , the order of computational complexity of computing the aforementioned Ω ̄ n is O(Kn log  n) with storage O( max{n, K}), where K is the number of random projections.

The proof of the above theorem is omitted because it is straightforward from Theorem 2.2. The statistical properties of the proposed estimator Ω ̄ n will be studied in the subsequent section (specifically in Section 4.4).

By a later result (cf. Theorem 4.19), we can apply Ω ̄ n in the independence tests. The corresponding asymptotic distribution of the test statistic Ω ̄ n can be approximated by a Gamma(α, β) distribution with α and β given in Eq. 4.7. We can compute the significance level of the test statistic by permutation and conduct the independence test accordingly. Recall that we have potentially multivariate X 1 , … , X n ∈ R p and Y 1 , … , Y n ∈ R q . Recall that K denotes the number of Monte Carlo iterations in our previous algorithm. Let α _s denote the prescribed significance level of the independence test. Let L denote the number of random permutations that we will adopt. We would like to test the null hypothesis H 0 —X and Y are independent—against its alternative. Recall Ω ̄ n is our proposed estimator in Eq. 3.1. The following algorithm describes a test of independence, which applies permutation to generate a threshold.

1) For each ℓ, 1 ≤ ℓ ≤ L, generate a random permutation of Y: Y ∗ , ℓ = ( Y 1 ∗ , … Y n ∗ ) ;

See Algorithm 2 in the Supplementary Appendix for a stand-alone description.

It is notified that one can use the approximate asymptotic distribution information to estimate a threshold in the independence test. The following describes such an approach. Recall that random vectors X 1 , … , X n ∈ R p and Y 1 , … , Y n ∈ R q , number of random projections K, and a prescribed significance level α _s have been mentioned earlier.

1) For each k (1 ≤ k ≤ K), randomly generate u _k and v _k from uniform ( S p − 1 ) and uniform ( S q − 1 ) , respectively.

The above procedure is motivated by the observation that the asymptotic distribution of the test statistic n Ω ̄ n can be approximated by a Gamma distribution, whose parameters can be estimated by Eq. 3.2 and Eq. 3.3. A stand-alone description of the above procedure can be found in Algorithm 3 in the Supplementary Appendix.

In this section, we will study some properties of distance covariances of randomly projected random vectors. We begin with a necessary and sufficient condition of independence.

Lemma 4.1

Suppose u and v are points on the hyper-spheres: u ∈ S p − 1 = { u ∈ R p : | u | = 1 } and v ∈ S q − 1 . We have random vectors X ∈ R p a n d Y ∈ R q a r e i n d e p e n d e n t if and only if V 2 ( u t X , v t Y ) = 0 , f o r a n y u ∈ S p − 1 , v ∈ S q − 1 .

The proof is relatively straightforward. We relegate a formal proof to the appendix. This lemma indicates that the independence is somewhat preserved under projections. The main contribution of the above result is to motivate us to think of using random projection, to reduce the multivariate random vectors into univariate random variables. As mentioned earlier, there exist fast algorithms of distance-based methods for univariate random variables.

The following result allows us to regard the distance covariance of random vectors of any dimension as an integral of distance covariance of univariate random variables, which are the projections of the aforementioned random vectors. The formulas in the following lemma provide the foundation for our proposed method: the distance covariances in the multivariate cases can be written as integrations of distance covariances in the univariate cases. our proposed method essentially adopts the principle of Monte Carlo to approximate such integrals. We again relegate the proof to the Supplementary Appendix.

Lemma 4.2

Suppose u and v are points on unit hyper-spheres: u ∈ S p − 1 = { u ∈ R p : | u | = 1 } and v ∈ S q − 1 . Let μ and ν denote the uniform probability measure on S p − 1 and S q − 1 , respectively. Then, we have for random vectors X ∈ R p and Y ∈ R q , V 2 ( X , Y ) = C p C q ∫ S p − 1 × S q − 1 V 2 ( u t X , v t Y ) d μ ( u ) d ν ( v ) , where C _p and C _q are two constants that are defined at the end of Section 1. Moreover, a similar result holds for the sample distance covariance: V n 2 ( X , Y ) = C p C q ∫ S p − 1 × S q − 1 V n 2 ( u t X , v t Y ) d μ ( u ) d ν ( v ) .

Besides the integral equations in the above lemma, we can also establish the following result for the unbiased estimator. Such a result provides the direct foundation of our proposed method. Recall that Ω _n , which is in Eq. 2.5, is an unbiased estimator of the distance covariance V 2 ( X , Y ) . A proof is provided in the Supplementary Appendix.

Lemma 4.3

Suppose u and v are points on the hyper-spheres: u ∈ S p − 1 = { u ∈ R p : | u | = 1 } and v ∈ S q − 1 . Let μ and ν denote the measure corresponding to the uniform densities on the surfaces S p − 1 and S q − 1 , respectively. Then, we have Ω n ( X , Y ) = C p C q ∫ S p − 1 × S q − 1 Ω n ( u t X , v t Y ) d μ ( u ) d ν ( v ) , where C _p and C _q are constants that were mentioned at the end of Section 1.

From the above lemma, recalling the design of our proposed estimator Ω ̄ n as in Eq. 3.1, it is straightforward to see that the proposed estimator Ω ̄ n is an unbiased estimator of Ω _n (X, Y). For completeness, we state the following without a proof.

Corollary 4.4

The proposed estimator Ω ̄ n in Eq. 3.1) is an unbiased estimator of the estimator Ω _n (X, Y) that was defined in Eq. 2.5 .

Note that the estimator Ω ̄ n in Eq. 3.1 evidently depends on the number of random projections K. Recall that to emphasize such a dependency, we sometimes use a notation Ω ̄ n , K ≜ Ω ̄ n . The following concentration inequality shows the speed that Ω ̄ n , K can converge to Ω _n as K → ∞.

Lemma 4.5

Suppose E[|X|²] < ∞ and E[|Y|²] < ∞. For any ϵ > 0, we have P | Ω ̄ n , K − Ω n | > ϵ ≤ 2 ⁡ exp − C K ϵ 2 T r [ Σ X ] T r [ Σ Y ] , where Σ _X and Σ _Y are the covariance matrices of X and Y, respectively, Tr[Σ _X ] and Tr[Σ _Y ] are their matrix traces, and C = 2 25 C p 2 C q 2 is a constant.

The proof is a relatively standard application of Hoeffding’s inequality [11], which has been relegated to the appendix. The above lemma essentially indicates that the quantity | Ω ̄ n , K − Ω n | converges to zero at a rate no worse than O ( 1 / K ) .

The asymptotic behavior of the sample distance covariance Ω _n in Eq. 2.5 of this paper, has been studied in many places, seeing [5,8,10,12]. We found that it is still worthwhile to present them here, as we will use them to establish the statistical properties of our proposed estimator. The asymptotic distributions of Ω _n will be studied under two situations: 1) a general case and 2) when X and Y are assumed to be independent. We will see that the asymptotic distributions are different in these two situations.

It has been showed in ([8], Theorem 3.2) that Ω _n is a U-statistic. In the following, we state the result without formal proof. We will need the following function, denoted by h ₄, which takes four pairs of input variables: h 4 ( ( X 1 , Y 1 ) , ( X 2 , Y 2 ) , ( X 3 , Y 3 ) , ( X 4 , Y 4 ) ) = 1 4 ∑ 1 ≤ i , j ≤ 4 , i ≠ j | X i − X j | | Y i − Y j | − 1 4 ∑ i = 1 4 ∑ 1 ≤ j ≤ 4 , j ≠ i | X i − X j | ∑ 1 ≤ j ≤ 4 , j ≠ i | Y i − Y j | + 1 24 ∑ 1 ≤ i , j ≤ 4 , i ≠ j | X i − X j | ∑ 1 ≤ i , j ≤ 4 , i ≠ j | Y i − Y j | . (4.1)

Note that the definition of h ₄ coincides with Ω _n when the number of observations n = 4.

Lemma 4.6

(U-statistics). Let Ψ₄ denote all distinct 4-subset of {1, …, n} and let us define X _ψ = {X _i |i ∈ ψ} and Y _ψ = {Y _i |i ∈ ψ}, then Ω _n is a U-statistic and can be expressed as Ω n = n 4 − 1 ∑ ψ ∈ Ψ 4 h 4 X ψ , Y ψ .

From the literature of the U-statistics, we know that the following quantities play critical roles. We state them here: h 1 ( ( X 1 , Y 1 ) ) = E 2,3,4 [ h 4 ( ( X 1 , Y 1 ) , ( X 2 , Y 2 ) , ( X 3 , Y 3 ) , ( X 4 , Y 4 ) ) ] , h 2 ( ( X 1 , Y 1 ) , ( X 2 , Y 2 ) ) = E 3,4 [ h 4 ( ( X 1 , Y 1 ) , ( X 2 , Y 2 ) , ( X 3 , Y 3 ) , ( X 4 , Y 4 ) ) ] , h 3 ( ( X 1 , Y 1 ) , ( X 2 , Y 2 ) , ( X 3 , Y 3 ) ) = E 4 [ h 4 ( ( X 1 , Y 1 ) , ( X 2 , Y 2 ) , ( X 3 , Y 3 ) , ( X 4 , Y 4 ) ) ] , where E 2,3,4 stands for taking expectation over (X ₂, Y ₂), (X ₃, Y ₃) and (X ₄, Y ₄); E 3,4 stands for taking expectation over (X ₃, Y ₃) and (X ₄, Y ₄); and E 4 stands for taking expectation over (X ₄, Y ₄); respectively.

One immediate application of the above notations is the following result, which quantifies the variance of Ω _n . Since the formula is a known result, seeing ([13], Chapter 5.2.1 Lemma A), we state it without proof.

Lemma 4.7

(Variance of the U-statistic). The variance of Ω _n could be written as Var ( Ω n ) = n 4 − 1 ∑ l = 1 4 4 l n − 4 4 − l Var ( h l ) = 16 n Var ( h 1 ) + 240 n 2 Var ( h 1 ) + 72 n 2 Var ( h 2 ) + O 1 n 3 , where O(⋅) is the standard big O notation in mathematics.

From the above lemma, we can see that Var(h ₁) and Var(h ₂) play indispensable roles in determining the variance of Ω _n . The following lemma shows that under some conditions, we can ensure that   Var(h ₁) and Var(h ₂) are bounded. A proof has been relegated to the appendix.

Lemma 4.8

If we have E [ | X | 2 ] < ∞ , E [ | Y | 2 ] < ∞ and E [ | X | 2 | Y | 2 ] < ∞ , then we have Var(h ₄) < ∞ . Consequently, we also have Var(h ₁) < ∞ and Var(h ₂) < ∞ .

Even though as indicated in Lemma 4.7, the quantities h ₁(X ₁, Y ₁) and h ₂((X ₁, Y ₁), (X ₂, Y ₂)) play important roles in determining the variance of Ω _n , in a generic case, they do not have a simple formula. The following lemma gives the generic formulas for h ₁(X ₁, Y ₁) and h ₂((X ₁, Y ₁), (X ₂, Y ₂)). Its calculation can be found in the Supplementary Appendix.

Lemma 4.9

(Generic h ₁ and h ₂). In the general case, assuming (X ₁, Y ₁), (X, Y), (X′, Y′), and (X ^″, Y ^″) are independent and identically distributed, we have h 1 ( ( X 1 , Y 1 ) ) = 1 2 E [ | X 1 − X ′ | | Y 1 − Y ′ | ] − 1 2 E [ | X 1 − X ′ | | Y 1 − Y ″ | ] + 1 2 E [ | X 1 − X ′ | | Y − Y ″ | ] − 1 2 E [ | X 1 − X ′ | | Y ′ − Y ″ | ] + 1 2 E [ | X − X ″ | | Y 1 − Y ′ | ] − 1 2 E [ | X ′ − X ″ | | Y 1 − Y ′ | ] + 1 2 E [ | X − X ′ | | Y − Y ′ | ] − 1 2 E [ | X − X ′ | | Y − Y ″ | ] . We have a similar formula for h ₂((X ₁, Y ₁), (X ₂, Y ₂)) in (B.7). Due to its length, we do not display it here.

If one assumes that X and Y are independent, we can have a simpler formula for h ₁, h ₂, as well as their corresponding variances. We list the results below, with detailed calculations relegated to the appendix. One can see that under independence, the corresponding formulas are much simpler.

Lemma 4.10

When X and Y are independent, we have the following. For (X, Y) and (X′, Y′) that are independent and identically distributed as (X ₁, Y ₁) and (X ₂, Y ₂), we have h 1 ( ( X 1 , Y 1 ) ) = 0 , (4.2) h 2 ( ( X 1 , Y 1 ) , ( X 2 , Y 2 ) ) = 1 6 | X 1 − X 2 | − E [ | X 1 − X | ] − E [ | X 2 − X | ] + E [ | X − X ′ | ] ( | Y 1 − Y 2 | − E [ | Y 1 − Y | ] − E [ | Y 2 − Y | ] + E [ | Y − Y ′ | ] ) , (4.3) Var ( h 2 ) = 1 36 V 2 ( X , X ) V 2 ( Y , Y ) , (4.4) where E stands for the expectation operators with respect to X, X and X′, Y, or Y and Y′, whenever appropriate, respectively.

If we have 0 <  Var(h ₁) < ∞, it is known that the asymptotic distribution of Ω _n is normal, as stated in the following. Note that based on Lemma 4.10, X and Y cannot be independent; otherwise one should have h ₁ = 0 almost surely. The following theorem is based on a known result on the convergence of U-statistics, seeing ([13], Chapter 5.5.1 Theorem A). We state it without a proof.

Theorem 4.11

Suppose 0 <  Var(h ₁) < ∞ and Var(h ₄) < ∞, then we have Ω n → P V 2 ( X , Y ) moreover, we have n ( Ω n − V 2 ( X , Y ) ) → D N ( 0,16 Var ( h 1 ) ) , a s n → ∞ . When X and Y are independent, the asymptotic distribution of n Ω n is no longer normal. In this case, from Lemma 4.10, we have h 1 ( ( X 1 , Y 1 ) ) = 0  almost surely ,  and  V a r [ h 1 ( ( X 1 , Y 1 ) ) ] = 0 .

The following theorem, which applies a result in ([13], Chapter 5.5.2), indicates that nΩ _n converges to a weighted sum of (possibly infinitely many) independent χ 1 2 random variables.

Theorem 4.12

If X and Y are independent, the asymptotic distribution of Ω _n is n Ω n → D ∑ i = 1 ∞ λ i ( Z i 2 − 1 ) = ∑ i = 1 ∞ λ i Z i 2 − ∑ i = 1 ∞ λ i , where Z i 2 ∼ χ 1 2 i.i.d, λ _i ’s are the eigenvalues of operator G that is defined as G g ( x 1 , y 1 ) = E x 2 , y 2 [ 6 h 2 ( ( x 1 , y 1 ) , ( x 2 , y 2 ) ) g ( x 2 , y 2 ) ] , where function h ₂((⋅, ⋅), (⋅, ⋅)) was defined in (4.3).

Proof

The asymptotic distribution of Ω _n is from the result in ([13], Chapter 5.5.2).

See Section 4.3 for more details on methods for computing the value of λ _i ’s. In particular, we will show that we have ∑ i = 1 ∞ λ i = E [ | X − X ′ | ] E [ | Y − Y ′ | ] (Corollary 4.15) and ∑ i = 1 ∞ λ i 2 = V 2 ( X , X ) V 2 ( Y , Y ) (which is essentially from Eq. 4.4 and Lemma 4.7).

From Theorem 4.12, we see that the eigenvalues λ _i ’s play important role in determining the asymptotic distribution of Ω _n . We study its properties here. Throughout this subsection, we assume that X and Y are independent. Let us recall that the asymptotic distribution of sample distance covariance Ω _n , n Ω n → D ∑ i = 1 ∞ λ i ( Z i 2 − 1 ) = ∑ i = 1 ∞ λ i Z i 2 − ∑ i = 1 ∞ λ i , where λ _i ’s are the eigenvalues of the operator G that is defined as G g ( x 1 , y 1 ) = E x 2 , y 2 [ 6 h 2 ( ( x 1 , y 1 ) , ( x 2 , y 2 ) ) g ( x 2 , y 2 ) ] , where function h₂((⋅, ⋅), (⋅, ⋅)) was defined in Eq. 4.3. By definition, eigenvalues λ ₁ , λ ₂ , … corresponding to distinct solutions of the following equation G g ( x 1 , y 1 ) = λ g ( x 1 , y 1 ) . (4.5)

We now study the properties of λ _i ’s. Utilizing Lemma 12 and Eq. 4.4 in [12], we can verify the following result. We give details of verifications in the Supplementary Appendix.

Lemma 4.13

Both of the following two functions are positive definite kernels: h X ( X 1 , X 2 ) = − | X 1 − X 2 | + E [ | X 1 − X | ] + E [ | X 2 − X | ] − E [ | X − X ′ | ] and h Y ( Y 1 , Y 2 ) = − | Y 1 − Y 2 | + E [ | Y 1 − Y | ] + E [ | Y 2 − Y | ] − E [ | Y − Y ′ | ] .

The above result gives us a foundation to apply the equivalence result that has been articulated thoroughly in [12]. Equipped with the above lemma, we have the following result, which characterizes a property of λ _i ’s. The detailed proof can be found in the Supplementary Appendix.

Lemma 4.14

Suppose {λ ₁ , λ ₂ , …} are the set of eigenvalues of kernel 6h₂((x₁, y₁), (x₂, y₂)), { λ 1 X , λ 2 X , … } and { λ 1 Y , λ 2 Y , … } are the sets of eigenvalues of the positive definite kernels h_X and h_Y, respectively. We have the following: { λ 1 , λ 2 , … } = { λ 1 X , λ 2 X , … } ⊗ { λ 1 Y , λ 2 Y , … } ; that is, each λ _i satisfying (4.5) can be written as, for some j, j′, λ i = λ j X ⋅ λ j ′ Y where λ j X and λ j ′ Y are the eigenvalues corresponding to kernel functions h _X (X ₁ , X ₂ ) and h _Y (Y ₁ , Y ₂ ), respectively.

Above lemma implies that eigenvalues of h ₂ could be obtained immediately after knowing the eigenvalues of h _X and h _Y . But, in practice, there usually does not exist analytic solution for even the eigenvalues of h _X or h _Y . Instead, given the observations (X ₁ , …, X _n ) and (Y ₁ , …, Y _n ) , we can compute the eigenvalues of matrices K ̃ X = ( h X ( X i , X j ) ) n × n and K ̃ Y = ( h Y ( Y i , Y j ) ) n × n and use those empirical eigenvalues to approximate λ 1 X , λ 2 X , … and λ 1 Y , λ 2 Y , … , and then consequently λ ₁ , λ ₂ , …

We end this subsection with the following corollary on the summations of eigenvalues, which is necessary for the proof of Theorem 4.12. The proof can be found in the Supplementary Appendix.

Corollary 4.15

The aforementioned eigenvalues λ 1 X , λ 2 X , … and λ 1 Y , λ 2 Y , … satisfy ∑ i = 1 ∞ λ i X = E [ | X − X ′ | ] , a n d ∑ i = 1 ∞ λ i Y = E [ | Y − Y ′ | ] . As a result, we have ∑ i = 1 ∞ λ i = E [ | X − X ′ | ] E [ | Y − Y ′ | ] , and ∑ i = 1 ∞ λ i 2 = V 2 ( X , X ) V 2 ( Y , Y ) .

We have reviewed the properties of the statistics Ω _n in a previous section (Section 4.2). The disadvantage of directly applying Ω _n (which is defined in Eq. 2.5) is that for multivariate X and Y , the implementation may require at least O(n ² ) operations. Recall that for univariate X and Y, an O(n log  n) algorithm exists, cf. Theorem 2.2. The proposed estimator ( Ω ̄ n in Eq. 3.1) is the averaged distance covariances, after randomly projecting X and Y to one-dimensional spaces, respectively. In this section, we will study the asymptotic behavior of Ω ̄ n . It turns out that the analysis will be similar to the works in Section 4.2. The asymptotic distribution of Ω ̄ n will differ in two cases: (1) the general case and (2) the case when X and Y are independent.