Guillaume Allibert

Estelle Courtial

What can Prediction Bring to Image-Based Visual Servoing ?

The purpose of this paper is to show what image prediction can bring to Image-Based Visual Servoing. The visual feature prediction is obtained thanks to the interaction matrix. Based on a Model Predictive Control strategy, the visual servoing task is formulated into an optimization problem. The error between the reference features and the predicted features is to be minimized over a receding prediction horizon. Numerous simulations highlight the interest of prediction, especially for difficult configurations such as large motion and rotation.

I. INTRODUCTION

Since the eighties, visual servoing has been largely developed from a theoretical viewpoint but also from a practical viewpoint. Many applications of visual servoing, among others things, for aerial, submarine or medical robots have been reported in the literature [START_REF] Bourquardez | Image-based visual servo control of the translation kinematics of a quadrotor aerial vehicle[END_REF], [START_REF] Dunbabin | Low-cost vision-based AUV guidance system for reef navigation[END_REF], [START_REF] Gangloff | Hight speed visual servoing of a 6 dof manipulator using MIMO predictive control[END_REF]. The fundamental classification of visual servoing distinguishes three approaches: image-based control (2D), position-based control (3D) and a hybrid approach (2 1/2 D). Further details about visual servoing can be found in [START_REF] Chaumette | Visual Servo Control, Part I: Basic Approaches[END_REF], [START_REF] Chaumette | Visual Servo Control, Part II: Advanced Approaches[END_REF]. Here, we focus our interest on Image-Based Visual Servoing (IBVS). In this context, the control task consists in determining the control input applied to the robotic system according to the error between the reference image and the current image from the camera. We consider a free-flying camera with six degreeof-freedom (6 dof). The relationship between the camera velocity screw τ and the time variation of the visual features ṡ is given by: ṡ(t) = L s τ(t)

where L s is the interaction matrix related to s. The aim of visual servoing is to minimize the error e(t) between the reference features s * , assumed to be constant, and the measured features s(t). It is defined by:

e(t) = s(t) -s * (2)
In order to satisfy an exponential decay of the error (2), mathematically expressed by: ė(t) = -λ e(t) with λ > 0 [START_REF] Danès | Towards an LMI approach to multicriteria visual servoing in robotics[END_REF] and considering the open-loop model [START_REF] Allibert | Visual Predictive Control for Manipulators with Catadioptric Camera[END_REF], we obtain the classical feedback control law written as:

τ(t) = -λ L s + e(t) (4)
where L s + is the pseudo-inverse of the approximated interaction matrix.

G. Allibert and E. Courtial are with the Institut PRISME EA 4229, Polytech'Orléans, 8 rue Léonard de Vinci, 45072 Orléans Cedex France, guillaume.allibert@univ-orleans.fr, estelle.courtial@univ-orleans.fr Another way to deal with IBVS task is to use advanced control laws such as optimal control [START_REF] Hashimoto | LQ optimal and nonlinear approaches to visual servoing[END_REF], Linear Matrix Inequalities [START_REF] Danès | Towards an LMI approach to multicriteria visual servoing in robotics[END_REF] and Model Predictive Control (MPC) [START_REF] Allibert | Visual Predictive Control for Manipulators with Catadioptric Camera[END_REF], [START_REF] Allibert | Real-time visual predictive controller for image-based trajectory tracking of mobile robot[END_REF], [START_REF] Gangloff | Hight speed visual servoing of a 6 dof manipulator using MIMO predictive control[END_REF], [START_REF] Sauvée | Image Based Visual Servoing through Nonlinear Model Predictive Control[END_REF]. In [START_REF] Gangloff | Hight speed visual servoing of a 6 dof manipulator using MIMO predictive control[END_REF], a predictive controller is used for motion compensation in target tracking applications. The prediction of the target motion is used to reject perturbation in order to cancel tracking errors. In [START_REF] Sauvée | Image Based Visual Servoing through Nonlinear Model Predictive Control[END_REF], the predictive controller is used from ultrasound images for a medical application. The strategy proposed in this paper exploits MPC for visual servoing tasks. The IBVS objective is formulated as solving on-line a nonlinear optimization problem expressed in the image plane [START_REF] Allibert | Visual Predictive Control for Manipulators with Catadioptric Camera[END_REF], [START_REF] Allibert | Real-time visual predictive controller for image-based trajectory tracking of mobile robot[END_REF]. A nonlinear global model combining the robotic model and the camera model allows to predict the evolution of the visual features over a finite receding horizon. The error between the reference features and the predicted model features is to be minimized with respect to the robotic control inputs. This strategy, named Visual Predictive Control (VPC), offers two advantages. First, VPC can easily take into account constraints such as visibility constraints or state/input constraints. The constraint handling, especially visibility constraints, can be very useful to deal with obstacle avoidance. Secondly, the visual prediction can play a crucial role for difficult configurations. The purpose of this paper is to show what image prediction brings to classical visual servoing. We propose an alternative approach based on the VPC strategy where the visual prediction is obtained no longer by the nonlinear global model but by the interaction matrix. In this case, no 3D data are required contrary to the global model [START_REF] Allibert | Visual Predictive Control for Manipulators with Catadioptric Camera[END_REF]. On the other hand, no 3D constraint would directly be taken into account. The interest of the prediction is pointed out through many simulations describing difficult configurations such as large motion and rotation on a free-flying perspective camera. The paper is organized as follows. In section II, the context of the study is stated and the principle of MPC is briefly recalled. Then, in section III, the proposed strategy and the interaction matrix is developed. In Section IV, simulation experiments highlight the efficiency of the proposed approach. The influence of the prediction horizon and the weighted matrix are addressed. Finally conclusions and future tasks are detailed in the last section.

II. VISUAL PREDICTION IN IBVS BASED ON OPTIMIZATION PROBLEM

A. Context of the study

Let us consider a point visual feature denoted s. For a 3D point with coordinates P = (X,Y, Z) T in the camera frame, which is projected in the image plane as a 2D point in normalized coordinates such that s = (u, v) T with u = X/Z and v = Y /Z, it has been shown that the interaction matrix, related to s, is given by:

L s 1 = -1 Z 0 u Z uv -(1 + u 2 ) v 0 -1 Z v Z 1 + v 2 -uv -u (5)
The value Z is the depth of the point relative to the camera frame. To avoid the estimation of this parameter at each iteration, the depth computed or measured at the reference position Z * is generally used in the control scheme (4). Consequently, the interaction matrix (5) depends only on the current time measure of the visual feature s.

To control a 6 dof free-flying camera, a minimum of three points is necessary. To a value of s corresponds a single camera pose w.r.t. the target if four points are considered.

The global interaction matrix is then obtained by stacking the L s i for i ∈ [1, ..., 4]:

L s = ⎛ ⎜ ⎜ ⎝ L s 1 L s 2 L s 3 L s 4 ⎞ ⎟ ⎟ ⎠ (6)
In this work, we propose to use (1), ( 5) and ( 6) to predict, at the current time, the evolution of the visual features over a finite prediction horizon.

B. Nonlinear Model Predictive Control

Nonlinear Model Predictive Control, also named Receding Horizon Control (RHC), has been largely developed for the control of constrained nonlinear processes [START_REF] Mayne | Constrained model predictive control: Stability and optimality[END_REF]. The control problem (setpoint or trajectory tracking) is formulated as solving on-line a nonlinear optimization problem. Based on the process model, the controller predicts the behavior of the system over a prediction horizon N p . The difference between the reference trajectory and the predicted model behavior defines the cost function J to be minimized with respect to a control sequence ũ. Due to disturbances and model mismatches, this procedure is repeated at each sampling instant. Only the first component of the optimal control sequence is really applied to the process. At the next sampling instant, when the measurements are updated, the finite horizon moves one step forward and the procedure starts again. The control structure considered is the well-known Internal Model Control (IMC) structure (see Fig. 1). According to this control scheme, we can write:

y ref (j) + _ _ + u(j) y m (j) Process Model (j) y d (j) y p (j) Controller Optimization Model + _ ym(j)
y d ( j) = y re f ( j) -ε( j) y d ( j) = y re f ( j) -(y p ( j) -y m ( j)) y d ( j) -y m ( j) = y re f ( j) -y p ( j) (7)
The tracking of the reference trajectory y re f by the process output y p is equivalent to the tracking of the desired trajectory y d by the model output y m . The cost function J is defined as a quadratic function of states and control inputs. Due to the IMC structure, J can be written in discrete-time as:

J(x, u) = k+N p ∑ j=k+1 [y d ( j) -y m ( j)] T Q [y d ( j) -y m ( j)] (8)
The mathematical formulation of NMPC strategy is then given by: min

u J(x, u) (9)
subject to:

x( j + 1) = f (x( j), u( j)) y m ( j) = h(x( j)) ( 10
)
x ∈ R n , u ∈ R m , y ∈ R p are respectively the states, the inputs and the outputs.

The first nonlinear equation ( 10) describes the dynamics of the system where x( j + 1) represents the predicted state at time j + 1 from the current time k. The predicted states are initialized with the system states at time k which guarantees the implicit feedback of the IMC structure when modeling errors and disturbances are equal to zero. The error ε( j), j ∈ [k + 1, k + N p ] is assumed to be constant over the prediction horizon and equal to the measurement ε(k). ũ = {u k , u k+1 , ..., u k+N c , ..., u k+N p -1 } is the optimal control sequence. From u(k + N c + 1) to u(k + N p -1), the control input is constant and equal to u(k + N c ) where N c is the control horizon. The weighted matrix Q is a symmetric definite positive matrix. One of the main advantages of NMPC is the capability to explicitly handle constraints on the states (x( j)

∈ X ⊂ R n , j ∈ [k + 1, k + N p ],
X is the set of admissible states) and on the control inputs (u( j

) ∈ U ⊂ R m , j ∈ [k, k + N p -1],
U is the set of admissible inputs). These constraints are easily added to the problem [START_REF] Folio | Dealing with visual features loss during a vision-based task for a mobile robot[END_REF].

