Empirical likelihood for quantile regression models with response data missing at random

1 Introduction

Since the seminal work of Koenker [1], quantile regression (QR) has been an indispensable and versatile tool for statistical research due to its promising performance and elegant mathematical properties, and attracted immediately considerable attention, resulting in numerous papers (e.g., see [1–7]) devoted to various theoretical extensions of this significant topic. Moreover, QR has been widely applied to a variety of fields such as economics, finance, biology, and medicine. Compared to the mean regression model which is commonly used by the traditional least square (LS) methods, QR is able to directly estimate the effects of the covariates at different quantiles of the response variable and therefore provide more information about the distribution of the response variable. Furthermore, QR is less sensitive to outliers due to its specific estimation.

Despite the significant theoretical advances and a rapidly growing literature on QR, only scant attention has been paid to QR when the data samples contain missing values which may lead to substantial distortions on the results. In fact, missing data are a commonplace in practice due to various reasons such as loss of information caused by uncontrollable factors, unwillingness of some sampled units to supply the desired information, failure on the part of investigators to gather correct information, and so forth. Dating back to the early 1970s, spurred by the advances in computer technology that made it possible to perform laborious numerical calculations, the literature on statistical analysis of real data with missing values has flourished in applied work, see [8–12]. Although missing data analysis has a long history in statistics, little work on QR has taken missing data into account. Recently, Yoon [13] proposed an imputation method where the imputed values are drawn from the conditional quantile function of the response with which data are incomplete, but his method is valid only under independent and identically distributed (i.i.d.) errors. Wei [7] developed an iterative imputation procedure for the covariates with missing values in a linear QR model that is valid under non-i.i.d. error terms. Lv [5] discussed smoothed empirical likelihood analysis with missing response in partially linear quantile regression. Sherwood [6] recently proposed an inverse probability weighting QR approach for analyzing health care cost data when the covariates are MAR. Sun [14] studied QR for competing risk data when the failure type was missing. Chen [4] discussed efficient QR analysis with missing observations. Shu [15] proposed some imputation methods for quantile estimation under missing at random.

On the other hand, empirical likelihood (EL) method, introduced by Owen [16, 17], has many advantages over normal approximation methods for constructing confidence intervals. For example, the EL method produces confidence intervals or regions whose shape and orientation are determined entirely by the data, and the empirical likelihood regions are range preserving and transformation respecting. Many authors have used this method for linear, nonparametric and semiparametric regression models. About the quantile regression model, Chen [18] constructed the EL confidence intervals for population quantiles. Tang [19]developed an EL approach for estimating equations with missing data. Wang [20] considered EL for quantile regression models with longitudinal data. Whang [21] proposed a smoothed EL and discussed its higher-order properties with cross sectional data. Otsu [22] studied the first-order approximation of a smoothed conditional EL approach. The EL method has also been used for the analysis of censored survival data, for example, Zhao [23].

In this paper, a empirical-likelihood-based method is proposed to study quantile linear regression models with response data missing at random. A class of quantile empirical log-likelihood (QEL) ratios of the regression parameters are defined firstly which include QEL ratio with complete-case data, weighted QEL ratio and imputed QEL ratio. Then the statistical inference on the mean of the response for a given quantile level is further studied to obtain a bias-corrected QEL ratio of the mean of the response for a given quantile level. It is proved that the QEL ratios of both the regression parameters and the mean of the response for a given quantile level are asymptotically χ² distribution. To compare the quantile empirical likelihood method with a normal approximation method, we also construct a class of estimators for the regression parameters and the mean of the response for a given quantile level. It is shown that this class of estimators are asymptotically normal. Furthermore, we derive consistent estimators of asymptotic variance, their confidence intervals (regions) can be constructed directly of the regressions parameters and the mean of the response for a given quantile level.

