Sophie Allain

Stéphane Méric
email: stephane.meric@insa-rennes.fr

Hongquan Wang

Eric Pottier

First Investigations on Influence of Resolution Cell Size for Monitoring Bare Agricultural Soil by using Ground Based SAR (PoSAR) System Measurements

The objective of this study is to evaluate the influence of resolution cell size of a SAR image in order to monitor bare agricultural soil. In order to achieve this objective, we propose to use Ground Based (GB) SAR experiments, in which the specific levels of surface roughness and soil moisture are set deliberately. For these first investigations, we propose to study the backscattering coefficients σHH for different surface roughness values regarding multi-resolution. These first results we obtained with our GB-SAR system (called PoSAR for Pocket SAR) are in agreement with previous experiments on indoor SAR measurements (JRC dataset).

I. INTRODUCTION

Polarimetric SAR provides an opportunity to monitor the spatial and temporal distributions of the soil parameters in both regional and continental scale. Especially, polarimetric SAR response from land cover can be significantly influenced by the following two parameter categories: [START_REF] Zribi | A new empirical model to retrieve soil moisture and roughness from c-band radar data[END_REF] Target parameters, such as roughness, soil moisture, soil texture and vegetation status; and (2) Sensor parameters, such as frequency, polarization and incidence angle. Thus, theoretical or empirical models are normally used to convert the measured backscattering coefficients into soil moisture and surface roughness. Moreover, representative algorithms are established which relate the backscattering coefficients and co/cross polarization ratio directly with surface roughness and soil moisture using semiempirical coefficients. Consequently, on the one hand, retrieval of surface parameters from multi polarimetric SAR measurement has been extensively investigated in the past decades [START_REF] Zribi | A new empirical model to retrieve soil moisture and roughness from c-band radar data[END_REF]. On the other hand, the multi-angular SAR images have also potentials to improve the estimation of soil parameters [START_REF] Zribi | A new empirical model to retrieve soil moisture and roughness from c-band radar data[END_REF]. In that way, multi-angular and multi-polarization configuration are compared in [START_REF] Molineaux | Using multi-angular radar data to discriminate the influence of rough surface scattering on soil moisture inversions over bare soils[END_REF], indicating the sensitivity of multi-angular configuration to surface roughness is 10 times than multipolarization, which encourages the continuous exploration of multi-angular models. For instance, it was proposed in [START_REF] Zribi | A new empirical model to retrieve soil moisture and roughness from c-band radar data[END_REF] an original surface roughness parameter to include horizontal and vertical roughness statistical information, and found that backscattering coefficient difference between two SAR acquisitions depends more on surface roughness than soil moisture, in accordance with the results of [START_REF] Molineaux | Using multi-angular radar data to discriminate the influence of rough surface scattering on soil moisture inversions over bare soils[END_REF].

Nowadays, operational spaceborne and airborne SAR provide major observational images to monitor soil status. Nevertheless, during the time interval of multi-angular image acquisitions, the surface status (soil moisture or roughness) may changes. Moreover, it is quite difficult to accomplish all the ground truth measurements over large areas at the same time. These are obvious constraints with currently orbiting radar systems which can not acquire multi-angular images at the same time. To overcome this drawback, the Ground Based (GB) SAR system provides an opportunity to gap the spaceborne and airborne SAR system limitations. It is capable of acquiring multi-incidence angle images in quasi-invariant soil and roughness conditions.

In this study, we aim to analyze the influence of resolution cell size on the sensitivity of polarimetric parameters to bare soil. The GB-SAR experiments we achieved are analyzed and firstly compared with the behaviors of indoor polarimetric SAR measurements. Section II describes the GB-SAR PoSAR experiments and the soil conditions we set up for the measurements. In section III, the indoor polarimetric SAR measurements (JRC dataset) we used for the comparison are introduced. Moreover, the influence of SAR image resolution is analysed as well as the soil parameter regarding several roughness conditions. Finally, main conclusions are presented in section IV.

II. EXPERIMENT DESCRIPTIONS A. GB-SAR system measurement

A very high resolution GB-SAR PoSAR system is developed in the IETR [START_REF] Ferro-Famil | PoSAR: A VHR tomographic GB-SAR system application to snow cover 3-D imaging at X and Ku bands[END_REF]. This system is comprised of vector network analyzer, accurate positioning system, antennas, and 3-m rail which undertakes the system. The signal is generated by vector network analyzer in stepped-frequency mode, and then transmitted during system movement which is controlled by the accurate positioning system. The GB-SAR system obtain each line of the image in a stationary state, thus the Doppler centroid is exactly zero. In addition, the value of Doppler rate can be derived by the geometry configuration between the system and target. These accurate estimations of Doppler parameters is essential for signal synthesis. Two GB-SAR During each GB-SAR PoSAR operation, in situ measurements of soil moisture and surface roughness are implemented. In the one hand, surface roughness was measured by a simple and fast chain approach proposed by [START_REF] Saleh | Soil roughness measurements:chain method[END_REF]. It is based on the justification that when a chain of given length L1 (146.5 cm with 2.2 cm linkage length in our case) is placed straightly along a surface, the covering horizontal distance L2 decreases as surface roughness increases. Therefore, the Saleh Roughness Factor (SRF) is defined as: SRF = 100(1 -L2/L1). This SRF roughness descriptor is used in our study. The SRF was measured 10 times uniformly in two perpendicular directions so as to represent the entire field Fig. 1(a). Moreover, considering the surface roughness is frequently characterized as random roughness (s:root mean square height of surface profile) measured using laser/pin roughness meter, a transformation relationship s(cm) = a SRF b is proposed in [START_REF] Gilley | Random roughness assessment by the pin and chain method[END_REF]. The coefficients a and b are in function of rainfall amount [START_REF] Gilley | Random roughness assessment by the pin and chain method[END_REF], and also affected a lot by the linkage length as reported [START_REF] Merrill | Use of chain set for scale-sensitive and erosion-relevant measurements of soil surface roughness[END_REF]. In our study, the coefficients for roughness scale transformation are derived as a = 0.5072, b = 0.7867 based on the regression between our measurements using laser and synchronous sampling data using chain approach Fig. 1(b). For this roughness measurement validation, four roughness conditions (ks = 1.9, 2.8, 3.0, 3.9, where k wave number and s standard deviation of surface profile) are set by the tillage operation Fig. 2(a). On the other hand, soil moisture at depth 3.8cm was measured using a calibrated Time Domain Reflectometry (TDR) with samples distributed homogeneously over the test site so as to obtain a representative mean soil moisture value. Different soil moisture status are set by using the Deltalab Microprocessor Controlled Spray System, EID 330 that simulates rain events in still air. Rainfall is simulated by a constant speed oscillating nozzle (Deltalab, Tec Jet SS 6560) at a height of 3.8 m as shown in Fig. 2(b).

In the multi-angular SAR measurement campaign on July, 2013, the GB-SAR system, as shown in Fig. 3 • and the beam width in range direction is 54 • (see Fig. 4). Moreover, a bare soil is selected to carry out this experiment for characterizing different roughness status in a same SAR image (see Fig. 5). The soil moisture can be considered as isotropic over this bare surface. Our measurement is carried out at frequency between 4.9 GHz and 5.9 GHz, with 801 frequency samples. On 16 October 2013, the system is mounted on a scaffold with height 6.3 m as shown in Fig. 6(a). The four antennas from left to right shown in Fig. 6(b) correspond H, V, H, V polarization respectively. The operating central frequency is set to 5.4 GHz with a 500 MHz bandwidth. The angle of Look Of Sight (LOS) is 38 • and the beam width in range direction is 54 • (see Fig. 7).

III. SAR DATA ANALYSIS

At this time, we only present the first results we obtained about the influence of resolution cell size by processing GB- SAR PoSAR data and regarding the soil properties (roughness under one moisture status of around 20%).

A. GB-SAR PoSAR data processing

The received signal is synthesized using back projection algorithm [START_REF] Cumming | Digital processing of synthetic aperture radar data[END_REF] and the basic radar scene which is exhibited in the Fig. 8 is used for polarimetric and radiometric calibrations. Firstly, an antenna correction is applied to compensate the antenna pattern. Then the polarimetric calibration is implemented on the synthesized images using approach [START_REF] Ainsworth | Orientation angle preserving a posteriori polarimetcalibration SAR calibration[END_REF] to remove the cross talk and channel imbalance. Finally, the radiometric calibration is accomplished by calculating the theoretical radar cross section of trihedral corner reflector [START_REF] Van Zyl | Calibration of polarimetric radar images using only image parameters and trihedral corner reflector responses[END_REF] and comparing with the integral energy of the trihedral within the image. As shown in Fig. 9, the trihedral corner reflector deployed within the scene is well synthesized with azimuth and range resolution around 15 cm in July and around 30 cm in October.

B. JRC dataset

In order to compare and validate GB-SAR measurements, multi-resolution indoor SAR measurements performed at the European Microwave Signature Laboratory (EMSL) are analyzed. SAR data were acquired for different cell resolution at 14 GHz for an incidence angle of 40 • . The target is a smooth isotropic Gaussian surface with the following roughness parameters values: correlation length L c of 6 cm and rms heights of 0.4 cm and 2.5 cm for smooth and rough artificial surface respectively. The resolution values R are equal to 4.2 cm, 4.8 cm, 5.4 cm, 6 cm, 6.6 cm, 24 cm. The backscattering power in HH polarization (see Fig. 10) is observed for different resolution sizes. As described in [START_REF] Nesti | Comparison of backscattered signal statistics as derived from indoor scatterometric and sar experiments[END_REF], [START_REF] Allain | Influence of resolution cell size for surface parameters retrieval from polarimetric sar data[END_REF], the power decreases with resolution with the same behavior for the two surfaces. Moreover, the observation of two SAR images with spatial

C. Comparisons with our GB-SAR measurements

The GB-SAR images obtained during the October campaigns are focused with different resolutions by varying the frequency bandwidth and azimuth aperture extent [START_REF] Wang | Monitoring Bare Agricultural Soil: Comparison between Ground Based SAR PoSAR System Measurements and Multi-angular RADARSAT-2 Datasets[END_REF]. Four different surface roughness conditions are set: rms height equal to 1.7 cm, 2.5 cm, 2.7 cm and 3.5 cm respectively. In accordance with the behaviors of JRC dataset, the HH polarization signature in our GB-SAR measurements decreases with resolution (see Fig. 12). However, the behaviors of co polarization signature HH is different from the JRC dataset and especially concerning the roughness sensitivity. Indeed, the two surface roughness conditions in our GB-SAR experiments are not so sensitive as in JRC dataset. Furthermore, two SAR images in HH polarization with resolution 30 cm and 120 cm are shown in Fig. 13, indicating the energy distributions are more homogeneous at a larger resolution cell. The SAR responses vary due to the resolution cell differences.

IV. CONCLUSION

In this study, multi-resolution polarimetric SAR images acquired by GB-SAR and in anechoid chamber are analyzed to characterize the influence of resolution cell size on radar response over bare surfaces. The main objective of this first study is to validate our new GB-SAR dataset regarding the JRC dataset. Consequently, we can say that both processing we apply to our radar measurements make possible to deeply study our SAR images from bare surface we can monitor. Also, this first study shows that the reflected power decreases with the resolution. For the further version of this proposal, we plan to study the roughness surface and periodic surface effects on the SAR image with respect the range resolution. Moreover, SAR polarimetric descriptors will be studied in order to retrieve surface parameters in case of agricultural bare soils.