III. 2D VISUAL SERVOING THROUGH NMPC

The IBVS task is formulated into the minimization of an image error over a prediction horizon. Due to the IMC structure combined with visual servoing task, the equations ( 7) become:

ε( j) = s( j) -s m ( j) s d ( j) = s * -ε( j) ( 11
)
where s * , s and s m are respectively the reference features, the current features and the predicted model features. The latter are obtained thanks to equation ( 1), discretized with a first order approximation (T e is the sampling period):

s m ( j + 1) = s m ( j) + T e L s ( j)τ( j) (12)
In [START_REF] Folio | Dealing with visual features loss during a vision-based task for a mobile robot[END_REF], this dynamic equation is solved to reconstruct the image data in case of occlusion. Consequently, the 2D visual predictive controller can be written as:

min τ {J(s, τ) = k+N p ∑ j=k+1 [s d ( j) -s m ( j)] T Q [s d ( j) -s m ( j)]}
(13) with:

• τ = {τ k , τ k+1 , ..., τ k+N c , ..., τ k+N p -1 } : the optimal sequence of the camera velocity screw; • Q : a symmetric definite positive matrix;

• N p , N c : the prediction and control horizons. The optimization problem ( 13) is solved at each sampling time T e . Only the first component of the optimal control sequence, τ(k), is really applied. Numerous optimization routines are available in software libraries to solve this kind of problem: conjugate gradient methods, quasi-Newton's methods, etc. In our case, a large-scale algorithm is used and more precisely, the function fminunc from Matlab optimization toolbox. To sum up, the different steps of the procedure are:

Step 1: Image acquisition and initializations; Test 1: s(t)s * ≤ β ( for small β > 0)

Yes: End ; No: Step 2;

Step 2: Model initialization with measurements;

Step 3: Solution of the optimization problem (13);

Step 4: Application of the optimal control τ(k);

Step 5: New image acquisition and go to Test 1;

The setting parameters of this approach are the prediction horizon (N p ), the control horizon (N c ) and the weighted matrix (Q):

• the prediction horizon is chosen in order to satisfy a compromise between scheme stability (long horizon) and numerical feasibility in term of computational time required (short horizon); • the control input is usually kept constant over the prediction horizon which corresponds to a control horizon equal to 1; • the matrix Q is often the identity matrix or can be a time-varying matrix for stabilization task. In the next section, the role of these parameters is discussed for different configurations.

IV. SIMULATION EXPERIMENTS

Different simulations illustrate what the image prediction can bring to 2D visual servoing. The results are compared with the classical approach introduced in section I (cf. Equations ( 1), ( 2), (3)). For all simulations, the sampling time T e is equal to 40 ms. This choice allows to consider a real time application with an usual camera (25 frames per second). The control task consists in positioning a perspective freeflying camera with respect to the perspective projection of four points. These four points form a square of 20 cm in length in Cartesian space. The reference image is obtained when the target pose expressed in the camera frame (R C ) is equal to P T /C = (0, 0, 0.5, 0, 0, 0) T (see Fig. 2). The first three components correspond to the translation expressed in meters and the last three components are the roll, the pitch and the yaw angles expressed in radians. The coordinates of the four points in the reference image are: s * = (u d1 , v d1 , u d2 , v d2 , u d3 , v d3 , u d4 , v d4 ) T = (-0.2, 0.2, 0.2, 0.2, 0.2, -0.2, -0.2, -0.2) T (see Fig. 3). To avoid the estimation of the depth parameter Z(t) at each iteration, its value Z * measured at the reference position is used (Z * = 0.5m). The discrete-time formulation of the interaction matrix ( 5) is used in the predictive approach. The VPC requires to set three parameters: the prediction horizon N p , the control horizon N c and the weighted matrix Q.

• The control horizon is kept constant and equal to 1 (N c = 1). Only one control is calculated over the prediction horizon. • The weighted matrix is either the identity matrix Q = I 8×8 , constant over the prediction horizon, or a timevarying matrix Q = Q(t), weighting the error at each sampling instant more and more over the prediction horizon. It gives importance to the control objective at the end of the prediction horizon. The time variation of Q is given by: Q( j) = 2 * Q( j -1) with Q(1) = I 8×8 . • The choice of the prediction horizon is discussed below. Even if VPC can handle constraints, no constraint will be considered here in order to compare the proposed approach with the classical visual servoing. In both cases, the control inputs are normalized if needed. This is accomplished by computing a scaling factor based on the knowledge of the velocity bounds (25 cm/s for the translation speed and 0.25 rad/s for the rotation speed). The scaling factor is similar to multiply all the control inputs by a gain coefficient in order to satisfy the bounds.

A. Pure rotation around the optical axis of the camera

The required camera motion is a pure rotation of π 2 radians around the optical axis Z. We compare the classical IBVS with three cases:

• N p = 1 (see Fig. 5);

• N p = 10 with Q = I and Q = Q(t) (see Fig. 6,7);

• N p = 20 with Q = I and Q = Q(t) (see Fig. 8,9); For N p = 1, the results are similar to the classical IBVS (see Fig. 4). The visual features try to follow shortest path in the image plane which is the straight line. Consequently, the camera moves and rotates respectively along and around the optical axis. The only difference is that the control law of Fig. 4 has an exponential decrease.

Remark: if the depth Z is calculated at each sampling time (Z=Z(t)), the trajectories in the image

Camera pose errors

Seconds When the prediction horizon increases (N p = 10 or N p = 20), the optimization routine gives a solution which decreases the translation motion around the optical axis whatever the choice of Q. Indeed, the only control which minimizes the cost function over the prediction horizon, consists in achieving only the rotation. Consequently,the image prediction allows to decouple the control in order to effectuate only the desired motion: the pure rotation around the optical axis.

Errors ε T x (m) ε T y (m) ε T z (m) ε W x (rad) ε W y (rad) ε W z (rad)
In Fig. 9, the visual feature trajectories are almost circles. The translation along the optical axis is near to zero. Only the rotation is done.

In Fig. 10, the evolution of

Δ Z (Δ Z = max(Z(t) T /C -Z * ) for t ∈ [0, 10s]) is represented versus the horizon prediction for Q = I and Q = Q(t). We can remark that Q = Q(t)
enhances the control decoupling. In the classical IBVS approach, a π radians rotation around the optical axis is known to lead to a failure of the control law [START_REF] Chaumette | Potential problems of stability and convergence in image-based and position-based visual servoing[END_REF]. As IBVS strategy chooses the shortest path in the Seconds Errors v (pixels)

ε v 1 ε v 2 ε v 3 ε v 4

Camera pose errors

Seconds

Errors ε T x (m) ε T y (m) ε T z (m) ε W x (rad) ε W y (rad) ε W z (rad)
Fig. 6. N p = 10 and

Q = I -0.3 -0.2 -0.1 0 0.1 0.2 0.3 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 Image plane u v 1 2 3 4 5 6 7 8 9 10 -0.4 -0.2 0 0.2 0

Camera pose errors

Seconds Camera pose errors Seconds Errors v (pixels) image plane which is the straight line, the camera performs an infinite retreat corresponding to a singularity. In contrast, the predictive controller with N p = 20 and Q = Q(t) generates a decoupling control: the π rotation around the optical axis is effectuated almost without translation (see Fig. 11).

Errors ε T x (m) ε T y (m) ε T z (m) ε W x (rad) ε W y (rad) ε W z (rad)
Seconds Errors ε T x (m) ε T y (m) ε T z (m) ε W x (rad) ε W y (rad) ε W z (rad)
ε v 1 ε v 2 ε v 3 ε v 4
ε T x (m) ε T y (m) ε T z (m) ε W x (rad) ε W y (rad) ε W z (rad)
-0.

ε T x (m) ε T y (m) ε T z (m) ε W x (rad) ε W y (rad) ε W z (rad)

B. Large displacements

The initial target pose expressed in the camera frame is given by: P T /C = (0.08, -0.13, 0.66, 0.46, 0.25, -1.1) T .

For N p = 1 (see Fig. 12) and whatever the choice of Q, the obtained results are once again approximately the same as the ones obtained with the classical IBVS. The evolution of the target position (resp. orientation) expressed in the camera frame shows the coupling of the dof involving useless movement. At least, for N p = 10 and Q = Q(t) (see Fig. 14), with approximatively the same feature trajectories in the image plane (the trajectories seem to be smoother), we can remark that the controls obtained are decoupled. This can be verified by the camera position and orientation.

Fig. 1 .

 1 Fig. 1. Internal Model Control Structure

Fig. 3 .

 3 Fig. 2. 3D desired posture

  Fig. 4. Classical IBVS

Fig. 5 .

 5 Fig. 5. N p = 1 and Q = I

Fig. 7 .

 7 Fig. 7. N p = 10 and Q = Q(t)

Fig. 8 .

 8 Fig. 8. N p = 20 and Q = I

Fig. 9 .ΔFig. 10 .

 910 Fig. 9. N p = 20 and Q = Q(t)

Fig. 11 .

 11 Fig. 11. π rotation: N p = 20 and Q = Q(t)

Fig. 12 .Fig. 13 .

 1213 Fig. 12. N p = 1 and Q = IFor N p = 5 and Q = I (see Fig.13), the results are clearly better. The feature trajectories in the image plane are more straightforward and the useless displacements in the Cartesian space are minimized.

VI. ACKNOWLEDGMENTS

The authors gratefully acknowledge Franc ¸ois Chaumette for his original idea and for the fruitful discussions.

Fig. 14. N p = 10 and

C. Robustness

The robustness w.r.t modeling errors (25% on the intrinsic camera parameters in the model) and disturbances (white noise added to the output s) is tested (see Fig. 15). Due to the IMC structure, the VPC converges to the desired image in spite of a large displacement to achieve.

Camera pose errors