The rest of this paper is organized as follows. In Section 2, a class of QEL ratios and estimators for the regression parameters are constructed with missing response data and their asymptotic distributions are derived. In Section 3, a bias-corrected QEL ratio and the maximum empirical QEL estimator for mean of the response at a given quantile level are proposed and their asymptotic properties are studied. A simulation study is conducted in Section 4 to demonstrate the finite-sample performance of the proposed method. A real-world data set is analyzed to illustrated the applications of the proposed method in Section 5. The proof of the main results are postponed to Section 6.

2 Quantile empirical likelihood (QEL) method with missing response

Consider the quantile linear regression model

where Y_i is the ith observation of the response Y, X_i is the ith observation of the covariates X and a d × 1 vector, τ ∈ (0,1) is the quantile level of interest, β_τ is a d × 1 vector of unknown quantile regression parameters and ε_i is the error satisfying P(ε_i < 0|X_i) = τ for i = 1,2, ⋯, n.

For the model (1), we focus on the situation where some observations of Y in a sample of size n may be missing while X is observed completely. As a consequence, we have an incomplete sample { X i , Y i , δ i } i = 1 n with δ_i = 0 if Y_i is missing and δ_i = 1, otherwise. Throughout this paper, we assume that the observations of Y are missing at random (MAR) which implies that δ and Y are conditionally independent given X. That is, P(δ = 1|X,Y) = P(δ = 1|X) = p(X). As pointed out in [24], MAR is a common assumption for statistical analysis with missing data and is reasonable in many practical situations. Hereafter, we will simply write β_τ as β whenever no confusion is made.

3 Quantile empirical likelihood for the mean of the response

Some methods are provided firstly in this section to conduct a inference on the mean of the response θ by using empirical likelihood. Then a weighted quantile regression imputation is used to construct a weighted-corrected quantile empirical likelihood ratio of θ such that this ratio has an asymptotic χ² distribution.

4 A simulation study

In this section a simulation study is carried out to investigate the finite-sample performance of the proposed approaches. We consider the following two models.

Model 1 (homoscedastic): Y_i = X_iβ + ε_i, i = 1, 2, ⋯, n;

Model 2 (heteroscedastic): Y_i = X_iβ + ξX_i ε_i, i = 1,2, ⋯, n,

where the observation X_i (i = 1,2,⋯, n) of the covariates X were drawn the N(0,1), ε i ( i = 1 , 2 , ⋯ , n ) i i d ∼ N(0,1), and ξ and β were set to be 0.5 and 1, respectively. In simulation study, we focus on τ = 0.5,0.8. We considered the following three selection probability functions proposed by Wang and Rao (see [8]).

Case 1 : p₁(x) = 0.8 + 0.2|x – 1| if |x – 1| ≤ 1, and 0.95 otherwise.

Case 2 : p₂(x) = 0.9 – 0.2|x – 1| if |x – 1| ≤ 4.5, and 0.1 otherwise.

Case 3 : p₃(x) = 0.6 for all x.

The average missing rates corresponding to the preceding three cases are approximately 0.09, 0.26 and 0.40, respectively. For each of the three cases, we generated 2000 Monte Carlo random samples of size n = 50, 100 and 150. The kernel function K(x) in (5) was taken to be K(x) = 0.75(1 – x²) if |x| ≤ 1; K(x) = 0, otherwise. We used the cross-validation method to select the optimal bandwidths h_opt. The simulations were implemented in the following two situations.

(1) The confidence intervals of β. For the two models, we used four methods, namely, the quantile empirical likelihood with complete-case data (QCEL), the quantile weighted empirical likelihood (QWEL), the quantile empirical likelihood based on imputed values (QIEL) and the normal approximation(NA) in Theorem 2.2. For convenience, in what follows N A ( β ^ Q E L ) a n d N A ( β ^ Q ) denote the corresponding normal approximation confidence intervals for β ^ Q E L a n d β ^ Q . The average lengths of the confidence intervals and their corresponding empirical coverage probabilities, with a nominal level 1 – α = 0.95 and τ = 0.5, 0.8, were computed with 2000 simulation runs. The results are reported in Tables 1–4.

