The Priestley-Chao Estimator of Conditional Density with Uniformly Distributed Random Design

283 98 (3) STATISTIKA 2018 The Priestley-Chao Estimator of Conditional Density with Uniformly Distributed Random Design Kateřina Konečná 1 |Brno University of Technology, Brno, Czech Republic 1 Institute of Mathematics and Descriptive Geometry, Faculty of Civil Engineering, Žižkova 17, 602 00 Brno, Czech Republic. E-mail: [email protected]. 2 This article is based on contribution at the conference Robust 2018. Abstract The present paper is focused on non-parametric estimation of conditional density. Conditional density can be regarded as a generalization of regression thus the kernel estimator of conditional density can be derived from the kernel estimator of the regression function. We concentrate on the Priestley-Chao estimator of conditional density with a random design presented by a uniformly distributed unconditional variable. The statistical properties of such an estimator are given. As the smoothing parameters have the most significant influence on the quality of the final estimate, the leave-one-out maximum likelihood method is proposed for their detection. Its performance is compared with the cross-validation method and with two alternatives of the reference rule method. The theoretical part is complemented by a simulation study.2 Keywords Priestley-Chao estimator of conditional density, random design, uniform marginal density, bandwidth selection, maximum likelihood method, reference rule method JEL code C14 INTRODUCTION Kernel smoothing is still a popular non-parametric procedure, in theory as well as in practice. There are numerous monographs concerned with the kernel smoothing approach, e.g., Wand and Jones (1994). Computational implementations in MATLAB were developed by Horová et al. (2012). The present paper focuses on the kernel conditional density estimation. Several estimator types can be found in the literature with the Nadaraya-Watson one being probably best known (see Rosenblatt, 1969). The local linear estimator of conditional density was suggested by Fan et al. (1996) for its better statistical properties and boundary effects. Conditional density can be regarded as a generalization of regression, which models the conditional mean while conditional density models the whole distribution. This is the reason why a kernel regression estimator can be generalized to a kernel conditional density estimator. The present paper extends the ANALYSES 284 Priestley-Chao regression estimator (for detailed information see Priestley and Chao, 1972) to estimate even conditional densities. Each kernel estimator depends on the smoothing parameters called bandwidths, values which significantly influence the final estimation. This is the reason why so much importance is given to their selection. There are many methods discussed in the literature, most of them suggested for the NadarayaWatson estimator, with only a few of them for the local linear estimator. Introduced by Fan and Yim (2004), Hansen (2004) and Hall et al. (2004) and based on minimizing the Integrated Squared Error, cross-validation is a method typical of bandwidth selection. Bashtannyk and Hyndman (2001) suggested a reference rule method for normal underlying conditional density and for two marginal density choices – normal and uniform. Some methods extend the methods suggested for kernel regression. The iterative method proposed by Konečná and Horová (2014), for one, is motivated by the iterative method developed for kernel density estimation and for kernel regression (for detailed information, see Horová and Zelinka (2007), Horová et al. (2012), Koláček and Horová (2012)). Other examples include the bootstrap method by Bashtannyk and Hyndman (2001) and Fan and Yim (2004) as well as the fast dual-tree based algorithms using a maximum likelihood criterion (see Holmes et al., 2012). Kernel conditional density estimation is still employed in practice: Takeuchi et al. (2009) show its application in medicine (the relative change in spinal bone mineral density is explored as a function of the age of adolescents), Jeon and Taylor (2012) are interested in 1-to-72-hours ahead wind-power prediction from which the management of wind farms and electricity systems can profit. Another application, forecasting electricity smart meter data, helping consumers to analyze and to minimize their electricity consumption and enabling new pricing strategies for suppliers, is introduced by Arora and Taylor (2016). As mentioned above, papers are focused primarily on the Nadaraya-Watson or the local linear estimator. The present paper suggests the Priestley-Chao estimator for the uniformly distributed design, based on the estimator suggested for the equally spaced design (see Konečná, 2017). The leave-one-out maximum likelihood method follows the one proposed by Konečná (2018). The paper is organized as follows: Section 1 deals with the Priestley-Chao estimator of conditional density and its statistical properties. The optimal values of the smoothing parameters are derived, and the leave-one-out maximum likelihood method for their practical estimation proposed in Section 2. This method is complemented by the cross-validation method and by two alternatives of the reference rule method. A simulation study in Section 3 then presents the performance of the methods by a simulation study. The proofs of the statistical properties can be found in the Appendix. 1 THE PRIESTLEY-CHAO ESTIMATOR OF CONDITIONAL DENSITY The conditional density f(y|x) models the probability of a random variable Y given a random variable X, represented by a fixed observation X = x. Let {(Xi,Yi), i = 1,…,n} be an observed data sample of a pair of real random variables (X,Y). The kernel estimate of conditional density generally takes the form: (1) where wi(x) is a weight function, and K is a real, symmetric, nonnegative kernel function satisfying: (2) The present paper uses the Gaussian kernel. The smoothing parameters hx > 0, hy > 0 control the smoothness of the estimate. The estimate of conditional density is also influenced by the estimator type (1). 285 98 (3) STATISTIKA 2018 Our focus is on the Priestley-Chao estimator, originally proposed for the kernel regression estimation (Priestley and Chao, 1972). Konečná (2017) dealt with the Priestley-Chao estimator for conditional density with the fixed design, i.e., the fixed values of the design variable X were assumed. Next, we are concerned with the estimator for the random design specified by a uniformly distributed variable X on the interval [0,1]. The estimator can easily be extended for the design variable X on the interval [a, b], a < b. The statistical properties of the estimator will be given and methods for bandwidth detection proposed. Let X be a uniformly distributed random variable with the marginal density function: As the focused weight function in Formula (1) is , the Priestley-Chao estimator takes the form: (3) The Priestley-Chao estimator of the regression function is expressed by the conditional mean of the Formula (3): Theorem 1 Let X be a uniformly distributed random variable on the interval [0,1], Y be a random variable with density f(y|x) being at least twice continuously differentiable, and K(x) be a real, symmetric, nonnegative kernel function satisfying (2). For , the asymptotic bias (AB) and the asymptotic variance (AV) are given by: where: R2(K) = ∫RK2 (u) du. Proof. The proof can be found in the Appendix. The local quality of the estimate at the point [x,y] is given by the mean squared error (MSE) which is the simple decomposition to the variance (V) and the squared bias (SB). Considering the main terms only, the asymptotic MSE (AMSE) is obtained as: (4) The statistical properties of the Formula (3), particularly the global quality measure expressed by the asymptotic mean integrated squared error (AMISE), are necessary for assessing the quality of the estimate and the theoretical values of the smoothing parameters. AMISE is obtained by integrating (4) weighted by the marginal density g(x) as: ANALYSES 286 The following form of the AMISE is more succinct for further processing: (5) where the constants c1, c2, c3, c4 are given by: (6) Remark. Note that all the integrals with respect to x are computed over the support of the X variable, i.e., over the interval [0,1]. The integrals with respect to y are considered over R. 2 METHODS FOR BANDWIDTH SELECTION The values of the smoothing parameters have an essential significance for the final estimate of conditional density. First, the optimal widths of the smoothing parameters are derived as the values minimizing the AMISE. As the optimal bandwidths depend on the true conditional and marginal density function, it is necessary to develop a data-driven method for their estimation. In this section, the leave-one-out maximum likelihood method, the cross-validation method, and two alternatives of the reference rule method are suggested for their detection. 2.1 Optimal values of the smoothing parameters The optimal values of the smoothing parameters are given as the values which minimize AMISE given by (5). By differentiating (5) with respect to hx and hy and setting the derivatives to 0, we obtain the following system of non-linear equations: (7) Solving system (7), the optimal bandwidths are given by: (8) Both hx * and hy * are of order n–1⁄6 while the order of AMISE is n–2⁄3. 2.2 The leave-one-out maximum likelihood method As mentioned above, with a real dataset, a data-driven method is needed for bandwidth selection. We will modify the maximum likelihood method, which is a standard statistical procedure for estimating 287 98 (3) STATISTIKA 2018 unknown parameters. This method was originally proposed for kernel density estimation by LeivaMurillo and Artes-Rodriguez (2012), and their approach is generalized to include the Priestley-Chao estimator of conditional density. Since the objective function: (9) is considered for all n observations, the optimization problem has a trivial solution. If i = j in (9), the objective function (9) increases to infinity for hx → 0 and hy → 0. Of course, this is not the desired behaviour because, with very small values of bandwidths, the final estimate tends to be undersmoothed. This problem can be solved by leaving out one observation and employing the modified objective function: (10) If the natural logarithm of the likelihood function L* given by (10) is taken into account, the values of the smoothing parameters maximize: and are developed as: 2.3 The leave-one-out cross-validation method The cross-validation method is a standard procedure for bandwidth selection in kernel smoothing. Introduced by Fan and Yim (2004), Hansen (2004) and Hall et al. (2004), the method is associated with the global quality measure of the estimator, with the integrated squared error (ISE). With being the estimate at the point (Xi,Yi) using the points {(Xj,Yj), j ≠ i}, the cross-validation function: is the proper estimator of the ISE. The values of the smoothing parameters are given by: 2.4 The reference rule method The reference rule method was originally proposed for the Nadaraya-Watson estimator by Bashtannyk and Hyndman (2001). They assumed a normally distributed random variable Y|(X = x) with linear conditional mean and constant or linear standard deviation. Additionally, they distinguished two possibilities for the marginal distribution, considering uniform and truncated normal marginal densities. Our approach via the Priestley-Chao estimator corresponds to the choice of a uniform marginal density. We assume that the conditional distribution is normal with the mean m(x) and standard deviation σ(x). Hence, the conditional density of Y|(X = x) is: ANALYSES 288 (11) Two different situations are considered: (a) According to Bashtannyk and Hyndman (2001), the model with the linear conditional mean m(x) = p0 + p1 x and the linear standard deviation σ(x) = q0 + q1 x is assumed. The values of the constants c1, … , c4 are given by the expressions: (12) where . The values of the smoothing parameters are obtained by substituting (12) into (8). (b) The model with the quadratic conditional mean m(x) = p0 + p1 x + p2 x2 and the constant standard deviation σ is suggested. The constants c1, … , c4 are given by: (13) The values of the smoothing parameters are obtained by substituting the terms (13) into (8). Remark. The expressions (13) are obtained by differentiating (11) twice and substituting them into (6). As the computations of the integrals in (6) include many auxiliary derivations, only a sketch of them is presented. A conditional random variable Y|(X = x) ~ N(m(x), σ2) with the density function f(y|x) is assumed. The following equality: where f1(y|x) is a density function of a conditional random variable Y1|(X = x) ~ N(m(x), σ2), was used. For evaluating the integrals in (6), the below auxiliary expression was derived: 289 98 (3) STATISTIKA 2018 where A = R · [0,1] is the domain of integration, and ∙ denotes the ceiling function. 3 SIMULATION STUDY In this section, a simulation study comparing four methods for bandwidth estimation is conducted. The considered methods are the maximum likelihood method (ML), the cross-validation method (CV), the reference rule method with linear conditional mean and linear standard deviation (REF1), and the reference rule method with quadratic conditional mean and constant standard deviation (REF2). Two models are involved in the simulation study. To demonstrate the adaptability of the methods to various shapes of the regression function or conditional density, a changing shape of the conditional mean is presented in the first model. In the second model, a bimodal and non-symmetric conditional distribution is chosen as a mixture of two normal densities. The models are defined as: In both simulation studies, one hundred observations (n = 100) were generated. The described methods for bandwidth selection are compared from several points of view – the estimates of the smoothing parameters and the quality measure of the estimate. The quality of the final estimate was measured by an estimate of the integrated squared error (ISE) given by: ANALYSES 290 where y = (y1, … , yN) is a vector of equally spaced values over the sample space of Y and ∆ is the distance between two consecutive values of y, i.e. ∆ = yj+1 – yj, j = 1, … , N – 1. In both simulation studies, the number of y values was set to N = 100. For both models, five hundred repetitions have been made to obtain the described characteristics. The results are displayed in boxplots and supplemented by numerical values in the text. First, the results for the model M1 are summarized. A scatterplot of one set of the sample values of model M1 is displayed in Figure 1. Boxplots of the estimates of the smoothing parameters hx and hy are displayed in Figure 2. It can be seen, that the ML and CV methods lead to good values for hx (medians are 0.0613 for ML and 0.0288 for CV), the REF1 and REF2 methods produce highly variable bandwidths exceeding the optimal value hx * = 0.0455. Considering the estimates of the smoothing parameter hy, the values estimated using REF1 resemble those with REF2. Their medians 0.2927 and 0.2863 (in this order) are close to the optimal value hy * = 0.3213, and the estimates are characterized by a low variability (their standard deviations are close to 0.02). The ML method gives slightly higher values than the optimum hy *, the median of the values being 0.3631. The CV method tends to produce values well under the optimal value which results in a much undersmoothed estimate of conditional density. Figure 1 A scatterplot of one hundred observations of model M1 Figure 2 Estimates of the smoothing parameters hx and hy along with the optimal values (horizontal lines) for the ML, CV, REF1, and REF2 methods in model M1 Source: The author's own construction Note: The log-scale of the vertical axis in the left-hand-side panel. Source: The author's own construction Figure 3 Estimates of the ISE values (expressed in the logscale) for the ML, CV, REF1, and REF2 method and for the optimal bandwidth choice (OPT) in model M1 Source: The author's own construction 291 98 (3) STATISTIKA 2018 The ML method results in an ISE that is smaller than any other methods considered as well as the optimal bandwidth choice (OPT). Both reference rule methods produce ISE values slightly higher than OPT while the values of the CV method are well above the optimal bandwidth choice. Now, we focus on the results of model M2. The scatterplot of a sample generated by model M2 is shown in Figure 4. Boxplots of bandwidth estimates are shown in Figure 5. Both reference rule methods produce values of hx and hy well above the optimal values, which results in an overvalued final estimate with worse capability to adapt to bimodal conditional density. On the other hand, the low values of both smoothing parameters obtained by the CV method lead to the undersmoothing of conditional density and abundance of useless information in the data. The ML method performs the best in this simulation study. The ML method is also suitable in terms of the ISE (see Figure 6). The medians of the ISE values obtained by the reference rule methods (0.0555 for REF1 and 0.0527 for REF2) do not reach the median (0.0742) of the ISE values for optimal bandwidths, but the REF1 method suffers from the large variability (standard deviation is 0.3143). The CV method provides highly variable ISE values exceeding the OPT values. Figure 4 A scatterplot of 100 observations by model M2 Figure 6 Estimates of the ISE values (expressed in the log-scale) for the ML, CV, REF1, and REF2 methods and for the optimal bandwidths (OPT) in model M2 Figure 5 Estimates of the smoothing parameters hx and hy and the optimal values (horizontal lines) for the ML, CV, REF1, and REF2 method in model M2 Source: The author's own construction Source: The author's own construction Source: The author's own construction