Seconds For all simulations presented in this section, the computational time needed to solve the optimization problem is about 30ms. It can be largely reduced in using C/C++.

V. CONCLUSIONS

In this paper, an alternative predictive control approach for 2D visual servoing has been proposed. The control objective is formulated into an optimization problem in the image plane. The difference between the reference features and the predicted features is to be minimized over the predicted horizon in regard to the camera velocity screw inputs. The predicted features are computed thanks to the interaction matrix. This allows to avoid 3D data necessary in the case of a global nonlinear model [START_REF] Allibert | Visual Predictive Control for Manipulators with Catadioptric Camera[END_REF]. We have shown, through simulations, what image prediction can bring to IBVS: a decoupled control. This decoupling control can be of great interest for large displacements or rotations around the optical axis. Due to the IMC structure, the VPC approach is also robust in regard to modeling errors and disturbances. The visual servoing task formulated into an optimization problem is suited to deal with constraint handling. For instance, visibility constraints which ensure that the visual features stay in the camera field of view, can easily be added to the optimization problem.


Mohamad Yehya
email: mohamad.yehya@3sr-grenoble.fr

Edward Ando
email: edward.ando@3sr-grenoble.fr

Frédéric Dufour

Alessandro Tengattini

Edward And Òa

Fluid

Keywords: fluid flow, neutron imaging, porous media, heavy water, normal water

come

INTRODUCTION

The investigation of fluid flow in porous media, and specifically full field measurements of local processes, is of significant interest for many industrial and environmental applications such as the disposal of hazardous waste, water evaporation and infiltration in soil. Local measurements of fluid flow in 3D present a convenient way to access a permeability field in different engineering porous media (such as concrete, bricks and ceramics, rocks, fuel cells etc.). The ability to deduce a permeability field has important repercussions on the development of better hydraulic material models (which can be thermo-, hydro-mechanically coupled).

In the case of an incompressible fluid (water for specific pressures for instance), under laminar flow condition, the permeability of the medium can be directly obtained by applying Darcy's law. The concept of permeability defined by Darcy's law, expresses, on a macroscopic scale, the physics of the flow of a viscous Newtonian fluid (the forces due to the viscosity dominate inertial forces) at the scale of pores. In the literature, many authors have made gas-or water-permeability measurements on different materials such as concrete, 1-3 mortars, [START_REF] Loosveldt | Experimental study of gas and liquid permeability of a mortar[END_REF] rocks, [START_REF] Huo | Experimental investigation of stress-dependency of relative permeability in rock fractures[END_REF] etc.. Standard laboratory experiments consist in measuring flow from one end of a test specimen to the other without full knowledge of how the fluid flows in the sample. The current measurements only allow for the identification of average characteristics where extrapolation to the engineering scale (power stations and dams for concrete, rock mass for CO 2 reservoir, oil pumping, shale gas, etc.) is necessary. In the case of homogeneous materials and an assumed homogeneous flow, this can be sufficient. However, a large number of complex heterogeneous materials are routinely used in engineering practice today (concrete and reinforced concrete are an excellent example), whose permeability is locally inhomogeneous. In the case, the understanding of the details of fluid † Institute of Engineering Univ. Grenoble Alpes flow are pivotal because of the implications on the preferential paths of damage and weathering which, in turn affects durability and containment.

The literature reports several works investigating the porosity, the moisture content and the sorptivity of a porous specimen such as sandstone, 6, 7 concrete, [START_REF] De Beer | Testing the durability of concrete with neutron radiography[END_REF][START_REF] Kanematsu | Quantification of water penetration into concrete through cracks by neutron radiography[END_REF][START_REF] Zhang | Neutron imaging of water penetration into cracked steel reinforced concrete[END_REF][START_REF] Zhang | Neutron radiography, a powerful method to determine time-dependent moisture distributions in concrete[END_REF] porous asphalt, 12 clay rock 19 using neutron radiography. Most of the studies, if not all, focus on the rising water front by imbibition so that a local velocity field can be obtained by assessing the displacement of the front in time which is a product of the competition between capillary suction driving and the medium's permeability. This competition makes the acquired data difficult to analyse. In order to get permeability-dominated flow directly, a saturated and fluid-pressure controlled test is recommended. De Beer and Middleton 6 were the first to use neutron imaging for two-phase fluid flow through an oil bearing rock inside a Hassler high-pressure apparatus in order to reproduce deep reservoir pressure conditions. In their experiment, they injected D 2 O through a sandstone sample under a pressure of 0.02 MPa causing oil to move at 3 cm 3 /min. This important step for measuring permeability-dominated flow directly is significantly extended in this work to higher pressures allowing lower permeability materials to be tackled.

Quantification and analysis of the spatial distribution of fluid flow in porous media is to this day a major scientific challenge. In order to investigate low permeability porous media (Concrete has an intrinsic permeability of around 10 -17 m 2 ) 15 in reasonable experimental time (hours not days), a large pressure gradient is required. For example, for a material like concrete cited above (assuming no dependence of permeability on mean stress), by applying Darcy's law with a pressure gradient of 40 MPa, a water front would cross the sample of 7 cm height in around 200 minutes. Consequently, a higher confining pressure should be applied to avoid the damage by the fluid driving pressure. Making a measurement of fluid flow in saturated conditions will allow the geometric complexities inherent with capillary rise to be avoided. This paper proposes a full-field imaging method, using neutron imaging, to quantify fluid flow within a porous specimen, and applies this technique to the very challenging case of concrete (low permeability) in saturated conditions (where no obvious flow front is normally observed using regular imaging techniques).

The aim of this paper is to show the design of the new experimental set-up that allows hydraulic flow through cylindrical samples to be monitored under 50 MPa confinement using neutron radiography and tomography.

NEUTRON IMAGING

Neutron radiography has shown in previous studies its sensitivity to fluids with high hydrogen content (especially water) [START_REF] Perfect | Neutron imaging of hydrogen-rich fluids in geomaterials and engineered porous media: A review[END_REF] and multiphase fluid flow 16 in a specimen. It is a powerful tool in the determination of the hydraulic behaviour in geomaterials and concrete compared to x-rays (which itself is better suited for studying the solid skeleton [START_REF] Poinard | Compression triaxial behavior of concrete: The role of the mesostructure by analysis of X-ray tomographic images[END_REF] ).

The essential difference between x-rays and neutrons is that x-rays interact with the electron clouds of the atom and their attenuation scales with the atomic number Z (heavier atoms attenuates more), whereas the neutrons interact with the nuclei and their attenuation (and scattering) does not follow an intuitive pattern across the periodic table, as hinted by Figure 1. The difference between the two interaction mechanisms gives a difference in contrast when both types of radiation are used for tomography. In this work, all neutron imaging is performed on NeXT-Grenoble, the imaging instrument on D50 at the ILL. [START_REF] Tengattini | NeXT-Grenoble, a novel facility for Neutron and X-ray Tomography in Grenoble[END_REF] Figure 2 shows two horizontal slices of a D 2 O saturated cylindrical concrete sample of 7cm in diameter using both techniques * . Although the image is sharper with x-rays since the resolution is higher, the x-ray attenuation of the cement matrix and the aggregates is very close meaning that they are difficult to distinguish. While with neutrons, we can visualize aggregates, cement paste and pores but with lower resolution. Some heavy elements (such as titanium and aluminium) are relatively transparent to neutrons, which allows testing apparatus capable of withstanding high pressures to be used to perform in-situ neutron imaging tests. Unlike x-rays, neutrons can be highly sensitive to some light elements such as hydrogen and lithium and are furthermore isotope-sensitive. An excellent example of isotope sensitivity is between hydrogen (H) which is highly attenuating and deuterium (D, i.e., an H atom containing a neutron) which attenuates significantly fewer Difference between neutron and x-ray attenuation for some elements. Adapted from: http: //centrin.net.id/ ~nslbatan/, https://www.psi.ch/niag/comparison-to-x-ray, http://nmi3.eu/news-and-media/ neutron-imaging-past-present-and-future.html neutrons. Consequently, the significant contrast in neutron imaging between water (H 2 O) and heavy water (D 2 O) (liquids which are otherwise hydraulically similar † ) allows a fluid front in saturated conditions to be visualised with ease.

EXPERIMENTAL SETUP AND MEASURES

Neutron Hassler Cell

Figure 3 and Table 1 show the original design and the principal characteristics of the neutron Hassler cell. It was optimised to accommodate cylindrical samples (easy to pour or core) which also simplifies the acquisition of in-situ tomographies and furthermore takes advantage of the large availability of tubular membranes.

Material

Fluid content liquid

Confining volume 325 ml Table 1. General characteristics of the neutron Hassler cell.

In order to percolate the water along the axis of the cylindrical sample, a fluid injection from one side of the specimen is required, together with a pressure gradient to drive the flow. The exit pressure used here is atmospheric pressure, even though with the proposed apparatus any back-pressure can be imposed. To prevent † The difference in terms of neutron attenuation is due to the neutron capture cross-section of deuterium σ=0.00046 barn, almost 1000 times lower than the one of ordinary water. http://akilia.alwaysdata.net/scf_old/ produit-du-jour/eau-lourde.html In terms of mechanical properties, D2O is 10.6% denser and 20% more viscous than H2O at a temperature of 20 • C. However, the other properties (surface tension, boiling point, etc.) are quite similar. https://en.wikipedia.org/wiki/ Heavy_water Note that these two can mix and often there are even hydrogen/deuterium exchanges that lead to a semiheavy water (HDO). Nonetheless this process is sufficiently slow compared to the time of the test so it can be neglected. leaking on the radial edges of the specimen and to prevent the internal damage of the material due to internal pressurization, a confining pressure greater than the injection pressure should be applied on the membrane that covers the sample. One standard [START_REF]Test Method for Water Permeability of Concrete Using Triaxial Cell[END_REF] suggests for example a confining pressure at least 1.25 times the injection pressure, as well as limiting the confinement to one-half of the uniaxiam compressive strength of the concrete.

