Density estimation is indeed the process of constructing an estimate of the probability density function (pdf) from an available data. This estimation not only represents the data distribution but also provides summary statistics such as the mean, median, variance, moments and quantiles. Furthermore, density estimates provide information about distribution characteristics, including skewness, kurtosis, and multimodality within the data. The estimation of pdf is a fundamental concept in statistics and a widely researched topic. There are two commonly methods for density estimation: parametric and nonparametric methods. The parametric method assumes that the data is drawn from a known distribution, whereas the nonparametric method aims to estimate the density function directly from the data. Several nonparametric density estimation techniques commonly include the histogram, naïve density estimator, nearest neighbor method, and orthogonal series estimator. Kernel density estimation (KDE) is a widely used nonparametric density estimation. The KDE relies on the kernel function which determines the weight assigned to each data point, and the bandwidth which controls the smoothness of the estimate. Hence, the selection of the bandwidth is the most crucial in the context of kernel density estimation.

The selection of the bandwidth is a critical issue that arises in the context of KDE, as the performance of the KDE depends on the chosen bandwidth. A small bandwidth value results in an undersmoothed density, while a large bandwidth value leads to an oversmoothed density (Gramacki, 2018). There are various methods to determine a bandwidth for KDE. The primary categories of bandwidths for KDE include rules-of-thumb (ROT), cross-validation (CV), and plug-in (PI). The plug-in method has been demonstrated to offer excellent performance in many cases. Due to its demonstrated excellent performance, the plug-in method is the common first choice in practical applications. However, there is still room for further improvement in its implementation (Wand and Jones, 1995). Plug-in bandwidths are operated on the straightforward concept of substituting estimated values of the unknown quantities into formulas for achieving the asymptotically optimal bandwidth.

In order to overcome optimal bandwidth, this study introduces the proposed plug-in bandwidth. This bandwidth leverages the first kind shifted Chebyshev polynomials, providing a solution to the problem. The effectiveness of the methods relies on the estimation of integrated squared density derivative functionals, a subject that has been explored by many researchers. Silverman (1986) provided a comprehensive overview of density estimation techniques, including discussions on bandwidth selection and the role of integrated squared density derivative functionals. Sheather and Jones (1991) discussed a data-driven method for bandwidth selection in kernel density estimation, which related to integrated squared density derivative functionals. Raykar and Duraiswami (2006) developed the algorithms for estimating density derivatives using the univariate Gaussian kernel. These algorithms are utilized to calculate the optimal bandwidth for kernel density estimation. Tenreiro (2011, 2020) proposed direct plug-in bandwidth for the KDE based on the Fourier series and the Hermite series. In a recent study, Dharmani (2022) introduced a bandwidth selection by employing the near Gaussian assumption. This assumption enables the use of the Gram-Charlier A series as an approximation to the function for the purpose of estimating its density derivative. The objective of this paper is to derive a bandwidth by using the first kind shifted Chebyshev polynomials as an approximation to the density function. This is aimed at estimating the integrated squared density derivative functionals.

The remaining sections of this article are organized as follows. Section 2 provides an overview of the fundamental properties of kernel density estimation. In Section 3, various methods for bandwidth selection are discussed. These methods include least squares cross-validation bandwidth, an improved version of rules of thumb bandwidth, and the Sheather and Jones plug-in bandwidth. Section 4 offers a brief definition of the first kind shifted Chebyshev polynomials, then utilizes them to approximate the underlying density function, and finally presents the proposed plug-in bandwidth based on this estimator. Section 5 presents a simulation study of the proposed bandwidth, examining its performance under different distributions and sample sizes using the R programming language. Additionally, Section 6 applies the proposed bandwidth to real dataset. Finally, Section 7 concludes the article with a summary of the findings.

The kernel density estimator for a random sample X₁, X₂, . . . , X_n drawn from a common and typically unknown density f (x), as defined by Rosenblatt (1956) and Parzen (1962), is expressed as

where K (x) is the kernel function and h is the bandwidth with positive value.

The kernel function K (x) plays as the weight function and satisfies the following properties: K (x) ≥ 0, K (x) = K (−x) , ∫ K (x) dx = 1, ∫ xK (x) dx = 0 and k₂ = ∫ x²K (x) dx ≠ 0. The bandwidth h determines the level of smoothness of the density estimate.

In practice, it is common to consider a global error criterion that measures the distance between the estimated density function f̂ (x; h) and the true density function f (x). One such error criterion is the integrated squared error (ISE) given by ISE (f^(x;h))=∫-∞∞[f^(x;h)-f(x)]2 dx. A more appropriate approach would involve analyzing the expected value of this quantity, known as the mean integrated square error (MISE), which is defined as

where k₀ = ∫ (K (x))²dx and k₂ = ∫ x²K (x) dx (Gramacki, 2018).

The assumptions are f (x) is assumed to be sufficiently smooth: Its second derivative f″(x) is bounded, continuous and square integrable. Also, if (nh)⁻¹ → 0 and h → 0 as n → ∞ then MISE f̂ (x; h) → 0. It obtains the asymptotic mean integrated square error (AMISE) as follows:

