Year: 2020

Abstract:
We propose a generalized probability kernel(GPK) on discrete distributions with finite support. This probability kernel, defined as kernel between distributions instead of samples, generalizes the existing discrepancy statistics such as maximum mean discrepancy(MMD) as well as probability product kernels, and extends to more general cases. For both existing and newly proposed statistics, we estimate them through empirical frequency and illustrate the strategy to analyze the resulting bias and convergence bounds. We further propose power-MMD, a natural extension of MMD in the framework of GPK, illustrating its usage for the task of two-sample test. Our work connects the fields of discrete distribution-property estimation and kernel-based hypothesis test, which might shed light on more new possibilities.

text: We propose a generalized probability kernel(GPK) on discrete distributions with finite support. This probability kernel, defined as kernel between distributions instead of samples, generalizes the existing discrepancy statistics such as maximum mean discrepancy(MMD) as well as probability product kernels, and extends to more general cases. For both existing and newly proposed statistics, we estimate them through empirical frequency and illustrate the strategy to analyze the resulting bias and convergence bounds. We further propose power-MMD, a natural extension of MMD in the framework of GPK, illustrating its usage for the task of two-sample test. Our work connects the fields of discrete distribution-property estimation and kernel-based hypothesis test, which might shed light on more new possibilities.

heading: 1 INTRODUCTION
text: We focus on the two-sample problem, which is given two i.i.d samples {x1, x2, ...xn} , {y1, y2, ..., yn}, could we infer the discrepancy between underlying distributions they are drawn from. For such a problem, the option of hypothesis test(two-sample test) is most popular, and a variety of statistics in estimating the discrepancy is proposed. In recent years, RKHS based method such as maximum mean discrepancy(MMD) has gained a lot of attention. (Gretton et al., 2012) has shown that in a universal-RKHS F , MMD(F ,p, q) = 0 if and only if p = q, thus could be used for the two-sample hypothesis test. (Gretton et al., 2012) further provides unbiased estimator of MMD with fast asymptotic convergence rate, illustrating its advantages.
On the other hand, estimating distribution properties with plugin(empirical) estimators on discrete setting is an active research area in recent years, where people focus on problem settings with large support size but not so large sample size. The Bernstein polynomial technique is introduced to analyze the bias of the plugin estimators in (Yi & Alon, 2020), which provides remarkable progress on bias-reduction methods of the plugin estimators. It is thus interesting to ask if the plugin estimators could motivate new results for the RKHS-based two-sample test.
Another interesting topic is about the probability kernel, defined as kernel function over probabilities, instead of over samples. As is easily seen, any discrepancy measure of distribution p and q could potentially be valid probability kernels, not so much work focuses on this. While (Jebara et al., 2004) introduced the so called probability product kernels which generalize a variety of discrepancy measures, its properties remain further study.
Motivated by above observations, our work focuses on a specialized probability kernel function which is a direct generalization of sample-based RKHS methods such as MMD. We focus on using plugin-estimator as the default estimator of the kernel function we defined, and illustrate that with the help of Bernstein polynomial techniques, we could analyze the bias and convergence bounds of these plugin-estimators. Our work thus connects the fields of discrete distribution-property estimation and kernel-based hypothesis test, which brings interesting possibilities.

heading: 2 NOTATION
text: We use bold symbol p, q ∈ Rk to represent a probability function over a discrete support with support size k, and pi, qi represents the ith entry of p and q. We use {v1, v2, ..., vk}, vi ∈ Rd to represent the support of p, q. [k] := {1, 2, 3..., k} represents the set of indices of elements in {v1, v2, ..., vk}. We use φ ◦ (p, q) to denote an element-wise function from Rk × Rk to Rk, where (φ ◦ (p, q))i = φ ◦ (pi, qi) and φ ◦ p to denote an element-wise function from Rk to Rk, where (φ ◦ p)i = φ ◦ pi. With a slight abuse of notation, we denote pρ, p − q as element-wise function defined above. We use kernel(p, q) to denote kernel function which maps from Rk × Rk to real value R. And kernel(x, y), x, y ∈ Rd represents a kernel function from Rd × Rd to real value R. We use K to denote the gram matrix generated from kernel(x, y) on finite support {v1, v2, ..., vk}, where Kij = kernel(vi, vj). We use {x1, x2, ..., xn} ∼ p and {y1, y2, ..., yn} ∼ q to denote the samples from distribution p and q, where n is the sample size.

heading: 3 GENERALIZED PROBABILITY KERNEL
text: Probability kernel function, defined as kernel function between distributions instead of samples, is a natural extension of the idea of kernel function in sample space. Definition 1. Given distribution p and q belongs to a family of discrete distribution with the same finite support {v1, v2, ..., vk}, vi ∈ Rd, where k is the support size, we define the probability kernel function as PK(p, q), which is a kernel function maps from Rk × Rk to real value R.
Many discrepancy measures, such as MMD, can serve as probability kernel functions, but people usually don’t use the term of probability kernel function when describing them. The reason is that for most of the time, we only consider a limited number of distributions, and do not need or have the resources to navigate through all the distributions within the family. For example, when looking into the two-sample problem, we usually assume two samples {x1, x2, ..., xn} ∈ Rd and {y1, y2, ..., yn} ∈ Rd are i.i.d drawn from two distributions p and q, and use the discrepancy measure MMD[F ,p, q] to determine if p and q are indistinguishable in the RKHS F . We do not consider all other distributions in F that is irrelevant to our samples! So far the idea of kernel function between distributions is in practice not so much useful, however, here in this paper, we propose, when considering the plugin-estimator of many of the existing discrepancy measures, it is beneficial to view them as probability kernel functions.

