Maik Schwarz
email: maik.schwarz@uclouvain.be

Sébastien Van Bellegem

Consistent density deconvolution under partially known error distribution

Keywords: deconvolution, measurement error, density estimation, minimum distance estimator, identification 2000 MSC: Primary 62G07, secondary 62G20

We estimate the distribution of a real-valued random variable from contaminated observations. The additive error is supposed to be normally distributed, but with unknown variance. The distribution is identifiable from the observations if we restrict the class of considered distributions by a simple condition in the time domain. A minimum distance estimator is shown to be consistent imposing only a slightly stronger assumption than the identification condition.

Introduction

A standard problem in nonparametric statistics is to consistently estimate the distribution of some real random variable X from a statistical sample that is subject to an additive measurement error ε. The usual setting is to assume independent and identically distributed (iid) observations from a random variable Y that is such that Y = X + ε. Knowing the cumulative distribution function (cdf) of ε, a vast literature focuses on the accuracy of estimation of the cdf of X (e.g. [START_REF] Carroll | Optimal rates of convergence for deconvolving a density[END_REF]; [START_REF] Fan | On the optimal rate of convergence for nonparametric deconvolution problems[END_REF]). The full knowledge of the cdf of ε is of course a strong assumption that is rarely encountered in real data analysis. However this assumption is also an identification condition: without the full information on the cdf of ε, the cdf of X cannot be identified from the observations of Y .

In order to circumvent that issue, some authors have worked under other sampling processes. In addition to the observation of a sample of Y some papers assume to observe another independent sample from the measurement error ε. That new sample allows to estimate the cdf of ε in a first step, and therefore to recover the cdf of X in a second step (e.g. [START_REF] Neumann | On the effect of estimating the error density in nonparametric deconvolution[END_REF], [START_REF] Cavalier | Adaptive estimation for inverse problems with noisy operators[END_REF], [START_REF] Johannes | Convergence rates for ill-posed inverse problems with an unknown operator[END_REF]). In other studies, it is instead assumed to observe longitudinal, or multilevel, versions of Y (e.g. [START_REF] Li | Nonparametric estimation of the measurement error model using multiple indicators[END_REF]; [START_REF] Meister | Deconvolution from non-standard error densities under replicated measurements[END_REF]).

ACCEPTED MANUSCRIPT

The present paper develops another strategy to identify the cdf of X. We assume that the measurement error, ε, is normally distributed with mean zero and unknown variance σ 2 . In this setting the cdf of X is of course not identified from the observation of Y . However if we restrict the class of cdfs of X we show below that the model becomes identified in spite of the partial knowledge on the cdf of ε. In Section 2 below we prove that the distribution of X is identifiable from the observation Y = X + ε if there is a set of positive Lebesgue measure on which it puts no mass.

The normality of ε is not a restrictive assumption in real data analysis, in which the main issue is often to evaluate the level of the noise σ 2 . In [START_REF] Hall | Estimating a changepoint, boundary, or frontier in the presence of observation error[END_REF] a similar setting is considered but under the additional assumption that σ 2 depends on the sample size n in such a way that σ 2 → 0 as n → ∞. [START_REF] Matias | Semiparametric deconvolution with unknown noise variance[END_REF] and [START_REF] Butucea | Minimax estimation of the noise level and of the deconvolution density in a semiparametric deconvolution model[END_REF] also consider the consistent estimation of σ 2 and of the cdf of X under strong restrictions on the characteristic function of X. In Section 3 we show that if we slightly restrict the class of cdfs of X, then we can prove the consistency of a minimum penalized contrast estimator of the cdf and σ 2 . The originality of our approach is to exploit the qualitative prior on X (which is vanishing on some interval) in order to separate it from the measurement error.

The next section addresses the identification issue and in Section 3 we prove the consistency of the minimum penalized contrast estimator. A technical lemma on measure theory is deferred to an appendix.

Identification

Suppose we want to recover the probability distribution P X of a random variable X that is observed with an additive and independent random contamination error ε. This measurement error is assumed to be normally distributed with mean zero and unknown variance σ 2 . The resulting observational model is

Y = X + ε. (1)
The distribution of Y is the convolution P Y = P X * N σ , where N σ denotes the probability distribution of ε.

Writing ϕ X , ϕ Y and ϕ σ for the characteristic functions of P X , P Y and N σ , respectively, the convolution equation is equivalent to ϕ Y = ϕ X ϕ σ by virtue of the convolution theorem. Because of the uncertainty about the variance of the measurement error, not all probability distributions can be recovered from the model. Define the set of distributions

P 0 = P ∈ P ∃ A ∈ B(R) : |A| > 0 ∧ P (A) = 0 ,
where B(R) denotes the set of Borel sets in R and P the set of all probability distributions, and |A| the Lebesgue measure of A. In the following theorem, we show that all distributions belonging to P 0 are identifiable from the observational model.

ACCEPTED MANUSCRIPT

Theorem 2.1 (Identification). The model defined by ( 1) is identifiable for the parameter space P 0 × (0, ∞), that is, for any two probability measures P 1 , P 2 ∈ P 0 and σ 1 , σ 2 > 0, we have that P 1 * N σ1 = P 2 * N σ2

implies P 1 = P 2 and σ 1 = σ 2 .

The proof of this theorem is based on the following lemma.

Lemma 2.2. Let P 1 and P 2 be probability distributions and 0 < σ 1 < σ 2 . Then,