To obtain a water flow speed adapted to the imaging requirements (i.e., fast enough to maintain scanning times reasonable), low permeability materials require a large pressure gradient, and therefore a large confinement pressure. This in turn requires strong materials with low neutron attenuation. The AG3 aluminium alloy is an obvious choice for neutron attenuation but a cell respecting safety standards could not be designed, since this would require very large cell walls which, by increasing the sample-detector distance would severely affect the resolution of the system. [START_REF] Perfect | Neutron imaging and A reference for the imaging community[END_REF] On the other hand, the titanium alloy (TA6V) that contains almost 90% of pure titanium, 6% of Aluminium and 4% of Vanadium, attenuates relatively few neutrons due to the small attenuation of its elements. It has also a good compromise between mechanical strength and toughness in addition to its excellent resistance to corrosion. For these reasons, the body of the cell is chosen to be made by TA6V. Given the very high strength of this alloy (100 MPa), this cell can withstand a high confinement pressure for a minimum thickness. In this case, a confinement up to 50 MPa can be achieved with a thickness of TA6V of 19.5 mm.

To impose such confining pressures an incompressible liquid is required for safety reasons. In the proposed design D 2 O is selected as a confining fluid to minimise the overall attenuation of the apparatus. To further reduce attenuation, the volume of confining fluid is minimised, which also reduces the cost of each test (at the time of writing D 2 O with a purity of 99.8% is sourced at around 1000 €/l). Total confining volume of 325 cm 3 will be needed for a sample of 200 mm length. A cylindrical sample of 70 mm in diameter can be tested in this cell. The length of this sample can vary between 100 and 200 mm thanks to a movable shutter -labelled 8 in Figure 3. Shorter samples can be tested also by adding cylindrical shims at the ends of the specimen. In addition, a bespoke circular base has been designed to set the cell on the existing rotation stage in the beamline and to have the possibility to rotate it and make tomographies in the most repeatable conditions.

The sample is encased within a special membrane to prevent flow along the edges (see section 3.2), after which it is placed and connected to one end of the cell, specifically the fixed shutter (3 in Figure 3), which is equipped with 1/8" nozzles that give access to both the confinement chamber as well as the sample. Thereafter, the cell walls (a filleted hollow cylinder) are screwed on this first base. Finally, the remaining base, which is equipped has the same fittings as the opposite end, is added to the system. It is at this point straightforward to add the confinement fluid, impose the desired pressure and start injection into the pore space of the specimen.

Membranes

In order to seal the sample during the experiment and to prevent penetration of the confining fluid, a membrane is necessary. The imposition of the high confinement pressure requires a membrane with good mechanical resistance, while being sufficiently flexible to adapt itself to the external shape of the specimen and as little neutron-attenuating as possible.

Three membranes were tested with a concrete sample under a pressure of 50 MPa in order to choose the most suitable (see Figure 4 & Table 2 ). From an imaging point of view, the pure aluminium membrane is transparent to neutrons but it has not been possible to increase its thickness in order to not lose flexibility for the implementation. The tested aluminium membrane did not resist more than 20 MPa confinement pressure. There was a rupture at 2 cm above the O-ring and significant deformations around the external asperities of the sample. On the other hand, the neoprene membrane resisted a confining pressure of 50 MPa but attenuates significantly more than 55% of the neutron beam on NeXT because of the high hydrogen content ( Note that this membrane attenuates more than a cylindrical dry concrete sample of 7 cm diameter!). Teflon membranes of 2-mm thickness were found to be the optimal compromise, given their comparatively low neutron attenuation, good mechanical resistance and flexibility. Furthermore, these Teflon membranes are reusable and do not require their internal diameter to be within one millimetre of the sample's external diameter, which is instead a strong limitation of the aluminium membrane. Especially when using concrete, the control of the sample dimensions with one millimetre accuracy is cumbersome to ensure since some shrinkage phenomena necessarily take place during hardening and thus the external dimensions can slightly change during time.

Hydraulic equipment

Two high-pressure fluid pumps are needed to operate the cell: one to impose the confining pressure and the other for the driving pressure. During the test, two syringe pumps with accurate flow rates were used. These pumps were pressure-controlled. The pump used to maintain the confining pressure is a Sanchez VPSSH 6/700 and can reach a pressure of 70 MPa and has a small reservoir of 6 cm 3 . The injection pump is an ultra-high pressure Cetoni neMESYS XL 7000N that offers a number of different reservoirs which can be remotely refilled with valve control. With the 10 cm 3 reservoir used in this work, a pressure of 89 MPa can be reached. High-pressure tubes are used in order to ensure the connection between the pumps and the cell. Each pump should be connected to a valve in order to stop the injection (to do a tomography for example) or disconnect the confinement when necessary. A flow meter is required at the outlet of the cell to measure the exit flow.

Outline of the neutron facility

Neutron flux (cm -2 s -1 mA -1 ) 1x10 Exposure time for one image in the experiment (s) 10

Pixel Size (µm/pixel) in the experiment 100

Table 3. Characteristics of the neutron beam facility used in this project.