heading: 3.1 DEFINATION OF GENERALIZED PROBABILITY KERNEL
text: Definition 2 (Generalized probability kernel). Given the family S of discrete distribution on support {v1, v2, .., vk} where vi ∈ Rd. Let F be a unit ball in a universal-RKHS H with associated continuous kernel RK(x, y), where for any x ∈ Rd and y ∈ Rd, RK(x, y) maps from Rd ×Rd to R. We denote gram matrix Kij = RK(vi, vj).
The generalized probability kernel function on distribution p, q ∈ S is GPKF,φ(p, q) = φ ◦ (p, q)Kφ ◦ (q,p)T = ∑ i∈[k] ∑ j∈[k] φ ◦ (pi, qi)Kijφ ◦ (qj , pj)
where φ ◦ (p, q) is an element-wise mapping function on discrete distribution p, q ∈ S, which maps from Rk × Rk to Rk,
Obviously, under this definition, the GPK is a symmetric probability kernel function where GPKF,φ(p, q) = GPKF,φ(q,p)
Mapping function φ represent a great amount of possibilities. For most cases, we need to narrow down the region and equipped it with some convenient properties so that the GPK measure could be useful.
One example is for the measurement of discrepancy, where we want GPKF,φ(p, q) = 0 if and only if p = q. Definition 3 (discrepancy probability kernel). Let S be a family of discrete distribution p ∈ S on support {v1, v2, ..., vk}. A discrepancy probability kernel is a kernel function PK(p, q) that PK(p, q) = 0 if and only if p = q
Theorem 1. GPKF,φ(p, q) with the mapping function φ that satisfies:
1. symmetric or antisymmetric with respect to p and q: φ ◦ (p, q) = φ ◦ (q,p) or φ ◦ (p, q) = −φ ◦ (q,p) 2. ‖φ ◦ (p, q)‖2 = ‖φ ◦ (q,p)‖2 = 0 if and only if p = q, where ‖ · ‖2 represents L2 norm. is a discrepancy probability kernel.
Proof. GPKF,φ(p, q) = ∑ i∈[k] ∑ j∈[k] φ ◦ (pi, qi)Kijφ ◦ (qj , pj)
= φ ◦ (p, q)Kφ ◦ (q,p)T = ±φ ◦ (p, q)Kφ ◦ (p, q)T
= ±vKvT
K is a semipositive definite matrix, thus by definition of positive definite matrix, vKvT ≥ 0, where equality holds if and only if v = 0, and since v = φ ◦ (p, q), this condition further means φ ◦ (p, q) = 0, which holds if and only if p = q.
Another example is the polynomial GPK, which is our main focus of this paper. Such a subclass of GPK is interesting since we can build unbiased estimators of it using techniques of Bernstein polynomial in (Qian et al., 2011). As we will show in section 5., we also have analyzable convergence bounds for the resulting unbiased estimators, illustrating its potential usage for applications such as two-sample test. Definition 4 (polynomial GPK). The polynomial GPK is the subset of GPK that equipped with the mapping function φ that is polynomial in p and q: φ ◦ (p, q) = ∑o l=0 ∑o s=0 αl,sp
lqs where o ∈ Z is the degree of the polynomial, and al,s ∈ R is the coefficient
Below we give some examples of polynomial GPK, which include MMD proposed in (Gretton et al., 2012), and the newly proposed power-MMD in this paper, which is a natural extension of MMD, from the view point of probability kernels.

heading: 3.1.1 EXAMPLE 1: MMD AS MEMBER OF POLYNOMIAL GPK
text: Given discrete distribution p, q with support {v1, v2, ..., vk}, we can rewrite MMD with distribution probability function pi, qi:
MMD2F (p, q) = ‖Ex∼pf(x)−Ex′∼qf(x′)‖ 2 H
= ∥∥∥∥∥∥ ∑ i∈[k] f(vi)pi − ∑ i∈[k] f(vi)qi ∥∥∥∥∥∥ 2
H
= ∥∥∥∥∥∥ ∑ i∈[k] f(vi)pi − f(vi)qi ∥∥∥∥∥∥ 2
H = ∑ i∈[k] ∑ j∈[k] (pi − qi)f(vi)f(vj)(pj − qj) = ∑ i∈[k] ∑ j∈[k] (pi − qi)Kij(pj − qj)
= −GPKF,φl(p, q)
Where φl ◦ (p, q) = p − q, H is the RKHS defined in MMD literature, and f is the function that maps vi toH. GPKF,φl(p, q) is a special case of polynomial GPK where α1,0 = 1, α0,1 = −1, and all other coefficients are 0.

heading: 3.1.2 EXAMPLE 2: PRODUCT GPK AS MEMBERS OF POLYNOMIAL GPK
text: Definition 5 (product GPK). Let p and q be probability distributions on support {v1, v2, ..., vk}, and l ∈ Z be nonnegative integer. The product GPK is a subset of polynomial GPK where αl,0 = 1, and all other coefficients are 0. the corresponding mapping function is: φ(p, q) = pl
The probability product kernel as in (Jebara et al., 2004) is a special case of product GPK where K is a identity matrix.