Fig. 1 :

 1 Fig. 1: (a) Surface roughness SRF measurements; (b) Synchronous samplings using laser and chain.

  (a) is positioned on a tower with a height of around 12 m. The four antennas from left to right shown in Fig. 3(b) correspond V, H, V, H polarization respectively. The angle of Look Of Sight (LOS) is

Fig. 2 :Fig. 3 :

 23 Fig. 2: Different (a) surface roughness; (b) soil moisture setting up in GB-SAR experiments.

Fig. 4 :Fig. 6 :

 46 Fig. 4: Geometric slant configuration in July

Fig. 7 :

 7 Fig. 7: Geometric slant configuration in October

Fig. 8 :Fig. 9 :

 89 Fig. 8: The trihedral corner reflectors within scene

Fig. 10 :

 10 Fig. 10: Backscattering power in JRC datasets.

Fig. 11 :Fig. 12 :

 1112 Fig. 11: SAR images with JRC datasets

Fig. 13 :

 13 Fig. 13: Backscattering power in our GB-SAR experiments

ACKNOWLEDGMENT

The authors would like to really thank Cécile Leconte and Frédéric Boutet for their help in ground truth measurements.


Eric Beutner
email: e.beutner@maastrichtuniversity.nl

Laurent Bordes
email: laurent.bordes@univ-pau.fr

Laurent Doyen
email: laurent.doyen@imag.fr

Keywords: Recurrent event data, Semi-parametric statistical model, Effective age process, Profile likelihood inference, virtual age process

We consider a semi-parametric model for recurrent events. The model consists of an unknown hazard rate function, the infinite-dimensional parameter of the model, and a parametrically specified effective age function. We will present a condition on the family of effective age functions under which the profile likelihood function evaluated at the parameter vector θ, say, exceeds the profile likelihood function evaluated at the parameter vector θ, say, with probability p. From this we derive a condition under which profile likelihood inference for the finite-dimensional parameter of the model leads to inconsistent estimates. Examples will be presented. In particular, we will provide an example where the profile likelihood function is monotone with probability one regardless of the true data generating process.

The failure of the profile likelihood method for semi-parametric effective age models

Introduction

Recurrent event data arise from the study of processes that generate events repeatedly over time. Such processes occur in many settings such as biomedicine, clinical trials and engineering to mention a few. For a list of references and some examples of recurrent event data see, for instance, [START_REF] Nelson | Recurrent Events Data Analysis for Product Repairs, Disease Recurrences, and Other Applications[END_REF], [START_REF] Cook | The Statistical Analysis of Recurrent Events[END_REF] and [START_REF] Aalen | Survival and Event History Analysis: A Process Point of View[END_REF]. In this article, our starting point is a semi-parametric model for recurrent events that was introduced by [START_REF] Peña | Models for recurrent events in reliability and survival analysis[END_REF]; see also Equation (1) below. Among other things, the model incorporates the effects of interventions after each event occurrence through an effective age process (or virtual age process). Probably the best known effective age process is the one arising from a renewal process where, after each event occurrence, the effective age is set back to zero. For further information on effective age processes see also [START_REF] Hollander | Nonparametric Methods for Repair models[END_REF], [START_REF] Last | Asymptotic and monotonicity properties of some repairable systems[END_REF], [START_REF] Lindqvist | On the statistical modeling and analysis of repairable systems[END_REF], and [START_REF] Peña | Dynamic Modeling and Statistical Analysis of Event Times[END_REF]. Statistical results for the model introduced in Peña and Hollander ( 2004) can be found in [START_REF] Peña | Semiparametric inference for a general class of models for recurrent events[END_REF] as well as in [START_REF] Dorado | Nonparametric estimation for a general repair model[END_REF] and [START_REF] Adekpedjou | A general class of semiparametric models for recurrent event data[END_REF] who consider sub-models for which they prove consistency and derive weak convergence results. See also [START_REF] Gärtner | Estimation in a general repair model based on left-truncated data[END_REF] who considers a slightly different data collection process. The most general results on consistency and weak convergence were obtained very recently by [START_REF] Peña | Asymptotics for a class of dynamic recurrent event models[END_REF] who restricts the general model given in Equation ( 1) below only by considering the case without frailties. In these articles it is assumed that the effective age function is entirely known. This implies that the way the interventions influence the effective age must be known by the statistician.

Here we question whether this assumption can be weakened in a semi-parametric model. More precisely, we analyse whether the profile likelihood function can be used to derive consistent estimators when the effective age process is not assumed to be known but parametrically specified.

Inference based on the likelihood function and its variants has a long history; for general accounts and a recent review see, for instance, [START_REF] Barndorff-Nielsen | Parametric Statistical Models and Inference[END_REF], [START_REF] Barndorff-Nielsen | Inference and Asymptotics[END_REF], [START_REF] Davison | Statistical models[END_REF], [START_REF] Severini | Likelihood Methods in Statistics[END_REF] and [START_REF] Reid | Aspects of likelihood inference[END_REF]. When the parameter is of the form (ζ, η) with η being a nuisance parameter, inference for ζ is often based on the profile likelihood function or modifications and adjustments to it. This approach has been applied in both parametric and semi-parametric problems; see, for example, [START_REF] Barndorff-Nielsen | Parametric Statistical Models and Inference[END_REF], [START_REF] Barndorff-Nielsen | Inference and Asymptotics[END_REF], [START_REF] Davison | Statistical models[END_REF], [START_REF] Fraser | Likelihood for Component Parameters[END_REF], [START_REF] Mccullagh | Adjustments to profile likelihood[END_REF], [START_REF] Reid | Aspects of likelihood inference[END_REF], [START_REF] Scott | Fitting regression models to case-control data by maximum likelihood[END_REF], [START_REF] Severini | Likelihood Methods in Statistics[END_REF], [START_REF] Severini | Profile likelihood and conditionally parametric models[END_REF], and [START_REF] Slud | Efficient semiparametric estimators via modified profile likelihood[END_REF]. For some semiparametric models like Cox's proportional hazards model, asymptotic normality of the profile (partial) likelihood estimator has been known for awhile; see, for instance, [START_REF] Andersen | Statistical Models Based on Counting Processes[END_REF] and [START_REF] Huang | A maximum pseudo-profile likelihood estimator for the Cox model under length-biased sampling[END_REF] for some recent extension. A general result in the semi-parametric context was proved by [START_REF] Murphy | On profile likelihood (with discussion)[END_REF] who showed that profile likelihood inference for the finite-dimensional parameter behaves like ordinary likelihood inference whenever some functional-analytic conditions are satisfied; see also [START_REF] Hirosi | Efficiency of profile likelihood in semi-parametric models[END_REF] who gave a weaker set of conditions. The result by [START_REF] Murphy | On profile likelihood (with discussion)[END_REF] implies, for example, asymptotic normality of the profile likelihood estimator, and it was successfully applied by many authors and in different settings; see, for example, [START_REF] Breslow | Large sample theory for semiparametric regression models with two-phase, outcome dependent sampling[END_REF], [START_REF] Braekers | Cox's regression model under partially informative censoring[END_REF], [START_REF] Claeskens | An asymptotic theory for model selection inference in general semiparametric problems[END_REF], [START_REF] Xu | Using profile likelihood for semiparametric model selection with application to proportional hazards mixed model[END_REF], and [START_REF] Zeng | A general asymptotic theory for maximum likelihood estimation in semiparametric regression models with censored data[END_REF]. However, it is worth recalling that the standard approach, i.e. profiling out the infinite-dimensional parameter through a right-continuous step function, may lead to an inconsistent estimator for the finite-dimensional parameter of the model. For a single event model with covariates (accelerated failure time model) this can be easily seen; see, for instance, [START_REF] Zeng | Efficient estimation for the accelerated failure time model[END_REF].

As mentioned above, we address the following question: Suppose we profile out the infinite-dimensional parameter by using a right-continuous step function. Can we use the resulting profile likelihood function of the above mentioned semi-parametric model for recurrent events if the effective age process is parametrically specified? Here infinitedimensional parameter refers to the integrated λ that is used in Model 2.1 below, i.e. by the infinite-dimensional parameter we mean the integrated hazard rate function. Denote the set in which the finite dimensional parameter lies by Θ, and let θ ∈ Θ and θ ∈ Θ.

We shall give a condition on the family of effective age processes under which the profile likelihood function at θ is not less than at θ with probability p, say. Additionally, we shall present an extension of this condition under which the profile likelihood function at θ exceeds the profile likelihood function at θ with probability p , say. From this one can easily derive a corollary providing conditions that rule out the possibility to obtain a consistent estimator based on the profile likelihood function. Examples will be presented to which the conditions given in our main results can be easily applied. In particular, we provide an example where p equals one regardless of the true probability measure and of the sample size. In the same example, we find a lower bound for p that does not depend on the true probability measure or on the sample size. Still for the same example, we will infer from our main results that the profile likelihood function is monotonically decreasing with probability one whatever the true probability measure and the sample size. Furthermore, we present statistical models containing the renewal process and the non-homogeneous Poisson process as special cases for which it will turn out that the profile likelihood function at the parameter corresponding to the non-homogeneous Poisson process is never less than at the parameter corresponding to the renewal process regardless of whether the data come from a renewal process or an non-homogeneous Poisson process. The rest of this article is organized as follows. In Section 2 we define the model for recurrent events that we consider, explain its relation to the model introduced by [START_REF] Peña | Models for recurrent events in reliability and survival analysis[END_REF], detail how the profile likelihood is derived, and present our main results as well as examples to which they apply. Following the standard procedure the derivation of the profile likelihood is based on a formula valid if the true model is continuous, whereas profiling out the infinite-dimensional parameter is done w.r.t. a jump function see Section 2.2 and in particular Remark 2.5. In the literature this technique to profile out the infinite-dimensional parameter is often compared to a technique that profiles out the infinite-dimensional parameter using a formula valid for a 'discrete model'. We study this technique for the model considered in this article through a simulation study in Section 3. The theoretical derivations of this technique for the model considered here are carried out in Appendix B. Additional simulation results illustrating the conditions imposed in our main result are also presented in Section 3. Moreover, the results shown there will illustrate the decrease in the above mentioned example. All proofs are given in Appendix A. Appendix C contains results on the identifiability of the examples presented below.

Main result

Throughout, we shall use the following conventions:

N := {1, 2, . . .}, N 0 := N ∪ {0}, R + := {x ∈ R|x ≥ 0}
, the subscript 0 indicates the true, but unknown parameter, x ∧ y stands for the minimum of x and y, and for a function f we denote by f (x-) and f (x+) the left-hand and right-hand limit of f at x, respectively. Convergence in probability is denoted by P → and P means that convergence in probability does not take place. For a simple counting process N we denote by 0 = S 0 < S 1 < S 2 < . . . the sequence of jump times. In the next section we define the model we analyse, then in Section 2.2 we explain the estimators and show that the estimator for the cumulative hazard rate can be considered to be a non-parametric maximum likelihood estimator (NPMLE), and in Section 2.3 we give our main results.

The model

We will consider the following statistical model:

Model 2.1 Let N = {N (s), 0 ≤ s ≤ s * }, s * ∈ R + ,
be a (simple) counting process on some measurable space (Ω, F) endowed with a filtration F = {F s , 0 ≤ s ≤ s * } satisfying the usual conditions. Let P (λ,θ) , λ ∈ Γ, θ ∈ Θ, where Γ denotes the set of all hazard rate functions and Θ ⊂ R d , be a set of probability measures on (Ω, F) such that under

