The nonparametric bootstrap for the current status model

TL;DR

This paper demonstrates that while direct bootstrap of the nonparametric MLE in the current status model is inconsistent, bootstrapping functionals of the MLE yields valid confidence intervals, supported by convergence results.

Contribution

It introduces a method for valid bootstrap inference in the current status model by focusing on functionals of the MLE, overcoming previous inconsistency issues.

Findings

Bootstrapped MLE converges at the correct rate in Lp-distance.

Bootstrapping functionals of the MLE produces valid confidence intervals.

Results extend to the current status regression model.

Abstract

It has been proved that direct bootstrapping of the nonparametric maximum likelihood estimator (MLE) of the distribution function in the current status model leads to inconsistent confidence intervals. We show that bootstrapping of functionals of the MLE can however be used to produce valid intervals. To this end, we prove that the bootstrapped MLE converges at the right rate in the $L_p$-distance. We also discuss applications of this result to the current status regression model.

Tables

Table 1. Table 1 : Average length of the SMLE-based CIs defined in ( 4.1.3 ) for different bandwidth choices ( h ∼ n − 1 / 5 similar-to ℎ superscript 𝑛 1 5 h\sim n^{-1/5} and h ∼ n − 1 / 4 similar-to ℎ superscript 𝑛 1 4 h\sim n^{-1/4} ) and average length of the MLE-based CIs proposed by [ 4 ] and [ 32 ] at timepoints t = 0.5 , 1 , 1.5 𝑡 0.5 1 1.5 t=0.5,1,1.5 .

	Uniform			Exponential
Method	$t=0.5$	$t=1$	$t=1.5$	$t=0.5$	$t=1$	$t=1.5$
SMLE ($h\sim n^{-1/5}$)	0.064819	0.077020	0.064976	0.085540	0.087565	0.057716
SMLE ($h\sim n^{-1/4}$)	0.079671	0.092096	0.079757	0.085540	0.087565	0.057716
MLE ([4])	0.164767	0.184590	0.165699	0.204079	0.161122	0.104002
MLE ([32])	0.183982	0.202430	0.186452	0.225882	0.176159	0.118541

Table 2. Table 2 : Simulation model 1: mean, n 𝑛 n times the variance and n 𝑛 n times MSE. CP: coverage proportion of 95% confidence intervals (Wald-type intervals based on a kernel variance estimate and classical bootstrap intervals) that contain the true parameter value β 0 = 0.5 subscript 𝛽 0 0.5 \beta_{0}=0.5 , AL: Average length of the CI, for different samples sizes n 𝑛 n based on N = 1 , 000 𝑁 1 000 N=1,000 simulation runs and B = 1 , 000 𝐵 1 000 B=1,000 bootstrap samples. ϵ = 0.001 italic-ϵ 0.001 \epsilon=0.001 . SSE = simple score estimator, MRCE = maximum rank correlation estimator and ESE = efficient score estimator.

Estimate	$n$	mean	$n\times$var	$n\times$MSE	Wald-type CI		Bootstrap CI
					CP	AL	CP	AL
SSE	100	0.498943	0.310723	0.310968	0.978	0.265883	0.824	0.204163
	500	0.499717	0.220885	0.220925	0.982	0.097457	0.897	0.080317
	1000	0.500720	0.217415	0.217933	0.977	0.065837	0.924	0.055648
	5000	0.499993	0.195111	0.195112	0.977	0.027159	0.945	0.024423
MRCE	100	0.497996	0.308180	0.308582	0.979	0.268731	0.821	0.205522
	500	0.499761	0.251232	0.251260	0.978	0.098028	0.862	0.089143
	1000	0.500553	0.246388	0.246693	0.973	0.063990	0.911	0.053129
	5000	0.499876	0.208386	0.208462	0.965	0.027197	0.922	0.026987