heading: 3.1.3 EXAMPLE 3: POWER-MMD AS MEMBERS OF POLYNOMIAL GPK
text: Another interesting subset of polynomial GPK is the one extends MMD case into a power form and we denote it as power-MMD:
Definition 6 (power-MMD). Let p and q be probability distributions on support {v1, v2, ..., vk} and ρ ∈ Z be a positive integer. then the power-MMD is a subset of polynomial GPK where αρ,0 = 1, α0,ρ = −1, and all other coefficients are 0. the corresponding mapping function is: φ(p, q) = pρ − qρ
Apparently, MMD is a special case of power-MMD where ρ = 1, and power-MMD satisfies the requirement in Theorem 1, thus has the potential usage of discrepancy measure. In section 5., we will show that power-MMD has unbiased estimator with analyzable convergence bounds thus could be used for two-sample test.

heading: 3.2 DISCUSSION OF GPK IN DISCRETE SETTING
text: As one may easily notice, the definition of GPK includes a gram matrix generated by the kernel function RK(vi, vj) which measures the discrepancy between vi, vj ∈ {v1, v2, .., vk}. While considering the cases of categorical distribution, values of discrete variables does not relate to any notion of distance, this raises the question: how the introduced gram matrix will be beneficial in any cases?
The answer is twofold: 1. Many natural processes produce discrete distributions where there possibly exists a similarity measure in values which imply the similarity in frequencies of occurrence(probability values). For example, in the field of natural language process(NLP), one may treat words as atomic units with no notion of similarity between words, as these are represented as indices in a vocabulary. However, given large number of training samples, similarity measure between words could be made possible using techniques such as words2vec(Mikolov et al., 2013). Such techniques generally result in better performance and have become the important preprocessing techniques for NLP tasks(Goodfellow et al., 2016). 2. As there are cases where the values of discrete variables are totally irrelevant, or people may use kernel function RK(vi, vj) which doesn’t correctly imply the similarity in probability values, the GPK framework may still capture the similarity between distributions. One example is the case of MMD, which is, as we discussed above, an element of GPK family. As proved in (Gretton et al., 2012), MMD is a distribution free measurement between two samples, which means no matter what kind of p, q and kernel(x, y), we have, the MMD2F (p, q) measure will be 0 if and only if p = q. However, the bad choice of kernel function does have a negative effect on convergence bounds of the empirical estimator proposed in (Gretton et al., 2012), and will influence the results of two-sample test. For this reason, we mainly focus on dataset with known relativity measures in our experiment section.

heading: 4 PLUGIN-ESTIMATOR FOR GPK
text: So far we have defined the GPK and discussed some subsets of GPK with potential usage of two-sample test. Next we discuss how to build an estimator, given a member of GPK. In this section, we propose the plugin-estimator, which based on the count of occurrence of each value vi ∈ {v1, v2, ..., vk} in samples {x1, x2, ..., xn} ∈ p or q. We illustrate that by doing so, the techniques of Bernstein polynomial in (Qian et al., 2011) could be used to help building unbiased estimators for any members of polynomial GPK. Furthermore, we provide analyzable convergence bounds of these estimators.
We begin with the definition of plugin-estimators:
Definition 7. Suppose we have i.i.d samples of distribution p as Xn1 := {x1, x2, ..., xn1} ∼ p and Xn2 := {xn1+1, xn1+2, ..., xn1+n2} ∼ p. And also the i.i.d samples of distribution q as Ym1 := {y1, y2, ..., ym1} ∼ q and Ym2 := {ym1+1, ym1+2, ..., ym1+m2} ∼ q.
Let N (n1)i denotes the number of occurrence of value vi ∈ {v1, v2, ..., vk} in sample Xn1 , and Si,n1 := (N (n1) i , n1) denotes the collection of N (n1) i and n1. The same follows for Xn2 , Ym1 and Ym2
We define the plugin-estimator of GPKF,φ(p, q) as GPKE [F , φ,X, Y ] = ∑ i∈[k] ∑ j∈[k] fφ (Si,n1 , Si,m1)Kijfφ (Sj,m2 , Sj,n2)
where fφ is a function related to function φ, and K is the gram matrix brought by F .
Here our setting is different from the unbiased estimator MMD2u of (Gretton et al., 2012), where in their setting Xn1 , Xn2 represent the same sample from p and so do for Ym1 , Ym2 from q. Instead, we are using the same setting as the linear time statistic MMD2l proposed in (Gretton et al., 2012). Another way of viewing this is that for our setting, given two samples {x1, x2, ..., xn}, {y1, y2, ..., yn} from p and q, we depart each sample of x and y into two parts, yielding 4 different samples with size n1, n2, m1, m2, and then calculate the empirical frequencies for plugin-estimator defined above.

heading: 4.1 POLYNOMIAL GPK WITH UNBIASED PLUGIN-ESTIMATORS
text: One of our main contributions of this paper is the proposal that we can always find an unbiased plugin-estimator for any members in polynomial GPK family. The basic idea is that we can analyze the expectation of plugin-estimators through Bernstein polynomial, and use the existing results of (Qian et al., 2011) to build the unbiased plugin estimators. Theorem 2. Denote
gj(k, n) :=  gj(k, n) = ( k j )( n j )−1 , for j ≤ k
0, for j > k
Then any member of polynomial GPK[F , φ,p, q] equipped with polynomial mapping function φ(p, q) = ∑o l=0 ∑o s=0 αl,sp
lqs of degree o ∈ Z, has an unbiased plugin-estimator with mapping function fφ to be:
fφ(Si,n1 , Si,m1) = o∑ l=0 o∑ s=0 αl,sgl(N (n1) i , n1)gs(N (m1) i ,m1)
Proof. The basic idea is directly using the result of Bernstein polynomial in (Qian et al., 2011) to build unbiased estimators. We put our formal proof in appendix
For notation simplicity, we define the plugin-estimator discussed above to be the default-pluginestimator for polynomial GPK: Definition 8 (default-plugin-estimator for polynomial GPK). The plugin-estimator defined in Theorem 2 is the default-plugin-estimator for polynomial GPK.
This plugin-estimator, according to Theorem 2, is an unbiased estimator