P (λ,θ) the F-compensator A = {A(s), 0 ≤ s ≤ s * } of N is given by A(s) = ˆs 0 Y (u)λ(ε θ (u)) du,
where the process Y = {Y (s), 0 ≤ s ≤ s * } is predictable, non-increasing and fulfils Y (s) ∈ {0, 1}, ∀s ∈ [0, s * ]. For every θ ∈ Θ we have that ε θ = {ε θ (s), 0 ≤ s ≤ s * } is a predictable process with the following additional properties:

(a) ε θ (0, ω) = c 0 P (λ,θ) -a.s. for some

c 0 ∈ R + ; (b) s → ε θ (s, ω) is P (λ,θ) -a.s. non-negative; (c) We have P (λ,θ) -a.s. that s → ε θ (s, ω) is continuous on (S k-1 (ω), S k (ω)], k ∈ N, and differentiable on (S k-1 (ω), S k (ω)), k ∈ N, with positive derivatives. The restriction of ε θ to the random time interval (S k-1 (ω), S k (ω)] is denoted by ε θ k-1 , k ∈ N.
Examples for the Model 2. Model 2.1 is in the same spirit as the model introduced by [START_REF] Peña | Models for recurrent events in reliability and survival analysis[END_REF]; see also [START_REF] Peña | Dynamic Modeling and Statistical Analysis of Event Times[END_REF] and [START_REF] Peña | Semiparametric inference for a general class of models for recurrent events[END_REF]. In their model the compensator is assumed to be of the form

A(s|Z, X(u), 0 ≤ u ≤ s) = ˆs 0 Z Y (u) ρ(N (u-), α) ψ(β T X(u)) λ(ε(u)) du, (1)
where Z is a frailty variable, ρ is a mapping from N 0 to R + of known functional form depending on some unknown parameter vector α ∈ A ⊂ R p with ρ(0; α) = 1, for all α ∈ A, ψ is a known mapping from R to R + with β ∈ B ⊂ R q an unknown parameter vector and X = {X(s), 0 ≤ s ≤ s * } is an R q -valued stochastic process interpreted as the possibly time-varying covariates. The superscript T denotes the transpose. The predictable process ε = {ε(s), 0 ≤ s ≤ s * } fulfils the properties (a)-(c) mentioned in Model 2.1. The parameter of interest for the model defined by Equation ( 1) is thus (α 0 , β 0 , λ 0 ). As mentioned in the introduction statistical results on this model or submodels of it can be found in [START_REF] Dorado | Nonparametric estimation for a general repair model[END_REF], [START_REF] Peña | Semiparametric inference for a general class of models for recurrent events[END_REF], [START_REF] Adekpedjou | A general class of semiparametric models for recurrent event data[END_REF], and [START_REF] Peña | Asymptotics for a class of dynamic recurrent event models[END_REF]. Thus, on one hand the statistical model introduced in [START_REF] Peña | Models for recurrent events in reliability and survival analysis[END_REF]

Z ≡ 1, ρ(N (s-), α) ≡ 1, X(s) = constant, ψ(β T X) = exp(β T X)
, and ε(s) = s) for which profile likelihood inference for β leads to a consistent and asymptotically normally distributed estimator; see, for instance, [START_REF] Andersen | Cox's Regression Model for Counting Processes: A Large Sample Study[END_REF]. We now present some examples for Model 2.1.

Example 2.2 (ARA 1 or Kijima I with non-random repair) For an ARA 1 model we have

ε θ k-1 (s, ω) = s -θ • S k-1 (ω) with Θ = [0, 1]
. This is the same model as a Kijima I model with non-random repair; see, for instance, [START_REF] Kijima | Periodical replacement problem without assuming minimal repair[END_REF], [START_REF] Kijima | Some results for repairable systems with general repair[END_REF] and [START_REF] Dorado | Nonparametric estimation for a general repair model[END_REF]. Notice that this model contains renewal processes with θ = 1 and non-homogeneous Poisson processes with θ = 0. 3

Example 2.3 (ARA ∞ or Kijima II with non-random repair) For an ARA ∞ model we have:

ε θ k-1 (s, ω) = s -θ k-1 l=1 (1 -θ) k-1-l S l (ω) with Θ = [0, 1]
. This is the same model as a Kijima II model with non-random repair; see again [START_REF] Kijima | Periodical replacement problem without assuming minimal repair[END_REF], [START_REF] Kijima | Some results for repairable systems with general repair[END_REF] and [START_REF] Dorado | Nonparametric estimation for a general repair model[END_REF]. Clearly, as in Example 2.2 we see that for θ = 1 we get renewal processes and taking θ = 0 results in non-homogeneous Poisson processes. 3 Example 2.4 Gonzáles et al. ( 2005) fitted a sub-model of the one given in Example 2.3 taking into account covariates via the equation as given in (1) to the data of 63 patients having a subtype of indolent non-Hodgkins lymphomas. They restricted Θ = [0, 1] in Example 2.3 to the discrete set {0, 0.5, 1}. Then a 0 stands for no response to the therapy/intervention, 0.5 means a partial remission, and 1 indicates a perfect intervention.

3

The notion of ARA (Arithmetic Reduction of Age) has been introduced in [START_REF] Doyen | Classes of imperfect repair models based on reduction of failure intensity or virtual age[END_REF]. Further examples for effective age processes can be found in [START_REF] Doyen | Classes of imperfect repair models based on reduction of failure intensity or virtual age[END_REF], [START_REF] Lindqvist | On the statistical modeling and analysis of repairable systems[END_REF][START_REF] Peña | Semiparametric inference for a general class of models for recurrent events[END_REF].

The estimators and NPMLE

We now begin by introducing a profile likelihood method for estimating the unknown parameter vector θ. Let N 1 , . . . , N m be m independent copies of a counting process as described in Model 2.1. For each counting process denote by Y i , 1 ≤ i ≤ m, the process as introduced in Model 2.1 and let 0 = S i,0 < S i,1 < . . . be the jump times of the process N i . By ε θ i we denote the predictable process arising in the definition of N i , 1 ≤ i ≤ m, and by ε θ i,j-1 its restriction to the time interval (S i,j-1 , S i,j ]. Following Jacod (1975) (see also Andersen et al. (1993, II.7) and [START_REF] Peña | Semiparametric inference for a general class of models for recurrent events[END_REF]) the full likelihood process equals

L m,F (s|λ, ε θ , D m (s)) = m i=1 s u=0 Y i (u)λ(ε θ i (u)) N i (∆u) × exp - m i=1 ˆs 0 Y i (u)λ(ε θ i (u)) du , (2)
where at time s the data

D m (s) equals D m (s) := {N 1 (u), . . . , N m (u), Y 1 (u), . . . , Y m (u), 0 ≤ u ≤ s}.
To obtain the profile likelihood function from (2) we first introduce an estimator for Λ 0 the cumulative hazard rate of λ 0 . In doing so, following a technique of [START_REF] Peña | Nonparametric estimation with recurrent event data[END_REF] who extended an idea of [START_REF] Gill | Testing with replacement and the product-limit estimator[END_REF] and Selke (1988), we define double indexed processes; see also [START_REF] Selke | Sequential Analysis of the Proportional Hazards Model[END_REF] who seem to be the first to consider double indexed processes in survival analysis. Below we demonstrate that the resulting estimator may be considered to be a NPMLE. Firstly, define the double indexed process

N θ i , 1 ≤ i ≤ m, by N θ i (s, t) := ˆs 0 Z θ i (u, t) dN i (u), 0 ≤ s ≤ s * , 0 ≤ t < ∞, with Z θ i (u, t) := 1 {ε θ i (u)≤t} , 1 ≤ i ≤ m. N θ i (s, t)
denotes the number of events over the period (0, s] for the ith unit whose effective age at time of occurrence was at most t. Thus, the first time variable s of N θ i stands for the observation time and the second time variable t for the effective age time. Notice that N θ i depends on the effective age parameter θ in contrast to N i . Secondly, we define what has been called the adjusted at risk process (or generalized at risk process)

Y θ i = {Y θ i (s, t), 0 ≤ s ≤ s * , 0 ≤ t < ∞} by Y θ i (s, t) := N i (s-) j=1 γ θ i,j-1 (t) • 1 (ε θ i,j-1 (S i,j-1 +),ε θ i,j-1 (S i,j )] (t) + γ θ i,N i (s-) (t) • 1 (ε θ i,N i (s-) (S i,N i (s-) +),ε θ i,N i (s-) (s∧τ i )] (t), (3)
where the functions γ θ i,j-1 are defined by

γ θ i,j-1 (t) := 1 ε θ i,j-1 ε θ i,j-1 -1 (t)
with (ε θ i,j-1 ) denoting the derivative of ε θ i,j-1 w.r.t. observational time s and (ε θ i,j-1 ) -1 denoting the inverse w.r.t. observational time. With the help of Y θ i one can rewrite the integral arising in the full likelihood (cf. Equation ( 2)) in terms of λ evaluated at the observational time s instead of at the effective age time ε i (s); see Equation ( 6) below. Notice that the jth event of the ith unit contributes to the risk set at the time pair (s, t) if it fulfils three conditions: Firstly, it occurred during the observation period [0,s), secondly the effective age of the jth event of the ith unit is larger than or equal to t, and thirdly the effective age of the ith unit immediately after the intervention succeeding the (j -1)th event is less than t. For further information on the adjusted at-risk process see [START_REF] Peña | Asymptotics for a class of dynamic recurrent event models[END_REF]. We define

S θ m (s, t) := m i=1 Y θ i (s, t). (4)
For fixed θ we define the following method-of-moments estimator Λ m for Λ 0 :

Λ m (s, t|θ) := ˆt 0 J θ m (s, u) S θ m (s, u) m i=1 N θ i (s, du) , (5)
where

J θ m (s, u) := 1 {S θ m (s,u)>0} .
A justification for calling Λ m a method-of-moment estimator can be found in [START_REF] Peña | Semiparametric inference for a general class of models for recurrent events[END_REF] after their Proposition 1. Moreover, for θ known and s fixed this estimator is consistent and converges, after being suitably normalized, to a Gaussian process (cf. [START_REF] Peña | Asymptotics for a class of dynamic recurrent event models[END_REF]). We now demonstrate that for observational time s and effective age parameter θ fixed the estimator Λ m (s, t|θ) can be seen to be a NPMLE. For this, notice first of all that the substitution rule

ˆb a f (x) dx = ˆφ(b) φ(a) f (φ -1 (x))(φ -1 (x)) dx
implies (with f = λ • ε θ and φ = ε θ ) that the full likelihood (2) can be written as

L m,F (s|λ, ε θ , D m (s)) = m i=1 s u=0 Y i (u)λ(ε θ i (u)) N i (∆u) × exp - ˆ∞ 0 S θ m (s, u) dΛ(u) . (6)
Now if we take Λ to be a jump function with jumps at

ε θ k, -1 (S k, ), 1 ≤ k ≤ m, 1 ≤ ≤ N k (s), i.e
. the effective age of the kth unit at the time of the th event, and if we denote these jumps by λ θ k, the log of the full likelihood becomes