ESE	100	0.500145	0.337755	0.337757	0.964	0.252687	0.824	0.223849
	500	0.499671	0.217428	0.217482	0.978	0.094390	0.896	0.080003
	1000	0.500742	0.207401	0.207953	0.973	0.063990	0.911	0.053129
	5000	0.500228	0.185614	0.185874	0.972	0.026396	0.904	0.022285

Table 3. Table 3 : Simulation model 2: mean, n 𝑛 n times the variance and n 𝑛 n times MSE. CP: coverage proportion of 95% confidence intervals (Wald-type intervals based on a kernel variance estimate and classical bootstrap intervals) that contain the true parameter value β 0 = 1 subscript 𝛽 0 1 \beta_{0}=1 , AL: Average length of the CI, for different samples sizes n 𝑛 n based on N = 1 , 000 𝑁 1 000 N=1,000 simulation runs and B = 1 , 000 𝐵 1 000 B=1,000 bootstrap samples. ϵ = 0.001 italic-ϵ 0.001 \epsilon=0.001 . SSE = simple score estimator, MRCE = maximum rank correlation estimator and ESE = efficient score estimator.

Estimate	$n$	mean	$n\times$var	$n\times$MSE	Wald-type CI		Bootstrap CI
					CP	AL	CP	AL
SSE	100	0.935732	4.525330	4.938096	0.922	1.000283	0.855	0.79952
	500	0.966217	4.676249	5.246881	0.926	0.399728	0.902	0.364210
	1000	0.977799	5.032432	5.525339	0.933	0.279928	0.914	0.262449
	5000	0.989466	4.580756	5.135616	0.945	0.124375	0.948	0.121388
MRCE	100	1.038510	8.500588	8.648890	0.925	1.125225	0.889	1.364034
	500	1.006050	6.443404	6.461690	0.932	0.429007	0.912	0.473787
	1000	1.002680	6.294143	6.301326	0.939	0.296537	0.903	0.320908
	5000	0.998502	5.160694	5.171915	0.962	0.129512	0.954	0.136487
ESE	100	0.974199	5.722576	5.789144	0.768	0.604649	0.827	0.910229
	500	0.998806	5.984291	5.985003	0.823	0.290297	0.902	0.430819
	1000	1.005545	6.032743	6.063495	0.841	0.214280	0.928	0.302124
	5000	1.002462	5.244373	5.274692	0.892	0.104281	0.951	0.131427

Equations

$$\mathbb{C}=\arg \underset{t}{\max }\left\{W(t)-t^2\right\},$$

$$n^{1/3}\{4F_0(t)(1-F_0(t))f_0(t)/g(t){\}}^{-1/3}\{{\hat{F}}_n(t)-F_n(t)\}$$

$$\stackrel{\mathcal{D}}{\to }\arg \underset{t}{\max }(W(t)+\hat{W}(t)-t^2)-\arg \underset{t}{\max }(W(t)-t^2),$$

$$p_{F_0}(t,\delta )=\left[\delta F_0(t)+(1-\delta )\{1-F_0(t)\}\right]g(t).$$

$${\ell }_n(F)=n^{-1}\sum _{i=1}^n\left[{\Delta }_i\log F(T_i)+(1-{\Delta }_i)\log \{1-F(T_i)\}\right],$$

$$\left(i,\sum _{j\le i}{\Delta }_{(j)}\right),$$

$$V_n(t)=n^{-1}\sum _{i=1}^n{\Delta }_i{1}_{\{T_i\le t\}},$$

$$U_n(a)=\text{mbox}argmin\{t\in \mathbb{R}:V_n(t)-a{\mathbb{G}}_n(t)\}.$$

$$F_n(t)\ge a⟺U_n(a)\le t,$$

$$\displaystyle E\left\{\left\|n^{1/3}\left\{\hat{F}_{n}-F_{0}\right\}\right\|_{p}\Bigm{|}Z_{1},\ldots,Z_{n}\right\}=O_{p}(1),$$

$$\displaystyle\sup_{t\in[0,R]}E\left\{\left.n^{1/3}\bigl{|}\hat{F}_{n}(t)-F_{0}(t)\bigr{|}\right|Z_{1},\dots,Z_{n}\right\}=O_{p}(1).$$