The preliminary test presented in this article was performed at the newly developed NeXT instrument at the Institut Laue Langevin (ILL) at Grenoble-France (https://next-grenoble.fr/). [START_REF] Tengattini | NeXT-Grenoble, a novel facility for Neutron and X-ray Tomography in Grenoble[END_REF] The ILL is the world's most powerful neutron source with 1.5x10 15 thermal neutrons per second per cm 2 (n/s/cm 2 ) and a mean neutron energy of 12.3 meV. This novel instrument, which at the time of this writing is still under development, currently allows for tomographies down to 30 microns resolution (and is expected to achieve soon a resolution of 10 microns) and has field of view up to 170x170 mm 2 . The collimation length is about 10 m and the pinhole (D) ranges between 1.5 mm and at the time of the test 15 mm, although the radiation protection has recently been upgraded to allow pinholes up to 30 mm. The neutron flux at the instrument, for a collimation ratio of L/D = 660, is about 1x10 8 n/s/cm 2 . Table 3 shows the characteristics of the neutron radiography facility used in this project.

EXPERIMENTAL METHOD

The experimental procedure is sketched in Figure 5. The following steps are required to prepare and execute a permeability test with in-situ neutron scanning:

• A neutron and x-ray tomography of the sample before starting the test is recommended in order to know its initial state and the distribution of aggregates and pores.

• The sample should be saturated with heavy water (if local permeability in the sense of Darcy's Law is to be calculated later 1). Heavy water D 2 O is selected as the ideal initial pore fluid because its significantly lower neutron attenuation (7 times less than the normal water H 2 O 15 ) simplifies the imaging of the initial state and microstructure of the sample.

• The sample is set on the base of the Hassler cell and covered with the membrane before closing the cell hermetically 2,3 (See section 3.1).

• The system of pumps, cables, valves and flow meter is then connected and set as indicated 4.

Figure 5. Experimental approach of the experiment.

• The confinement chamber is filled with heavy water 5 and air is expelled from the system by pushing the water until bubbles stop appearing from the exit nozzles.

• Confining pressure is slowly increased until the desired confinement is achieved.

• The pressure in the injection pump is slowly increased until the desired hydraulic head is achieved, injecting normal water under constant pressure while monitoring the outflow using the flow meter 6,7.

• While injecting, neutron radiographies can be acquired in order to measure the advancement of the normal water front due to the contrast between H 2 O and D 2 O. Since this technique only provides a projection of the sample, it allows for fast acquisitions, so that it is possible to follow the front "live".

• It is possible to stop the injection, which allows for the stable conditions necessary to acquire tomographies while maintaining the same confining pressure. These tomographies provide extra informations than the radiographies about the shape of the 3D front in the sample.

• By assessing the rising front in the sample, a speed map can be obtained. A speed map can be used as input to advanced inverse analysis techniques to yield a permeability field.

It should be noted that neutron images are normalised by the "dark" and "flat" calibration images (standard practice), and artefacts from gamma radiation are detected and removed from radiography by using a "remove outliers" procedure.

PRELIMINARY TEST

In this section we present a representative test performed with the proposed setup to clarify the potential of the cell and highlight some key findings enabled by the use of the cell. A standard cylindrical concrete sample of 7 cm in diameter and 7 cm in height was selected to contain sufficient aggregates to be representative, and was cast and kept saturated with heavy water. A tomography of the sample in its membrane was taken before mounting it in the cell. This tomography was acquired by taking 5034 radiographies (3 projections per position) of 2 sec exposure with no binning ‡ (images of 2048 x 2048 pixels) over 360 degrees. Figure 7 shows the attenuation of the saturated D 2 O concrete sample inside the membrane. In the image it is clear that there is a highly neutron-attenuation region, corresponding to a film of heavy water interfacing the two, which was added to avoid evaporation and to keep the sample saturated while the test is performed. It is presumed that this water is pushed outside the sample when the confining pressure forces the membrane in contact with the sample in order to seal it well.

The porosity is greater at the surface because of the wall effect in the molding of a granular medium. This porosity is probably full of heavy water that attenuates more. As we can see the sample was not fully saturated (the light zones correspond in fact with air-filled pores) since it was kept under water for one month before the test without performing vacuum water saturation methods.

The flow experiment set-up detailed in section 4 was prepared, as shown in Figure 6, at which point the test was performed. During the test, the confining pump was filled with D 2 O and the pressure was set at 30 MPa. On the other hand, the injection pump was saturated with H 2 O and refilled remotely during the experiment thanks to an H 2 O reservoir. After a stable confinement was achieved, H 2 O was injected at 20 MPa. This pressure was optimised based on the speed of the imaging instrument and the available beamtime. The neutron beam was ‡ Binning is the procedure of merging a cluster of pixels into a single pixel. As such, in a binning 2, an array of 4 pixels becomes a single larger pixel which reduce the overall number of pixels. In our case, using binning 2 images leads to images of 1024 x 1024 pixels. Binning reduces noise (at the cost of a doubled pixel size). Reading the camera directly in a binned mode may also reduce the impact of readout noise. oriented by a collimator before it interacts through the cell and the sample. The transmitted neutrons are then converted into visible light by a scintillator screen. The light emitted is deflected by the mirror and recorded by the camera.

During the water percolation process, a series of radiographies with an exposure time of 10 seconds were acquired through the injection phase, which lasted around 1hr 40mins. The neutron exposure time has been increased to compensate for the attenuation of the cell and confining fluid. Periodically the injection was paused to acquire longer neutron-tomographies in order to characterize the H 2 O front in 3D. 2-hour tomographies were launched during the experiment. 720 projections were acquired, each the result of the average of 3 projections, over an 180 degrees rotation (see section 6). The exposure time was also 10 secs and the images were acquired with binning 2. It should be noted that the chosen number of projections is an optimisation between the time of the test and the resolution/quality of the image.

RESULTS

Figure 8 shows the neutron radiographies of the normal water flow (blue -high attenuation) into the concrete sample saturated with heavy water (green -low attenuation) during the injection test.

A series of horizontal sections of the 3D tomography acquired after 39 min of H 2 O injection is shown in Figure 9, the distance between two slices being of about 1 cm in the sample. It can therein be observed that the normal water is flowing more pronouncedly in the left part of the sample. Figure 8 shows that the water percolates directly in the higher porosity zones of the concrete, notably in the large pore visible at 1cm from the bottom which was observed to rapidly fill.

CONCLUSION

This contribution presents a new experimental method and preliminary results from neutron imaging under the high pressures required to impose significant fluid flows in low permeability media, taking advantage of the high contrast between the normal water flushed in the sample and the heavy water used to saturate it. The key aim of the work was to follow a rising front in the saturated specimen under high pressure, which provides the possibility to acquire a velocity field and therefore to estimate a local permeability field by inverse analysis in the future. In this paper we presented the design of the titanium Hassler cell that allows for studying large samples with a very low permeability since high pressure can be used. The methodology and the experimental set-up used to study the driving fluid flow within the porous media by neutron imaging was presented. The test performed on concrete showed that the front is not homogeneous contrary to expectations, which opens the door to many questions concerning the transfer of properties in concrete and the key role played by the porosity. This study is a proof of concept opening new venues due to the combination of high neutron flux instruments allowing tests to be performed on low-permeability materials in reasonable time.

Figure 1 .

 1 Figure 1.Difference between neutron and x-ray attenuation for some elements. Adapted from: http: //centrin.net.id/ ~nslbatan/, https://www.psi.ch/niag/comparison-to-x-ray, http://nmi3.eu/news-and-media/ neutron-imaging-past-present-and-future.html

10 •

 10 C to 150 • C Internal diameter of the sample 70 mm Length of the sample variable, up to 200 mm

Figure 2 .

 2 Figure 2. A tomography of a cylindrical concrete sample of 7 cm in diameter using neutrons (left) and x-rays (right).

Figure 4 .

 4 Figure 4. Neutron attenuation of the Neoprene and Teflon membranes

Figure 6 .

 6 Figure 6. Experimental set-up at the NeXT-Grenoble beamline on D50 at the ILL.

Figure 7 .

 7 Figure 7. Horizontal slices of the reconstructed 3D tomography of the heavy water saturated concrete sample within its Teflon membrane before the beginning of the normal water injection in different positions of the sample (1 cm, 3 cm and 5 cm).

Figure 8 .Figure 9 .

 89 Figure 8. Time-lapse neutron radiography images of normal water injection into heavy water concrete sample. Radiographies are averaged over 11 minutes. Neutron attenuation is the raw attenuation measured net of cell, confining fluid and membrane.

Table 2 .

 2 Mechanical resistance and attenuation of the membranes.

ACKNOWLEDGMENTS

We would like to thank Erika Tudisco for the suggestion of Teflon membranes. This work was supported by EDF-SEPTEN/DTG/CIH within the Chair PERENITI agreement with the Grenoble INP Partnership Foundation. The Chair PERENITI partners shall not in any circumstances be deemed liable for the content of this publication which is only binding its author. The Laboratoire 3SR is part of the LabEx Tec 21 (Investissements d'Avenir-Grant agreement n • ANR -11-LABX-0030).


P Gentile

A Valero

D Gaboriau

M Boniface

D Aldakov

S Sadki

D Aradilla

ALD Alumina Passivated Silicon Nanotrees electrodes for new Ultrastable Microsupercapacitors

The current trend towards miniaturized and autonomous electronic devices requires innovative energy storage solutions. For instance, autonomous micro-sensor networks or implantable medical devices would need a robust power source with high cyclability and a large power density, which might be out of the scope of conventional battery technologies. For such applications, microsupercapacitors (µSCs) are promising alternatives, and their integration "on-chip" could allow significant innovations. 1 However, finding a suitable "on-chip" µSCs technology implies addressing key challenges, such as temperature resistance, silicon industry compatibility and good electrochemical performances on a small footprint. Following this trend, our work focuses on µSCs using highly doped silicon nanowires (SiNWs) and nanotrees 2,3 (SiNTs) as current collector. The fine morphological tuning of the nanostructure allowed by the bottom-up approach permits a careful design of the electrodes architectures, with a considerable liberty compared to other techniques. Such latitude allows optimizing porosity and ionic and electronic pathways while keeping robust mechanical performances, depending on the target application or other parameters like surface modification, functionalization by pseudo-capacitive material, electrolyte… Nanostructures such as SiNWs and SiNTrs demonstrated excellent cyclability with more than 1 million cycles of galvanostatic charge/discharge under a 4 V wide electrochemical windows in EMI-TFSI ionic liquid, with large power densities and good capacitance values. 3,4 Moreover, the use of silicon for electrode material allows extremely interesting developments towards "on-chip" integration and potential scaleup production using standard silicon industry processes for small micro-sized energy storage devices. Furthermore, we have also investigated the impact of the addition of a high-k dielectric layer, such as Al2O3 as protective films on silicon nanotrees. The electrochemical performances was enhanced, allowing symmetric 2 electrodes device to reach an unprecedented cell voltage of 5.5 V, improving energy and maximum power densities compared to unmodified nanostructured silicon. The cyclability was also largely enhanced, with only 3% capacitance fade after 10 6 galvanostatic charge/discharge cycles at 4 V, and no degradation even after several 10 5 cycles over 5 V 5 . In addition, the protective alumina layer makes it possible to use aqueous electrolytes, not usable with crude silicon, in order to significantly increase the capacities of the µSCs from 1mF/cm 2 to 10mF/cm 2 and open the door to use metal oxides. This project has received funding from


Bouhadjera Feriel

Ould Saïd

Elias Corresponding

Mohamed Riad Remita

Nonparametric relative error estimation of the regression function for censored data

Keywords: Asymptotic normality, Censored data, Kernel estimate, Relative regression error, Uniform almost sure convergence, V-C classes

Let (T i ) i be a sequence of independent identically distributed (i.i.d.) random variables (r.v.) of interest distributed as T and (X i ) i be a corresponding vector of covariates taking values on R d . In censorship models the r.v. T is subject to random censoring by another r.v. C. In this paper we built a new kernel estimator based on the so-called synthetic data of the mean squared relative error for the regression function.

We establish the uniform almost sure convergence with rate over a compact set and its asymptotic normality. The asymptotic variance is explicitly given and as product we give a confidence bands. A simulation study has been conducted to comfort our theoretical results.

Introduction

Let (X i , T i ) i be a R d × R * + , (d ≥ 1) valued sequence of random vectors that we assume drawn from the pair (X, T ) which is defined on a probability space (Ω, F, P). The purpose of this work is to study the effect of a random covariable X on a r.v. T which is subject to right censoring by another r.v. C. This relation of regression is modeled by:

T = r(X) + , (1.1)
where r(•) is the regression function and a sequence of error independent to X. Usually, r(•) = E[T |X] is estimated by minimizing the mean squared loss function E[(T -r(X)) 2 X]. However, this loss function is based on some restrictive conditions that is the variance of the residual is the same for all the observations, which is inadequate when the data contains some outliers. Therefore, in order to overcome this drawback we consider an alternative approach allow to construct an efficient predictor even if the data is affected by the presence of outliers. So, in this paper the limitations of the classical regression are counteracted by estimating the regression function with respect to the minimization of the following mean squared relative error, for T > 0,

E T -r(X) T 2 X . (1.2)
The latter is a more meaningful measure of performance of a predictor than the usual error in the presence of outliers. It is easy to see that the solution of the minimization problem of (1.2) is given by

r(X) = E T -1 |X E T -2 |X . (1.3)
Park and Stefanski (1998) have shown that the solution given by (1.3) satisfies

E T -1 |X E T -2 |X
≤ E T X a.s. (1.4) provided that the first two conditional inverse moments are finite. The authors consider parametric approaches to estimate the regression function r(•) which focused on estimating the mean and variance functions modeling methods [START_REF] Carroll | Transformation and Weighting in Regression[END_REF] of the inverse response T -1 as function of X. Without claiming to be exhaustive, we can quote [START_REF] Narula | Prediction, linear regression and the minimum sum of relative errors[END_REF] who studied an estimation method for minimizing the sum of absolute relatives residuals. [START_REF] Farum | Improving the relative error of estimation[END_REF] developed an estimation method designed to reduce absolute relative error. Khoshgoftaar et al.

(1992) studied the asymptotic properties of the estimators minimizing the sum of squared relative errors. In this contribution, we focus on nonparametric approach. To the best of our knowledge, only the paper of [START_REF] Park | Relative error prediction via kernel regression smoothers[END_REF] study the nonparametric regression using the relative error as loss function. They studied the asymptotic properties of an estimator minimizing the sum of squared relative errors by applying local linear approach.

In many estimation problems, it is not always possible, to make complete measurements when the available sample data is incomplete in the sense that measures are not available for all members of a random sample. For example, in medical follow-up studies, it often happens for various reasons, that the duration of interest can not be observed. This may be due to the loss of view of the patient at the beginning or end of the study period. These values are censored. The censored values, although unknown, must be taken into account to obtain a correct estimate and precise conclusions. For such practical observations, conventional statistical procedures are no longer valid and more elaborate techniques are used to model such observations.

One of the classical cases for incomplete data is the right-censored data. In this case, we observe another r.v. C with continuous distribution function (d.f.) G, we can only observe a sample account the effect of censorship on the distribution. For that we consider the sample (Y i , δ i ) i , for 1 ≤ i ≤ n and we put :

(X i , Y i , δ i ) where Y i = T i ∧ C i and δ i = 1 {T i ≤C i } , for i = 1, . . . ,
T * ,- i = δ i Y - i G(Y i ) , = 1, 2 (1.5)
where G is the survival function of the censoring rv C.

All along this paper, we suppose that:

(T i , X i ) i and (C i ) i are independent for 1 ≤ i ≤ n. (1.6)
Then from the equation (1.5) and the hypothesis (1.6), we get,

E[T * ,- 1 |X 1 ] = E δ 1 Y - 1 G(Y 1 ) |X 1 = E E δ 1 Y - 1 G(Y 1 ) |T 1 , X 1 |X 1 = E T - 1 G(T 1 ) E 1 {T 1 ≤C 1 } |T 1 |X 1 = E[T - 1 |X 1 ]
. This paper offers then an alternative approach to traditional estimation models by considering the minimization of the least relative error for regressions models when the data are randomly right censored. We establish the strong and uniform consistencie (with rate) of the constructed estimator and then the asymptotic normality has been shown. At the best of our acknowledge there is no result concerning the nonparametric regression function for censoring data using the relative error.

The rest of the paper is organized as follows: Section 2 is devoted to the presentation of the new estimator of the mean squared relative error of the regression function. The assumptions and main results are given in Section 3. Simulations are drawn in Section 4. Finally, the proofs are relegated to Section 5 with some auxiliary results.

Definition of the new estimator

Let (T i ) 1≤i≤n be an i.i.d. n-sample of r.v. of interest with commun unknown continuous d.f. F and let (X i ) 1≤i≤n be a corresponding vector of covariates with joint density function f (•). As mentioned before the solution of (1.2) is given by

r(x) = E T -1 |X = x E T -2 |X = x := r1 (x) r2 (x) , (2.1) with r (•) = r (•) f (•) where r (x) = R * + t -f T,X (t, x)dt for = 1, 2.
Recall that, in the case of complete data, a well-known Nadaraya Watson (N-W) estimator of r(•)

is given by

r * (x) = n i=1 T i W i,n (x)
with

W i,n (x) =            K x-X i hn n i=1 K x-X i hn , if n i=1 K x-X i hn = 0; 1/n, otherwise,
where h n is a sequence of positive real numbers (bandwidth) that decreases to zero when n goes to infinity and K is a kernel function defined in R d . Thus, a natural estimator of (2.1) is given by

r n (x) = n i=1 T -1 i K x -X i h n n i=1 T -2 i K x -X i h n , (2.2)
this is the analogous N-W estimator which is nothing other than a special case of the censored case.

As mentioned before, when the r.v. T is subject to right censoring by another r.v. C, we define rn (x) as a "pseudo-estimator" of r(x) that is, for any x ∈ R, we have,

rn (x) = n i=1 δ i Y -1 i G(Y i ) K x -X i h n n i=1 δ i Y -2 i G(Y i ) K x -X i h n =: r1,n (x) r2,n (x) . (2.
3)

The latter can not be calculated as G is unknown. Then to define a genuine estimator of r(•), we replace G by its [START_REF] Kaplan | Nonparametric estimation from incomplete observations[END_REF] estimator which is defined by

G n (t) =        n i=1 1 - 1 -δ i n -i + 1 1 {Y i ≤t} if t < Y (n) , 0 otherwise, (2.4)
where

Y (1) ≤ Y (2) ≤ • • • ≤ Y (n)
are the order statistics of the Y i and δ i is the indicator of noncensoring. The properties of Ḡn (t) have been studied by many authors. So a calculable estimator of r(•) is given by

r n (x) = n i=1 δ i Y -1 i G n (Y i ) K x -X i h n n i=1 δ i Y -2 i G n (Y i ) K x -X i h n =: r 1,n (x) r 2,n (x) (2.5)
where

r ,n (x) = r ,n (x) f n (x) = n i=1 δ i Y - i G n (Y i ) K x -X i h n n i=1 K x -X i h n , (2.6)
for = 1, 2 and f n (•) is the well-known kernel estimator of the joint density function f (•).

Hypotheses and main results

In order to state our results, we introduce some notations. For any d.f. L, let τ L = sup{y, L(y) < 1} be a upper endpoint of L. Assume that τ F < ∞, G(τ F ) > 0. All along the paper, when no confusion is possible, we denote by M any generic strictly positive constant such that M ≥ T -and by

r (•) = E[T -|X = •] the conditional -inverse moments of T given X and = 1, 2. Furthermore r (•) = r (•)f (•), with f is the density of X. On the other hand, log 2 (•) = log log(•) denotes the iterated logarithm function. Finally denote C 0 = {x ∈ R/f (x) > 0} the open set and C be a compact subset of C 0 .
We will make use of the following hypotheses.

H. The bandwidth h n satisfies:

i) lim n→∞ h n = 0, lim n→∞ nh n = +∞, lim n→∞ log n nhn = 0 ii) h n log 2 n = o(1). iii) lim n→∞ nh 5 n = 0. K. The kernel K(•) is: i) Continuously differentiable compactly supported density function. ii) R |t|K(t)dt = 0 and R t 2 K(t)dt < ∞ holds. iii) R tK 2 (t)dt < ∞ D. i) The function r (•), for = 1, 2
, is twice continuously differentiable and