heading: 4.2 DEVIATION BOUND OF PLUGIN-ESTIMATORS
text: Another topic about plugin-estimator is its deviation bound. We directly use the McDiamid’s inequality to derive the default-plugin-estimator for polynomial GPK: Theorem 3. The default-plugin-estimator of GPK[F , φ,p, q] equipped with polynomial mapping function φ(p, q) = ∑o l=0 ∑o s=0 αl,sp
lqs of degree o ∈ Z has the convergence bound:
∀a > 0,Pr(|GPKE [F , φ,X, Y ]− E[GPKE [F , φ,X, Y ]]| ≥ a) ≤ 2e− 2a2 Z
where
Z = (( n1 ( τ (1)n1,m1 )2 +m1 ( τ (2)n1,m1 )2) Φ2m2,n2 + ( m2 ( τ (1)m2,n2 )2 + n2 ( τ (2)m2,n2 )2) Φ2n1,m1 ) K2max
Φn,m = ∑ i∈[k] ∣∣∣∣∣ o∑ l=0 o∑ s=0 αl,sgl(N (n) i , n)gs(N (m) i ,m) ∣∣∣∣∣ τ (1)n,m = sup
i∈[k] ( o∑ l=0 o∑ s=0 l N (n) i · |αl,s| · gl(N (n)i , n)gs(N (m) i ,m) )
τ (2)n,m = sup i∈[k] ( o∑ l=0 o∑ s=0 s N (m) i · |αl,s| · gl(N (n)i , n)gs(N (m) i ,m) )
Kmax is the largest value of entries in K
Proof. The basic idea is to use the McDiamid’s inequality, and we put our formal proof into the appendix.

heading: 5 EXAMPLE: POWER-MMD AS A NATURAL EXTENSION TO MMD FROM
text: GPK VIEWPOINT
In this section, we mainly discuss power-MMD as defined in 3.1.3. We analyze the bias and convergence bound of its plugin-estimators using the techniques we introduced so far, illustrating that such a natural extension to MMD from GPK viewpoint could be beneficial for two-sample test.

heading: 5.1 PLUGIN-ESTIMATORS OF POWER-MMD
text: As we already discussed in section 4.1.3, power-MMD is a subset of polynomial GPK. According to Theorem 2, any member GPKF,φρ(p, q) in power-MMD has a default-plugin-estimator with the mapping function fφ(Si,n1 , Si,m1) = gρ(Ni, n1)− gρ(Mi,m1) Remark 3.1. When ρ = 1, the power-MMD return to the original MMD case. Remarkably, the default-plugin-estimator of this case is equivalent to the linear time statistic MMD2l proposed in (Gretton et al., 2012): GPKE [F , φl, X, Y ] = MMD2l [F , X, Y ] For details of the derivation, see appendix

heading: 5.2 DEVIATION BOUND OF PLUGIN-ESTIMATORS OF POWER-MMD
text: Corollary 3.1. Denote τn = supi∈[k] ( ρ
N (n) i
gρ(N (n) i , n) ) The default-plugin-estimator of power-MMD GPKF,φρ(p, q) has uniform convergence bound defined in Theorem 3. with τ (1)n,m = τn and τ (2) n,m = τm
Corollary 3.2. Consider the case where n1 = n2 = m1 = m2 = n The default-plugin-estimator of power-MMD GPKF,φρ(p, q) has uniform convergence bound:
Pr(|GPKE [F , φρ, X, Y ]− E[GPKE [F , φρ, X, Y ]]| ≥ a) ≤ 2e −na2 (ρ2Φ2n1,m1 +ρ2Φ2m2,n2 )K2max ≤ 2e −na2 8ρ2K2max
Proof. The first inequality above comes from:
ρ
Ni gρ(Ni, n) = ρ(Ni − 1)(Ni − 2)...(Ni − ρ+ 1) n(n− 1)...(n− ρ+ 1) ≤ ρ n Nρ−1i nρ−1 ≤ ρ n
where supi∈[k] ( ρ Ni gρ(Ni, n) ) = ρn only stands for extreme case such that there exist Ni = n, i.e. all the samples belongs to the same value vi ∈ {v1, v2, ..., vk}. And the second inequality above comes from:
Φn,m = ∑ i∈[k] ∣∣∣gρ(N (n)i , n)− gρ(N (m)i ,m)∣∣∣ ≤∑ i∈[k] ∣∣∣gρ(N (n)i , n)∣∣∣+ ∑ i∈[k] ∣∣∣gρ(N (m)i ,m)∣∣∣ ≤ ∑ i∈[k] ( N (n) i n )ρ + ∑ i∈[k] ( N (m) i m )ρ ≤ ∑ i∈[k] N (n) i n ρ + ∑ i∈[k] N (m) i m ρ = 2
Remark 3.2. recall in (Gretton et al., 2012), the deviation bound for linear time estimator MMD2l [F , X, Y ] is
∀a > 0,Pr( ∣∣MMD2l [F , X, Y ]− E[MMD2l [F , X, Y ]]∣∣ ≥ a) ≤ 2e −na28K2max
Interestingly, this bound is the same as the case ρ = 1 in Corollary 3.2. Note that according to section 6.1, the default-plugin-estimator of power-MMD with ρ = 1 is actually in equivalent to MMD2l case in (Gretton et al., 2012). Our bound generalize the bound in (Gretton et al., 2012) and provide a tighter version. Note that the bounds for special case of ρ = 1 has simpler derivation, and the reader may refer to appendix for more details.