$${\hat{\mathbb{P}}}_n=n^{-1}\sum _{i=1}^n{M}_{ni}{1}_{Z_i},$$

$$M_n=(M_{n1},\dots ,M_{nn})\sim \text{multinomial}(n,n^{-1},\dots ,n^{-1}),$$

$$\left(\sum _{j=1}^i{M}_{n(j)},\sum _{j=1}^i{M}_{n(j)}{\Delta }_{(j)}\right),$$

$$P^{\ast }\left(P_{M|Z}\left\{|{\Gamma }_n|>ϵ\right\}>\eta \right)\to 0\text{ as }n\to \infty ,$$

$$0<\underset{t\in [0,R]}{\inf }{f}_0(t)<\underset{t\in [0,R]}{\sup }{f}_0(t)<\infty .$$

$$U(a)=F_0^{-1}(a)0<a<1,$$

$${\hat{U}}_n(a)=\text{argmin}\{t\in [0,R]:{\hat{V}}_n(t)-a{\hat{\mathbb{G}}}_n(t)\}0<a<1,$$

$${\hat{V}}_n(t)={\int }_{u\in [0,t]}\delta d{\hat{\mathbb{P}}}_n(u,\delta )\text{ and }{\hat{\mathbb{G}}}_n(t)={\int }_{u\in [0,t]}d{\hat{\mathbb{P}}}_n(u,\delta )t\in [0,R].$$

$$\left\{\exists x\in [0,R]:P_{M|Z}\left\{n^{1/3}\left|{\hat{U}}_n(a)-U(a)\right|\ge x\right\}>K_1{e}^{-K_2{x}^{3/2}}\right\}=o_p(1),$$

$$P_{M|Z}\left\{n^{1/3}\left|{\hat{U}}_n(a)-U(a)\right|\ge x\right\}\le K_1{e}^{-K_2{x}^{3/2}}$$

$$E_{M|Z}{\left[n^{1/3}\{{\hat{F}}_n(t)-F_0(t){\}}_{+}\right]}^p={\int }_0^{\infty }{P}_{M|Z}\left\{n^{1/3}\{{\hat{F}}_n(t)-F_0(t)\}\ge x\right\}p{x}^{p-1}dx,$$

$$P_{M|Z}\left\{{\hat{U}}_n\left(a+n^{-1/3}x\right)\le t\right\}$$

$$\displaystyle\quad=P_{M|Z}\Biggl{[}n^{1/3}\left\{\hat{U}_{n}\left(a+n^{-1/3}x\right)-U\left(a+n^{-1/3}x\right)\right\}$$

$$\displaystyle\qquad\qquad\qquad\qquad\qquad\qquad\qquad\leq n^{1/3}\left\{t-U\left(a+n^{-1/3}x\right)\right\}\Biggr{]},$$

$$P_{M|Z}\left\{n^{1/3}\{{\hat{F}}_n(t)-F_0(t)\}\ge x\right\}=P_{M|Z}\left\{{\hat{U}}_n\left(a+n^{-1/3}x\right)\le t\right\},$$

$$\left\{\exists t\in [0,R]:E_{M|Z}{\left|{\hat{F}}_n(t)-F_0(t)\right|}^p>K{n}^{-p/3}\right\}=o_p(1).$$

$$P\left\{\underset{t\in [0,R]}{\sup }{E}_{M|Z}\left|{\hat{F}}_n(t)-F_0(t)\right|>K_1{n}^{-1/3}\right\}⟶0,n\to \infty ,$$

$$\displaystyle P\left\{E_{M|Z}\bigl{\|}\hat{F}_{n}-F_{0}\bigr{\|}_{2}>K_{2}n^{-1/3}\right\}\longrightarrow 0,\qquad n\to\infty.$$

$${\tilde{F}}_{nh}(t)=\int \mathbb{K}\left((t-x)/h\right)d{F}_n(x),$$