sup x∈C |r (x)| < +∞.
ii) The function Υ (•), for = 2, 3, 4, is continuously differentiable and

sup x∈C Υ (x) < +∞.
iii) There exists Γ > 0 such that r2 (x) > Γ for all x ∈ C.

Discussions on the hypotheses

1. The independence assumption between (C i ) i and (T i , X i ) i may seem to be strong and one can think of replacing it by a classical conditional independence assumption between (C i ) i and (T i ) i given (X i ) i . However in the conditionally hypothesis we propose the following estimator for the regression function r(x) where r ,n for = 1, 2 are given by

r ,n (x) = n i=1 δ i Y - i G n (Y i X i ) K x -X i h n (3.1)
where G n (Y i X i ) is Beran's estimator of the survival conditional distribution of the censored r.v. C given X. Then we get an analogous estimator as in (2.5) using (3.1). As mentioned before and as far as we know there is no rate of convergence for this estimate as in the unconditional case (see [START_REF] Deheuvels | Functional limit laws for the increments of Kaplan-Meier product limit processes and applications[END_REF]. We think that this issue has to be addressed if we aim to get rates of convergence. Moreover our framework is classical and was considered by [START_REF] Carbonez | Partitioning estimates of a regression function under random censoring[END_REF]

Results

We can now present our results. The proofs of these are established in Section 5. We first state a uniform consistency result with rate for r n (•).

Theorem 1. Under hypotheses H i), K and D i), iii), we have:

sup x∈C | r n (x) -r(x) |= O a.s. max log n nh 2 n 1/2 , h 2 n as n -→ ∞.
Remark 1. It is clear that we can give the same result in R d , d > 1 without difficulties. The proofs are analogous. Therefore, the Theorem 1 becomes:

sup x∈C | r n (x) -r(x) |= O a.s. max log n nh d+1 n 1/2 , h 2 n as n -→ ∞.
In what follows we will state the asymptotic normality result. For this, let

Σ(x) =    Υ 2 (x) Υ 3 (x) Υ 3 (x) Υ 4 (x)   
be the covariance matrix, with

Υ 2 (x) = t -2 G(t) f T,X (t, x)dt and Υ 3 (x) = t -3 G(t) f T,X (t, x)dt, for = 1, 2.
Now we are in position to give our asymptotic normality result.

Theorem 2. Suppose that hypotheses H , K and Di), ii) holds. Let A = x ∈ C and r (x) = 0, = 1, 2 and Υ j (x) = 0, j = 2, 3, 4 , we have

nh n (r n (x) -r(x)) D --→ N (0, σ 2 (x)) as n -→ ∞ where σ 2 (x) = κ Υ 2 (x)r 2 2 (x) -2Υ 3 (x)r 1 (x)r 2 (x) + Υ 4 (x)r 2 1 (x) r4 2 (x)
for κ = K 2 (t)dt and D --→ denotes the convergence in distribution.

Confidence interval

The determination of confidence interval requires the estimation of the unknown quantity σ n (x). A plug-in estimate and using the following estimate of Υ 2 (x), for = 1, 2, and Υ 3 (x) given by

Υ 2 (x) = 1 nh n n i=1 Y -2 i G n (Y i ) K x -X i h n and Υ 3 (x) = 1 nh n n i=1 Y -3 i G n (Y i ) K x -X i h n
respectively and (2.5) we get a consistent estimate of σ 2 (x). This yields a confidence interval of asymptotic level 1 -ζ for r(x) given by

r n (x) -t 1-ζ/2 σ n (x) √ nh n ; r n (x) + t 1-ζ/2 σ n (x) √ nh n
where t 1-ζ/2 denotes the 1 -ζ/2 quantile of the standard normal distribution.

Comeback to complete data

At the best of our knowledge there are no analogous results for the complete.The analogous results can be state by putting C = +∞ and therefore G(•) = 1.

To give an overview of the performance of our estimator, we graph it in the next section.

Simulation Study

The main objective of this part is to evaluate the good behavior of our estimator for different censoring rates and sample sizes and to show the efficiency of this approach compared to the classical one.

Consistency

Simulations settings