P 1 * N σ1 = P 2 * N σ2 =⇒ P 1 = P 2 * N σ3 , where σ 3 = σ 2 2 -σ 2 1 .
Proof. First, apply the convolution theorem on both sides of the equation, then divide by ϕ σ1 which is positive everywhere. To conclude, it suffices to remark that ϕ σ3 = (ϕ σ2 /ϕ σ1 ).

Proof of Theorem 2.1. Suppose that (P 1 , σ 1 ), (P 2 , σ 2 ) ∈ P 0 × (0, ∞) are chosen in such a way that

P 1 * N σ1 = P 2 * N σ2
. It has to be shown that this implies (P 1 , σ 1 ) = (P 2 , σ 2 ). First, we prove by contradiction that σ 1 = σ 2 . Suppose that σ 1 = σ 2 . Without loss of generality, say σ 1 < σ 2 . By virtue of Lemma 2.2, this implies P 1 = P 2 * N σ3 .

We show now that this is only possible if P 1 is not in P 0 which contradicts the assumption. Indeed, denotes the product measure. We have that a

1 -b 1 < a 2 -b 2 because of |B| < |A|. It is easily verified that B × [a 1 -b 1 , a 2 -b 2 ] ⊂ S A and hence P 1 (A) = (P 2 * N σ3 )(A) ≥ P 2 (B) N σ3 ([a 1 -b 1 , a 2 -b 2 ]) > 0.
This contradicts the assumption that P 1 ∈ P 0 , showing that σ 1 = σ 2 .

The characteristic function of the normal distribution being positive everywhere, an application of the convolution theorem completes the proof.

Remark 2.3. The identification theorem assumes that the measurement error ε is normally distributed with an unknown variance σ 2 . Although this is a most natural assumption from a practical point of view, it should be noticed that the proof essentially exploits the infiniteness of the support of ε. Therefore, the identification result may be extended to general (but known) scale families of error distributions, provided some counterpart of Lemma 2.2 holds. This is the case for Cauchy distributions with location parameter µ = 0 or, more generally, for stable distributions with fixed exponent α ∈ (0, 2], skewness parameter β = 0, and location µ = 0, for example.

It is worth noticing that if in Theorem 2.1 we do not suppose both P 1 and P 2 to belong to P 0 , the conclusion is false in general as the following counterexample illustrates. Let P 1 be the uniform distribution on the unit interval and ϕ 1 its characteristic function. Clearly, P 1 ∈ P 0 . If we let P 2 be the probability distribution with characteristic function ϕ 2 := ϕ 1 ϕ σ /ϕ (σ/2) , then, in view of the convolution theorem, we have P 1 * N σ = P 2 * N (σ/2) , but P 1 = P 2 .

ACCEPTED MANUSCRIPT

We conclude this section by remarking that the probability distributions in question are not required to have densities. For those having one, the identification condition can be equivalently expressed by requiring that the density has to vanish on a set of positive Lebesgue measure.

Estimation

Now suppose we have an i.i. -→ P Y , which is hence a necessary condition of consistency. We call an estimator satisfying this condition admissible.

Minimum distance estimation

Let ϕ Y n (t) = 1 n n k=1 exp(itY k ) be the empirical characteristic function of the observations. For characteristic functions ϕ X , ϕ σ , and ϕ Y let us define a distance ρ,

ρ( ϕ X , ϕ σ ; ϕ Y ) := R | ϕ X (t) ϕ σ (t) -ϕ Y (t)| h(t) dt, (2)
where h is some continuous and strictly positive probability density ensuring the existence of the integral.

With a slight abuse of notation we do not make the dependence of the distance on h explicit, as it does not have any influence on the results derived in this work. The estimation consists in choosing P X n and N b σn such that their characteristic functions ϕ X n and ϕ b σn minimize ρ(• , • ; ϕ Y n ). Since this minimum is not necessarily attained, we give the following definition.

Definition 3.1. Let (δ n ) n∈N be a vanishing sequence and C a set of probability distributions. In the context of the deconvolution model, we call a random sequence ( P X n , σ n ) depending on the observations {Y 1 , . . . , Y n } a minimum distance estimator on C if it is such that the corresponding characteristic functions yield

ρ( ϕ X n , ϕ b σn ; ϕ Y n ) ≤ inf e ϕ X ∈ΦC e σ≥0 ρ( ϕ X , ϕ e σ ; ϕ Y n ) + δ n , (3)
where we denote by Φ C the set of all characteristic functions of some class of distributions C. Let further

P X+ε n := P X n * N b σn , the characteristic function of which is ϕ X n ϕ b σn .
Our aim is to prove the consistency of this estimator. Obviously, this requires further assumptions on the class C. In the first instance, we show that the minimum distance estimator is always admissible.

Lemma 3.2 (Admissibility). Any minimum distance estimator ( P X n , σ n ) on the set P of all probability distributions is admissible.

ACCEPTED MANUSCRIPT

Lemma 3.4. For any R, η > 0, weakly convergent sequences in P η R have their limit in P 0 .

Proof. Let (P n ) n∈N be a sequence in P η R . Then, we have

∀n ∈ N ∃ interval A n ⊂ [-R, R] : P n (A n ) = 0 ∧ |A n | ≥ η. (6)
Suppose further that (P n ) n∈N converges in distribution to some P ∞ . We have to show that there is some