log L m,F (s|λ θ , ε θ , D m (s)) = m k=1 N k (s) =1 log λ θ k, - m k=1 N k (s) =1 λ θ k,    (i,j)∈I θ k, γ θ i,j-1 (ε θ k, -1 (S k, ))    - m k=1 N k (s) =1 λ θ k,    i∈I θ,τ k k γ θ i,N i (s-) ε θ k, -1 (S k, )    , (7)
where for every pair (k, ), 1 ≤ k ≤ m, 1 ≤ ≤ N k (s), the sets I θ k, are defined by

I θ k, := {(i, j), 1 ≤ i ≤ m, 1 ≤ j ≤ N i (s-)|ε θ i,j-1 (S i,j-1 +) < ε θ k, -1 (S k, ) ≤ ε θ i,j-1 (S i,j )}
and for every k, 1 ≤ k ≤ m, the sets I θ,τ k k are defined by

I θ,τ k k := {i, 1 ≤ i ≤ m|ε θ i,N i (s-) (S i,N i (s-) +) < ε θ k, -1 (S k, ) ≤ ε θ i,N i (s-) (s ∧ τ i )}.
Clearly, maximizing (7) with respect to λ θ k, has solution given by

λ θ k, = 1 (i,j)∈I θ k, γ θ i,j-1 (ε θ k, -1 (S k, )) + i∈I θ,τ k k γ θ i,N i (s-) ε θ k, -1 (S k, ) . (8)
Hence, upon substituting the NPMLE Λ m for Λ in the full likelihood we obtain from Equation ( 6) that for every fixed θ the resulting log profile likelihood function m,P , up to a constant, equals

m,P (s|θ, Λ m , D m (s)) = - ˆs 0 m i=1 log S θ m (s, ε θ i (w)) dN i (w),
because replacing Λ by Λ m in the argument of the exponential function in the full likelihood (cf. Equation 2), we obtain

exp - ˆ∞ 0 S θ m (s, u) d Λ(u) = exp - m i=1 N i (s) .
It is worth mentioning that the argument of exp being free of the finite-dimensional parameter is not a peculiarity of the model we consider here. For instance, for Cox's proportional hazards model, when plugging in the NPMLE into the full likelihood the argument of the exponential function is free of the regression parameter; see, for instance, [START_REF] Johansen | An Extension of Cox's Regression Model[END_REF].

Remark 2.5 As mentioned in the introduction we followed the standard procedure to derive a NPMLE to profile out the infinite-dimensional component that is to say our NPMLE was taken to be a jump function whereas formula ( 2) is valid only in the continuous case. It is, therefore, common to additionally study the profile likelihood function if instead of formula (2) the corresponding formula for the discrete case is used. We shall do the same here. In Appendix B we derive all relevant formulas as well as the resulting profile likelihood function and in Section 3 this profile likelihood function is studied through a simulation study. 3

Remark 2.6 Notice that (8) boils down to well-known NPMLEs in special cases. Recall that the model in Example 2.2 contains, for instance, the renewal process (corresponding to θ = 1) with effective age process s-S k-1 (ω) and the non-homogeneous Poisson process (corresponding to θ = 0) whose effective age process equals s during the observation period. In these cases noting that γ i,j-1 ≡ 1 we obtain from (8) with T i,j := S i,j -S i,j-1 and τ (i) the ith smallest value among the τ 1 , . . . , τ m the well-known NPMLEs

λ θ=1 k, = 1 m i=1 N i (s-) j=1 1 {T k, ≤T i,j } + m i=1 1 {T k, ≤(s∧τ i )-S i,N i (s-) }
, and

λ θ=0 k, = 1 m -i + 1 , whener τ (i-1) < s k, ≤ τ (i)
with τ (0) := 0; see, for instance, [START_REF] Peña | Nonparametric estimation with recurrent event data[END_REF] and [START_REF] Lawless | The analysis of recurrent events for multiple subjects[END_REF], respectively. 3

Main result

In this section we present our main result. We first state an assumption that is needed in the theorem and provide examples when the assumption is satisfied. We further analyse this assumption in the simulation study in Section 3. The other assumptions made in our main result are illustrated below Theorem 2.10.

Assumption 2.7 For Y θ i , 1 ≤ i ≤ m, as defined in (3) we have for every θ ∈ Θ with probability 1

Y θ i (s * , t) = J i (s * ) j=1 γ θ i,j-1 (t) • 1 (ε θ i,j-1 (S i,j-1 +),ε θ i,j-1 (S i,j )] (t), 0 ≤ t < ∞,
where

J i (s * ), 1 ≤ i ≤ m, are random variables taking values in N 0 . Assumption 2.7 means that γ θ i,j-1 (t) • 1 (ε θ iN i (s-) (S iN i (s-) +),ε θ iN i (s-) (s * ∧τ i )] (t) is either of no relevance for Y θ i (s * , t) or of the form γ θ i,j-1 (t) • 1 (ε θ i,j-1 (S i,j-1 +),ε θ i,j-1 (S i,j )] (t)
. We now give two examples that fulfil Assumption 2.7.

Example 2.8 (Type-II censoring) Let τ i = S i,n i with n i ∈ N, 1 ≤ i ≤ m, and s * ≥ max i τ i . In this case we have s * ∧ τ i = s * ∧ S i,n i = S i,n i . Moreover, if s * > S i,n i , then (ε θ i,N i (s * -) (S i,N i (s * -) +), ε θ i,N i (s * -) (s * ∧ τ i )]
equals the empty set and the representation in Assumption 2.7 holds with

J i (s * ) = N i (s * -) = N i (s * ). Finally, if s * = S i,n i , then we have (ε θ i,N i (s * -) (S i,N i (s * -) +), ε θ i,N i (s * -) (s * ∧ τ i )] = (ε θ i,N i (s * -) (S i,N i (s * -) +), ε θ i,N i (s * -) (S i,N i (s * ) )],
and the representation in Assumption 2.7 holds with

J i (s * ) = N i (s * -) + 1 = N i (s * ). 3
Example 2.9 (Compact support and finite number of interventions) Let λ 0 be such that ´v 0 λ 0 (u) du = ∞ for some v ∈ R + . Additionally, suppose that we consider the model of Example 2.2 with τ i = s * ∧ S i2 , where s * > 2v. By definition of τ i we observe at most two events for the ith unit. Moreover, in the model of Example 2.2, the largest value we can observe for

S i2 equals v + [v -v(1 -θ)] = v + θv which is maximal for θ = 1. Hence, since s * > 2v and ´v 0 λ 0 (u) du = ∞, we have τ i = S i2 , 1 ≤ i ≤ m, so that (ε θ i,N i (s * -) (S i,N i (s * -) +), ε θ i,N i (s * -) (s * ∧ τ i )] = (ε θ i,2 (S i,2 +), ε θ i,2 (S i,2 )]
equals the empty set, and the representation in Assumption 2.7 holds with

J i (s * ) = N i (s * -). 3
We now state our main result whose proof, as mentioned in the introduction, is given in Appendix A.

Theorem 2.10 Let Assumption 2.7 be satisfied, N i , 1 ≤ i ≤ m, etc. be as above, and denote by (P λ 0 ,θ 0 ) m the m-fold product measure of P λ 0 ,θ 0 . Moreover, let θ and θ be such that there exists a c > 0 with the following property: For every t ≥ 0 we have that

(P λ 0 ,θ 0 ) m (γ θ i,j-1 (t) ≤ c, γ θ i,j-1 (t) ≥ c, 1 ≤ i ≤ m, 1 ≤ j ≤ J i (s * ), 0 ≤ t < ∞) = 1.
Then (a) Denote by A m,θ, θ the set of all ω's such that for all pairs (i, j),

1 ≤ i ≤ m, 1 ≤ j ≤ J i (s * ), and all pairs (k, ) 1 ≤ k ≤ m, 1 ≤ ≤ J k (s * ), we have that ε θ i,j-1 (S i,j-1 (ω)) < ε θ k, -1 (S k, (ω))
implies that

ε θ i,j-1 (S i,j-1 (ω)) < ε θ k, -1 (S k, (ω)).
Then we have

(P λ 0 ,θ 0 ) m P,m s * |θ, Λ m , D m (s) ≥ P,m s * | θ, Λ m , D m (s) ≥ (P λ 0 ,θ 0 ) m A m,θ, θ .
(b) Denote by B m,θ, θ the set of all ω ∈ A m,θ, θ for which we additionally have that there are at least two pairs (i, j),

1 ≤ i ≤ m, 1 ≤ j ≤ J i (s * ), and (k, ), 1 ≤ k ≤ m, 1 ≤ ≤ J k (s * ), such that ε θ i,j-1 (S i,j-1 (ω)) < ε θ k, -1 (S k, (ω)) but ε θ i,j-1 (S i,j-1 (ω)) ≥ ε θ k, -1 (S k, (ω)).
Then we have

(P λ 0 ,θ 0 ) m P,m s * |θ, Λ m , D m (s) > P,m s * | θ, Λ m , D m (s) ≥ (P λ 0 ,θ 0 ) m B m,θ, θ .
From Theorem 2.10 one can easily derive a criterion for inconsistency of the maximizer of the log-likelihood function denoted by θm . We state the result as a corollary whose proof is omitted.

Corollary 2.11 Denote by B(θ 0 , ) an -ball around θ 0 and assume that θ is such that for some m ∈ N we have for all m ≥ m that

(P λ 0 ,θ 0 ) m B m,θ, θ ≥ c, c > 0, ∀ θ ∈ B(θ 0 , ). Then θm P θ 0 , as m → ∞.
Before giving examples to which Theorem 2.10 and Corollary 2.11 apply we present the following example that makes Theorem 2.10 and Corollary 2.11 plausible in a simplified setting.

Example 2.12 Consider the model of Example 2.2 (the same reasoning applies to the model of Example 2.3), take m = 2, let both samples be Type-II censored and take arbitrary event times s 1,1 , . . . , s 1,n 1 and s 2,1 , . . . , s 2,n 2 . Then with s = max{s 1,n 1 , s 2,n 2 } the function S 0 2 (s, •) (see Equation ( 4)) equals:

1 (0,s 1,1 ] (•) + 1 (s 1,1 ,s 1,2 ] (•) + . . . + 1 (s 1,n 1 -1 ,s 1,n 1 ] (•) + . . . + 1 (s 2,n 2 -1 ,s 2,n 2 ] (•),
whereas the function S 1 2 (s, •) equals

1 (0,s 1,1 ] (•) + 1 (0,s 1,2 ] (•) + . . . + 1 (0,s 1,n 1 ] (•) + . . . + 1 (0,s 2,n 2 ] (•).
Clearly, plugging in an arbitrary t ≥ 0 into S 0 2 there are at most two intervals that contain t. This is because the intervals coming from the first and the second sample, respectively, do not overlap. On the other hand, looking at the function S 1 2 we see that for small values of t there are n 1 + n 2 overlapping intervals. As t increases the number of overlapping intervals decreases from n 1 + n 2 to n 1 + n 2 -1 to n 1 + n 2 -2 and so on. Apparently, this behaviour does not depend on how the samples s 1,1 , . . . , s 1,n 1 and s 2,1 , . . . , s 2,n 2 were generated. Of course, such a behaviour rules out the possibility to get a consistent estimator based on the profile likelihood function.

3

We now provide examples to which Theorem 2.10 and Corollary 2.11 can be applied.