where θ₂ = ∫( f″ (x))²dx, k₀ = ∫ (K (x))²dx and k₂ = ∫ x²K (x) dx.

Several methods are available for determining appropriate bandwidth for KDE. The three main types of bandwidths are as follows: Cross-validation (CV), rules-of-thumb (ROT) and plug-in (PI) (Gramacki, 2018). Cross-validation involves techniques like least squares cross-validation (LSCV), biased cross-validation (BCV), and smoothed cross-validation (SCV). Rules-of-thumb includes approaches such as Silverman’s rule of thumb and its improved version. Plug-in methods have been explored by various authors, including Park and Marron (1990), Sheather and Jones (1991), and Hall et al. (1991), among others. This article provides concise explanations of the following methods chosen for study: least squares cross-validation, the improved version of Silverman’s rule of thumb, and Sheather and Jones plug-in bandwidth.

3.1. Least squares cross validation

A well-known method for selecting the bandwidth is least squares cross-validation (LSCV), proposed by Rudemo (1982) and Bowman (1984). The main objective is to find the optimal bandwidth h that minimizes the ISE using the estimator f̂ (x; h) for density f (x). The integrated squared error of f̂ (x; h) is represented as

(Silverman, 1986; Wand and Jones, 1995).

The first term ∫ f̂ (x; h)²dx of (3.1) can be calculated from the data which was proved by Härdle (1991) as

where K * K (x) is the convolution of K (x).

The second term ∫ f̂ (x; h) f (x) dx of (3.1), which depends on h and involves the unknown density f (x), has to be estimated. Notice that, ∫ f̂ (x; h) f (x) dx is the expected value of f̂ (x; h) which can be estimated by

where f^-i(Xi;h)=(1/((n-1)h)) Σj≠inK((x-Xj)/h) is the estimator based on the sample with X_i deleted.

Finally, the last term ∫ f (x)²dx of (3.1) is independent of the bandwidth h. Therefore, the last term can be moved to the left side of the equation and can be written as

As a result, Equations (3.2) and (3.3) are inserted into the Equation (3.4), leading to the least squares cross-validation function as

(Härdle et al., 2004). The bandwidth that minimizes the function LSCV (h) is denoted by h_LSCV.

3.2. Rules of thumb

The optimal bandwidth minimizes AMISE f̂ (x; h) with respect to h, and solving the first partial derivative with respect to h yields

where θ₂ = ∫( f ″(x))²dx, k₀ = ∫ (K (x))²dx and k₂ = x²K (x) dx.

The rule-of-thumb bandwidth is determined by replacing the density function f (x) with the normal distribution having zero mean and variance σ² and using the Gaussian kernel function in Equation (3.6). The rule-of-thumb bandwidth, denoted by h_ROT1, is calculated using the formula

The rule-of-thumb bandwidth is sensitive to outliers, which cause an overestimation of σ and lead to a larger bandwidth. To make the estimator more robust, the interquartile range (IQR) is used. The improved version of the rule-of-thumb bandwidth, denoted as h_ROT2, is defined as

(Härdle et al., 2004).

3.3. Plug-in

The concept of plug-in bandwidth was originally introduced by Woodroofe (1970). This concept is based on the idea of using an optimal bandwidth that minimizes AMISE( f̂ (x; h_SCBS)). Plug-in bandwidth is based on the substitution of the unknown quantity θ₂ = ∫ ( f″ (x))²dx. Sheather and Jones (1991) also proposed the plug-in bandwidth, denoted as h_SJDP. They provided a solution for this bandwidth selection method as

The pilot bandwidth for the estimation of ψ₄ is a function γ of h. The choice of γ is defined by

where ψ̂₄ (g₁) and ψ̂₆ (g₂) are kernel estimates of ψ₄ and ψ₆. The choice of g₁ and g₂ are formulated by

where ψ^6=-15/(16 πσ^7), ψ^8=105/(32 πσ^9), K(4)(0)=3/ 2π and K(6) (0)=-15/ 2π (Wand and Jones, 1995).

Bandwidth selection in KDE using the AMISE criteria involves estimating the second-order derivative of the unknown density being estimated. The first kind shifted Chebyshev series expansion can be used as an approximation method for an unknown density function. This section will cover the necessary background on the first kind shifted Chebyshev polynomials and derive the bandwidth.

4.1. The shifted Chebyshev polynomials

The first kind Chebyshev polynomials of degree m, where m ∈ {0, 1, 2, . . .}, are denoted as T_m (x) and defined on the interval [−1, 1] . A more general recurrence relation is

In order to use the first kind Chebyshev polynomials on a finite range [a, b], a transformation can be applied to generate the so-called the first kind shifted Chebyshev polynomials. This transformation involves using the equation

Afterward, the first kind shifted Chebyshev polynomials are generated by

In the context of the interval [a, b], the approximating function can be approximated using the first kind shifted Chebyshev series as

where the coefficients are defined via the formula

where

is the Chebyshev zero nodes.

The integration of the squared function of the second-order derivative of the first kind shifted Chebyshev series expansion θ₂ = ∫ ( f″ (y))²dy can be expressed as