$$\mathbb{K}(u)={\int }_{-\infty }^{u}K(x)dx,$$

$$K(u)=\frac{35}{32}{\left(1-u^2\right)}^3{1}_{[-1,1]}(u).$$

Full text

62G09,

62N01,

bootstrap,

current status,

MLE,

smooth functionals,

arXiv:1701.07359 \startlocaldefs

\setattributejournalname \endlocaldefs

1 Introduction

In the current status model, the variable of interest is a survival variable $X$ with distribution function $F_{0}$. However, instead of observing the exact survival time $X$, a censoring variable $T\sim G$ is observed together with the indicator $\Delta=1_{X\leq T}$. Such data arise naturally in clinical trials when a patient can only be checked at one measurement due to destructive testing. A lot of research has been published on the behavior of the maximum likelihood estimator (MLE) $F_{n}$ of the distribution function $F_{0}$. The limiting distribution of $n^{1/3}(F_{n}(t)-F_{0}(t)$) is after scaling by the constant $\kappa=\{4F_{0}(t)(1-F_{0}(t))f_{0}(t)/g(t)\}^{1/3}$ given by

[TABLE]

where $W$ is a two-sided Brownian motion with $W(0)=0$ (see [19]). Other estimators with similar asymptotic properties are Chernoff’s estimator of the mode ([6]), the Grenander estimator ([10]) of a nonincreasing density, Manski’s maximum score estimator ([27]) and Rouseeuw’s least median of squares estimator ([29]). A general framework for cube-root $n$ asymptotics is given in [25].

In this paper we investigate the behavior of Efron’s nonparametric bootstrap method ([9]) for constructing confidence intervals for smooth functionals of the MLE. It is known that the nonparametric bootstrap is inconsistent for generating the limit distribution of the MLE. The authors of [2] prove that (conditional on the data),

[TABLE]

where $\hat{F}_{n}$ is the bootstrap MLE and $W$ and $\hat{W}$ are two independent two-sided Brownian motions originating at zero. A similar result is obtained in [26] and in [31] for the Grenander estimator. The maximum score estimator of [27] is another example of a cube-root $n$ statistic with asymptotic distribution derived in [25], where inconsistency of the nonparametric bootstrap for this estimator is shown in [2].

Constructing asymptotic confidence intervals for the distribution function in the current status model based on Chernoff’s distribution and the normalizing constant $\kappa$ is complicated by the need to compute the critical values of $\mathbb{C}$ and to estimate the density $f_{0}$ consistently. Since this turns out to be a rather difficult task several alternative bootstrap methods have been proposed based on resampling from a smooth estimate. [32] consider a smooth kernel estimate $\tilde{F}$ of $F_{0}$ and resample the $\Delta_{i}$ from a Bernoulli distribution with probability $\tilde{F}(T_{i})$, while keeping the censoring variables $T_{i}$ fixed and center the values of the bootstrap samples by subtracting the smooth estimate of the distribution function. [26] and [31] propose similar smooth respampling schemes for the Grenander estimator and a model-based smoothed bootstrap procedure for making inference on the maximum score estimator is developed in [28]. All methods result in consistent estimation of the (suitably standardized) distribution $\mathbb{C}$ conditional on the original data.

A drawback of this approach is that smoothness conditions of $F_{0}$ are used which allow faster than cube-root $n$ estimation of $F_{0}$. This raises the question if one should really use confidence intervals based on the MLE instead of on a faster converging estimate.