For this purpose, simulation data are generated from model (1.1) where covariates X have normal distribution on N (5, 2) and random effect have standard normal distribution. For the rest, we proceed in the following way:

• Generate the censoring variable C according to the normal law with (µ = 11, σ = 1).

• Calculate the response variable T = αX + β + c with (α = 2, β = 1 and c = 0.2).

• The censored data are calculated as Y = T ∧ C and δ = 1 {T ≤C} .The observed data therefore becomes (X, Y, δ).

• The Kaplan-Meier estimator is calculated for the distribution function of censorship variable C in (2.4).

• The choice of K is not decisive, we choose then the standard Gaussian kernel (K(u)

= 1 √ 2π exp (-u 2 2
)). In contrast, the choice of bandwidth is crucial that's why we take the optimal one h n = 0.55 log(n) n 0.33 .

• Finally, we calculate the expression of our estimator obtained from (2.5) for a compact set

C = [1, 4].
Under each simulation setting, 100,300 and 500 replications are conducted.

Simulation results

Effect of sample size with fixed censorship rate. From Figure 1, we can see that the quality of fit increases with n when censoring rate (CR) and bandwidth kept unchanged. Nonlinear functions We consider the case of nonlinear regression by choosing this three kinds of model:

(1) Parabolic T = x 2 + 1 + ε,

(2) Sinusoidal T = sin 1 2 x 2 + 1 + ε, ( 3
) Exponential T = exp 1 2 x + ε.
The curves are shown in Figure 3. Note that the quality of fit deteriorates when the period is very small.

Classical regression versus Relative error regression with respect to the censorship rate

In order to highlight the efficiency of relative error estimation, we draw up a comparative study.

For that, we simulated the classical regression estimator for randomly right censored data defined in Guessoum and Ould Saïd (2008) by

r(x) = n i=1 δ i Y i Ḡn (Y i ) K x -X i h n n i=1 K x -X i h n ,
for the same parameters listed below. From Figure 4 below, it is clear that the classical regression estimator deteriorates when the censorship rate increases considerably. Effect of outliers for the two methods with a fixed sample size and censorship rate.

To show the robustness of our approach, we generate the case where the data contains outliers.

For that we set both that sample size and censorship rate (n = 500 and T C ≈ 50). To create this outlier effect, 20 values of this sample are multiplied by a factor called M F . From Figure 5, we can see that our estimator is very close to the theoretical curve with respect to the classical one. Then, it is very clear that our approach is widely better than the classical one in the presence of outliers.

M SE 1 = 1 n n i=1 (T i -r(X i )) 2 and M SE 2 = 1 n n i=1 (T i -r n (X i )) 2
for three sample sizes and censoring level.

It can be seen from Table 1

Asymptotic normality

The purpose of this part is to highlight the theoretical results obtained in Theorem 2, by studying by simulation the asymptotic normality. To do this, we compare the shape of the estimated density to that of the standard normal density in the case of a linear regression model:

T = 2X + 1 + ,
we reproduce the same steps as in the previous subsection for X exp(1.5) and C exp [START_REF] Chow | Probability theory. Independence, interchangeability, martingales[END_REF].

Throughout this subsection, we fix x = 0 and replicate m independent n-sample size. Then, we calculate the asymptotic variance. For that we replace the Υ k (.) for k = 2, 3, 4 by their estimators in (3.3) and r ,n (.) for = 1, 2 by their estimators in (2.6). A calculable estimator of the normalized deviation is given by:

σ 2 n (0) = κ Υ2 (0)r 2 2,n (0) -2 Υ3 (0)r 1,n (0)r 2,n (0) + Υ4 (0)r 2 1,n (0) r4 2,n (0)
, we consider now the sequence:

A j = nh n σ 2 n,j (0) 1/2
(r n,j (0) -1), which under Theorem 2, A j follows asymptotically to N (0, 1). Then, we build a kernel density estimator for the A k that we compare with the standard normal law for different values of n and

h * = c log(m) m 0.2
where the constant c is chosen appropriately. Finally, for a sample size m = 200 and a censorship rate (CR ≈ 66%), we conduct n = 100, 300 and 500 replications. The figure 6 show the quality of goodness of fit.

Auxiliary results and proofs

Proof. Using (2.1), (2.3), (2.5), (2.6) and for x ∈ R, we consider the following decomposition :

r n (x) -r(x) = 1 r2,n (x) r1,n (x) -r1,n (x) + ( r1,n (x) -E[ r1,n (x)]) + (E[ r1,n (x)] -r1 (x)) + r(x) ( r2,n (x) -r2,n (x)) + (E[ r2,n (x)] -r2,n (x)) + (r 2 (x) -E[ r2,n (x)])
which by triangle inequality, we have

sup x∈C |r n (x) -r(x)| ≤ 1 inf x∈C |r 2,n (x)| {sup x∈C {|r 1,n (x) -r1,n (x)| + | r1,n (x) -E[ r1,n (x)]| + |E[ r1,n (x)] -r1 (x)|} + sup x∈C |r(x)|{| r2,n (x) -r2,n (x)| + |E[ r2,n (x)] -r2,n (x)| + |r 2 (x) -E[ r2,n (x)]|}}.
In the sequel, we give a sequence of lemmas that are helpful in proving our results.

Lemma 1. Under hypotheses H i) and K , we have, for = 1, 2, that

sup x∈C |r ,n (x) -r ,n (x)| = O a.s. log 2 n n 1/2
as n -→ ∞.

(5.1)

Proof. For = 1, 2, we have

|r ,n (x) -r ,n (x)| = 1 nh n n i=1 δ i Y - i G n (Y i ) K x -X i h n - δ i Y - i G(Y i ) K x -X i h n = 1 nh n n i=1 T - i G n (T i ) 1 {T i ≤C i } K x -X i h n - T - i G n (T i ) 1 {T i ≤C i } K x -X i h n ≤ 1 nh n n i=1 T - i K x -X i h n 1 G n (T i ) - 1 G(T i ) ≤ 1 nh n n i=1 T - i K x -X i h n G n (T i ) -G(T i ) G n (T i )G(T i ) ≤ sup t<τ F | G n (t) -G(t) | G n (τ F )G(τ F ) M nh n n i=1 K x -X i h n ,
then by using the strong law of large numbers (SLLN) and law of iterated logarithm (LIL) on the censoring law (see formula (4.28) in [START_REF] Deheuvels | Functional limit laws for the increments of Kaplan-Meier product limit processes and applications[END_REF], we get,

sup x∈C | r ,n (x) -r ,n (x) |≤ M G 2 (τ F ) E 1 h n K x -X 1 h n log 2 n n .
Then hypotheses H i) and K complete the proof of the lemma.

Lemma 2. Under hypotheses H i), K and Di), we have, for = 1, 2, that:

sup x∈C | E[ r ,n (x)] -r (x) |= O(h 2 n ) as n -→ ∞. (5.2)
Proof. Using the conditional expectation properties, we get,

E [ r ,n (x)] = 1 h n E δ 1 Y - 1 G(Y 1 ) K x -X 1 h n = 1 h n E K x -X 1 h n E δ 1 Y - 1 G(Y 1 ) |X 1 = 1 h n K x -u h n E δ 1 Y - 1 G(Y 1 ) |X 1 = u f (u)du,
and as

1 {T 1 ≤C 1 } Y - 1 = 1 {T 1 ≤C 1 } T - 1 , we get E T - G(T 1 ) E 1 {T 1 ≤C 1 } |T 1 |X 1 = u = r (u). (5.3)
By a change of variable and using r (•) = r (•)f (•), we get

E 1 h n δ 1 Y - G(Y ) K x -X h n -r (x) = 1 h n K x -u h n r (u)f (u)du -r (x) = 1 h n K x -u h n r (u)du -r (x) = K(t)[r (x -h n t) -r (x)]dt .
We use a Taylor expansion to r (•) for ζ ∈]x -h n t, x[, we get

sup x∈C |(E( r ,n (x)) -r (x)| = sup x∈C K(t)[-h n tr (x) - h 2 n t 2 2 r (ζ)]dt ≤ h n sup x∈C tK(t)r (x)dt + h 2 n sup x∈C t 2 2 K(t)r (ζ)dt ≤ h n sup x∈C tK(t)r (x)dt + h 2 n 2 sup x∈C r (ζ)t 2 K(t)dt .
Under Hypotheses H i) and K ii), the first term is equal to zero. The second term goes to zero for n → ∞ from hypotheses Di) and K ii). The last result complete the proof of the lemma. (5.4)

Proof. Let us consider the i.i.d sequence (X 1 , Y 1 , δ 1 ), . . . , (X n , Y n , δ n ) and define

Φ n = θ x : R × R * + × {0, 1} → R + /θ x (u, y, δ) = δy - nh n G(y) K x -u h n , x ∈ R .
By Lemma (3b) in Giné and Guillou (1999), Φ n is Vapnik-Cervonenkis (V-C) class of no-negative measurable functions. These are uniformly bounded with respective envelopes

Θ = M K ∞ nh n G(τ F ) .
Moreover,

sup x∈C E [θ x (X 1 , Y 1 , δ 1 )] ≤ M K ∞ nh n G(τ F ) =: U n
In the same way, we get

sup x∈C V ar [θ x (X 1 , Y 1 , δ 1 )] ≤ sup x∈C E θ 2 x (X 1 , Y 1 , δ 1 ) ≤ M 2 K 2 2 f ∞ n 2 h n G 2 (τ F ) =: σ 2 n with σ n ≤ U n for
P sup θx∈Φn n i=1 (θ x (X i , Y i , δ i ) -E[θ x (X 1 , Y 1 , δ 1 )]) > A log n nh n ≤ L exp        - A log n nhn L M K ∞ nhn Ḡ(τ F ) log        1 + A √ log n K ∞ (nhn) 3/2 Ḡ(τ F ) L M K 2 √ f ∞ √ nhnG(τ F ) + M K ∞ nhn Ḡ(τ F ) log B K ∞ √ hn K 2 2 √ f ∞ 2              
, and using log(1 + x) ≈ x( for x → 0), the right-hand of the last equation becomes an order of

