A Novel Nonparametric Distance Estimator for Densities with Error Bounds

entropy 2013, 15, 1609-1623; doi:10.3390/e15051609 entropy ISSN 1099-4300 www.mdpi.com/journal/entropy Article A Novel Nonparametric Distance Estimator for Densities with Error Bounds Alexandre R.F. Carvalho 1,*, João Manuel R. S. Tavares 1 and Jose C. Principe 2 1 Instituto de Engenharia Mecânica e Gestão Industrial, Faculdade de Engenharia, Universidade do Porto; Rua Dr. Roberto Frias, s/n 4200-465 Porto, Portugal; E-Mail: [email protected] 2 Computational Neuro Engineering Laboratory, University of Florida; EB451 Engineering Building, University of Florida, Gainesville, FL 32611, USA; E-Mail: [email protected] * Author to whom correspondence should be addressed; E-Mail: [email protected]; Tel.: +351-225-081-719; Fax: +351225-081-440. Received: 19 December 2012; in revised form: 25 April 2013 / Accepted: 28 April 2013 / Published: 6 May 2013 Abstract: The use of a metric to assess distance between probability densities is an important practical problem. In this work, a particular metric induced by an α-divergence is studied. The Hellinger metric can be interpreted as a particular case within the framework of generalized Tsallis divergences and entropies. The nonparametric Parzen’s density estimator emerges as a natural candidate to estimate the underlying probability density function, since it may account for data from different groups, or experiments with distinct instrumental precisions, i.e., non-independent and identically distributed (non-i.i.d.) data. However, the information theoretic derived metric of the nonparametric Parzen’s density estimator displays infinite variance, limiting the direct use of resampling estimators. Based on measure theory, we present a change of measure to build a finite variance density allowing the use of resampling estimators. In order to counteract the poor scaling with dimension, we propose a new nonparametric two-stage robust resampling estimator of Hellinger’s metric error bounds for heterocedastic data. The approach presents very promising results allowing the use of different covariances for different clusters with impact on the distance evaluation. Keywords: generalized differential entropies; generalized differential divergences; Tsallis entropy; Hellinger metric; nonparametric estimators; heterocedastic data PACS Codes: 02.50.-r, 02.50.Cw, 89.70.-a, 89.70.Cf OPEN ACCESS Entropy 2013, 15 1610 1. Introduction Distances measures between two probability densities have been extensively studied in the last century [1]. These measures address two important main objectives: how difficult it is to distinguish between one pair of densities in the context of others and to assess the closeness of two densities, compared to others [2]. In learning scenarios essentially associated with the test of a single hypothesis, the use of a divergence to represent the notion of distance is efficient. However, in scenarios involving multiple hypotheses, such as clustering, image retrieval, or pattern recognition and signal detection, for instance, the non-symmetric and non-metric nature of divergences becomes problematic [3]. When deciding the closest or the farthest among three or more clusters, the use of a metric is important. In this work, a novel nonparametric metric estimator for densities with error bounds is presented. Shannon’s entropy has a central role in information-theoretic studies. However, the concept of information is so rich that perhaps there is no single definition that will be able to quantify information properly [4]. The idea of using information theory functional, such as entropies or divergences, in statistical inference is not new. In fact, the so-called statistical information theory has been the subject of much research over the last half century [5]. How to measure the distance between two densities is an open problem with several proposals since the work of Hellinger in 1909 with Hellinger’s distance [1], Kullback and Leibler (1951), with Kullback-Leibler’s divergence [6], Bregman (1967) with the Bregman’s divergence [7], Jeffreys (1974) with J-distance [8], RAO (1985) and Jianhua Lin (1991) with Jensen-Shannon’s divergence [9,10], Menéndez et al. (1997) with (h, Φ)-entropy differential metric [11], Seth and Principe (2008) with correntropy [12], among others. This work looks into Hellinger’s metric that is the preferred [13,14] or natural model metric [15]. In 2007 Puga found that Hellinger’s metric is one particular α-divergence [16]. Here, we propose a new measure change to solve the nonparametric metric estimation and a two stage robust estimator with error bounds. 2. Theory Background Following Hartley’s (1928) and Shannon’s (1948) works [17,18], Alfred Rényi introduced in 1960 the generalized α-entropy [19] of probability density function   xf :     dxxffR    ln 1 1, α (1) The corresponding generalized differential divergence between two densities  xf1 and   xf2 is:    1 12 1 2 1 ,ln 1 Rfx D ff dx fx      (2) Gell-Mann and Tsallis considered another family of α-entropies [20]:        dxxffT   1 1 1 (3) being the corresponding α-divergences given as:    1 12 1 2 1 ,1 1 Tfx Dff dx fx                (4) Entropy 2013, 15 1611 Making 1  , one easily can conclude that:      fHfTfR S        11 limlim (5) and:      12 12 12 11 lim , lim , , RT KL D ff D ff D ff     (6) where:        dxxfxffHSln (7) is Shannon’s differential entropy and:      dx xf xf xfffDKL 2 1 121 ln, (8) is Kulback-Leibler’s divergence. Another member of these families is the Rényi’s quadratic entropy (α=2) that is defined as:      dxxffR 2 2ln (9) while the respective divergence is:    2 1 212 2 ,ln Rfx Dff dx fx  (10) Rényi’s quadratic entropy, given by Equation (9), is particularly interesting because it accepts a close form nonparametric estimator, saving computational time compared to numerical integration or resampling [21,22]. α-Entropy families given by Equations (1) and (3) are monotonically coupled (Ramshaw [23]) through:     111 R Te       (11) Therefore, an optimization in one family has equivalence in the other. 2.1. Square-Root Entropy Let us consider 21  in Equations (1)–(4). Then, the square-root entropy in the form of Tsallis is:   12 22Tf fxdx    (12) with the corresponding divergence given as:     12 1 2 1 2 ,22 T Dff fxfxdx  (13) In Rényi’s form one finds, respectively, the entropy:   12 2ln , R ffxdx  (14) and the divergence:   12 1 2 1 2 ,2ln R Dff fxfxdx  (15) It should be noted that, from Equation (13), one obtains: Entropy 2013, 15 1612       21 2 212121 ,, ffMdxxfxfffDT  (16) where  21,ffM is a information theoretic derived metric that, among other properties, verifies the triangular inequality. This particular -divergence, by means of a monotonous transformation, induces the Hellinger’s distance, which is a metric [13,14,24]:  12 12 ,22(,) M ff Iff (17) On the other hand, information theoretic derived metrics given by Equations (15) and (16) are also related with Hellinger’s affinity or Bhattacharya’s coefficient (   1,0 21   ffI ):   12 1 2 ,. I ff fxfxdx  (18) Considering the expected cross-value of two probability density functions  21,ffC :         12 12 2 1 1 2 ,ff C f f E f E f f xf xdx  (19) the Hellinger’s affinity given by Equation (18) can be then written as:     12 12 12 ,, , I ff Cff f xdx CffHf     (20) where  xf is the normalized product density:       12 12 , f xf x fx Cff  (21) and   fH the corresponding entropy of the information theoretic derived metric. This metric has bounds that can be directly computed from the samples as shown by Puga [16]. These bounds often present overlapping hypothesis intervals, and resampling estimation is a necessary tool to remove ambiguities and access distances between densities. 2.2. Nonparametric Hellinger’s Affinity Estimation Let us focus on the application of the previous measures on two Parzen’s nonparametric densities [25] from two data clusters       11 2 1 1 )1( 1 ,...,, N xxxCl  and         22 2 2 1 )2( 2 ,...,, N xxxCl :       1 1 1 1 1 1,, 1N j j xxG N xf  (22) and:    22 22 1 2 1,, N j j fx Gx x N    (23) where    ,,xG is the Parzen’s Gaussian kernel, also known as kernel bandwidth, with the approximation of covariance I 2  , and mean  given as:    2 2 () 2 2 1 1 ,, 2 xi i i Gx e         (24) Entropy 2013, 15 1613 where  is the dimension. Notice that the two clusters in Equations (22) and (23) may have two different Gaussian kernel covariances. The kernel covariance may be obtained directly from the a priori knowledge of the instruments used to produce the data; for instance, two different instruments with different precisions may produce the same data, but the densities should reflect the measurements error through the bandwidth, (covariances). To estimate the bandwidth without instrumental apriori knowledge it is possible to estimate the kernel bandwidth with a suitable method, such as k-Nearest Neighbor (k-NN), Silverman [25] or Scott [26]. Now, let us adopt the summing convention  12 11, N i N jji and define the following auxiliary variables:   222 12 2   (25) 2222 *12 2     (26)     12 22 2 ,2 1 2 ij i j sxx    (27)       2 21 ,jiji xxd  (28) and:  22 ,, 22 , , ij kl dd Dij kl Fd e e    (29) The nonparametric estimator    f ˆ results in:      2 ,*, , ˆ,, D ij ij ij f Fd G s   (30) 2.3. The Resampling Estimator The bootstrap resampling is reached through the distribution of probability given by Equation (30) combined with the random generation of samples   k  from nonparametric Parzen’s density with diagonal covariance, which is a well-established as well as a computationally efficient procedure [27]. Then, the synthesized samples are directly usable in the estimator:         1 111 lim lim ˆ K fK KK kk f Hf f d d f EHf K ff                         (31) with  f dii k...~  . However, the use of Equation (31) is associated with serious practical difficulties because the second moment:  2 dHf      (32) Entropy 2013, 15 1614 has infinite variance, which is a condition where the central limit theorem is not valid. In this work, we use measure theory and propose the following change of measure:   zf    (33) with the associated density )(zfZ:                      max 0 )( ˆ 1 lim 11 k Z z kz k K ff dzzf zz E f E  with Z dii kfz ...~ (34) This new density presents a finite second moment:  max 2 0 1() z z f zdz H f z   (35) having z f a limited support between 0 (zero) and max z. This is a density with an abrupt jump in max z end of the density. However, the approximation properties of a histogram are not affected by a simple jump at the end of the density [26], hence the histogram estimator was used to estimate )( ˆkZ zf with . The product probability density function   z f must be estimated from the random variable    fz , but it ensures finite variance, which is a requisite of the central limit theorem and the t-student confidence interval may be used Equation (39). Figure 1. A logarithmic scale for the 21 /   coefficient variation and metric measure change between the two respective densities   21,ff .  k k fz   ˆ Entropy 2013, 15 1615 To test the algorithm, we consider the simplest case of Hellinger’s metric (17) associated with the nonparametric densities of Equations (22) and (23). In this particular case, we have access to the analytical value of Hellinger’s metric:  2 1,2 2 2 12 12 22 12 2 ,22 d Mff e            (36) Using Equations (34) and (36), we can quantify the computational behavior of the resampling estimator. Let us first consider the behavior of the Parzen’s density estimator with two distinct kernel sizes: 12 ,   . In the simplest case with only two kernels, located at the same coordinates, despite the same location, different Parzen’s windows in Equation (36) provide different distances, as can be observed in Figure 1. It is possible to verify the symmetric behavior of the distance estimator and realize that the bandwidth of the Parzen’s kernel is important to access the distance between clusters. This is a relevant characteristic, especially when the experimental data have different instrumental origins with different measurement precisions; the use of different bandwidths in the Parzen’s kernels may reflect this importance feature of the density, and this implies that the data is heterocedastic. To quantify the error bounds estimation performance, we propose the generation of N1 distance samples m M ~ from resampling the density of Equation (34), and to estimate Z f ~ we use a discrete histogram with N1 bins, obtaining the ordered k z and k Z f ~ . As such, the metric estimator becomes:  1 ~ 12 1... 1 22 , () mk Z kmN MCfffzz z         (37) which can be written as:  max ~ 12 0 1 22 , () k z k Z k M Cff f z z z     (38) To assess the error bounds estimation, we use the t-student 95% confidence interval (39), which is a maximum entropy distribution [28,29] and provides a parametric approach to robust statistics [30], and allows the following calculation of the confidence limits:   1 1 2 1,0.5 0.95/ 2 1 11 1 ,(1) N Nm k LU M t M M NN       (39) We calculate the 95% confidence limits, the upper (U) and the lower ( L ) for the respective density resampling. The variance of this new estimator is well controlled in one dimension. The unexpected drawback of this estimator is its poor scaling performance with increased dimension, as depicted in Figure 2. The new variable    fz may be seen as a projection of the multidimensional Parzen’s kernels into a 1-Dimensional function. This insight allowed the design of a two-stage estimator for     f that circumvents both problems: infinite variance and poor scalability with dimensionality. M ~ Entropy 2013, 15 1616 Figure 2. Illustration of the metric estimator behavior for dimensions 1 (one) to 5. 2.4. The Two Stage Resampling Estimator We propose the generation of N1 distance samples )( ~ n k M from resampling the density )( k f  , which constitutes one trial )(n:   1 12 () () 1... , 22 n kn kkN Cff M f          (40) It is possible to estimate )( ~ n M as:   112 () () 1 1 , 2 2 N n n kk Cff MNf     (41) For each trial )(n, the 95% confidence limits, the upper )(n Uand the lower )(n L for the respective density resampling, can be calculated:  1 1 2 () () () () () 1,0.5 0.95/2 1 11 1 ,(1) N nn n n n Nk k LU M t M M NN         (42) It may seem that this step is enough to estimate the metric  ),( ~ 22, ~ 2121 ffIffM  , but the theoretically predicted undesired behavior associated to Equation (32), with large confidence intervals is present in this estimator. To demonstrate this drawback, we have simulated 100 trials of the simplest case of nonparametric Hellinger’s metric, Equation (36), with Euclidean distance 1 21  xxd . As can be observed in Figure 3, the large confidence intervals are present, hence the motivation for the two-stage error bound estimator. Entropy 2013, 15 1617 Figure 3. t-Student 95% confidence intervals for Hellinger’s metric defined by dots; the exact value is represented by a continuous line; the predicted large intervals are marked with triangles; and the miss-estimated intervals are marked with circles. To achieve a robust error bound estimator   RR UL ~ , ~ for the expected value of ),( 21 ffM with similar results of )(zfZ in one dimension, we propose a new two-stage method. Comparing the results of the two densities resampling, we found that 31 selected trials out of 33 from )( ˆk f   was in good agreement with )( ~zfZ. With 33 trails )(n generated with N1 random samples each as:   1 12 () 1... 1...33 , 22 ˆn kkN n Cff f                (43) sorting the amplitude )()( nn LU  and keeping the 31 smallest intervals with the correspondent estimated affinities ( )( ~ n s M), we obtain the estimator for the second stage with: 31 () 12 1 1 (, ) 31 n ss n Mff M    (44) Then, we calculated the respective t-student 95% confidence interval   ss UL ~ , ~ with the selected trials )( ˆn s M. To overcome the miss-estimated intervals, we have defined a second estimator for the lower limit of the interval ( 2 ~ L) and a second estimator for the upper limit of the interval ( 2 ~ U): 31 31 () () 22 11 11 ,, 31 31 nn ss nn LU L U             (45) which is a potentially asymmetric interval, guided by the selected first-stage interval limits. The robust estimator for the lower limit of the interval ( R L ~ ) and the robust estimator for the upper limit of the interval ( R U ~ ) were defined as:   22 ,min,,max, RR s s LU LL UU         (46) ),( ~ 21 ffMs