Example 2.13 (ARA 1 or Kijima I with non-random repair) Consider again the model of Example 2.2 and notice first of all that γ θ i,j-1 ≡ 1 for every θ ∈ [0, 1]. Let 0 ≤ θ < θ ≤ 1. We first discuss part (a) of Theorem 2.10 for this model. For arbitrary positive real numbers x, y and z with y < z we have that

x -θx < z -θy ⇒ x -θx < z -θy.
Indeed, the functions f 1 , f 2 : [0, 1] → R, defined by f 1 ( θ) := x -θx and f 2 ( θ) := z -θy, respectively, are both monotonically decreasing. By assumption f 1 (θ) < f 2 (θ). Moreover, f 1 (1) = 0 and f 2 (1) > 0. Hence, f 1 (θ) < f 2 (θ), ∀θ > θ, because f 1 and f 2 are linear. Hence, for this model (P λ 0 ,θ 0 ) m (A m,θ, θ) = 1 for every m ∈ N whenever θ < θ regardless of (λ 0 , θ 0 ) (recall that we consider a simple counting process so that we have S k, > S k, -1 with probability one and notice that S k, corresponds to z and S k, -1 to y).

We now turn to part (b) of Theorem 2.10. For x, y and z as before and 0

≤ θ ≤ θ < 1 the condition x -θx < z -θy, but x -θx ≥ z -θy is equivalent to z -y 1 -θ + y ≤ x < z -y 1 - θ + y. (9)
Now, if we consider the data generating process of Example 2.8 with n i ≥ 2 and m ≥ 2, then we have for θ < θ from ( 9)

(P λ 0 ,θ 0 ) m B m,θ, θ ≥ (P λ 0 ,θ 0 ) m S 2,2 -S 2,1 1 -θ + S 2,1 ≤ S 1,1 < S 2,2 -S 2,1 1 - θ + S 2,1 = ˆR2 (P λ 0 ,θ 0 ) m s 2,2 -s 2,1 1 -θ + s 2,1 ≤ S 1,1 < s 2,2 -s 2,1 1 - θ + s 2,1 dF S 2,2 ,S 2,1 λ 0 ,θ 0 (s 2,2 , s 2,1 ) = ˆR2 F S 1,1 λ 0 ,θ 0 s 2,2 -s 2,1 1 - θ + s 2,1 -F S 1,1 λ 0 ,θ 0 s 2,2 -s 2,1 1 -θ + s 2,1 dF S 2,2 ,S 2,1 λ 0 ,θ 0 (s 2,2 , s 2,1 ), ( 10
)
where the first equality follows from the independence of the units, and F S 2,2 ,S 2,1 λ 0 ,θ 0 denotes the joint distribution function of (S 2,1 , S 2,2 ) under (P λ 0 ,θ 0 ) m and similar F S 1,1 λ 0 ,θ 0 denotes the distribution function of S 1,1 under (P λ 0 ,θ 0 ) m . Clearly, if λ 0 is such that the corresponding cumulative distribution function F 0 , which equals here F S 1,1 λ 0 ,θ 0 , is strictly increasing on R, the integrand in ( 10) is positive whenever θ < θ. With slightly more effort other cases as, for instance, an F 0 which is constant on some intervals can be discussed. Let us finally consider the condition of Corollary 2.11. Assume that θ 0 ∈ [0, 1] is not equal to 0. Consider [θ 0 -, θ 0 + ] with [θ 0 -, θ 0 + ] [0, 1] and let 0 ≤ θ < θ 0 -. Then, from Equation ( 10) we have a lower bound for (P λ 0 ,θ 0 ) m (B m,θ,θ 0 -) that does not depend on m, m ≥ 2. Clearly, this lower bound also holds for (P λ 0 ,θ 0 ) m B m,θ, θ , θ ∈ [θ 0 -, θ 0 + ], as the integrand in Equation ( 10) is non-decreasing in θ for θ fixed.

3

Example 2.14 (ARA ∞ or Kijima II with non-random repair) Consider again the model of Example 2.3 and notice that we again have that γ θ i,j-1 ≡ 1 for every θ ∈ [0, 1]. Let 0 ≤ θ < θ ≤ 1. Then for arbitrary s 1 < ... < s i-1 and s1 < ... < sk (positive) real numbers the condition in part (a) of Theorem 2.10 reads as

s i-1 -θ i-1 =1 (1 -θ) i-1-s < sk -θ k-1 =1 (1 -θ) k-1-s ⇒s i-1 - θ i-1 =1 (1 -θ) i-1-s < sk - θ k-1 =1 (1 -θ) k-1-s . ( 11
)
This implication may not hold for every pair (θ, θ) with θ < θ regardless of s 1 < ... < s i-1 and s1 < ... < sk ; see Section 3 for more details on that. However, we see that it holds for 0 ≤ θ < 1 and θ = 1 so that (P λ 0 ,θ 0 ) m (A m,θ,1 ) = 1 and Theorem 2.10 now implies

(P λ 0 ,θ 0 ) m (l P,m (s * |θ) ≥ l P,m (s * |1)) = 1, 0 ≤ θ < 1.
Furthermore, for i = 2 and k = 2 we find that Equation (11) boils down to Equation (9) (with x = s 1 , y = s1 and z = s2 ). It therefore follows from Example 2.13 that (P λ 0 ,θ 0 ) m (B m,θ,1 ) > 0, 0 ≤ θ < 1. 3

Simulation results

In this section, as mentioned in the introduction, we provide some simulation results illustrating Assumption 2.7 imposed in Theorem 2.10. We have seen in Example 2.8 that under Type-II censoring Assumption 2.7 holds. It is also clear from Example 2.9 that, in general, under Type-I censoring Assumption 2.7 is not fulfilled. Moreover, in Example 2.14 we left some questions open. Here, we will analyse these questions further by Monte Carlo simulations. Furthermore, as mentioned in the introduction and in Remark 2.5 the behaviour of the "discrete log profile likelihood" is studied by Monte Carlo simulations as well. In all the simulations the hazard rate function used in the definition of Model 2.1 was taken to come from a right truncated Weibull distribution with reliability function S d defined for d > 0 by

S d (t) := e -t 2 -e -d 2 1 -e -d 2 1 [0,d) (t), t > 0.
Hence the corresponding hazard rate function λ d equals

λ d (t) = 2te -t 2 e -t 2 -e -d 2 1 [0,d) (t), t > 0.
We simulated data following an ARA 1 model (see Example 2.2) and an ARA ∞ model (see Example 2.3), respectively, for various values of (θ, d). In addition we introduced two types of censoring: Type-I and Type-II censoring. The sample size m was always taken to be equal to 100. Under Type-I censoring the data observed were given by S i,1 < • • • < S i,r i for (random) r i , i = 1, . . . , 100, and S i,r i +1 , i = 1, . . . , 100, was right censored by τ i = τ , i = 1, . . . , 100, for some non-random τ > 0. In case of Type-I censoring Assumption 2.7 is not fulfilled. If the censoring time is random and equal to S i,r , i = 1, . . . , 100, for some non-random r the data are Type-II censored (see Example 2.8) and Assumption 2.7 holds. Combining the two models, ARA 1 and ARA ∞ , and the two types of censoring one obtains four possible combinations. The results, i.e. ten realizations of the function P,100 , for Type-II censoring and the following values of (θ 0 , d, r) ∈ {0, 0.5, 1} × {5} × {2, 5} are given in Fig. 1 and3. Here θ 0 denotes as in Examples 2.13 and 2.14 the true parameter. Type-I censored data were simulated for (θ 0 , d, τ ) ∈ {0.1, 0.5, 1} × {2} × {2.2} and (θ 0 , d, τ ) ∈ {0, 0.5, 1} × {6} × {5.9}. The reason for taking (θ 0 , d, τ ) = (0.1, 2, 2.2) instead of (θ 0 , d, τ ) = (0, 2, 2.2) is that θ 0 corresponds to a non-homogeneous Poisson process for which with d = 2 the observations are not censored. The results, i.e. again ten realizations of the function P,100 , are given in Fig. 2 and4.

Examples 2.8 and 2.13 together imply that conditions (a) and (b) of Theorem 2.10 are met for an ARA 1 model under Type-II censoring with probability one and a positive probability, respectively. It is also clear that Assumption 2.7 is not fulfilled for Type-I censoring. The results for Type-II censoring are given in Fig. 1. Clearly, the ten realizations of the function P,100 are decreasing in θ as it was proved in Example 2.13 by verifying condition (a) of Theorem 2.10 for all pairs (θ, θ) with θ < θ. We also see that all realizations of P,100 in Fig. 1 seem to be strictly decreasing as a function of θ so that the probabilities with which condition (b) of Theorem 2.10 hold might be quite large for the above λ d s. The picture is slightly different for the combination ARA 1 and Type-I censoring. Of course assumptions (a) and (b) of Theorem 2.10 on the statistical model are not affected by considering Type-I censored data instead of Type-II censored data, but assumption (Assumption 2.7) on the sampling procedure is not met. Comparing the left-and right-hand side of Fig. 2 we see that on the left-hand side not all realizations of P,100 are monotonically decreasing whereas on the right-hand side all realizations seem to lead to a monotonically decreasing P,100 . This is probably a result of the fact that due to the increased observation period (τ = 2.2 on the left-hand side and τ = 5.9 on the right-hand side) the part of P,100 stemming from Type-I censored observations (exactly these observations seem to prevent P,100 from being monotonically decreasing in case of Type-I censoring) gets outweighed by the number of observed failure times. Now we briefly discuss the simulation results for the ARA ∞ model and the two types of censoring. In Example 2.14 we did not prove that conditions (a) and (b) of Theorem 2.10 hold for all pairs (θ, θ) with θ < θ. The simulation results for an ARA ∞ model under Type-II censoring (Assumption 2.7 is then fulfilled) suggest that condition (a) (and maybe even (b)) of Theorem 2.10 may also hold with probability one for all pairs (θ, θ) with θ < θ as for the ARA 1 model. However, this is not the case as can be seen from Fig. 5 where we plotted the difference between the right-and left-hand side in the first displayed equation of Example 2.14 for i = 9, s 1 = 5.0, s 2 = 7.1, s 3 = 12.2, s 4 = 16.3, s 5 = 17.0, s 6 = 20.5, s 7 = 22.5, s 8 = 27, and k = 9, s1 = 3.4, s2 = 7.9, s3 = 10, s4 = 14.0, s5 = 19.6, s6 = 22.6, s7 = 23.3, s8 = 26.0, s9 = 27.1 as a function of θ. However, the realizations shown in Fig. 3 suggest that the probabilities of the events in condition (a) and (b) of Theorem 2.10, respectively, are relatively large. In case of an ARA ∞ model and Type-I censoring Assumption 2.7 is not met and condition (a) of Theorem 2.10 does not hold for all pairs (θ, θ) with θ < θ. Nevertheless, the simulation results shown in Fig. 4 suggest that the profile likelihood estimator remains inconsistent. Here θ denotes the true parameter. 20 q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q 0.0 0.2 0.4 0.6 0.8 1.0

0.0 0.2 0.4 0.6 0.8 1.0 θ g -f Figure 5: Difference between g(θ) := sk -θ k-1 =1 (1 -θ) k-1-s and f (θ) := s i-1 - θ i-1 =1 (1 -θ) i-1-s
, where i = 9, s 1 = 5.0, s 2 = 7.1, s 3 = 12.2, s 4 = 16.3, s 5 = 17.0, s 6 = 20.5, s 7 = 22.5, s 8 = 27, k = 9 and s1 = 3.4, s2 = 7.9, s3 = 10, s4 = 14.0, s5 = 19.6, s6 = 22.6, s7 = 23.3, s8 = 26.0, s9 = 27.1, see Example 2.14.