This latter procedure is followed in [14], where the authors consider constructing confidence intervals around the smoothed maximum likelihood estimator (SMLE) of $F_{0}$ in the current status model. The SMLE is a kernel estimate based on the MLE with an asymptotic normal distribution, instead of Chernoff’s limiting distribution ([16]). The bootstrap method proposed in [14] is however still based on the smooth bootstrap procedure described in [32] and not on Efron’s nonparametric bootstrap. We show in this paper that the construction of confidence intervals around the SMLE based on the nonparametric bootstrap can also be proved to be valid, where one does not resample from a smooth estimate of $F_{0}$, but just resamples with replacement from the pairs $(T_{i},\Delta_{i})$ in the original sample. This method already has been used without proof in [17] and also in [18] and the present manuscript intends to fill the gap of the missing proofs here. An important difference with the smooth bootstrap in [14] is that for the centering of the estimates in the nonparametric bootstrap samples the SMLE of the original sample is used, whereas this will not work for the resampling as proposed in [14]; in the latter case one needs to center the estimates in the bootstrap samples by a kernel convolution of the SMLE in the original sample. It is not clear which method is better, and the most striking fact is the similarity of the results of the two methods in our simulations. An advantage of the purely nonparametric bootstrap, discussed in the present paper, might be its conceptual simplicity and the absence of the need to center with a convolution of the SMLE in the centering of the bootstrap samples instead of the SMLE itself. An advantage of the smooth bootstrap, discussed in [14] might be the fact that only the indicators $\Delta_{i}$ are being resampled, and that in this sense one stays closest to the sample distribution of the observation times $T_{i}$, which stay fixed in this procedure.

Although it is argued in [8] that the naive bootstrap will not work for their goodness-of-fit test for monotone functions, based on the Grenander estimator, no theoretical justification for this conjecture is given. Other examples where a smooth bootstrap procedure is used, are the likelihood ratio type two-sample test for current status data proposed by [11] and the test for equality of functions under monotonicity constraints proposed by [7]. Both tests establish asymptotic normality for the test statistic considered.

The paper is organized as follows: In Section 2 we introduce the current status model and review some interesting properties of the MLE. The validity of the nonparametric bootstrap is discussed in Section 3. In Section 4 we provide two examples to illustrate the applicability of our result. In the first example we construct pointwise confidence intervals based on the smoothed MLE in the current status model. The second example deals with doing inferences for a finite dimensional regression parameter in the current status linear regression model. For both examples, the theoretical and finite sample behavior of the nonparametric bootstrap is discussed. Section 5 presents some concluding remarks. The proofs of our results are given in Section 6.

2 The current status model and the MLE

Let $Z_{1}=(T_{1},\Delta_{1}),\ldots,Z_{n}=(T_{n},\Delta_{n})$ be an i.i.d. sample from the probability space $([0,R]\times\{0,1\},{\cal A},P)$, where $\Delta_{i}=1_{X_{i}\leq T_{i}}$ and $R>0$. The $X_{i}$ are interpreted as (nonnegative) survival times with distribution function $F_{0}$. Instead of observing $X$, a censoring variable $T\sim G$ is observed (with density $g$) independent of $X$. One could say that in the current status model, each observation $Z_{i}$ represents the current status of the item $i$ at time $T_{i}$. The density of $Z_{i}$ with respect to the product of Lebesgue measure and counting measure on $[0,R]\times\{0,1\}$ is given by

[TABLE]

The maximum likelihood estimator $F_{n}$ is defined as the maximizer of the log likelihood given by (up to a constant not depending on $F$),

[TABLE]

over all distribution functions $F:[0,\infty]\mapsto[0,1]$. [19] show that the MLE can be characterized as the left-continuous slope of the greatest convex minorant of a cumulative sum diagram consisting of the points (0,0) and

[TABLE]

where we let $T_{(j)}$ denote the $j$th order statistic of the $T_{i}$ and $\Delta_{(j)}$ be the $\Delta_{i}$ corresponding to it (assuming no ties are present in the data). An important property of the MLE is the so-called switch relation, see [17] p. 69. Let ${\mathbb{G}}_{n}$ be the empirical distribution function of $T_{1},\ldots,T_{n}$ and define the process $V_{n}$ by

[TABLE]

and the process (in $a$) $U_{n}$ by

[TABLE]

Then, taking $a=F_{0}(t)$, we get the switch relation:

[TABLE]

see also Figure 2.