heading: 5.3 TWO-SAMPLE TEST USING POWER-MMD
text: Corollary 3.3. A hypothesis test of level α for the null hypothesis p = q has the acceptance region∣∣∣GPKE [F,φρ,X,Y ]√ Z ∣∣∣ <√ 12 log((α2 )−1) Where Z is defined in Corollary 3.1 The two-sample test for power-MMD then follows this procedure: 1. calculate v =∣∣∣GPKE [F,φρ,X,Y ]√
Z ∣∣∣. 2. check if v < √ 12 log((α2 )−1), if so, accept the null hypothesis, otherwise reject the null hypothesis.
Next we analyze the performance of our proposed two-sample test under two cases: ρ = 1 and ρ > 1
5.3.1 ρ = 1 CASE
For ρ = 1 case, since GPKE [F , φ1, X, Y ] is equivalent to MMD2l [F , φ1, X, Y ], the only difference between our proposal and that of Gretton et al. (2012) is the convergence bound. According to Remark 3.2, we provide a tighter bound for the test statistic, thus we will certainly have a better performance using power-MMD.
5.3.2 ρ > 1 CASE
We need to answer two questions for the ρ > 1 case: 1. when applying power-MMD in practice, is the proposed statistics numerical stable? 2. will the performance of two-sample test gets better when ρ gets larger?
For the question of numerical stability, since gρ(Ni, n) ≤ ( Ni n )ρ , the term will exponentially decrease with the increase of ρ. This effect will cause numerical problem when Ni n and ρ is large. One solution is to find an upper-bound of ∣∣∣GPKE [F,φρ,X,Y ]√ Z
∣∣∣ which is numerical stable. Corollary 3.4. Consider the simplest case where n1 = n2 = m1 = m2 = n. Define
CN := {N (n1)1 , N (n1) 2 , ..., N (n1) k , N (n2) 1 , N (n2) 2 , ..., N (n2) k , N (m1) 1 , N (m1) 2 , ...,
N (m1) k , N (m2) 1 , N (m2) 2 , ..., N (m2) k }
to be the set of all counts of occurrence in the four samples Xn1 , Xn2 , Xm1 , Xm2 . Denote SN = supNi∈CN (Ni) to be the maximium value in the set CN We have:∣∣∣∣GPKE [F , φρ, X, Y ]√Z ∣∣∣∣ ≤ ∣∣∣∣ SN ·GPK′EKmaxρ√2n · Φ′ ∣∣∣∣
where GPK′E := ∑ i,j∈[k] ( gρ(N (n1) i , n)− gρ(N (m1) i , SN ) ) Kij ( gρ(N (n2) j , SN )− gρ(N (m2) j , SN ) ) and
Φ′ = ∑ i∈[k] ∣∣∣gρ(N (n)i , SN )− gρ(N (m)i , SN )∣∣∣ For cases when Ni are not far less from SN , GPK′E will be much more numerical stable than GPKE [F , φρ, X, Y ]. To answer the question related to the performance of two-sample test when ρ get larger, we need to analyze the case when p 6= q, if ∣∣∣GPKE [F,φρ,X,Y ]√ Z
∣∣∣ increase with the increase of ρ. Unfortunately, there is no clear answer to this.

heading: 6 SUMMARY
text: To summarize, we introduce the framework of generalized probability kernel(GPK). While GPK represents a large family of probability kernels, we focus on polynomial GPK since all members of such subset of GPK have unbiased plugin-estimators. Remarkably, a natural extension of MMD from the viewpoint of polynomial GPK, which we call power-MMD, could be used for two-sample test. Theoretical study shows that for ρ = 1 case, power-MMD outperforms linear time MMD proposed in Gretton et al. (2012), and the performance of ρ > 1 case is left for future work. For members of GPK which do not belong to polynomial GPK, it is not easy to design an unbiased estimators. However, bais reduction techniqes proposed in (Yi et al., 2018) and (Yi & Alon, 2020) could be used, and we still have the chance to apply two-sample test with the resulting estimators. Such a possibility is also left for future work.