Tables 1–4 show the following results. Firstly, for Case 1, QIEL yields lightly longer interval lengths but higher coverage probabilities than the other three methods. For Cases 2 and 3, QIEL performs better than the other three methods in the sense that its confidence intervals have uniformly shorter average lengths and higher coverage probabilities, which indicates that quantile regression imputation is necessary when the missing rate is large. Secondly, both QCEL and QWEL result in slightly longer interval lengths but higher coverage probabilities than N A ( β ^ Q E L ) a n d N A ( β ^ Q ) . In addition, the confidence inervals obtained by N A ( β ^ Q E L ) a n d N A ( β ^ Q ) show nearly equal lengths and coverage accuracies in the same case. Thirdly, as expected all the interval lengths decrease and the empirical coverage probabilities increase as n increases for every given missing rate. Observably, the missing rate also affects the interval length and coverage probability. Generally, the interval length increases and the coverage probability decreases as the missing rate increases for every fixed sample size. However, the two values fail to change by a large amount for the QIEL method because the quantile regression imputation is used in QIEL. Furthermore, it is also seen that than the other methods for the heteroscedastic model QIEL still performs much better.

(2) The confidence intervals of θ. The weighted-corrected empirical likelihood(WCQEL) based on l ^ ( θ ) and NA were considered. θ ^ Q W I defined in (9) is used to estimat θ. The empirical coverage probabilities and average lengths of the confidence intervals, with a nominal level 1 – α = 0.95 were computed with 2000 simulation runs. The results are reported in Table 5.

It is seen from Table 5 that WQCEL produces slightly longer interval lengths, but higher coverage probabilities than NA does. All the coverage probabilities increases and the average lengths decrease as n increase. In addition, the coverage probabilities and average lengths depend on the selection probability function p(x) and the quantile level τ.

5 A real-data example

The data originally presented by [25] is investigated in this section to support the proposition that food expenditure constitutes a declining share of personal income. This data that has not any missing data consists of 235 budget surveys of 19th century European working class households. More details of the discussion on this data can be found in [26]. We consider the following linear QR model: Y_i = β₀(τ) + β₁(τ)X_i,i = 1,2,⋯, 235, where Y is the centered annual household food expenditure and X is the centered annual household income in Belgian francs. In order to use the data set to illustrate our method, artificial missing data was created by deleting some of the response values at random. Assume that 25% of the response values in this data are missed. The missing indicator δ is generated from the probability function p(x) = 0.9 – 0.2|x – 1| if |x – 1| ≤ 4.5, and 0.1 otherwise.

We now present the estimator and the 95% confidence interval of β based on the proposed QILE method and the normal approximation method (NA) based on Theorem 2.2 with τ = 0.4 and 0.7. The results are shown in Table 6. From Table 6, we can see that the confidence interval obtained by the QIEL method has much shorter confidence interval than that obtained by the NA method, which shows that the former method is superior to the latter one.

6 Proofs of the main results

Let r > 2 be an integer. g(x), f(·|x) and F(·|x) are used to denote the density function of X,the density and distribution functions of ε conditional on X_i = x, respectively. Let c be a positive constant which is independent of n and may take a different value in different place. The following conditions will be used in this section.

(C1) {(Y_i, X_i): i = 1,2,⋯, n} are independent and identically distributed random vectors.

(C2) Both p(x) (the selection probability function) and g(x) have bounded derivatives up to order r almost surely and inf_x p(x) > 0:

(C3) K(·) is a kernel function of order r and is bounded and compactly supported on [– 1, 1]. Furthermore, there exist positive constants C₁, C₂ and ρ such that C 1 I [ | | u | | ≤ ρ ] ≤ K ( u ) ≤ C 2 I [ | | u | | ≤ ρ ] .

(C4) P ( ∥ X ∥ > M n ) = o ( n − 1 / 2 ) , where 0 < M_n →∞ as n → ∞.

(C5) The positive bandwidth parameter h satisfies nh^2r →0 when n →∞.

(C6) The matrices A and B defined in Theorem 2.2 are both nonsingular. Firstly, some lemmas are introduced to derive the main results.