A ∞ ∈ B([-R, R]) of positive Lebesgue measure such that P ∞ (A ∞ ) = 0, that is P ∞ ∈ P 0 .
Firstly, we deduce from ( 6) that there exists an x 0 ∈ [-R, R] which lies in infinitely many A n , or in other words,

∃x 0 ∈ [-R, R] ∃(n k ) k∈N ∀k ∈ N : x 0 ∈ A n k .
As all A n k are intervals of length at least η, there is an interval containing x 0 which is a null set for infinitely many measures of the sequence P n k . More precisely, there is a subsequence n

′ k of n k such that (x 0 -η/2, x 0 ] ⊂ k∈N A n ′ k or [x 0 , x 0 + η/2) ⊂ k∈N A n ′ k . Hence, we can choose A ∞ = (x 0 -η/2, x 0 ) or A ∞ = (x 0 , x 0 + η/2) such that |A ∞ | = η/2 > 0 and P n ′ k (A ∞ ) = 0 for all k ∈ N.
The latter assertion implies that lim inf n→∞ P n (A ∞ ) = 0. Recall that the P n converge weakly to P ∞ and that A ∞ is an open set. Therefore, the Portmanteau theorem allows us to conclude that P ∞ (A ∞ ) = 0.

Before proving consistency, recall the definition of the Lévy distance: For probability distributions P 1 , P 2 with cumulative distribution functions F 1 , F 2 , define

d(P 1 , P 2 ) := inf{δ > 0 | F 1 (x -δ) -δ ≤ F 2 (x) ≤ F 1 (x + δ) + δ ∀x ∈ R}.
For a sequence P n of probability distributions, one has that P n D -→ P if and only if d(F n , F ) → 0 as n → ∞.

In other words, d metrizes the weak convergence. Now define, for probability distributions P X and real numbers σ, ∆( P X , σ; P X , σ) = d( P X , P X ) + | σ -σ|.

Remark that ∆(P X n , σ n ; P X , σ) → 0 if and only if P X n D -→ P X and σ n → σ (and hence N σn D -→ N σ ).

Theorem 3.5 (Consistency). Let R, η > 0 and suppose that in the deconvolution model (1), we have P X ∈ P η R . Then, any minimum distance estimator ( P X n , σ n ) on P η R is consistent, that is, we have ∆( P X n , σ n ; P X , σ) → 0 almost surely.

Proof. We have seen in Lemma 3.2 that the considered estimator is admissible. Now we show that under the assumptions of this theorem, this implies ∆( P X n , σ n ; P X , σ) → 0. The proof is by contradiction. Assume there is a δ > 0 and an increasing sequence (n k ) k∈N in N such that ∆( P X n k , σ n k ; P X , σ) ≥ δ ∀ k ∈ N. Lemma A.1 furnishes a subsequence (n ′ k ) k∈N of (n k ) k∈N , a probability measure P X ∞ and a constant σ ∞ ≥ 0 such that