We finish this section by the comparison of the "discrete" and the "continuous log profile likelihood" as explained in Remark 2.5. For the above given right-truncated Weibull distribution and various values of (θ 0 , d, r) and a sample size of m = 100 we simulated five samples from an ARA 1 model. For each simulated sample the two versions (continuous and discrete) of the profile likelihood function were calculated for θ ∈ [0, 1]. The continuous profile likelihood curves are plotted with solid black lines whereas the discrete profile likelihood curves are plotted with dotted red lines. The results are given in Fig. 6. The obtained results are consistent with Theorem 2.10 where we use the continuous version of the profile log-likelihood function (solid lines). In addition we observe that using a version of the profile log-likelihood accounting the fact that the NPMLE of the baseline cumulative hazard is discrete (Appendix B) still leads to monotone profile log-likelihood functions. Finally notice that we see in Fig. 6 that the differences of the realizations of the discrete and continuous version are for all θ ∈ [0, 1] close to log exp (-m i=1 N i (s)) = -m•r. This term is just a product integral (the second double product in (44)). That the differences of the realizations are almost constant as a function of θ is due to the fact that m is large compared to r. For m small and r large the differences are far from being constant as a function of θ.

Conclusion

We have seen in Example 2.13 a model for recurrent events for which the profile likelihood approach with a right-continuous step function as an estimator for Λ fails in all respects as the profile likelihood function in that example is monotonically decreasing with probability one regardless of the sample size and the true underlying probability measure. As mentioned in the introduction a similar behaviour may occur in a single event model with covariates. It is also clear from the two sentences preceding Remark 2.5 that this behaviour cannot only be attributed to the fact that the argument of the exponential in the profile likelihood function is free of the virtual age parameter θ. This is simply because the same goes for the profile likelihood function in Cox's proportional hazards model for which the profile likelihood approach works. Moreover, Appendix C rules out the possibility that the failure of the profile likelihood approach results from an identifiability issue. To conclude: It seems to be a fine line that divides semi-parametric models for which the profile likelihood approach works from those semi-parametric models for which it fails. An exact description, i.e. an 'if and only if' statement, might be difficult or even impossible to obtain. Nevertheless, a few general features are worth summarizing. When the NPMLE Λ(•; θ) is a right-continuous step function the profile likelihood method fails for the semi-parametric accelerated failure time model (with a single event). This is because for this model θ influences only the locations of the jumps of Λ(•; θ) but not the jump heights. This is in contrast to Cox's proportional hazards model where the jump heights of Λ(•; θ) do depend on θ whereas the locations of the jumps are the same for all θ. Here we have proved a result that gives conditions under which profile likelihood inference does not work for the effective age parameter in a semi-parametric recurrent event model. In the examples presented the jump heights of the NPMLE Λ(•; θ) depend on θ. Nevertheless, the profile likelihood method fails. This might be a result of the fact that θ does not only affect the jump heights, but also the locations. One may wonder whether the locations being independent of θ (together with some regularity conditions) is a sufficient condition for the profile likelihood method to work.

A. Proofs and auxiliary lemmas

Before we give the proof of Theorem 2.10 we state two lemmas. The first lemma will be used in the proof of part (a) and part (b) whereas the second lemma will only be used in the proof of part (b).

Lemma A.1 Let I 1 and I 2 be two subsets of a finite set J with

I 1 = J and |I 1 | = |I 2 |,
where for any set K we denote by |K| its cardinality. Moreover, assume that for at least one element i 1 of I 1 we have i 1 / ∈ I 2 . Then there is at least one element i 2 ∈ J such that i 2 ∈ I 2 , but i 2 / ∈ I 1 .

Proof The claim is immediate from the facts that I 1 J and I 2 J have the same cardinality and that i

1 ∈ I 1 , but i 1 / ∈ I 2 . 2
The following lemma is obvious and its proof is therefore omitted.

Lemma A.2 Let x i , y i , 1 ≤ i ≤ I, and xi , ỹi , 1 ≤ i ≤ I, be (non-negative) real numbers.

(a) Define the functions G and G both with domain {1, . . . , I} by

G(j) := I i=1 1 (x i ,∞) (y j )
and

G(j) := I i=1 1 (x i ,∞) (ỹ j ).
Then:

(i) If for a given j ∈ {1, . . . , I} we have that for all i ∈ {1, . . . , I} the relation x i < y j implies the following relation xi < ỹj , then G(j) ≤ G(j).

(ii) If additionally to the assumption in part (i) we have that there is an i ∈ {1, . . . , I} such that xi < ỹj but x i ≥ y j then G(j) < G(j).

(b) Denote by y (i) , 1 ≤ i ≤ I, and ỹ(i) , 1 ≤ i ≤ I, the increasingly ordered values of the y i , 1 ≤ i ≤ I, and of the ỹi , 1 ≤ i ≤ I, respectively. Define the functions G ord and Gord both with domain {1, . . . , I} by

G ord (j) := I i=1 1 (x i ,∞) (y (j) )
and

Gord (j) := I i=1 1 (x i ,∞) (ỹ (j) ).
Then we have for j ≤ k G ord (j) ≤ G ord (k) and Gord (j) ≤ Gord (k).

Proof of Theorem 2.10 Throughout the proof whenever appropriate we suppress the dependence on ω, otherwise it is made explicit. We first prove part (a). Notice first that under Assumption 2.7 we have for every θ that the log profile likelihood

P,m (s * | θ, Λ m , D m (s)) equals - m k=1 J k (s * ) =1 log S θ m (s * , ε θ k (S k )) = - m k=1 J k (s * ) =1 log S θ m (s * , ε θ k, -1 (S k, )) . (12)
Now, let ω be arbitrary. Put U (ω) := m k=1 J k (s * , ω). Denoting the increasingly ordered values of

ε θ k, -1 (S k, (ω)), 1 ≤ k ≤ m, 1 ≤ ≤ J k (s * , ω), by ε θ k, -1 (S k, (ω)) (p)
, p = 1, . . . , U (ω), the right-hand side of Equation ( 12) can be rewritten as

- U (ω) p=1 log S θ m s * , (ε θ k, -1 (S k, (ω))) (p) . (13)
From Equation ( 4) we now obtain that under Assumption 2.7 the quantity in (13) equals

- U (ω) p=1 log m i=1 J i (s * ,ω) j=1 γ θ i,j-1 ((ε θ k, -1 (S k, (ω))) (p) ) × 1 (ε θ i,j-1 (S i,j-1 (ω)),ε θ i,j-1 (S i,j (ω))] ((ε θ k, -1 (S k, (ω))) (p) ) . (14)
We see that under the assumption (P

λ 0 ,θ 0 ) m (γ θ i,j-1 (t) ≤ c, γ θ i,j-1 (t) ≥ c, 1 ≤ i ≤ m, 1 ≤ j ≤ J i (s * ), 0 ≤ t < ∞) = 1 our claim (P λ 0 ,θ 0 ) m P,m s * |θ, Λ m , D m (s) ≥ P,m s * | θ, Λ m , D m (s) ≥ (P λ 0 ,θ 0 ) m A m,θ, θ would follow if - U (ω) p=1 log m i=1 J i (s * ,ω) j=1 1 (ε θ i,j-1 (S i,j-1 (ω)),ε θ i,j-1 (S i,j (ω))] ((ε θ k, -1 (S k, (ω))) (p) ) ≥ - U (ω) p=1 log m i=1 J i (s * ,ω) j=1 1 (ε θ i,j-1 (S i,j-1 (ω)),ε θ i,j-1 (S i,j (ω))] ((ε θ k, -1 (S k, (ω))) (p) ) , (15)
were true for (almost) all ω ∈ A m,θ, θ, because we have

- U (ω) p=1 log m i=1 J i (s * ,ω) j=1 γ θ i,j-1 ((ε θ k, -1 (S k, (ω))) (p) ) × 1 (ε θ i,j-1 (S i,j-1 (ω)),ε θ i,j-1 (S i,j (ω))] ((ε θ k, -1 (S k, (ω))) (p) ) ≥ - U (ω) p=1 log c m i=1 J i (s * ,ω) j=1 1 (ε θ i,j-1 (S i,j-1 (ω)),ε θ i,j-1 (S i,j (ω))] ((ε θ k, -1 (S k, (ω))) (p) ) (16)
and

- U (ω) p=1 log c m i=1 J i (s * ,ω) j=1 1 (ε θ i,j-1 (S i,j-1 (ω)),ε θ i,j-1 (S i,j (ω))] ((ε θ k, -1 (S k, (ω))) (p) ) ≥ - U (ω) p=1 log m i=1 J i (s * ,ω) j=1 γ θ i,j-1 ((ε θ k, -1 (S k, (ω))) (p) ) × 1 (ε θ i,j-1 (S i,j-1 (ω)),ε θ i,j-1 (S i,j (ω))] ((ε θ k, -1 (S k, (ω))) (p) ) . (17)
Now notice that for arbitrary θ we can rewrite

- U (ω) p=1 log m i=1 J i (s * ,ω) j=1 1 (ε θ i,j-1 (S i,j-1 (ω)),ε θ i,j-1 (S i,j (ω))] ((ε θ k, -1 (S k, (ω))) (p) ) as - U (ω) p=1   log   m i=1 J i (s * ,ω) j=1 (ε θ k, -1 (S k, (ω))) (p) -(p -1)     , (18)
because the condition

(ε θ k, -1 (S k, (ω))) (p) ≤ ε θ i,j-1 (S i,j (ω))
fails for exactly (p -1) pairs (i, j), 1 ≤ i ≤ m, 1 ≤ j ≤ J i (s * ). As the representation (18) holds for every element of Θ, we see (15) would follow if for every p and (almost) every ω ∈ A m,θ, θ it were true that for every pair (i, j)

ε θ i,j-1 (S i,j-1 (ω)) < ε θ g θ (p) ,h θ (p) -1 (S g θ (p) ,h θ (p) (ω)) ⇒ ε θ i,j-1 (S i,j-1 (ω)) < ε θ k θ (p) , θ (p) -1 (S k θ (p) , θ (p)
(ω)).

(

) 19
Here ε θ

g θ (p) ,h θ (p) -1 (S g θ (p) ,h θ (p)
(ω)) denotes the p smallest value among the ε θ g,h-1 (S g,h (ω)),

1 ≤ g ≤ m, 1 ≤ h ≤ J g (s * ), and ε θ k θ (p) , θ (p) -1 (S k θ (p) , θ (p)
(ω)) denotes the p smallest value

among the ε θ k, -1 (S k, (ω)), 1 ≤ k ≤ m, 1 ≤ ≤ J k (s * ). Consequently, (g θ (p) , h θ (p)
) denotes the index pair of the p smallest observation among the ε θ g,h-1 (S g,h (ω)), 1 ≤ g ≤ m, 1 ≤ h ≤ J g (s * ) and similar for the pair (k θ (p) , θ (p) ). Before continuing we need some more notation. For two pairs (i, j) and (k, ) of natural numbers (i, j) = (k, ) means i = k and j = . Moreover, for any pair (i, j), 1 ≤ i ≤ m, 1 ≤ j ≤ J i (s * ), and θ ∈ {θ, θ} we denote by rk θ(i, j) the rank of ε θ i,j-1 (S i,j (ω)) among the