By substituting the values of θ₂ from (4.7) into the expressions for the optimal bandwidth derived in (3.6), the resulting bandwidth can be described as the first kind shifted Chebyshev series-based bandwidth:

4.2. The optimal value of m

The bandwidth h_SCBS is influenced by the number of terms m in the first kind shifted Chebyshev expansion used in estimating θ_2,m. The performance and smoothness of this expansion are also affected by m, acting as a smoothing parameter and representing the number of terms in the series expansion. The best choice for the number of terms is the smallest m that results in the lowest ISE. Then, the error in function expansion using the first kind shifted Chebyshev expansion is determined by MISE ( f̂ (x; h_SCBS)) = E ∫ [ f̂ (x; h_SCBS) − f (x)]²dx.

In this section, the aim is to evaluate the performance of the proposed bandwidth h_SCBS and compare with three other bandwidths used for density estimation: least squares cross-validation bandwidth (h_LSCV), the improved version of the rules of thumb bandwidth (h_ROT2), and the Sheather and Jones plug-in bandwidth (h_SJDP). This study compares different bandwidths for density estimation using fifteen normal mixture densities constructed by Marron and Wand (1992). These densities include various shapes, such as unimodal, bimodal, trimodal, and multimodal, each defined and visualized in Marron and Wand’s work. For each distribution, sample sizes of n = 25, 50, 100, 150, and 200 are considered, and the MISE of the estimator is computed over 500 replications. The Gaussian kernel function is used in all cases.

The main idea is to find the optimal number of terms in the expansion (m) that allows a good approximation of f (x) using the first kind shifted Chebyshev series expansion, in the sense of the mean integrated squared error (MISE). The results have been presented through box plots showing the estimated MISE as a function of the number of terms in the expansion, as shown in Figure 1. This graph illustrates the influence of the terms m on MISE across three density distributions (#4, #6, and #12). The x-axis labels represent the sequential terms in the expansion, and the medians of MISE are displayed in the box plots. Furthermore, a solid circle is used to indicate the optimal number of terms in the expansion, which corresponds to the smallest MISE. Figure 1(a) displays the MISE of Density #4 with sample size of n = 50, while Figure 1(b) shows the MISE of Density #6 with sample size n = 50. Figure 1(c) and Figure 1(d) show MISE for Density #12 with two different sample sizes n = 25 and n = 200, respectively. The performance is influenced by the combination of density and sample size, and specific combinations yield better results.

The simulation results evaluate the performance of the plug-in bandwidth by finding the bandwidth that minimizes the mean integrated squared error, MISE ( f̂ (x; h_SCBS)). Table 1 shows the MISE values for different bandwidths, and the bold text value indicates the smallest MISE associated with the bandwidths h_LSCV, h_ROT2, h_SJDP and h_SCBS. Overall, the proposed bandwidth h_SCBS demonstrates good performance compared to the other bandwidths, except for Density #6 (n = 50) and #11 (n = 50), where the improved version of rules of thumb bandwidth (h_ROT2) shows better performance. The simulation results indicate that the suggested bandwidth is a good choice for various scenarios. It provides excellent performance compared to all the other bandwidths under consideration, even though the calculation is quite complex.

In this section, kernel density estimation is applied to real datasets. The performance of the proposed bandwidth h_SCBS is verified against other bandwidths. All calculations are performed using the R programming language.

A real dataset named “flywheels” from Anderson-Cook (1999) and comprising 60 observations on flywheel imbalance angles, will be utilized. This analysis focuses on how different bandwidth choices influence kernel density estimates and histograms, serving as methods to understand data distribution. Figure 2 displays a histogram with 14 bins for this dataset. The density seems to exhibit asymmetric bimodal behavior. The kernel density estimate, using different bandwidth options such as h_LSCV, h_ROT2, h_SJDP and the proposed h_SCBS, is also overlaid on Figure 2. Different bandwidth options will be compared to find the best approach. The kernel density estimate using the suggested bandwidth h_SCBS fits well across the dataset, as confirmed in Table 2 by the lowest mean square error (MSE) value for h_SCBS.

When selecting bandwidth for kernel density estimation, the direct plug-in method is the common initial approach, but there is room for enhancement. Estimating the bandwidth involves finding the integration of the squared function of the second-order derivative of the unknown density to be estimated. This article introduces a bandwidth selection technique by incorporating the estimation of θ_2,m = ∫ ( f″ (x))²dx through the first kind shifted Chebyshev polynomials as an approximation to the function f (x).

The simulation studies revealed that the first kind shifted Chebyshev series-based plug-in bandwidth (h_SCBS) performs well among other bandwidths, such as least squares cross-validation bandwidth, improved rule of thumb bandwidth, and Sheather and Jones plug-in bandwidth. When applying kernel density estimation to estimate the density of the “flywheels” dataset, the results indicate that the proposed bandwidth (h_SCBS), with the lowest MSE, offers the best performance compared to other bandwidth methods. Even though obtaining the proposed bandwidth might be complex, it is still a favorable choice for practical application.

The authors would like to thank Kasetsart University and Rajamangala University of Technology Rattanakosin for the support.