L exp      - A Ḡ(τ F ) log n nhn LM K ∞ nh n A √ log n K ∞ (nhn) 3/2 Ḡ(τ F ) L M K 2 √ f ∞ √ nhnG(τ F ) 2      = Ln - G 2 (τ F ) M 2 K 2 2 f ∞ ( A L ) 2 ,
which by an appropriate choice of the constant A, can be made O(n -3/2 ). The latter being a general term of summable series and by Borel-Cantelli's lemma we conclude the proof.

Then Lemma 1-Lemma 3 permit to conclude the proof of the Theorem 1.

Next we proceed to the proof of the Theorem 2.

Proof. Our goal is to show

nh n (r n (x) -r(x)) D --→ N (0, σ 2 (x)) as n -→ ∞.
Note that, for = 1, 2,

nh n (r ,n (x) -r (x)) = nh n (r ,n (x) -r ,n (x)) + nh n ( r ,n (x) -E[ r ,n (x)]) + nh n (E[ r ,n (x)] -r (x)) =: Λ ,n (x) + Γ ,n (x) + Σ ,n (x).
First, we consider the negligible terms Λ ,n and Σ ,n . Proof. From Lemma 1, under H ii), we get,

Λ ,n (x) = nh n (r ,n (x) -r ,n (x)) = O a.s. h n log 2 n = o a.s. ( 1
).

(5.5)

In the same way, from Lemma 2, under H iii), we have,

Σ ,n (x) = nh n (E[ r ,n (x)] -r (x)) = O nh 5 n = o(1). (5.6)
Now we consider the dominant terms Γ ,n (x) for ∈ {1, 2} and prove Lemma 5.

Lemma 5. Under hypotheses H i),K and Di), ii), we have

(Γ 1,n (x), Γ 2,n (x)) T D --→ N (0, κΣ(x))
as n → ∞.

Proof. We first estimate the asymptotic variance, for = 1, 2, we get

V ar [Γ ,n (x)] = nh n E r2 ,n (x) -E 2 [ r ,n (x)] = h -1 n E δ 1 Y -2 1 G 2 (Y 1 ) K 2 x -X 1 h n -E 2 δ 1 Y - 1 G(Y 1 ) K x -X 1 h n = V 1 -V 2 .
For V 2 proceeding as in Lemma 2 and under H i),K and Di), we have

V 2 = h -1 n h n K(s)r (x -h n s)ds 2 = o(1). (5.7) Furthermore for V 1 , we get, V 1 = h -1 n E δ 1 Y -2 1 G 2 (Y 1 ) K 2 x -X 1 h n = h -1 n E K 2 x -X 1 h n E T -2 1 G(T 1 ) |X 1 = h -1 n K 2 x -u h n t -2 G(t) f T 1 |X 1 (t|u)dtf (u)du = h -1 n K 2 x -u h n t -2 G(t) f T 1 ,X 1 (t, u)dtdu = h -1 n K 2 x -u h n Υ 2 (u)du
by a change of variable and Taylor expansion, by K iii), we get,

= h -1 n K 2 (t)Υ 2 (x -h n t)h n dt = Υ 2 (x) K 2 (t)dt + o(h n ) = Υ 2 (x)κ + o(h n ) (5.8)
which together with (5.7) and (5.8) gives under Dii), for = 1, 2

V ar [Γ ,n (x)] -→ Υ 2 (x) K 2 (t)dt as n → ∞.

In addition under Dii) and K iii), we get easily

Cov(Γ 1,n (x), Γ 2,n (x)) = E[Γ 1,n (x)Γ 2,n (x)] -E[Γ 1,n (x)]E[Γ 2,n (x)] = h -1 n E δ 1 Y -3 1 G 2 (Y 1 ) K 2 x -X i h n -E δ 1 Y -1 1 G(Y 1 ) K x -X 1 h n E δ 1 Y -2 1 G(Y 1 ) K x -X 1 h n = Υ 3 (x)κ + o(1).
Next, we will show that any linear combinations are asymptotically gaussian. For let (z 1 , z 2 ) t be a real numbers, we put

∆ n (x) = 2 =1
z Γ ,n (x) =:

n i=1 z 1 ∆ 1 i,n (x) + z 2 ∆ 2 i,n (x) (5.9)
where

∆ i,n (x) := (nh n ) -1/2 δ i Y - i G(Y i ) K x -X i h n -E δ i Y - i G(Y i ) K x -X i h n for = 1, 2.
Now in order to show that (5.9) is asymptotically normal we verify the Berry-Essèen condition (Chow and Teicher (1997), p. 322). For that we need to prove :

ρ 3 n := n i=1
E ∆ i,n (x) 3 -→ 0 (5.10) with

E |∆ i,n (x)| 3 = h n n 3/2 E δ i Y - i h n G(Y i ) K x -X i h n -E δ i Y - i h n G(Y i ) K x -X i h n 3 .
Applying the C r inequality (see [START_REF] Loève | Probability theory[END_REF], p. 155), we get

E ∆ i,n (x) 3 ≤ 4 h n n 3/2 E δ i |Y i | -3 h 3 n G 3 (Y i ) K 3 x -X i h n + E δ i Y - i h n G(Y i ) K x -X i h n 3 ≤ 4 h n n 3/2 {M 1 + M 2 } .

n, with ∧ denotes the minimum and 1 A

 1 is the indicator function of the event A. When we talk about censored data, several authors like Carbonez et al. (1995), Kohler et al. (2002), Delecroix et al. (2008) and Guessoum and Ould Said (2008) uses the synthetic data that take into

Figure 1 :

 1 Figure 1: r(x) = 2x + 1, CR ≈ 50% and n = 100, 300 and 500, respectively.

Figure 2

 2 is obtained by varying the censoring rate for a fixed sample size (n = 300) and for that we push the variable of interest on the right by increasing the average of the normal distribution to observe more censorship variable (the number of complete observation decreases). It can be seen that the forecasting quality decreases when the CR rises in particular in the border.

Figure 2 :

 2 Figure 2: r(x) = 2x + 1, n = 300 and CR ≈ 15%, 50% and 80%, respectively.

Figure 3 :

 3 Figure 3: CR ≈ 50%, n = 500 for parabolic, sinus and exponential, respectively.

Figure 4 :

 4 Figure 4: r(x) = 2x + 1, n = 100 and CR ≈ 15%, 50% and 75%, respectively.

Figure 5 :

 5 Figure 5: r(x) = 2x + 1, n = 500, CR ≈ 50% and M F = 10, 50 and 100 respectively.

Figure 6 :

 6 Figure 6: CR ≈ 66% and n = 100, 300 and 500, respectively.

Lemma 3 .

 3 Under hypotheses H i) and K i), we have, for = 1, 2sup x∈C | r ,n (x) -E[ r ,n (x)] |= O a.s.log n nh n as n -→ ∞.

Lemma 4 .

 4 Under H ii),iii) and by Lemma 1, Lemma 2, both √ nh n Λ ,n (x) and √ nh n Σ ,n (x) are o(1) as n → ∞.

and as M 1 2 =

 12 and M 2 are bounded under K which gives that ρ 3 n o(1). Now, under H i) the property (5.10) is satisfied which proves the asymptotic normality of Γ ,n (x) and, together with (5.5) and (5.6), complete the proof of Lemma 5. Now to complete the proof of Theorem 2, consider the mapping θ from R × R * + to R defined by θ(x, y) = x/y. We deduce from Mann-Wald's Theorem (see Rao 1965, p. 321) that:nh n (r n (x) -r(x)) D --→ N (0, κ∇θ T Σ(x)∇θ)where the gradient ∇θ T = ∂θ ∂x , ∂θ ∂y is evaluated at (r 1 (x), r2 (x)). Simple algebra gives then the varianceσ 2 (x) = κ Υ 2 (x)r 2 2 (x) -2Υ 3 (x)r 1 (x)r 2 (x) + Υ 4 (x)r 2 1 (x) r4 2 (x), which completes the proof.

  and Kohler et al. (2002) among others. Note finally that this assumption implies the independence between (C i ) i and (T i ) i which ensures the identifiability of the model. 2. The hypothesis τ F < τ G is classical for asymptotic normality results in the censorship framework. It implies that G(T ) ≥ G(τ F ) > 0 a.s.3.The Hypotheses H i) and K concern the smoothing parameter h n and the kernel K(.) and are standard in nonparametric regression estimation for complete or incomplete data. Moreover,

D i) is needed to study the bias term. On another side, hypothesis D iii) is used to state the uniform consistency of the constructed estimator. Finally, hypotheses H ii), iii) and D ii) are needed for get asymptotic normality.

  Table of mean squared error. To close this part, we compare the values of the mean squared error (M SE) for the classical regression (M SE 1 ) and the relative error regression (M SE 2 ) methods

	represented by

Table 1 :

 1 that the variability of the mean squared error (M SE) of the two methods for a low censoring rate is not significantly considerable, i.e. the performance is the same for both methods. However, when the data is affected by the presence of censoring the M SE of relative error regression becomes smaller than the classical regression. It means that the relative error regression model is more stable than the classical regression in the presence of censorship. The M SE errors according to the censoring rates.

	Sample	CR	M SE 1	M SE 2
	size	(≈ %) (Classical Regression) (Relative Error Regression)
		20	0.0150	0.0027
	n=100	50	0.0802	0.0339
		80	0.4009	0.0885
		20	0.0103	0.0057
	n=300	50	0.0275	0.0138
		80	0.1284	0.0425
		20	0.0078	0.0008
	n=500	50	0.0138	0.0032
		80	0.1359	0.0136

  n large enough. Now applying Talagrand's inequality [see Proposition 2.2 in Giné and Guillou (2001)], with t ≥ A log n nhn , there exist two positives constants L and B such that