ε θ k, -1 (S k, (ω)), 1 ≤ k ≤ m, 1 ≤ ≤ J k (s * ), i.e. rk θ(i, j) is equal to 1 if ε θ i,j-1 (S i,j (ω)) is the smallest among the ε θ k, -1 (S k, (ω)), 1 ≤ k ≤ m, 1 ≤ ≤ J k (s * ).
For the rest of the proof of part (a) let ω ∈ A m,,θ, θ.

Case 1: Suppose that g θ (p) = k θ (p) and h θ (p) = θ (p) . Then we immediately see that ( 19) is just the assumption stated in part (a).

Case 2: Suppose that (g θ (p) , h θ (p) ) = (k θ (p) , θ (p) ) and that we additionally have rk θ (k θ (p) , θ (p) ) > p. Then (19) holds, because we have

ε θ k θ (p) , θ (p) -1 (S k θ (p) , θ (p) (ω)) > ε θ g θ (p) ,h θ (p) -1 (S g θ (p) ,h θ (p) (ω))
due to the fact that rk θ (k θ (p) , θ (p) ) > p and the implication

ε θ i,j-1 (S i,j-1 (ω)) < ε θ k θ (p) , θ (p) -1 (S k θ (p) , θ (p) (ω)) ⇒ ε θ i,j-1 (S i,j-1 (ω)) < ε θ k θ (p) , θ (p) -1 (S k θ (p) , θ (p)
(ω))

is just the assumption made in part (a).

Case 3:

Let (g θ (p) , h θ (p) ) = (k θ (p) , θ (p)
) and rk θ (k θ (p) , θ (p) ) < p and assume additionally that rk θ(g θ (p) , h θ (p) ) < p. Then ( 19) is true, because we have

ε θ g θ (p) ,h θ (p) -1 (S g θ (p) ,h θ (p) (ω)) < ε θ k θ (p) , θ (p) -1 (S k θ (p) , θ (p)
(ω))

due to the fact that rk θ(g θ (p) , h θ (p) ) < p and the implication

ε θ i,j-1 (S i,j-1 (ω)) < ε θ g θ (p) ,h θ (p) -1 (S g θ (p) ,h θ (p) (ω)) ⇒ ε θ i,j-1 (S i,j-1 (ω)) < ε θ g θ (p) ,h θ (p) -1 (S g θ (p) ,h θ (p) (ω))
is again just the assumption stated in part (a).

Case 4: Let (g ) > p together imply that there is at least one pair (e, f ), 1 ≤ e ≤ m, 1 ≤ f ≤ J e (s * ) such that rk θ (e, f ) > p and rk θ(e, f ) < p. Indeed, let

I θ := {(v, w), 1 ≤ v ≤ m, 1 ≤ w ≤ J v (s * )|rk θ (v, w) < rk θ (g, h) = p} and I θ := {(v, w), 1 ≤ v ≤ m, 1 ≤ w ≤ J v (s * )|rk θ(v, w) < rk θ (g, h) = p}.
Because This proves rk θ (e, f ) > p and therefore there is indeed at least one pair with the stated properties. Now we get from the assumption in part (a)

ε θ ij-1 (S i,j-1 (ω)) < ε θ e,f -1 (S e,f (ω)) ⇒ ε θ i,j-1 (S i,j-1 (ω)) < ε θ e,f -1 (S e,f (ω)).
This clearly implies (19), because

ε θ e,f -1 (S e,f (ω)) > ε θ g θ (p) ,h θ (p) -1 (S g θ (p) ,h θ (p) (ω)) and ε θ e,f -1 (S e,f (ω)) < ε θ k θ (p) , θ (p) -1 (S k θ (p) , θ (p)
(ω)).

This finishes the proof of part (a).

We now begin with the proof of the statement in part (b). Consider an arbitrary ω ∈ B m,θ, θ. To shorten the notation we introduce the functions

F θ ω , θ ∈ {θ, θ}, defined on {(g, h)|1 ≤ g ≤ m, 1 ≤ h ≤ J g (s * , ω)} by F θ ω (g, h) := m i=1 J i (s * ,ω) j=1 1 (ε θ i,j-1 (S i,j-1 (ω)),∞) (ε θ g,h-1 (S g,h (ω))), θ ∈ {θ, θ}.
Before continuing with the actual proof, we state the following facts about the functions F θ ω and F θ ω . Firstly, under the assumption made in part (a) of Theorem 2.10 we obtain from Lemma A.2 part (a), (i) that for every pair (g, h)

F θ ω (g, h) ≤ F θ ω (g, h). (20)
Secondly, the assumptions in part (a) and part (b) of Theorem 2.10 together imply in view of Lemma A.2 (a), (ii) that

F θ ω (k, ) < F θ ω (k, ). (21)
Thirdly, notice that part (b), (i) of Lemma A.2 implies for θ ∈ {θ, θ}

F θ ω (g, h) ≤ F θ ω (k, ), if rk θ(g, h) ≤ rk θ(k, ). (22)
Fourthly, in the proof of part (a) of Theorem 2.10 we showed ( 19) so that we have

F θ ω (g, h) ≤ F θ ω (k, ), if rk θ (g, h) = rk θ(k, ). (23)
We now start considering all possible cases. In each case considered it is sufficient to show that the inequality in ( 23) is strict for two pairs (g, h) and (k, ). This follows by

Case 1:

Case 3b:

Case 2:

Case 3a: combining ( 16) and ( 17). Of course, the two pairs might be equal. To ease the reading of the proof each case is illustrated by a graphic; see Fig. 7.

Case 1: rk θ(k, ) < rk θ (k, ). Then let the pair (u, v) be such that rk θ(u, v) = rk θ (k, ) and notice that we consequently have rk θ(u, v) > rk θ(k, ). We now obtain

F θ ω (u, v) (22) ≥ F θ ω (k, ) (21)
> F θ ω (k, ).

Hence, the inequality in ( 23) is strict for the pairs (u, v) and (k, ) that fulfil p = rk θ(u, v) = rk θ (k, ).

Case 2: rk θ(k, ) = rk θ (k, ). Then we immediately see that the inequality in ( 23) is strict for p = rk θ(k, ) = rk θ (k, ) by ( 21).

Case 3: rk θ(k, ) > rk θ (k, ). Before continuing with Case 3 we remark that we will assume that

F θ ω (g, h) = F θ ω (k, ), if rk θ (g, h) = rk θ(k, ). (24)
The reason why we can assume Equation ( 24) is that if (24) did not hold, we would have found two pairs for which the inequality in ( 23) would be strict. However, this would finish the proof. To make the remainder of the proof easily accessible we state the following lemma.

Lemma A.3 Under the assumption stated in (24) we have that

rk θ (i, j) ≥ rk θ (a, b) = rk θ(α, β) ≥ rk θ(i, j) implies F θ ω (i, j) = F θ ω (a, b) = F θ ω (α, β) = F θ ω (i, j).
Proof of Lemma A.3

F θ ω (i, j) (22) ≥ F θ ω (a, b) (24) = F θ ω (α, β) (22) ≥ F θ ω (i, j) (20) ≥ F θ ω (i, j).
Hence, we can replace all ≥ by = which finishes the proof. 2

Proof of Theorem 2.10 continued. Let the pairs (u, v) and (c, d) be such that rk θ (u, v) = rk θ(k, ) and rk θ (k, ) = rk θ(c, d).

Here R d + := {x ∈ R d |x i ≥ 0, i = 1, . . . , d}. For every d ∈ N we now define a probability measure Q d on R d by defining for j = 1, . . . , d transition kernels Q j s j |s 1 , . . . , s j-1 , ε θ j-1 := P d {ε θ j-1;s 1 ,...,s j-1 (s j )} P d (ε θ j-1;s 1 ,...,s j-1 (s j-1 ), ∞) , = λ d ε θ j-1;s 1 ,...,s j-1 (s j )

P d [ε θ j-1;s 1 ,...,s j-1 (s j ), ∞) P d (ε θ j-1;s 1 ,...,s j-1 (s j-1 ), ∞)

= λ d ε θ j-1;s 1 ,...,s j-1 (s j )

S d ε θ j-1;s 1 ,...,s j-1 (s j )-S d ε θ j-1;s 1 ,...,s j-1 (s j-1 ) for s j > s j-1 ,

where we use the convention that 0/0 := 0, (ε θ j-1;s 1 ,...,s j-1 ) -1 denotes the inverse of ε θ j-1;s 1 ,...,s j-1 , s 0 = 0, and λ d and S d are the hazard rate and the survival function corresponding to P d , respectively. In Equation ( 41) we have implicitly assumed that P d and the functions ε θ j-1;s 1 ,...,s j-1 match properly, i.e. P d (ε θ j-1;s 1 ,...,s j-1 (s j-1 ), ∞) > 0, j = 1, 2, . . .. Remark B.1 Constructing effective age models through Equation ( 41) is similar to the approaches in [START_REF] Dorado | Nonparametric estimation for a general repair model[END_REF] and [START_REF] Last | Asymptotic and monotonicity properties of some repairable systems[END_REF].

Example B.2 For the model discussed in Example 2.2 we have for θ ∈ [0, 1] that ε θ 0 (s) := s, s ∈ R + . Moreover, we have for θ ∈ [0, 1] that ε θ d-1 : R j + → R + , j ≥ 2, with ε θ j-1 (s) := s j -θs j-1 , s = (s 1 , . . . , s j-1 , s j ) ∈ R j + . Furthermore, for the model discussed in Example 2.3 we have for θ ∈ [0, 1] that ε θ 0 (s) := s, s ∈ R + . Moreover, we have for θ ∈ [0, 1] that ε θ j-1 : R j + → R + , j ≥ 2, with ε θ j-1 (s) := s j -θ j-1 =1 (1 -θ) j-1-s , s = (s 1 , . . . , s j-1 , s j ) ∈ R d + .

Now, let [0, s], D m (s), N i (s) and τ i be as in Section 2 and denote the event occurrence times by sij , 1 ≤ i ≤ m, 1 ≤ j ≤ N i (s). Then the full likelihood corresponding to the model defined by ( 41) is given by with d = max{N 1 (s), . . . , N m (s)} L F (s|P d , (ε θ i,j-1 ) 1≤i≤m,j∈N , D m (s)) = m i=1 N i (s) j=1 λ d ε θ i,j-1;s i,1 ,...,s i,j-1 (s i,j )