P X n ′ k D -→ P X ∞ and N b σ n ′ k D -→ N σ∞ . (7)

  let A = [a 1 , a 2 ] be some interval of positive length |A| = a 2a 1 and B = [b 1 , b 2 ] another interval with |B| < |A| and P 2 (B) > 0. By definition of the convolution and in view of the independence of X and ε in our model, we can write (P 2 * N σ3 )(A) = (P 2 ⊗ N σ3 )(S A ), where S A = {(x, y) ∈ R 2 | x + y ∈ A} and ⊗

  d. sample {Y 1 , . . . , Y n } from the model (1). Let

	D -→ denote convergence n , σ n ) of (P X , σ) is called consistent if, almost surely, P X in distribution. An estimator ( P X n D -→ P X and σ n → σ as n → ∞. For a consistent estimator, we always have P X n * N b σn D

Acknowledgements

The authors would like to thank the two reviewers for their comments that helped to improve the manuscript. This work was supported by the IAP research network grant no. P6/03 of the Belgian government (Belgian Science Policy) and by the "Fonds Spéciaux de Recherche" from the Université catholique de Louvain.

ACCEPTED MANUSCRIPT

Proof. The empirical characteristic function ϕ Y n converges almost surely pointwise to ϕ Y . By Lebesgue's Theorem, this implies ρ(ϕ X , ϕ σ ; ϕ Y n ) → 0 almost surely. Applying the triangle inequality and using ϕ X ϕ σ = ϕ Y , we obtain ρ( ϕ X n , ϕ b σn ; ϕ Y ) ≤ ρ( ϕ X n , ϕ b σn ; ϕ Y n )+ρ(ϕ X , ϕ σ ; ϕ Y n ). Because of (3), we can write ρ( ϕ X n , ϕ b σn ; ϕ Y n ) ≤ ρ(ϕ X , ϕ σ ; ϕ Y n ) + δ n , so that we conclude that, almost surely, ρ( ϕ X n , ϕ b σn ; ϕ Y ) ≤ 2ρ(ϕ X , ϕ σ ; ϕ Y n ) + δ n → 0. We choose an element ω ∈ Ω of the underlying probability space such that this convergence holds. As the integrand in ( 2) is non-negative, it follows that a 0 ϕ X n (t)ϕ b σn (t)dt -→ a 0 ϕ Y (t)dt for all a ∈ R as n → ∞. We have shown that the integrated characteristic functions of the measures P X+ε n converge to the integrated characteristic function of the probability measure P Y , so that, applying Theorem 6.3.3 from [START_REF] Chung | A course in probability theory[END_REF] ) → 0 almost surely by Lebegue's dominated convergence theorem. Thus, the sequence ( P X n , σ n ) is a candidate for a minimum distance estimator when δ n decreases sufficiently slow. It is easily verified that this estimator is admissible but not consistent. The following consideration shows in which way we have to restrict the class P 0 in order to obtain consistency. Assume that P X ∈ P 0 and let ( P X n , σ n ) be an admissible estimator. By virtue of Lemma A.1 below, admissibility implies the existence of an increasing sequence (n k ) k∈N , some probability measure P X ∞ , and σ ∞ ≥ 0 such that

and hence, due to admissibility and by uniqueness of the weak limit,

It follows from (4) that a necessary condition for P X n to be consistent is P X ∞ = P X . In view of ( 5) and Theorem 2.1, this is equivalent to P X ∞ ∈ P 0 . But this may not be the case in spite of all P X n k lying in P 0 as this class is not closed under convergence in distribution as the above counterexample shows.

Consistency

In Remark 3.3 we have seen that in order to show consistency of the minimum distance estimator, we need to restrict the set of considered distributions. For R, η > 0, let

Indeed, this choice avoids the problem of possibly obtaining a sequence of estimators the weak limit of which lies outside the identified class P 0 , as the following lemma shows. Note that for a positive random variable X, we have P X ∈ P η R for any choice of η and R.

ACCEPTED MANUSCRIPT

Denote the characteristic functions of P X ∞ and N σ∞ by ϕ X ∞ and ϕ σ∞ , respectively. Since weak convergence implies pointwise convergence of the corresponding characteristic functions, we obtain by Fatou's Lemma

As h is strictly positive and characteristic functions are uniformly continuous on R, we conclude that

Lemma 3.4 allows us to deduce from (7) that P X ∞ ∈ P 0 . Consequently, Theorem 2.1 ensures that P X ∞ = P X and σ ∞ = σ. Together with (7), this implies ∆(

This contradicts the assumption and completes the proof.

Conclusion

We have considered the problem of density deconvolution from one single contaminated sample with uncertainty in the error distribution and we have shown a minimum distance estimator to be consistent in this model. The estimation procedure presented here is inspired by a similar estimator suggested by [START_REF] Neumann | Deconvolution from panel data with unknown error distribution[END_REF] in the context of panel data. Neumann proposes an identification assumption which also implies consistency. This condition is expressed in terms of characteristic functions.

Unlike this, the focus of the present note is on a weak identification condition in the time domain which reflects the properties of the involved distributions in a more natural way. In Section 2 we have proposed such an assumption. An additional difficulty arises from the fact that this condition is too weak to imply consistency, which motivates the definition of admissibility. However, a slight restriction of the considered class of distributions is sufficient to conclude, as the theorem shows.

As far as convergence rates under additional assumptions on the characteristic functions are concerned, the reader may consult the work of [START_REF] Butucea | Minimax estimation of the noise level and of the deconvolution density in a semiparametric deconvolution model[END_REF] in which rates are developed and shown to be minimax-optimal.

As for practical computability, the infinite-dimensional minimization problem (3) could be reduced to a finite-dimensional one by considering σ ∈ Σ n and ϕ X ∈ Φ Pn only, where Σ n becomes dense in R as n → ∞ and P n = { Nn j=1 α j δ xj | α j ∈ R} is a collection of purely atomic probability distributions which grows with n. Such a procedure obviously involves the choice of a number of parameters. A further development of such a procedure as well as an investigation of its behavior and sensitivity to these choices goes beyond the scope of this work and appears to be a relevant and intersting topic for further research.

A. Technical lemma

Lemma A.1. Let Q n be a sequence of probability distributions and σ n a sequence of positive real numbers.

Suppose further that (Q n * N σn ) n∈N converges weakly to some probability distribution. Then, there exist an

as n → ∞, where N 0 := δ 0 denotes the Dirac measure by convention.

Proof. By Helly's selection theorem, there is a subsequence (n k ) k∈N and a subprobability measure [START_REF] Chung | A course in probability theory[END_REF]). We show below that the σ n k are bounded from above such that they have a convergent subsequence; without loss of generality, say σ n k → σ ∞ for some σ ∞ ≥ 0. Proposition 3.1 from [START_REF] Jain | Vague convergence of sums of independent random variables[END_REF]

By assumption, the same sequence converges weakly, and hence vaguely, to some distribution, so the uniqueness of the vague limit of measures on locally compact spaces implies Q ∞ (R) = 1 because (µ * ν)(R) = µ(R) ν(R) for any two finite measures µ and ν on (R, B(R)). The Portmanteau Theorem then implies Q n k D -→ Q ∞ , which was our claim.

It remains to prove that σ n k is bounded from above. We show that otherwise the sequence (Q n k * N σn k ) k∈N would not be tight, which contradicts its weak convergence. Random variable notation is more convenient for this argument, so let

Fix δ = (1/12) and J > 0, and let J = [-J, J]. Put I + j = [3jJ, (3j + 1)J] and I - j = [-(3j + 2)J, -(3j -1)J], and let I := j≥0 I + j be the disjoint union of the I + j . Because of the monotony of the normal density on [0, ∞), we have P[V k ∈ I + j ] > (1/3)P V k ∈ [3jJ, 3(j + 1)J) . The disjoint union over j ≥ 0 of the intervals on the right hand side of this inequality is [0, ∞), and P[V k ≥ 0] = (1/2). Thus, we have P[V k ∈ I] > (1/6).

We can now write P[W k ∈ J ] < (5/6) + (1/6)P[W k ∈ J | V k ∈ I], and it is sufficient to prove that the conditional probability appearing in this inequality is less than (1/2) for some k.

It is easy to check

. By construction, V k ∈ I + j and W k ∈ J together imply U k ∈ I - j . Using further the monotony of the normal density on [0, ∞), we deduce that

As the I - j are pairwise disjoint, the sum is bounded from above by 1, and hence

. If σ n k is unbounded, k can be chosen in such a way that the right hand side of this inequality is less than (1/2), which completes the proof.


Piet Groeneboom
email: p.groeneboom@tudelft.nl

Geurt Jongbloed
email: g.jongbloed@tudelft.nl

Generalized continuous isotonic regression

Keywords: Maximum smoothed likelihood estimation, greatest convex minorant, smooth monotone hazard function

The standard isotonic regression of a vector in IR n is the solution to a least squares projection problem on the cone C of vectors with 'increasing components' in IR n . Generalized isotonic regression problems are isotonic optimization problems that seem to be quite different from isotonic regression problems, but in fact have the same solution. In problems of maximum smoothed likelihood estimation, often continuous versions of generalized isotonic regression problems emerge. In this paper we show that the solutions to such problems are often given by the relatively simple continuous isotonic regression of a function. We illustrate the result using the problem of maximum smoothed likelihood estimation of an increasing hazard function.

Introduction

Consider a vector y ∈ IR n and define the closed convex cone C ⊂ IR n by C = {x = (x 1 , . . . , x n ) :

x 1 ≤ x 2 ≤ • • • ≤ x n }.
Let furthermore w = (w 1 , w 2 , . . . , w n ) be a weight vector in (0, ∞) n . The isotonic regression of y, with weights w, is defined as the minimizer over C of the following quadratic function:

φ q (x) = 1 2 n i=1 (x i -y i ) 2 w i . (1)
Early references on this problem include [START_REF] Ayer | An empirical distribution function for sampling with incomplete information[END_REF] and [START_REF] Grenander | On the theory of mortality measurement, part II[END_REF]. Comprehensive books on isotonic regression are [START_REF] Barlow | Statistical inference under order restrictions[END_REF] and Robertson et al

A C C E P T E D M A N U S C R I P T

ACCEPTED MANUSCRIPT (1988). The most basic result from these references is that the isotonic regression of y with weights w exists, is unique and can be constructed graphically, in two dimensions. Indeed, defining the set of points P i , i = 0, 1, . . . , n, by P 0 = (0, 0) and

P i = i j=1
w j , i j=1 w j y j , i = 1, 2, . . . , n, the k-th component of the isotonic regression of y is given by the left derivative of the greatest convex minorant of the set of points {P i : 0 ≤ i ≤ n}, evaluated at P k . See e.g. [START_REF] Robertson | Order restricted statistical inference[END_REF], Theorem 1.2.1.

Generalized isotonic regression problems are also concerned with minimizing an objective function over the cone C. The objective function is different from ( 1), but of a specific structure. Due to that structure, the solution to this seemingly more complex optimization problem turns out to be exactly the isotonic regression of y on C with some choice of weight vector w. Examples of these problems are the Poisson extremum problem and the Gamma extremum problem that occur in order constrained maximum likelihood estimation. See section 1.5 in [START_REF] Robertson | Order restricted statistical inference[END_REF].

A straightforward continuous version of isotonic regression is concerned with L 2 projections of functions on a compact interval [a, b]. Let w > 0 be an integrable weight function on [a, b] and g ∈ L 2 (w). Moreover, denote by K the subset of L 2 (w) consisting of the nondecreasing functions on [a, b]. The continuous isotonic regression of g on K is defined as the minimizer of the function

ψ q (f ) = 1 2 (f (x) -g(x)) 2 w(x) dx. (2)
In [START_REF] Mammen | Estimating a smooth monotone regression function[END_REF], it is mentioned that in case w ≡ 1, the solution to this problem is given by the right continuous slope of the convex minorant of the parameterized graph in IR 2 given by

t a w(x) dx, t a g(x)w(x) dx : t ∈ [a, b] (3)
In [START_REF] Anevski | Monotone spectral density estimation[END_REF], Lemma 2, this statement is proved for continuous functions g and weight function w ≡ 1. In fact, the proof of that result reveals that it also holds if w ≡ 1 is a continuous strictly positive weight function on [a, b]. This characterization, together with continuity of g, immediately shows that the right-continuous version of the solution is unique on [a, b). In [START_REF] Tantiyaswasdikul | Isotonic smoothing splines under sequential designs[END_REF], a penalized version of the continuous isotonic regression problem is considered, leading to isotonic smoothing splines.

In this paper we show that for a class of specific continuous projection problems, which we call generalized continuous isotonic regression problems, the solution is given by the continuous isotonic regression in a similar way as in the discrete situation. This implies that so-called maximum smoothed likelihood estimators in various models can be characterized and constructed via a simple convex minorant interpretation.

In section 2, we state the generalized continuous isotonic regression problem, give two examples from the literature concerning its use and prove that its solution is precisely that of a related continuous isotonic regression problem. In section 3, we introduce a maximum smoothed likelihood estimator for an increasing hazard rate, characterize it in terms of a convex minorant of a smooth graph in IR 2 and use that characterization to compute it based on a data set of ball bearings from [START_REF] Lieblein | Statistical investigation of the fatugue life of deep-groove ball bearings[END_REF] and also used in [START_REF] Gijbels | Nonparametric testing of a monotone hazard function via normalized spacings[END_REF].

Main result and some examples

In this section we describe the generalized continuous isotonic regression problem. We will follow closely the notation and ideas given for the discrete situation in [START_REF] Robertson | Order restricted statistical inference[END_REF], section 1.5.

Let g be a fixed real-valued function on [a, b] and I ⊂ IR be an interval such that g(x) ∈ I for all x ∈ [a, b]. Moreover, let Φ be a convex function on I with (right-) derivative φ and define for u, v ∈ I

∆ Φ (u, v) = Φ(u) -Φ(v) -(u -v)φ(v).
By the subgradient inequality for convex functions, ∆ Φ (u, v) ≥ 0 for all u, v ∈ I. If Φ is strictly convex, this inequality is strict for u = v. Now, taking a strictly positive continuous weight function w on [a, b], we can define the objective function

ψ(f ) = ∆ Φ (g(x), f (x)) w(x) dx (4)
and minimize this function over the class K. The main result in this paper is that the minimizer of (4) over K is exactly the minimizer of (2) over K.

Theorem 1. Suppose that g is a continuous function on [a, b] and w is a strictly positive, continuous function on [a, b]. Then the minimizer g * of ψ q defined in (2) over K is continuous and unique. Its value at x ∈ (a, b) is given by the slope of the greatest convex minorant of the graph defined by ( 3), evaluated at t = x. Moreover, for all f ∈ K,

ψ(f ) ≥ ψ(g * ),
where ψ is defined in (4). If Φ in the definition of ( 4) is strictly convex, g * is the unique continuous minimizer of ψ over K.

Before proving this theorem, we give two examples of situations where this theorem can be applied. A third example will follow in section 3. All examples are concerned with maximum smoothed likelihood estimation in the spirit of [START_REF] Eggermont | Maximum penalized likelihood estimation[END_REF].

Denoting the empirical distribution function of a sample by IF n and the sampling density (w.r.t. some dominating measure) parameterized by a parameter f by g f , maximum likelihood estimation is usually concerned with maximizing a function of the form

f → log g f (x) dIF n (x)
over f in the set of parameters. This function is the log likelihood function. The smoothed log likelihood function is a similar function, only the empirical distribution function IF n is replaced by some smooth alternative estimator. This could be a kernel estimator, spline estimator based on the data or a linear interpolation of IF n .

Example 1. Maximum Smoothed Likelihood estimation of a decreasing density. Let g be some smooth density estimate on [0, M ] (e.g. a kernel estimate based on a sample). The MSLE of the underlying density based on this g is then defined as the maximizer of

f → log f (x)g(x) dx -f (x) dx
over the class K d , the class of non-increasing nonnegative functions on [0, M ]. Note that the second term in this expression is a relaxation term, guaranteeing the solution to have integral equal to 1. The isotonized kernel estimator fn discussed in Van der [START_REF] Vaart | Smooth estimation of a monotone density[END_REF] is in fact an MSLE, where the initial density estimate is the standard kernel estimator as studied in [START_REF] Silverman | Density estimation for statistics and data analysis[END_REF]. To see this, note that Φ(u) = u log u yields

∆ Φ (g(x), f (x)) = g(x) log g(x) -f (x) log f (x) -(g(x) -f (x))(1 + log f (x)) = c g (x) -g(x) log f (x) + f (x)
and observe that the MSLE based on g coincides with the continuous isotonic regression of g, with weight function w ≡ 1.

Example 2. Maximum Smoothed Likelihood estimation in the current status model. This model is studied in Groeneboom et al. (2009). The maximizer of their smoothed log likelihood (3.1),

F → g n,1 (x) g n (x) log F (x) + 1 - g n,1 (x) g n (x) log(1 -F (x)) g n (x) dx
is characterized in their Theorem 3.1 via a direct argument. Viewing their optimization problem as generalized continuous isotonic regression problem with Φ(u) = u log u + (1u) log(1u), the characterization follows immediately. Indeed, this choice leads to

∆ Φ g n,1 (x) g n (x) , F (x) = c(x) - g n,1 (x) g n (x) log F (x) -1 - g n,1 (x) g n (x) log(1 -F (x))
where c does not depend on F . Using their notation, the solution of this optimization problem is the continuous isotonic regression of g n,1 /g n with weight function g n , i.e. the minimizer of

F → 1 2 g n,1 (x) g n (x) -F (x) 2 g n (x) dx
over all functions F in K.

We now prove Theorem 1 along the lines of the proof given in [START_REF] Robertson | Order restricted statistical inference[END_REF] in the discrete situation. In order to do so, we use the following two lemmas. The proofs of these results can be found in the appendix. First we prove that the continuous isotonic regression g * of a continuous function g is continuous.

Lemma 1. Let g be a continuous function on . Its value at x is given by the slope of the greatest convex minorant of the graph defined by ( 3), evaluated at t = x.

The next lemma is a generalization of Theorem 1.3.5 in [START_REF] Robertson | Order restricted statistical inference[END_REF]. (5)

Proof of Theorem 1. Let g * be the continuous isotonic regression of g with weight function w. Its continuity, uniqueness and representation as derivative of the convex minorant of the graph defined in (3) follow from Lemma 1. Before showing that this g * also minimizes ψ over K, we first establish a useful equality.

The number of 'level intervals' of g * of positive length is at most countable. Continuity of g * implies that these intervals are closed. Denote the intervals by [σ n , τ n ], n = 1, 2, . . .. From Lemma 2, it immediately follows that for any real-valued function Ψ,

a (g(x) -g * (x)) Ψ(g * (x))w(x) dx = {x : g * (x) =g(x)} (g(x) -g * (x)) Ψ(g * (x))w(x) dx = n τn σn (g(x) -g * (x)) Ψ(g * (x))w(x) dx = n Ψ(g * (τ n )) τn σn (g(x) -g * (x)) w(x) dx = 0.
(6)

A C C E P T E D M A N U S C R I P T ACCEPTED MANUSCRIPT

This means that for each f ∈ K

ψ(f ) -ψ(g * ) = (∆ Φ (g(x), f (x)) -∆ Φ (g(x), g * (x))) w(x) dx = Φ(g * (x)) -Φ(f (x)) -(g(x) -f (x))φ(f (x)) +(g(x) -g * (x))φ(g * (x)) w(x) dx (i) = (∆ Φ (g * (x), f (x)) -(g(x) -g * (x))φ(f (x))) w(x) dx (ii) ≥ 0.
In (i) we use ( 6), in (ii) that ∆ Φ (•, •) ≥ 0, the fact that x → φ(f (x)) is monotone and the fact that for each nondecreasing function m on [a, b],

(g(x) -g * (x))m(x)w(x) dx = lim ↓0 -1 (ψ q (g * ) -ψ q (g * + m)) ≤ 0.
In case Φ is strictly convex, inequality (ii) is strict if f = g * on a set of positive Lebesgue measure. 2

Application: smooth increasing hazard

Consider an i.i.d. sample X 1 , . . . , X n from an unknown density g on [0, ∞) with distribution function G and increasing hazard function f . The functions g, G and f satisfy the following relations:

G(x) = x 0 g(y) dy = 1 -exp - x 0 f (y) dy , f (x) = - d dx log(1 -G(x)) = g(x) 1 -G(x)
.

One approach to estimate f would be to use the maximum likelihood estimator. This estimator is well defined and described in [START_REF] Robertson | Order restricted statistical inference[END_REF], section 7.4. It is a step function, piecewise constant between successive observed data points. In some situations it is reasonable to assume that the underlying hazard rate is smooth and monotone. See, e.g., [START_REF] Hall | Nonparametric estimation of hazard rate under the constraint of monotonicity[END_REF] for examples and a "biased bootstrap" procedure to estimate a hazard rate under such an assumption. Alternatively, one could study a maximum smoothed likelihood estimator in the spirit of [START_REF] Eggermont | Maximum penalized likelihood estimation[END_REF] where a smoothed log likelihood function is maximized. Like in the examples described in section 2, the difference between this function and the log likelihood function is that the empirical distribution function of the data that occurs in the log likelihood, is replaced by a smooth estimator of the distribution function. More precisely, let x 1 < x 2 < • • • < x n be the ordered realized sample. Let G n be a smooth estimator of the distribution function G, with continuous density g n . The smoothed log likelihood is then defined by Maximizing this function over K is equivalent to computing the continuous isotonic regression of the function g n /(1 -G n ) with weight function w = 1 -G n . This follows from Theorem 1, using the same function Φ as in Example 1 in the previous section and taking weight function w = 1 -G. Note that g n /(1 -G n ) is the plug-in estimator of the hazard rate based on g n . Having computed the MSLE f of the hazard rate, the plug-in principle can be used to compute the MSLE of G and g within this class: Ĝ

(f ) = g n (x) log g f (x) dx = g n (x) log f (x) dx -g n (x) x y=0 f (y) dy dx = g n (x) log f (x) dx -f (x)(1 -G n (x)) dx.
(x) = 1 -exp - x 0 f (y) dy , ĝ(x) = (1 -Ĝ(x)) f (x)
In Groeneboom and Jongbloed (2009), the MSLE of f is used to construct a test on monotonicity of a hazard function. This approach is an alternative to approaches to this testing problem recently adopted by [START_REF] Gijbels | Nonparametric testing of a monotone hazard function via normalized spacings[END_REF] and [START_REF] Hall | Testing for monotone increasing hazard rate[END_REF].

To illustrate the MSLE, we use the ball bearing dataset from [START_REF] Lieblein | Statistical investigation of the fatugue life of deep-groove ball bearings[END_REF], also used in [START_REF] Gijbels | Nonparametric testing of a monotone hazard function via normalized spacings[END_REF]. This dataset consists of n = 23 failure times (although more could be said on this, see [START_REF] Catoni | The correct "ball bearings" data[END_REF]). The left panel of Figure 1 shows the MLE of [START_REF] Robertson | Order restricted statistical inference[END_REF] on the interval [0, 150], a step function, as well as the MSLE based on a kernel density estimator with bandwidth h = 100n -1/4 ≈ 46 (using a boundary kernel near zero). The middle panel of Figure 1 shows the empirical distribution function of the data as well as the estimates based on the MLE and the MSLE. The right panel of Figure 1 shows the density estimates corresponding to the MLE and the MSLE.

Appendix

In this section we prove Lemma 1 and Lemma 2.

Proof of Lemma 1: Existence of g * and its representation as derivative of the convex minorant of the graph defined in (3) follows from Lemma 1 and Lemma 2 in Anevski Suppose the right continuous function g * has a point of discontinuity x ∈ (a, b), i.e. g * (x)g * (x -) = δ > 0. Then x belongs to the support of the Stieltjes measure d(g * ) and hence, by Lemma 1 in [START_REF] Anevski | Monotone spectral density estimation[END_REF], t a (g(y)g * (y)) w(y) dy = 0 if t = x.

Since for each t ∈ [a, b] this integral is greater than or equal to zero and it is (as an integral) differentiable on (a, b), g * (x) = g(x) necessarily. Now define, for small > 0,

g * (y) = g * (y) 1 [a,x-) (y) + 1 [x,b] (y) + g(x)1 [x-,x)) (y).
Note that g * ∈ K. Using continuity of g and monotonicity of g * , we can choose sufficiently small to guarantee that |g(y)g(x)| < δ/4 and |g * (y)g(y)| > 3δ/4 for y ∈ [x -, x). This yields for > 0 sufficiently small

ψ q (g * ) -ψ q (g * ) = 1 2 x x- (g(x) -g(y)) 2 -(g * (y) -g(y)) 2 w(y) dy ≤ δ 2 32 - 9δ 2 32 x x- w(y) dy = - 1 4 δ 2 x x- w(y) dy < 0 contradicting ψ q (g * ) ≥ ψ q (g * ). 2
Proof of Lemma 2. If L c = τ is a singleton, τ is a point of increase of g * , and hence belongs to the support of the Stieltjes measure d(g * ). By Lemma 1 in [START_REF] Anevski | Monotone spectral density estimation[END_REF], this implies that τ a g * (y)w(y) dy = τ a g(y)w(y) dy and similar to in the proof of our Lemma 1, it follows that g * (τ ) = g(τ ).

If L c has positive Lebesgue measure, it is necessarily a closed interval. This follows from monotonicity of g * and continuity of both g and g * . Consider such an interval [σ, τ ], with a < σ < τ < b and define c as in (5). Now suppose that the value of g * on the interval [σ, τ ] is strictly smaller than c. For > 0, define g * (x) = g * (x) if x ∈ (σ, τ + ) g * (τ + ) if x ∈ (σ, τ + ) Then ψ q (g * )-ψ q (g * ) = 1 2 τ σ + τ + τ (g(x)g * (x)) 2 -(g(x)g * (x)) 2 w(x) dx = I (1) +I (2) where (using that g * (x) = g * (τ ) for x ∈ [σ, τ ]) Define δ = cg * (τ ) > 0. Using continuity of g * , we can make sure by taking > 0 sufficiently small, that

I (1) = 1 2 (g * (τ + ) -g * (τ ))
I (1) = 1 2 (g * (τ + ) -g * (τ ))) τ σ
(2cg * (τ + )g * (τ ))) w(x) dx > 1 4 (g * (τ + )g * (τ ))δ.