heading: A APPENDIX
text: A.1 BERNSTEIN POLYNOMIAL
Drawing i.i.d. samples Y m from any distribution p the expected value of the empirical estimator for a distribution property is
E [ ĤE (Y m) ] = ∑ i∈[k] E Mi∼bin(m,pi) [ h ( Mi m )]
Note that for any function f,m ∈ N, and x ∈ [0, 1], the degree-m Bernstein polynomial of f is
Bm(f, x) := m∑ j=0 f ( j m )( m j ) xj(1− x)m−j
Therefore, we can express the expectation of the empirical property estimator as
E Ym∼p
[ ĤE (Y m) ] = ∑ i∈[k] Bm (h, pi)
A.2 PROOF OF THEOREM 2
Proof. Recall the definition of polynomial GPK: GPK(F , φ,p, q) = ∑ i∈[k] ∑ j∈[k] φ(pi, qi)Ki,jφ(qj , pj) = ∑ i∈[k] ∑ j∈[k] Ki,j o∑ l,s,r,t=0 αl,sαr,tp l iq s i p t jq r j
recall in (Qian et al., 2011)
xj = n∑ k=j gj(k, n)bk,n(pi) = Ek∼bin(pi,n)gj(k, n)
where gj(k, n) is defined in the beginning of the theorem
GPK(F , φ,p, q) = ∑ i∈[k] ∑ j∈[k] Ki,j o∑ l,s,r,t=0 αl,sαr,tp l iq s i p t jq r j
= ∑ i∈[k] ∑ j∈[k] Ki,j o∑ l,s,r,t=0 αl,sαr,tEN(n1)i ∼bin(pi,n1)gl(Ni, n1)
E N
(m1) i ∼bin(qi,m1) gs(N (m1) i , n2)EN(n2)j ∼bin(pj ,n2)gt(N (n2) j , n2)EN(m2)j ∼bin(qj ,m2)gr(N (m2) j ,m2)
= E ∑ i∈[k] ∑ j∈[k] Ki,j o∑ l,s,r,t=0 αl,sαr,tgl(N (n1) i , n1)gs(N (m1) i ,m1)gt(N (n2) j , n2)gr(N (m2) j ,m2)  = E [GPKE(F , φ,X, Y )]
A.3 PROOF OF THEOREM 3.
Lemma 4. Let SNi,n1 := (N (n1) i , n1) denotes the collection of N (n1) i and n1. The same follows for Xn2 , Ym1 and Ym2 . Also for notation simplicity, let Si,n1 := SNi,n1 For the plugin-estimator
GPKE [F , φ,X, Y ] = ∑ i∈[k] ∑ j∈[k] fφ (SNi,n1 , SNi,m1)Kijfφ ( SNj ,m2 , SNj ,n2 ) with mapping function fφ(SNi,n1 , SNi,m1) having the following properties:
• fφ(SNi,n1 , SNi,m1) is a monotonic function related to N (n1) i and N (m1) i :
∀Ni > N ′i or ∀Ni < N ′i : fφ(SNi,n1 , SNi,m1) > fφ(SN ′i ,n1 , SNi,m1)
The same follows for N (m1)i
• |fφ(SNi±1,n1 , SNi,m1)− fφ(SNi,n1 , SNi,m1)| ≤ τn1
and
|fφ(SNi,n1 , SNi±1,m1)− fφ(SNi,n1 , SNi,m1)| ≤ τm1 where τn1 is a constant related to sample size n1, τm1 is a constant related to sample size m1, the same follows for τn2 and τm2
We have:
∀a > 0,Pr(|GPKE [F , φ,X, Y ]− E[GPKE [F , φ,X, Y ]]| ≥ a) ≤ 2e −2a2 Z
where
Z = (( n1τ 2 n1 + n2τ 2 n2 ) Φ22 + ( m1τ 2 m1 +m2τ 2 m2 ) Φ21 ) K2max
Φ1 = ∑ j∈[k] ∣∣fφ(SNj ,n1 , SNj ,m1)∣∣ ,Φ2 = ∑ j∈[k] ∣∣fφ(SNj ,n2 , SNj ,m2)∣∣ Kmax is the largest value of entries in K
Proof. recall McDiamid’s inequality
Theorem 5. Let Y1, . . . , Ym be independent random variables taking values in ranges R1, . . . , Rm and let F : R1 × . . . × Rm → C with the property that if one freezes all but the wth coordinate of F (y1, . . . , ym) for some 1 ≤ w ≤ m, then F only fluctuates by most cw > 0, thus | F (y1, . . . , yw−1, yw, yw+1, . . . , ym) − F (y1, . . . , yw−1, y′w, yw+1, . . . , ym) |≤ cw for all yj ∈ Rj and y′w ∈ Rw for 1 ≤ j ≤ m
Then for any a > 0, one has Pr(|F (Y )− E[F (Y )]| ≥ a) ≤ 2e − 2a2∑n i=1 c2 i
considering the plugin-estimator of GPK family:
GPKE [F , φ,X, Y ] = ∑ i∈[k] ∑ j∈[k] fφ (Si,n1 , Si,m1)Kijfφ (Sj,m2 , Sj,n2)
Without loss of generality, we rewrite the function fφ as:
fφ (SNi,n1 , SNi,m1) = F (x1, x2, ..., xs, ..., xn1 , N (m1) i ,m1) = FNi
Assume we freeze all but one element in Xn1 := {x1, x2, ..., xs, ..., xn1}, and only xs is allowed to change its value.
obviously, no matter how this element change, it always lies in the finite set of support {v1, v2, ..., vk}, without loss of generality, we assume xs changes its value from vi to vii, thus the corresponding count of occurrence Ni changes to Ni − 1 and Nii changes to Nii + 1
we have for xs ∈ Xn1 cs = supxs ∣∣GPKE(x1, x2, ..., xn)−GPKE(x1, x2, ..., xs, ..., xn)∣∣ = supxs ∣∣∣∣∣∣ ∑ i,j∈[k] ( F (x1, x2, ..., xs, ..., xn1 , N (m1) i ,m1)− F (x1, x2, ..., xs, ..., xn1 , N (m1) i ,m1) ) Kijfφ (Sj,m2 , Sj,n2)
∣∣∣∣∣∣ = supi,ii∈[k]
∣∣∣∣∣∣ ∑ j∈[k] ((FNi−1 − FNi)Ki,j + (FNii+1 − FNii,1)Kii,j) fφ(SNj ,n2 , SNj ,m2) ∣∣∣∣∣∣ ≤
∣∣∣∣∣∣ ∑ j∈[k] τn1(−Ki,j +Kii,j)fφ(SNj ,n2 , SNj ,m2) ∣∣∣∣∣∣ ≤ |τn1Kmax| ∑ j∈[k]
∣∣fφ(SNj ,n2 , SNj ,m2)∣∣ = τn1KmaxΦ2
where Φ2 = ∑ j∈[k] ∣∣fφ(SNj ,n2 , SNj ,m2)∣∣ Note that Ki,j := Kij
Similarly, for xs ∈ Xn2 cs ≤ τn2KmaxΦ1
where Φ1 = ∑ j∈[k] |φ (p1j , q1j)|
for xs ∈ Ym1 cs ≤ τm1KmaxΦ2
for xs ∈ Ym2 cs ≤ τm2KmaxΦ1
thus according to McDiamid’s inequality, we have
Pr(|GPKE [F , φ,X, Y ]− E [GPKE [F , φ,X, Y ]]| ≥ a) ≤ 2e −2a2∑n1+n2+m1+m2 i=1 c2 i
We set ci = τn1KmaxΦ2 for xi ∈ {Xn1},
ci = τm1KmaxΦ2
for xi ∈ {Ym1},
ci = τn2KmaxΦ1
for xi ∈ {Xn2},
ci = τm2KmaxΦ1
for xi ∈ {Ym2} and get
n1+n2+m1+m2∑ i=1 c2i = n1 · (τn1KmaxΦ2) 2 + n2 · (τn2KmaxΦ1) 2
+m1 · (τm1KmaxΦ2) 2 +m2 · (τm2(KmaxΦ1) 2 = (( n1τ 2 n1 +m1τ 2 m1 ) Φ22 + ( n2τ 2 n2 +m2τ 2 m2 ) Φ21 ) K2max
We are ready to proof theorem 3:
Proof. Define
Div (1) n,l,s = |gl(N (n) i + 1, n)gs(N (m) i ,m)− gl(N (n) i , n)gs(N (m) i ,m)|
= ∣∣∣∣ (Ni + 1)Ni(Ni − 1)...(Ni − ρ+ 2)n1(n1 − 1)...(n1 − ρ+ 1) − Ni(Ni − 1)...(Ni − ρ+ 1)n1(n1 − 1)...(n1 − ρ+ 1) ∣∣∣∣ · ∣∣∣gs(N (m)i ,m)∣∣∣
= ∣∣∣∣∣ lN (n)i − l + 1gl(N (n)i , n) ∣∣∣∣∣ · ∣∣∣gs(N (m)i ,m)∣∣∣
= l
N (n) i − l + 1
gl(N (n) i , n)gs(N (m) i ,m)
Div (2) n,l,s = |gl(N (n) i − 1, n)gs(N (m) i ,m)− gl(N (n) i , n)gs(N (m) i ,m)|
= ∣∣∣∣ (Ni − 1)(Ni − 2)...(Ni − ρ)n1(n1 − 1)...(n1 − ρ+ 1) − Ni(Ni − 1)...(Ni − ρ+ 1)n1(n1 − 1)...(n1 − ρ+ 1) ∣∣∣∣ · ∣∣∣gs(N (m)i ,m)∣∣∣
= ∣∣∣∣∣ lN (n)i gl(N (n)i , n) ∣∣∣∣∣ · ∣∣∣gs(N (m)i ,m)∣∣∣
= l
N (n) i
gl(N (n) i , n)gs(N (m) i ,m)
Divn,l,s = max(Div (1) n,l,s, Div (2) n,l,s) =
l
N (n) i
gl(N (n) i , n)gs(N (m) i ,m)
Recall the mapping function of default-plugin-estimator of polynomial GPK fφ (Si,n, Si,m) :=∑o l=0 ∑o s=0 αl,sgl(N (n) i , n)gs(N (m) i ,m)
Apparently fφ is a monotonic function with respect to Ni and Mi, thus condition 1. for Theorem 3. is satisfied
Since we also have:
|fφ(SN(n)i ±1,n, SN(m)i ,m)− fφ(SNi,n1 , SMi,m1)| ≤ o∑ l=0 o∑ s=0 |αl,sDivn,l,s|
= o∑ l=0 o∑ s=0 |αl,s| l N (n) i gl(N (n) i , n)gs(N (m) i ,m)
thus condition 2. for Theorem 3. is satisfied
A.4 PROOF OF COROLLERY 3.4
Proof. Define τ = supNi∈CN ρ Ni gρ(Ni, n) = ρ SN gρ(SN , n), Φn,m = max{Φn1,m1 ,Φm2,n2}. Since we have gρ(N,n)gρ(M,n) = gρ(N,M), we could get∣∣∣∣GPKE [F , φρ, X, Y ]√Z ∣∣∣∣ ≤ ∣∣∣∣∣∣ ∑ i,j∈[k] ( gρ(N (n1) i , n)− gρ(N (m1) i , n) ) Kij ( gρ(N (n2) j , n)− gρ(N (m2) j , n) ) KmaxτΦn,m √ 2n
∣∣∣∣∣∣ = ∣∣∣∣∣∣ ∑ i,j∈[k] ( gρ(N (n1) i , n)− gρ(N (m1) i , n) ) Kij ( gρ(N (n2) j , n)− gρ(N (m2) j , n) ) Kmax
ρ SN gρ(SN , n) ∑ i∈[k] ∣∣∣gρ(N (n)i , n)− gρ(N (m)i , n)∣∣∣√2n ∣∣∣∣∣∣
= ∣∣∣∣∣∣ SN ∑ i,j∈[k] ( gρ(N (n1) i , SN )− gρ(N (m1) i , SN ) ) Kij ( gρ(N (n2) j , SN )− gρ(N (m2) j , SN ) ) Kmaxρ √ 2n · ∑ i∈[k] ∣∣∣gρ(N (n)i , SN )− gρ(N (m)i , SN )∣∣∣ ∣∣∣∣∣∣
A.5 DETAILS OF REMARK 3.1 GPKE [F , φl, X, Y ] = ∑ i∈[k] ∑ j∈[k] ( g1(N (n1) i , n1)− g1(N (m1) i ,m1) ) Kij ( g1(N (n2) j , n2)− g1(N (m2) j ,m2) )
= ∑ i∈[k] ∑ j∈[k]
( N
(n1) i
n1 − N
(m1) i
m1
) Kij ( N (n2) j
n2 − N
(m2) j
m2
)
= 1
m1m2 m1∑ i=1 m2∑ j=1 k (xi, xj) + 1 n1n2 n1∑ i=1 n2∑ j=1 k (yi, yj)− 1 m1n2 m1∑ i=1 n2∑ j=1 k (xi, yj)− 1 m2n1 m2∑ i=1 n1∑ j=1 k (xi, yj)
= MMD2l [F , X, Y ]
A.6 THE CONVERGENCE BOUNDS OF POWER-MMD WITH ρ = 1
When ρ = 1, τn1 = ρ Ni gρ(N,n1) = 1 Ni · Nin1 = 1 n1 For simplicity consider the case where n1 = n2 = m1 = m2 = n, from Corollary 3.1. we have
Pr(|GPKE [F , φl, X, Y ]− E[GPKE [F , φl, X, Y ]]| ≥ a) ≤ 2e −na2 K2max2Φ 2
Since Φ = ∑ j∈[k] |fφ(Sj,n, Sj,m)| = ∑ j∈[k] ∣∣∣∣N(n)jn − N(m)jm ∣∣∣∣ ≤∑j∈[k] |N(n)jn | +∑j∈[k] |N(m)jm | = 2
we have
Pr(|GPKE [F , φl, X, Y ]− E[GPKE [F , φl, X, Y ]]| ≥ a) ≤ 2e −na2 K2max2Φ 2 ≤ 2e −na2 8K2max
recall in (Gretton et al., 2012), the deviation bound for linear time estimator MMD2l [F ,p, q] is
Pr( ∣∣MMD2l [F , X, Y ]− E[MMD2l [F , X, Y ]]∣∣ ≥ a) ≤ 2e −na28K2max
Thus our bound generalize the bound in (Gretton et al., 2012) and provide a tighter version.
A.7 IS BOUNDS GET TIGHTER WHEN ρ GETTING LARGER?
As we’ve already known, ρ = 1 case is equivalent to MMDl in (Gretton et al., 2012), one question rises: Would the performance of cases ρ > 1 better than widely used ρ = 1 case?
According to Corollary 3.2., since e −na2 8K2max ≤ e −na2
8ρ2K2max , the convergence bounds for ρ > 1 cases seem looser than ρ = 1 case, and this may give a negative answer to the question above.
However, the bound above is based on the worst cases where supi(Ni) = n, such that τn ≤ ρn and Φ ≤ 2. In practice, we are less likely to come across such a phenomena, instead, we may assume the supi(Ni) to be far smaller.
Without loss of generality, assume we have max( supi(N (n) i )
n , supi(N
(m) i )
m ) ≤ 1 α , where α ≥ 1, it is
easily seen:
τn ≤ ρ
n Nρ−1i nρ−1 ≤ ρ αρ−1n
and
Φn,m ≤ ∑ i∈[k]
( N
(n) i
n )ρ + ∑ i∈[k] ( N (m) i m )ρ ≤ 2 (( 1 α )ρ + ( 1− 1 α )ρ) define τ := max (τn1 , τn2 , τm1 , τm2) and Φ := max (Φn1,m1 ,Φm2,n2)
We have
Z ≤ 4nτ2Φ2K2max ≤ 16 ρ2
α2ρ−2n
(( 1
α
)ρ + ( 1− 1
α
)ρ)2 K2max = Zb
Plotting the Zb value with respect to variety value of ρ and α in Fig. 1, we can see that for α = 1, the bound will be looser given larger ρ. However, for α larger than around 1.25, which means the supi(Ni) is slightly smaller than the sample size, the bound will become tighter when ρ is large. This illustrate the benefit of using power-MMD with larger ρ in practice.
We could also get a tighter bound according to Corollary 3.1. Practically, it will be
much more beneficial to calculate the τn = supi∈[k]
( ρ
N (n) i
gρ(N (n) i , n) ) and Φn,m =∑
i∈[k] ∣∣∣gρ(N (n)i , n)− gρ(N (m)i ,m)∣∣∣ on the fly. That is to say, we do not estimate the convergence bounds before we receive the samples, instead, the calculation of the bounds is carried out together with the calculation of default-plugin-estimators. Remarkably this is still a distribution-free bounds, since we make no assumptions on the probability functions we apply hypothesis test upon.
However, the issue is although Z decreases when ρ increases, GPKE [F , φρ, X, Y ] also decreases when ρ increases. It is not clear how ρ will influence the value of ∣∣∣GPKE [F,φρ,X,Y ]√ Z ∣∣∣.

Review summaries:
The paper under review proposes to generalize MMD for discrete random variables whose labels take value in $\mathbb{R}^k$. They propose to estimate these generalized probability kernel distance using empirical estimators. Their properties are studied for two particular examples, namely a kernelized Stein discrenpancy and polynomials versions. Consistency and bias of both estimators are studied and bias corrected.