× m i=1 N i (s) j=1
S d ε θ i,j-1;s i,1 ,...,s i,j-1 (s i,j )-S d ε θ i,j-1;s i,1 ,...,s i,j-1 (s λ d ε θ i,j-1;s i,1 ,...,s i,j-1 (s i,j )

×    m i=1 N i (s) j=1 k∈I θ i,j (1 -λ d (z k ))    ×    m i=1 k∈I θ,τ i i (1 -λ d (z k ))    , (42)
where for every pair (i, j), 1 ≤ i ≤ m, 1 ≤ j ≤ N i (s), the sets I θ i,j are defined by I θ i,j := {k ∈ {1, . . . , K}|ε θ i,j-1;s i,1 ,...,s i,j-1 (s i,j-1 ) < z k < ε θ i,j-1;s i,1 ,...,s i,j-1 (s i,j )} and for every i, 1 ≤ i ≤ m, the sets I θ,τ i i are defined by I θ,τ i i := {k ∈ {1, . . . , K}|ε θ i,N i (s);s i,1 ,...,s i,N i (s) (s i,N i (s) ) < z k ≤ ε θ i,N i (s);s i,1 ,...,s i,N i (s) (τ i ∧ s)}.

To see that Equation ( 43 1 -λ d ε θ i,j-1;s i,1 ,...,s i,j-1 (u)

   ×    m i=1 u∈(I θ,τ i i ) -1
1 -λ d ε θ i,N i (s);s i,1 ,...,s i,N i (s) (u)

   , (43)
where for every pair (i, j), 1 ≤ i ≤ m, 1 ≤ j ≤ N i (s), the intervals (I θ i,j ) -1 are defined by (I θ i,j ) -1 := {u ∈ R + |s i,j-1 < u < si,j } and for every i, 1 ≤ i ≤ m, the intervals (I θ,τ i i ) -1 are defined by (I θ,τ i i ) -1 := {u ∈ R + |s i,N i (s) < u ≤ τ i ∧ s}, and the product w.r.t. u stands for the product integral. Now we are going to maximise (42) for θ fixed. As in Section 2.2, maximisation is w.r.t. all discrete probability measures that put (positive) mass at the points ε θ k, -1;s k,1 ,...,s k, -1 (s k, ), 1 ≤ k ≤ m, 1 ≤ ≤ N k (s), which we assume to be different. Also as in Section 2.2 we denote the corresponding hazard rates at these points by λ θ k, . Then the full likelihood (42) after rearranging can be written as

m k=1 N k (s) =1 λ θ k, m k=1 N k (s) =1 (1 -λ θ k, ) |I θ,B k, |+|I θ,τ k ,B k | , (44)
where for every pair (k, ), 1 ≤ k ≤ m, 1 ≤ ≤ N k (s), the sets I θ,B k, are defined (cf. Section 2.2) as

I θ,B
k, := {(i, j), 1 ≤ i ≤ m, 1 ≤ j ≤ N i (s)| ε θ i,j-1;s i,1 ,...,s i,j-1 (s i,j-1 ) < ε θ k, -1;s k,1 ,...,s k, -1 (s k, ) < ε θ i,j-1;s i,1 ,...,s i,j-1 (s i,j )} and for every k, 1 ≤ k ≤ m, the sets I θ,τ i ,B k are defined (cf. Section 2.2) by I θ,τ i ,B k := {i, 1 ≤ i ≤ m| ε θ i,N i (s);s i,1 ,...,s i,N i (s) (s i,N i (s) ) < ε θ k, -1;s k,1 ,...,s k, -1 (s k, ) ≤ ε θ i,N i (s);s i,1 ,...,s i,N i (s) (τ i ∧ s)}.

Moreover, as in Appendix A for any set I we denote by |I| its cardinality. As the function x → x(1-x) k , k ∈ N 0 with x ∈ [0, 1], is maximised at x = 1/(k +1) we see that Equation ( 44) is maximised at

λ θ k, = 1

Figure 1 :Figure 2 :Figure 3 :Figure 4 :

 1234 Figure 1: Each of the 6 figures contains ten graphs of θ → m,P (s|θ, Λ m , D m (s)) for θ ∈ [0, 1] and m = 100 obtained from simulated data for which λ = λ d , the effective age follows the ARA 1 assumption of Example 2.2, and for τ i = S i,r for 1 ≤ i ≤ m (Type II censoring, see Example 2.8). Here θ denotes the true parameter.

Figure 6 :

 6 Figure 6: Each of the 6 figures contains five graphs of θ → m,P (s|θ, Λ m , D m (s)) (solid lines) for θ ∈ [0, 1] and m = 100 obtained from simulated data for which λ = λ d , the effective age follows the ARA 1 assumption of Example 2.2, and for τ i = S i,r for 1 ≤ i ≤ m (Type II censoring, see Example 2.8). The dotted red lines are the realizations of the discrete profile likelihood curves. Here θ denotes the true parameter.

Figure 7 :

 7 Figure 7: Illustration of the different cases considered in the proof of part (b) of Theorem 2.10. The ranks of the pairs under θ are given on the upper axis in ascending order and under θ on the lower axis. For Case 3b the three graphics on the right-hand side illustrate the process of extracting n with the property as given in (32).

  ) is indeed the discrete time analogue of Equation (2

.

  Notice that this NPMLE coincides with the one given in Section 2.2 if the derivatives of the effective functions w.r.t. observational time are equal to one, because we have the relations|I θ,B k, | + 1 = |I θ k, | and |I θ,τ k ,B k | = |I θ,τ k k |.

  1 will be given below; see Examples 2.2, 2.3, and 2.4. Notice

	that we do not require
	that ε θ (s, ω) ≤ s. This means that the effective age might be larger than the current
	time s and that we do not exclude harmful interventions, i.e. interventions that increase
	the effective age instead of reducing it. To clarify the meaning of the effective age process
	let us consider two well-known processes:
	1. Renewal process: Replacing the observational unit by a new one results in a renewal
	process with effective age process that equals s -S k-1 (ω) on (S k-1 (ω), S k (ω)] at time s;
	2. Non-homogeneous Poisson process: Somehow on the other side of the spectrum is
	the non-homogeneous Poisson process with intensity function λ, because its effective age
	process at time s equals just s.

that, as usual, N (s) denotes the number of events over the period (0, s] for an observable unit, and the time of the ith recurrent event is denoted by S i . The interpretation of the predictable process Y is as follows: the unit is still under observation, i.e. at risk, if and only if Y (s) = 1. In the following we assume that Y is of the form Y = {Y (s), 0 ≤ s ≤ s * } with Y (s) = 1 {τ ≥s} , where τ is some positive (random) variable so that [0, τ ] is the (random) observation interval. Here and in the following 1 denotes the indicator function. Further, we refer to ε θ , θ ∈ Θ, as the effective age process and to θ as the effective age parameter. It describes the effect of interventions applied to the observational unit after experiencing a recurrent event. Notice

  is more general than Model 2.1, because Model 2.1 takes the functions ρ and ψ to be identically equal to one. On the other hand Model 2.1 is more general than the one introduced by Peña and Hollander, because it allows for a class of predictable processes ε θ , θ ∈ Θ, whereas Peña and Hollander and the above references assume that the process ε = {ε(s), 0 ≤ s ≤ s * } is known.

	Clearly, Model 2.1 could be
	extended to contain the model defined by Equation (1). However, as our focus here is
	on the extent to which profile likelihood inference for effective age models is possible we
	restrict ourselves to Model 2.1. At this point, it is worth mentioning that the model
	given by Equation (1) contains, for instance, Cox's proportional hazards model (with

  i,j-1 ) (s);s i,1 ,...,s i,N i (s) (τ i ∧ s) S d ε θ i,N i (s);s i,1 ,...,s i,N i (s) (s i,N i (s) )

	× i,N (40) m i=1 S d ε θ = m N i (s)
	i=1	j=1

i

* This work was partially supported by the French National Agency of Research (ANR), project AMMSI, number ANR-2011-BS01-021.

From (24) we then get

because rk θ (u, v) = rk θ(k, ). Define now the sets

1 and (k, ) / ∈ I θ 1 , it follows from Lemma A.1 that there must be a pair (c 1 , d 1 ) such that (c 1 , d 1 ) ∈ I θ 1 and (c 1 , d 1 ) / ∈ I θ 1 , i.e.

rk θ(c 1 , d 1 ) ≤ rk θ (k, )

= rk θ(c, d) and rk θ (c 1 , d 1 ) > rk θ (k, ). ( 27)

Now make the following assumption that we call Case 3a:

Under this assumption combining the left-hand side of ( 25) with the left-hand side of ( 27) we obtain from Lemma A.3

which implies that

Combining the left and right-hand side of (27) we have

and we obtain from Lemma A.3

However, in view of (29) the first equality in Equation (30) gives F θ ω (k, ) = F θ ω (k, ) which this is a contradiction to (21). Hence the assumption made in Equation ( 26) does not hold and we must have

which finishes the proof if we are in Case 3a, i.e. the assumption stated in (28) holds. Now if (28) does not hold, we consider Case 3b:

In this case continue extracting pairs (c n , d n ), n ≥ 2, with

Clearly, the pairs (c n , d n ) and the sets as given in (31) exist by Lemma A.1 at least as long as rk

is strictly increasing for every ω (with the property that we have to consider Case 3 and that we are not in Case 3a, i.e. the assumption stated in (28) does not hold) there is an n ≥ 2 such that

Notice that Equation (32) together with (31) implies

Let the pair (e n-1 , f n-1 ) be such that

Then

Applying Lemma A.3 to (35) we now obtain for every n, 2 ≤ n ≤ n

Furthermore, combining (32) and (33) we have

Applying once again Lemma A.3 we obtain

In particular,

and we know from Equation (30) which is valid here as well as it was derived from Equation ( 27)

Together with (36) Equations ( 38) and ( 39) now imply that

However, this is a contradiction to (21) showing that also in Case 3b the assumed equality in (24) cannot hold. This finishes the proof of part (b). 2

B. Discrete hazard rates

Recall that for a probability measure P on R with a finite support on {z 1 , . . . , z K }, say, the (discrete) hazard rate λ and the survival function S are given by the formulas

, and P((t, ∞)) =:

Recall also that the relation

Now we shall define an effective age model corresponding to a sequence of discrete probability measures P d with supports {z 1 , . . . , z K 1 , z K 1 +1 , . . . , z K d }, d ∈ N (the method could also be used for a continuous probability measure). This model has exactly the same transition kernels as Model 2.1. To do so, assume that for every θ ∈ Θ we have given a sequence of functions ε θ 0 , ε θ 1 , . . ., θ ∈ Θ, that are interpreted as the effective age functions for the time interval between the jth and the (j -1)th event, with the following properties: (c) The functions ε θ j-1;s 1 ,...,s j-1 (•) := ε θ j-1 (s 1 , . . . , s j-1 , •) are strictly increasing for every fixed s 1 , . . . , s j-1 .

C. Identifiability in Examples 2.2 and 2.3

First note that the distribution of (S 1 , S 2 ) (or (T 1 , T 2 )) is the same under Kijima I or Kijima II models; see Examples 2.2 and 2.3. Indeed, the joint distribution of (T 1 , T 1 ) is defined by its density function

Here f is an unknown probability density function with support [0, ∞) (S is the corresponding survival function) and θ is an unknown Euclidean parameter in [0, 1]. Let µ k be the Lebesgue measure on R k (k ≥ 1). Proving the identifiability of (θ, f ) requires to show the following one-to-one property: g

T 1 ,T 2 (t 1 , t 2 ) with respect to t 2 on R. Thus identifiability reduces to proving that θ = θ results from

Integrating ( 45) with respect to t 2 on [s, ∞) leads to

Using (45) (with t 2 = s) and ( 46) we obtain

The latter holds if and only if

Then for µ 2 almost all (t 1 , s) ∈ {(x, y) ∈ R 2 + |x < (b -a)/| θ -θ|, a -(1 -θ ∨ θ)t 1 < y < b -(1 -θ ∧ θ)t 1 } we must have s + (1 -θ)t 1 = s + (1 -θ)t 1 . Hence θ = θ which proves the semi-parametric identifiability of Kijima I and Kijima II models whenever the first two failures can be observed.