On the other hand,

|I (2) | ≤ ( g ∞ + g * ∞ )(g * (τ + ) -g * (τ )) w ∞
Combining these inequalities, we obtain that, for > 0 sufficiently small, ψ q (g * )ψ q (g * ) > 1 8 τ σ w(x) dx • (g * (τ + )g * (τ ))δ, contradicting that g * = argmin f ∈K ψ q (f ). If g * (τ ) > c, a similar argument, using the point σinstead of τ + , can be used. If σ = a and / or τ = b, the argument is a bit simpler. Instead of using σor τ + to define g * , this function can be defined directly on the first (in case of σ = a) or last (in case τ = b) interval.

2

  [a, b], w a continuous strictly positive function on[a, b]. Then the continuous isotonic regresion g * of g with weight function w exists, is unique and continuous on[a, b]

Lemma 2 .

 2 Let g be a continuous function on [a, b] and w a continuous strictly positive function on [a, b]. Let g * be the continuous isotonic regression of g with weight function w. Define, for each c ∈ IR, the level set L c = {x ∈ [a, b] : g * (x) = c}. If L c consists of exactly one point τ , necessarily g * (τ ) = g(τ ). If L c has positive Lebesgue measure, then L c = [σ, τ ], where a ≤ σ < τ ≤ b. Moreover,

  b

Figure 1 :

 1 Figure 1: The MLE's and MSLE's based on the ball bearing data on the interval [0, 150], h = 46. From left to right: the hazard functions; the empirical distribution function, the estimate corresponding to the MLE (red) and the MSLE (blue); the density estimates.

A

  C C E P T E D M A N U S C R I P TACCEPTED MANUSCRIPTand Soulier (2009), with obvious modification to incorporate strictly positive continuous w ≡ 1. Uniqueness of the right continuous minimizer follows from strict convexity of the function ψ q .

  x)g * (τ + )g * (τ )) w(x) dx A C C E P T E D M A N U S C R I P T (τ + )g * (x))(2g(x)g * (x)g * (τ + ))w(x) dx.

Acknowledgement We thank two anonymous referees for their suggestions and remarks.