Mohamed Elobaid
email: mohamed.elobaid@uniroma1.it

Mattia Mattioni
email: mattia.mattioni@uniroma1.it

Salvatore Monaco
email: salvatore.monaco@uniroma1.it

Dorothée Normand-Cyrot
email: dorothee.normand-cyrot@centralesupelec.fr

Digital path-following for a car-like robot

Keywords: Robotic systems, Novel control theory and techniques, Sampled-data nonlinear control

This paper deals with path following under digital control for a car-like robot via TFL. Assuming the existence of a continuous-time feedback, a multi-rate sampled-data control strategy is proposed and its effectiveness validated through simulations.

INTRODUCTION

Path following control deals with the confinement of the evolutions of a nonholonomic mechanical system to a prescribed geometric path. The assignment of time requirements on the way the system moves along the path represent an extra constraint which can be handled separately from the geometric requirement, allowing for better performances compared to trajectory tracking design [START_REF] Aguiar | Performance limitations in reference tracking and path following for nonlinear systems[END_REF]). Several approaches solving path following problems are discussed in the continuous-time literature. Among others, let us mention the work relying on energy shaping through Immersion and Invariance by [START_REF] Yi | Path following of a class of underactuated mechanical systems via immersion and invariance-based orbital stabilization[END_REF], sliding mode path following by [START_REF] Dagci | Path following controller design using sliding mode control theory[END_REF], smooth timevarying feedback and input scaling by De [START_REF] De Luca | Control of wheeled mobile robots: An experimental overview[END_REF] and, more relevant to this work, Transverse Feedback Linearization (TFL) by [START_REF] Banaszuk | Feedback linearization of transverse dynamics for periodic orbits[END_REF]; [START_REF] Altafini | Following a of curvature as an output regulation problem[END_REF]; [START_REF] Nielsen | On local transverse feedback linearization[END_REF]. Among these, the latter one aims to stabilize a general set which is made the zero dynamics sub-manifold associated with a suitable output function with well-defined relative degree. Accordingly, static (or possibly dynamic) feedback linearization is applied to stretch the trajectories onto the target set [START_REF] Nielsen | Output stabilization and maneuver regulation: A geometric approach[END_REF]). In this sense, path following admits a solution via TFL when the path can be implicitly represented by the equations which specify the zero-dynamics sub-manifold [START_REF] Akhtar | Path following using dynamic transverse feedback linearization for car-like[END_REF]). When approaching these problems in a digital context, one has to face well-known limitations due to the sampling process as, for instance, the loss of the relative degree and the rise of the unstable sampling zero-dynamics (see [START_REF] Aström | Zeros of sampled systems[END_REF]; [START_REF] Monaco | On the linearizing feedback in nonlinear sampled data control schemes[END_REF]) so directly affecting TFL in general. For those reasons, a first single rate solution preserving TFL in an approximate sense was proposed in [START_REF] Elobaid | Approximate transverse feedback linearization under digital control[END_REF] providing, under suitable assumptions, a solution to the path following problem as well. For instance, when dealing with the car-like robot approximate TFL under sampling is preserved for the kinematic model under a preliminary continuous-time dynamic extension [START_REF] Akhtar | Path following using dynamic transverse feedback linearization for car-like[END_REF]). Beyond TFL, a multirate technique has been proposed in Di [START_REF] Giamberardino | Piecewise continuous control for a car-like robot: implementation and experimental results[END_REF] to solve steering problems for mobile robots, when assuming a preliminary continuous-time feedback making the dynamics finitely discretizable. In the method we propose, both these demands are weakened. In this work, following [START_REF] Monaco | An introduction to motion planning under multirate digital control[END_REF], we propose a solution which circumvents the need for the dynamic extension, making use of a multi-rate sampled-data control strategy. Simulation results validate the proposed design approach in a comparative way with respect to the continuous-time solution and its direct implementation through Zero Order Holding (ZOH) referred to as emulation as well as the design approach previously recalled [START_REF] Elobaid | Approximate transverse feedback linearization under digital control[END_REF]). The proposed multi-rate solution provides efficient and improving results with respect to the preliminary ones.

The paper is organized as follows: Section II provides some background material and states the problem. The proposed control solution is developed in Section III in a constructive way. Simulations are discussed in Section IV. Concluding remarks end the manuscript.

Notations: S 1 denotes the unit disk, i.e. S 1 = {z ∈ C : |z| ≤ 1}. Z ≥0 denotes the set of non-negative integers. Given a pair of matrices (A, B), col(A, B) = A B while diag(A, B) is the block diagonal matrix with blocks A, B. The couple (A n , B n ) with

A n = 0 (n-1)×1 I n-1 0 0 1×(n-1) , B n = 0 (n-1)×1
1 is said to be in Brunovský's canonical form with I n-1 being the identity matrix of dimension n -1 and 0 m×n the zero matrix of dimension m × n. Given a manifold M and a closed connected set N ⊂ M , N is said to be invariant under the dynamics q = f (q) + g(q)u if for all q(0) ∈ N and any control u(•), q(t) ∈ N , ∀t. N is controlled invariant if there exists a feedback u making N invariant for the closed loop system. The point-to-set distance is denoted q 0 M = inf q -q 0 , q ∈ M . L f denotes the operator

L f = n i=1 f i (•) ∂ ∂qi , L f L g their composition.
Given a real valued function h(•) on R n , e L f h(q)| q(k) denotes the application of the Lie series operator e L f to the function h(q) evaluated at the state q(k). A continuous function R(x, δ) is of order O(δ p ) with p ≥ 1 if, whenever it is defined, it can be written as R(x, δ) = δ p-1 R(x, δ) and there exists a function β(δ) of class κ ∞ and δ > 0 such that ∀δ ≤ δ , | R(x, δ)| ≤ β(δ).

PRELIMINARIES AND PROBLEM STATEMENT

Path following for a car-like robot in continuous time

Consider the kinematic model of a car-like robot [START_REF] Siciliano | Robotics: modelling, planning and control[END_REF]) q = g 1 (q)v + g 2 (q)ω p = (x y)

with (see Figure 1)

g 1 (q) =      cos θ sin θ 1 tan φ 0      , g 2 (q) =    0 0 0 1    q = (q 1 q 2 q 3 q 4 ) = (x y θ φ) ∈ R 4 , v ∈ R
, the forward linear velocity, ω ∈ R, the angular velocity and the distance between the wheels. The position on the plane p = h(q) := (q 1 q 2 ) ∈ R 2 is the output of the system.

The desired path is given by a regular parameterized curve : D → R 2 with no self-intersections. Following [START_REF] Nielsen | Maneuver regulation via transverse feedback linearization: Theory and examples[END_REF]; [START_REF] Akhtar | Path following using dynamic transverse feedback linearization for car-like[END_REF], we refer to for a precise characterization, let the path (D) be an embedded sub-manifold of R 2 of dimension 1; i.e., there exists a function s(•) : R 2 → R such that 0 is a regular value of s and (D) = {w ∈ R 2 s.t. s(w) = 0}.

The path following problem for a car-like robot asks for the design of a feedback control law maneuvering the output of system (1) to approach and move along a given curve in a desired way. The problem is formally stated below.

Problem 2.1. Given a regular parameterized curve C = Im{ (D)}, find if possible, a smooth feedback law for system (1) such that for a set of some initial conditions X , with C ⊂ X the following holds true.

(1)

Invariance: if p(0) = h(q(0)) ∈ C then ∀t ≥ 0 p(t) C = 0. (2) Attractivity: system (1) under feedback is s.t ∀t ≥ 0, p C → 0 as t → ∞.
(3) Motion on the curve: system (1) traverses the curve C with a given desired velocity or acceleration profile ( πref (t), πref (t)).

It is known [START_REF] Altafini | Following a of curvature as an output regulation problem[END_REF]; [START_REF] Nielsen | Output stabilization and maneuver regulation: A geometric approach[END_REF]) that TFL is a natural tool to handle path following problems. To see this, denote the n = 3 dimensional sub-manifold Γ ⊆{q ∈ R 4 : q = (s • h) -1 (0)} as the control invariant subset for (1); i.e., the set of all initial Fig. 1. Kinematics of a car-like robot conditions q 0 ∈ R 4 under which the car-like robot is forced to remain on the curve (i.e., p(t) ∈ C for all t ≥ 0) under a suitably designed feedback. Consequently, Γ is referred to as the path following sub-manifold associated with the curve. With this in mind, Problem 2.1 is equivalent in the context of TFL to the problems of stabilizing Γ and zero dynamics assignment. The problem has been solved in [START_REF] Akhtar | Path following using dynamic transverse feedback linearization for car-like[END_REF] making use of continuoustime TFL through dynamic control defining the function

β : R 4 → R 2 β(q) = (α(q) π(q)) (2) with α(q) = (s • h)(q) ∈ R the transverse output function (i.e., the curve implicit function) and π(q) = arctan( q2 q1 ) is the tangent output function. Under the control dynamic extension η1 = η 2 , η2 = u 1 v = η 1 , ω = u 2 (3)
the feedback system q = f (q) + g1 (q)u 1 + g2 (q)u 2 f (q) = η 1 cos θ η 1 sin θ η 1 tan φ 0 η 2 0 g1 (q) = (0 0 0 0 0 1) , g2 (q) = (0 0 0 1 0 0) (4) with q = (q η 1 η 2 ) ∈ R ñ, ñ = 6 possesses strong vector relative degree r = (ñ -n n ) = (3 3) with respect to the output function (2). Accordingly, the following result is recalled.

Theorem 2.1. [START_REF] Akhtar | Path following using dynamic transverse feedback linearization for car-like[END_REF]). Given a regular, three times differentiable curve C in the plane with Γ the path following sub-manifold associated with C. Consider the dummy output function (2), then the Path-Following Problem 2.1 is solved by the feedback (3) with

u = γ(q, ν) = A -1 (q)(ν -B(q)) (5) with u = (u 1 , u 2 ) , A(q) = L g1 L 2 f α(q) L g2 L 2 f α(q) L g1 L 2 f π(q) L g2 L 2 f π(q) , B(q) = L 3 f α(q) L 3
f π(q) and ν = (ν 1 , ν 2 ) are respectively the transverse and tangential external stabilizing controls.

In the continuous-time case, dynamic extension is necessary for guaranteeing a well-defined relative degree ñ -n and, thus, TFL with respect to Γ via the input ω ∈ R n .

Remark 2.1. Setting the feedback (5), and the coordinates change φ(q) = (φ 1 (q) φ 2 (q)) with ξ = φ 1 (q) = (α(q) L f α(q) L 2 f α(q)) , z = φ 2 (q) = (π(q) L f π(q) L 2 f π(q)) the extended system (4) takes the form ξ

= A 3 ξ + B 3 ν 1 (6a) ż = A 3 z + B 3 ν 2 (6b)
where

A 3 ∈ R 3×3 , B 3 ∈ R 3 . From the result above, ξ ∈ R 3
is the transverse component to Γ while z describes the motion when the dynamics is restricted to Γ i.e. Γ = {(ξ z) : ξ = 0}. It is intuitively understood that Γ is precisely the zero dynamics sub-manifold of the system (4) with output α(q). Remark 2.2. Theorem 2.1 states that, under TFL with dynamics extension, ( 1) and ( 2) in Problem 2.1 are solved by ν 1 = -Kξ so to make (6a) asymptotically stable.

Denoting z = (z 1 z 2 z 3 ) , requirement (3) is solved setting ν 2 = -f 2 (z 3 -πref (t)) -f 1 (z 2 -πref (t)) + ... π ref (t) with f 1 , f 2 > 0 over the tangential dynamics (6b).

Path following for the car-like robot under sampling

In the sequel, we provide a solution to the path following problem for the car-like robot under sampled-data control and digital transverse feedback linearization as set in the following problem.

Problem 2.2. Design a digital control (v(k) ω(k)) = γ 3 δ (q(k), ν(k))
with external inputs ν solving the Path Following Problem 2.1 for the car-like robot (1) at all sampling instants t = kδ, k ≥ 0 and δ > 0 the sampling period.

In the present context, digital control design refers to design over the sampled-data equivalent model for which measures of the states are available at periodic sampling instants t = kδ, k ∈ Z ≥0 , where δ is the sampling period, and controls kept constant over δ. A first solution to Problem 2.2 can be carried out over the extended model (4) with inputs (u 1 , u 2 ) by directly applying the results in [START_REF] Elobaid | Approximate transverse feedback linearization under digital control[END_REF]. More in detail, let u i (t), i = 1, 2 be constant over the sampling period, i.e. u i (t) = u i (kδ) = u i (k) for t ∈ [kδ, (k + 1)δ[. Then the sampled-data model equivalent to (4) takes he form

q(k + 1) = F δ (q(k), u(k)) (7) with F δ (q, u) =e δ(L f +u1Lg 1 +u2Lg 2 ) q =q + j≥1 δ j j! (L f + u 1 L g1 + u 2 L g2 ) j q
with the function F δ (•, u) defined by its series expansion in powers of δ, [START_REF] Monaco | On nonlinear digital control[END_REF].

At this point, it is easily verified that the relative degree of the sampled-data model ( 7) with the output (2) falls to (1 1) for all δ > 0 [START_REF] Monaco | Minimumphase nonlinear discrete-time systems and feedback stabilization[END_REF]).

As a matter of facts, one gets that ∂α(q(k + 1))

∂u j (k) = δ ñ-n (ñ -n )! L gj L ñ-n -1 f α(q) q(k) + O(δ ñ-n +1 ) ∂π(q(k + 1)) ∂u j (k) = δ ñ n ! L gj L ñ -1 f π(q) q(k) + O(δ ñ +1 )
for n-n = 3 j = 1, 2 are non-zero (at least not simultaneously) by definition of the continuous-time relative degree.

Thus, the path following sub-manifold Γ is no longer the zero dynamics sub-manifold for the sampled-data model (7) with output α(q). In fact, the (typically unstable) zero dynamics of the sampled-data equivalent model evolves over a sub-manifold containing Γ ⊂ R 3 that is given by

Z SD = {q : α(q) = 0} ⊂ R 5 Z SD ⊃ Γ ={q : α(q) = L f α(q) = L 2 f α(q) = 0}
. Accordingly, transverse feedback linearization is lost for the sampled-data model ( 7) with respect to the output function (2).

In [START_REF] Elobaid | Approximate transverse feedback linearization under digital control[END_REF] TFL under sampling is achieved, up to a prescribed approximation order, by means of a δ-dependent dummy output computed from (2) as

α δ (q) = α(q) + δL f α(q) - δ 2 3 L 2 f α(q) + O(δ 3 ) (8a) π δ (q) = π(q) + δL f π(q) - δ 2 3 L 2 f π(q) + O(δ 3 ). ( 8b
)
The vector relative degree (3 3) is preserved under the dynamics (7) together with the zero dynamics sub-manifold Γ in O(δ 4 ), so that a digital solution can be computed to approximately solve the TFL problem under singlerate sampling. In fact, consider the coordinates change q → (ξ δ z δ ), with ξ δ = T 3 (δ)φ 1 (q), z δ = T 3 (δ)φ 2 (q), where φ 1 (•), φ 2 (•) are as in Remark 2.1, and T 3 (δ) as in [START_REF] Elobaid | Approximate transverse feedback linearization under digital control[END_REF]. This coordinates change, together with the piecewise continuous feedback for t ∈ [kδ, (k+1)δ[

η(t) = A 2 η(t) + B 2 u 1 (k) (9a) u(k) = A -1 (q(k))(ν δ (k) -B(q(k))) (9b
) recovers the TFL normal form in an approximate sense, with (v, ω) piecewise continuous. The control ν δ = -K δ col(ξ δ , z δ ), K δ : σ(I 3 +δ(A 3 +B 3 K δ )) ⊂ S 1 stabilizes the system to the path following sub-manifold and allows for solving Problem 2.2 in an approximate sense.

In the following we present an exact (and fully digital) solution for Problem 2.2 based on multi-rate sampling with no need of a preliminary continuous-time dynamic extension. Namely, we design a digital control (v(k) ω(k)) = γ 3 δ (q(k), ν(k)) for solving path following based on the sampled-data equivalent model to the kinematic car-like robot in (1); i.e., one gets for ω

(t) = ω(k) and v(t) = v(k) for t ∈ [kδ, (k + 1)δ[ the sampled-data equivalent kinemat- ics q(k + 1) = F δ (q(k), v(k), ω(k)) (10) with F δ (q, v, ω) =e δ(vLg 1 +ωLg 2 ) q = q + j≥1 δ j j!
(vL g1 + ωL g2 ) j q.

EXACT PATH FOLLOWING FOR A CAR-LIKE ROBOT UNDER MULTI-RATE SAMPLING

As shown in the previous section, dynamic extension is used over system (1) to guarantee v 2 L g2 L 2 g1 α(q) = 0, that is relative degree 3 with respect to the dummy output component in (2) associated with Γ (i.e., α(q)). However, the relative degree can also be guaranteed for the kinematic model (1) under multi-rate sampling without dynamic extension, as shown in [START_REF] Monaco | An introduction to motion planning under multirate digital control[END_REF]. To this end, we set in (10) a multi-rate of order n -n = 3 over the input ω that is δ = δ 3 and ω(t) = ω i (k) for t ∈ [kδ + (i -1) δ, kδ + i δ[ with i = 1, 2, 3 and v(t) = v(k) for t ∈ [kδ, (k + 1)δ[. Accordingly, the multi-rate model of (1) gets the form

q(k + 1) =F δ 3 (q(k), v(k), ω(k)) (11) with F δ 3 (q, v, ω) =F δ (•, v, ω 3 ) • F δ (•, v, ω 2 ) • F δ (q, v, ω 1 ) = j1,j2,j3≥0 δj1+j2+j3
j 1 !j 2 !j 3 ! (vL g1 + ω 1 L g2 ) j1 • (vL g1 + ω 2 L g2 ) j2 • (vL g1 + ω 3 L g2 ) j3 q.
In this respect, with reference to Problem 2.2, we seek for a digital piecewise constant control that preserves TFL with respect to Γ for the multi-rate equivalent model ( 11) based on the dummy output (2).

The problem is set for (11) based on the augmented output vector [START_REF] Monaco | An introduction to motion planning under multirate digital control[END_REF])

H(q) = (α(q) α(q) α(q)) ( 12) ensuring that the multi-rate model (11) possesses vector relative degree (1, 1, 1) so that transverse feedback linearizability under sampling is guaranteed. Accordingly, requirements (1),( 2) of Problem 2.1 at all sampling instants are satisfied by a digital feedback ω = ω δ (q) solution to the equality

H(F δ 3 (q, v, ω)) = A 3 δ H(q(k)) + B 3 δ ν 1 (k) ( 13
)
with A 3 δ = e 3 δA3 , B 3 δ = 3 δ 0 e τ A3 B 3 dτ , A 3 , B 3 as in (6a) and ν 1 the external stabilizing transverse control. More in detail, the TFL feedback solution to (13) is the one making the H(q) -ν 1 link in (11) linear. The following result asserts the existence of such feedback. Proposition 3.1. Consider the kinematic model of the carlike robot (1), and a regular parameterized curve : D → R 2 under the hypotheses of Theorem 2.1. Then, requirement (1) of Problem 2.1 is guaranteed, at all t = kδ with k ≥ 0, by the feedback ω = ω δ (q, v, ν 1 ) of the form

ω δ (q, v, ν 1 ) = ω 0 (q, v, ν 1 ) + i>0 δi (i + 1)! ω i (q, v, ν 1 ) (14)
defined as the unique solution to (13). In addition, requirement (2) of Problem 2.1 is guaranteed, at all t = kδ with k ≥ 0, setting

ν 1 (k) = -K δ H(q(k)) (15)
with

K δ such that σ(A 3 δ -B 3 δ K δ ) ⊂ S 1 .
Proof: Rewriting (13) as a formal series equality in powers of δ so getting

  e δ(v(k)Lg 1 +ω1(k)Lg 2 ) . . . e δ(v(k)Lg 1 +ω3(k)Lg 2 ) α(q) e δ(v(k)Lg 1 +ω1(k)Lg 2 ) . . . e δ(v(k)Lg 1 +ω3(k)Lg 2 ) α(q) e δ(v(k)Lg 1 +ω1(k)Lg 2 ) . . . e δ(v(k)Lg 1 +ω3(k)Lg 2 ) α(q)   = A 3 δ H(q(k)) + B 3 δ ν 1 (k) (16)
The equations ( 16) rewrite as S δ (q, ω, v, ν 1 ) = 0 with

S δ (q, ω, v, ν 1 ) (17) = δ3 S δ 1 (q, ω, v, ν 1 ) δ2 S δ 2 (q, ω, v, ν 1 ) δS δ 3 (q, ω, v, ν 1 ) and δ4-i S δ i = e δ(vLg 1 +ω1Lg 2 ) • . . . • e δ(vLg 1 +ω3Lg 2 ) α (i-1) (q) -α (i-1) (q) - 2 =i 3 δ ! α ( ) (q) - 3 4-i δ4-i (4 -i)! ν 1
Accordingly, because v(k) constant over the sampling interval, denoting f (q) = g 1 (q)v, one looks for ω satisfying ( 17), where each term can be written as S i (q, ω, v, ν 1 ) = j≥0 δj S i,j (q, ω, v, ν 1 ) with i = 1, 2, 3 and

S i,0 (q, ω, v, ν 1 ) = ∆ 4-i L g2 L 2 f α(q)ω + 1L 3 f α(q) -1ν 1 (k) where 1 = (1 1 1) , ∆ = col(∆ 3 ∆ 2 ∆ 1 ), and ∆ j = 1 j! j 4-j -(j -1) 4-j , (j -1) 4-j -(j -2) 4-j , 1
with j = 1, 2, 3. Following [START_REF] Monaco | On nonlinear digital control[END_REF] , it results that the matrix

∂ ∂ω S δ (q, v, ω)|δ →0 → ∆L g2 L 2 f α(q)
is full rank because L g2 L 2 f α(q) = 0 and ∆ is invertible. Hence, by the Implicit Function Theorem, the existence of ω δ unique solution to (13) of the form ( 14) can be deduced (Mattioni et al., 2017, Proposition 4.1). Under the coordinates transformation (ξ z) := q → φ(q) = (H(q) π(q)) , the controlled dynamics reads

ξ(k + 1) = A 3 δ ξ(k) + B 3 δ ν 1 (k) z(k + 1) = ψ(ξ(k), z(k), v(k), ν 1 (k))
where ψ(ξ, z, v, ν 1 ) = arctan q2(k+1) q1(k+1) | (ξ z) =φ -1 (q) bounded on the curve by definition and hence the result follows.

From the statement above, denoting πref (k) = πref (kδ) and πref (k) = πref (kδ), requirement (3) at all t = kδ is fulfilled by the discrete-time feedback

v(k + 1) =v(k) + δa(k) + δ 2 2 ν 2 (k) (18a)
a(k + 1) =a(k) + δν 2 (k) (18b) over the tangential component with v(k), a(k) are the linear velocity and acceleration over the path respectively and

ν 2 (k) = 1 δ 2 (a δ 0 -a δ 1 -1) 1 2δ (a δ 0 + a δ 1 -3) e(k) (19)
+ πref (k + 1) -πref (k) so to guarantee, asymptotically,

e(k) = v(k) -πref (k) a(k) -πref (k) → 0
for a δ 0 , a δ 1 given by a δ 0 = e -λ1δ e -λ2δ a δ 1 = -(e -λ1δ + e -λ2δ ), λ 1 , λ 2 > 0. Those arguments, together with Proposition 3.1, constitute the proof of the result below. Theorem 3.1. Consider the kinematic model of the car-like robot (1), and a regular parameterized curve : D → R 2 under the hypotheses of Theorem 2.1. Then, Problem 2.2 admits a solution under multi-rate control that is given by the feedback (v, ω) with: (i) ω = ω δ (q, v, ν 1 ) defined as in Proposition 3.1 with transverse control (15); (ii) v generated by the discrete dynamics (18) with tangential control (19).

Remark 3.1. It is important to stress the fact that, in the digital context, TFL is achieved working directly on the kinematic model of the car-like robot with no preliminary dynamic extension, contrarily to the continuous-time case.

The above discussion ascertains the intuitive expectation that whenever the path following problem is solvable in continuous time using TFL, a digital multi-rate feedback solution under sampling exists. As noted, the feedback component ω = ω δ (q, v, ν 1 ) soltion to (13) comes in the form of a series expansion in powers of δ. All terms of such an expansion ( 14) are computable through an iterative and constructive procedure solving, at each step, a linear equality in the corresponding unknown. For the first terms, one gets

ω 0 (q, v, ν 1 ) = 1ω(q, v, ν 1 ), ω(q, v, ν 1 ) = L 3 f α(q) -ν 1 L g2 L 2 f α(q) (20a) ω 1 (q, v, ν 1 ) = 15 4 - 15 2 87 4 ω(q, v, ν 1 ) (20b)
when denoting f (q) = g 1 (q)v and ω(q, v, ν 1 ) = (vL g1 + ω(q, v, ν 1 )L g2 )ω(q, v, ν 1 ).

Remark 3.2. From the expression below, it is clear that the digital feedback ω = ω δ (q, v, ν 1 ) solution to ( 13) is an expansion around the continuous-time (static) TFL solution when no velocity and acceleration requirements on the path are enforced (i.e., when only ( 1) and ( 2) in Problem 2.1 are given).

Because a closed-form of ( 14) cannot be computed, only approximations can be implemented in practice as the truncation of the series at a finite order ρ > 0, i.e., ω δ, [ρ] (q, v, ν 1 ) = ω 0 (q, v, ν 1 )+ ρ i=1 δi (i + 1)! ω i (q, v, ν 1 ). ( 21) Such approximate controllers ( 21) with ( 18)-( 19) solve Problem 2.2 in a practical sense with trajectories of the closed-loop system converging to a neighbourhood of target set Γ in O(δ p+1 ), [START_REF] Mattioni | Immersion and invariance stabilization of strictfeedback dynamics under sampling[END_REF].

SIMULATIONS

Consider the case where the car-like robot is required to follow a circle of radius r = 2. In this case; : D → R 2 , λ → (r cos q 3 r sin q 3 ) , satisfies the hypotheses of Theorem 2.1. Consequently, one has α(q) = q 2 1 + q 2 2 -r 2 . Suppose we are given a constant reference on the linear velocity πref (q) = 2 with zero acceleration (i.e., πref (q) = 0). To follow the path, let v(k) be as in ( 18)-( 19), with λ i = 5, i = 1, 2, and the first order approximate feedback be as in ( 21) with ρ = 1. In addition, for the transverse dynamics, let ν 1 (k) be as in ( 15) with K δ placing the poles of (A 3 δ + B 3 δ K δ ) in (e -0.317δ , e (-1.34+1.16i)δ) , e (-1.34-1.16i)δ) ) for δ = 3 δ. In Figures 2 and3, the initial condition is fixed as q 0 = (1 1 π 0) with increasing values of δ. The former one highlights that, albeit all controllers ensure convergence to the circle, the results under multi-rate control are slightly better compared to the other controllers. The latter one shows that as δ increases one notices deteriorated performances by both the approximate feedback in [START_REF] Elobaid | Approximate transverse feedback linearization under digital control[END_REF] as well as emulation based control with higher peak values for the control effort. In all cases, the multi-rate control ensures satisfactory performances with improved control effort even with respect to the continuous time control law.

In Figure 4, the initial condition is fixed as q 0 = (3 2 π 0) and δ = 0.3235s. Whereas similar comments hold true in this case for the emulation and the approximate single-rate controllers (which strongly rely on the continuous-time design), for such initialization (outside the circle) both the aforementioned controllers are significantly sensible to variations of δ, with no guarantee for convergence to the path under emulation. In this case, the multi-rate solution provides still notable performances with a desired velocity profile and limited control effort. Existence of a multi-rate digital feedback solution to the path following problem for a car-like robot via transverse feedback linearization was studied assuming a continuous time solution to the problem exists. In addition, this multirate solution was shown, through simulations, to provide better performances, both in terms of path-following position error and the required magnitude of input velocities to the robot. Perspectives concern the extension to the case of general nonlinear systems being transverse feedback linearizable in continuous time and the consequent application to further case studies.

  Simulation compare, for different initial conditions and increasing values of δ, the proposed multi-rate controller (21) with ρ = 1 (MR Sampling) with both the emulationbased (ZOH emulation) and the one proposed in Elobaid et al. (2020) (Approx Sampling), based on the continuoustime preliminary dynamic extension. For completeness, the results under continuous-time control are reported as well.

Fig. 2 .

 2 Fig. 2. From top to bottom: position on the plane, pathfollowing error, control effort for δ = 0.1s

Fig. 3 .

 3 Fig. 3. From top to bottom: position on the plane, pathfollowing error, control effort for δ = 0.3s

Partially funded by Université Franco-Italienne/Università Italo-Francese (UFI/UIF) through the Vinci program.


Mohamed Khalil Jabri
email: mohamed-khalil.jabri@imt-atlantique.fr

Robot Manipulation Learning Using Generative Adversarial Imitation Learning *

Imitation learning allows learning complex behaviors given demonstrations. Early approaches belonging to either Behavior Cloning or Inverse Reinforcement Learning were however of limited scalability to complex environments. A more promising approach termed as Generative Adversarial Imitation Learning tackles the imitation learning problem by drawing a connection with Generative Adversarial Networks. In this work, we advocate the use of this class of methods and investigate possible extensions by endowing them with global temporal consistency, in particular through a contrastive learning based approach.

Introduction

In [START_REF] Ho | Generative adversarial imitation learning[END_REF], Generative Adversarial Imitation Learning (GAIL) was proposed to tackle the Imitation Learning (IL) problem by drawing a connection with generative adversarial networks (GANs) [START_REF] Goodfellow | [END_REF]. GAIL allows to recover the demonstrator's behavior without explicitly recovering the reward function by training a policy to produce trajectories resembling the expert's demonstrations and a discriminator to distinguish them.

In this project, we focus on this class of IL algorithms and their application on an assistive robot manipulator to assist frail people. While GAIL provides a natural and direct way for this goal, it gives rise to a number of challenges. We believe we can take use of the advantages brought by this approach, and investigate possible ways to improve it. In the following sections, we discuss a specific limitation of the majority of works on adversarial imitation learning and present how we envision to mitigate it.

Research Questions

The sequential nature of IL requires learning a good representation that incorporates the underlying temporal structure governing the task. It is important to be able to identify and capture the relevant high-level temporal features that govern the demonstrator behaviour.

Indicatively, consider the task of picking an object. When the end-effector is too far from the object, the robot does not have to exactly match the expert trajectory all the way to the goal. Instead, it should only reason about features that shape the long-term outcome of the actions it takes. If the agent is too focused on temporally local features that govern the immediate outcome of its actions, it would be harder for it to recover when it deviates from the demonstrator's trajectory.

Model-based approaches allows effective long-term planning by learning a model of the agent's environment and using it to plan. However learning such a model is a challenging and data-intensive task. On the other hand, most works on model-free imitation learning do not take into account this aspect in the learning process, and this is even more true with adversarial imitation learning. The idea we are currently working on is a GAIL-based attempt to incorporate learning long-term temporal features within the model-free adversarial imitation learning setup. We envision to achieve that by using a a more complex architecture, while, at the same time, learning better representations.

Proposed Approach

Instead of the baseline GAIL, we are building our work on GoalGAIL [START_REF] Ding | [END_REF]. In GoalGAIL, the authors extend GAIL by conditioning the policy on the goal and using "relabeling" as introduced in Hindsight experience Replay (HER) [Andrychowicz et al., 2017] This allows to leverage failed attempts to improve the efficiency and speed up the learning. The learning is done by alternating two things: training the discriminator to distinguish the generated trajectories from real ones, and updating the policy using an offpolicy RL algorithm and the discriminator output as the reward.

We seek to extend GoalGAIL with two interrelated components that simultaneously tackle the discussed limitation. These extensions are inspired from recent works on GANs architectures, as well as recent works on contrastive selfsupervised learning that has proven effective in learning good data representations.

At this point, we finished reproducing the baseline on which we build our work as well as the first extension, and we are working on the second one and evaluating the performance on a set of simulated tasks. Video generation is a task that shares a particular aspect with the imitation learning problem in the sense they both need the inter-temporal dependencies to be taken into consideration, and in which GANs have proven effective. In [START_REF] Vougioukas | [END_REF], the authors introduce a GAN architecture capable of generating realistic facial animation by using two separate discriminators: one for the individual frames and another for the generated sequences. In our setting, we similarly consider an architecture composed of a policy network and two discriminators: a transition and a sequence discriminator, denoted D trans and D seq respectively in figure 1. The former impels the policy behaviour to be locally similar to the demonstrator's behavior while the latter encourages the policy trajectories to be globally similar to those of the expert.

For this architecture to be effective, the transitions and sequences representations need to be disentangled, meaning they should capture distinct features. This would allow the discriminators to distinguish between the immediate and the long-term effects of the actions.

Extension 2: Contrastive learning of temporally local and global representations

InfoMaxGAN [START_REF] Lee | [END_REF] is a GAN framework that mitigates GANs instability using contrastive learning of local and global visual features. On one hand, maximizing the mutual information between the local and global features helps the discriminator learn long-term representation in such a way as to reduce catastrophic forgetting within the non-stationary training environment. On the other hand, distinguishing local and global features by contrastive learning helps the generator produce more diverse images, which helps the generator avoid mode collapse.

We take inspiration from this work by trying to augment the proposed architecture with a contrastive learning task where the projected embeddings of sequences and transitions lying within them are pushed closer, while the projected representations of transitions lying within distinct sequences are pulled apart. As shown in figure 1, the discriminators D trans and D seq are fed with embeddings produced by two encoders: E trans and E seq respectively. Those embedding are passed through two projectors P trans and P seq that project them into a latent space, in which we compute the contrastive loss.

The contrastive loss has a double aim. First we want the temporally local and distant features learned by the transition and sequence encoders to be disentangled in such a way as to allow each discriminator to focus only on the features that are relevant to its task. Second, in analogy with InfoMax-GAN, we want the contrastive loss to stabilize the adversarial learning.

Figure 1 :

 1 Figure1: The trajectories produced by both the expert and the policy are stored in a memory from which transtions and sequences are sampled and relabeled to train two discriminators, while their embeddings are projected into a latent space to compute the contrastive loss

project ROGAN funded by the region of Brittany,


Olivier Marti
email: olivier.marti@lsce.ipsl.fr

Sébastien Nguyen

Pascale Braconnot

Sophie Valcke

Florian Lemarié

Eric Blayo

Diagnosing the ocean-atmosphere coupling schemes by using a mathematically consistent Schwarz iterative method

Introduction

Coupling algorithms implemented in coupled general circulation models (CGCMs) are driven by the necessity to conserve energy and mass at the air-sea interface. The discretization of the coupling problem leads to inconsistencies in time. It splits the simulation into small time intervals (coupling periods) over which averaged-in-time boundary data are exchanged. Atmosphere computes the fluxes at the interface (heat, water and momentum), and ocean computes the sea surface properties (water and sea ice temperatures, sea ice fraction, albedos, surface current). Two main algorithms are used, the parallel and the sequential atmospherefirst. In the parallel algorithm, both models run a coupling period with the boundary conditions of the other model from the previous coupling period. In the sequential atmosphere-first algorithm, atmosphere runs with oceanic boundary conditions from the previous coupling time. Then ocean runs with atmospheric boundary conditions of the current coupling period. Both algorithms are lagged: there is a time lag (of one coupling period) between the model and its boundary conditions, for one (sequential algorithm) or both (parallel algorithm) models.

Schwarz algorithms are attractive iterative coupling methods to cure the temporal inconsistencies and provide tightly and mathematically consistent coupled solutions. As discussed in [START_REF] Lemarié | Algorithmes de Schwarz et couplage océan-atmosphère[END_REF], the standard coupling methods correspond to one single iteration of a global-in-time iterative Schwarz method.

The Schwarz iterative method

A Schwarz iterative method has been implemented in the IPSLCM5A Earth System model [START_REF] Sepulchre | IPSL-CM5A2 -an Earth system model designed for multimillennial climate simulations[END_REF][START_REF] Marti | Key features of the IPSL ocean atmosphere model and its sensitivity to atmospheric resolution[END_REF][START_REF] Marti | A Schwarz iterative method to evaluate ocean-atmosphere coupling schemes. Implementation and diagnostics in IPSL-CM6-SW-VLR[END_REF]. The model repeats each coupling period, with the same initial state (the result from the previous coupling period). The boundary conditions are updated at each iteration, using the result of the previous iteration, until convergence. The solution is mathematically consistent: during a coupling period, each model uses the boundary conditions from the other model for the same coupling period. There is no lag.

The method has a huge computational cost, as each coupling period is iterated. It is not affordable for production runs, but allows us to compare the legacy coupling algorithms used in state-of-the art GCM to a mathematically consistent and synchronous algorithm, which is used as a reference to compute the error made with legacy algorithms.

We run 5-day experiments with the Schwarz iterative method, using three algorithms (parallel, sequential atmosphere-first and sequential ocean-first), and two coupling periods ∆t = 1h and ∆t = 4h. The Schwarz method can be used over all these algorithms. If the method is run until perfect convergence, the result is not dependant of the underlying coupling algorithm. In practice, the iterative method will be stopped before perfect convergence with an ad hoc criterion, and the different algorithms will give slightly different results.

Results

Fig. 1 shows the relative error in the change of sea surface temperature (SST) during one coupling period when the Schwarz method is not used. At each Schwarz iteration, the model computes an occurrence of the SST trend. At the first iteration, the trend is the one that the model would calculate with the legacy lagged coupling. It is compared with the trend obtained after convergence. The comparison is done on a unique trajectory of the model. The error is computed as the ratio between the correction due to the iterative procedure and the solution change with no Schwarz iteration. We consider each case of the iterative process, i.e. each point of the atmosphere grid, and for each coupling period.

In the parallel -∆t = 1h experiment, the relative error is negligible (less than 0.01) in about 15 % of the cases. It is small (less than 0.1) in almost 50 % of the cases. But it is larger than 0.1 for the other half. The relative error is even larger than 0.5 in 25 % of the cases. The atmosphere-first shows strongly improved results, with a negligible error for 97 % of the points. The results of the ocean-first experiment are close to the atmosphere-first experiments, with a slight deterioration. For the ∆t = 4h experiments, the errors are larger than in the ∆t = 1h case, but with the same hierarchy between the algorithms.

The error shows a strong diurnal cycle (not shown) with the lowest errors during the night. Errors are larger at noon than at midnight. The error is maximum after sunset and before sunrise, when the change of the insolation forcing evolves at the fastest pace.

We propose two hypotheses to explain the atmosphere-first algorithm performance compared to ocean-first. First, the atmosphere has shorter characteristic time scales than the ocean, with a more marked diurnal cycle. The atmospheric lower boundary conditions evolves slowly, and the atmospheric solution after the first half-iteration is then already quite close to its converged value, and provides a relevant and synchronized forcing to compute the oceanic solution in the second half-iteration. Second, the better performance of the atmosphere-first case can also be linked to the phasing of the solar radiation, which is the only external forcing and constrains the diurnal cycle. In the ocean-first case, the ocean is forced by fluxes, including solar radiation, calculated by the atmosphere at the previous coupling period. For atmosphere-first, the solar forcing is correctly phased.

Most current GCMs use the parallel algorithm which appears to have very large errors (see [START_REF] Marti | A Schwarz iterative method to evaluate ocean-atmosphere coupling schemes. Implementation and diagnostics in IPSL-CM6-SW-VLR[END_REF] for a short review). Our analysis shows that implementing sequential algorithms are simple ways to strongly reduce the error, with the atmosphere-first algorithm showing the best performance. The sequential algorithms, however, have a major drawback. The models do not run concurrently as, while one model is running, the other model waits for its boundary conditions. This eliminates a level of parallelism, and increases the time to solution of the coupled model.

Acknowledgements

This study is part of the ANR project COCOA (https: //anr.fr/Projet-ANR-16-CE01-0007). This work was granted access to the HPC resources of TGCC under an allocation made by GENCI (Grand équipement National de Calcul Intensif, grant 2019-A0040100239). It benefits from the development of the common modelling IPSL infrastructure coordinated by the IPSL climate modeling center (https://cmc.ipsl.fr).

Figure 1 :

 1 Figure1: Relative error of the change of sea surface temperature during a coupling period. The error is computed as the ratio between the correction due to the iterative procedure and the solution change during a coupling period with no Schwarz iteration. The ordinates show the number of cases in percentage of the total number of cases in time × space.


Patrice Bertail
email: bertail@ensae.fr

On constructing accurate confidence bands for ROC curves through smooth resampling

Keywords: ROC curve, confidence band, functional estimation, bootstrap, smoothing, asymptotic accuracy

This paper is devoted to thoroughly investigating how to bootstrap the ROC curve, a widely used visual tool for evaluating the accuracy of test/scoring statistics s(X) in the bipartite setup. The issue of confidence bands for the ROC curve is considered and a resampling procedure based on a smooth version of the empirical distribution called the "smoothed bootstrap" is introduced. Theoretical arguments and simulation results are presented to show that the "smoothed bootstrap" is preferable to a "naive" bootstrap in order to construct accurate confidence bands.

Introduction

Since the seminal contribution of [START_REF] Green | Signal detection theory and psychophysics[END_REF], so-called ROC curves (ROC standing for Receiving Operator Characteristic) have been extensively used in a wide variety of applications (anomaly detection in signal analysis, medical diagnosis, search engines, credit-risk screening) as a visual tool for evaluating the performance of a test statistic regarding its capacity of discrimination between two populations [START_REF] Fawcett | ROC graphs: Notes and practical considerations for data mining researchers[END_REF]. Whereas the statistical properties of their empirical counterparts have been only lately studied from the asymptotic angle [START_REF] Hsieh | Nonparametric and semi-parametric statistical estimation of the ROC curve[END_REF][START_REF] Girling | ROC confidence bands: An empirical evaluation[END_REF][START_REF] Ghosal | Bayesian ROC curve estimation under binormality using a partial likelihood based on ranks[END_REF][START_REF] Hall | Improved methods for bandwidth selection when estimating ROC curves[END_REF], ROC curves also have recently received much attention in the machine-learning literature through the development of statistical learning procedures tailored for the ranking problem [START_REF] Freund | An efficient boosting algorithm for combining preferences[END_REF][START_REF] Clémençon | Ranking and scoring using empirical risk minimization[END_REF]. The latter consists of determining, based on training data, a test statistic s(X) (also called a scoring function) with a ROC curve "as high as possible" at all points of the ROC space. Given a candidate s(X), it is thus of prime importance to assess its performance by computing a confidence band for the corresponding ROC curve, in a data-driven fashion preferably. Indeed, in such a functional setup, resampling-based procedures should naturally be preferred to those relying on computing/simulating the (gaussian) limiting distribution, as first observed in [START_REF] Macskassy | Confidence bands for ROC curves: methods and an empirical study[END_REF] and [START_REF] Macskassy | ROC confidence bands: An empirical evaluation[END_REF]2005a), where the use of the bootstrap is promoted for building confidence bands in the ROC space.

By building on recent works [START_REF] Hall | Nonparametric confidence intervals for receiver operating characteristic curves[END_REF][START_REF] Ghosal | Strong approximations for resample quantile process and application to ROC methodology[END_REF], it is the purpose of this paper to investigate how the bootstrap approach should be practically implemented based on a thorough analysis of the asymptotic properties of empirical ROC curves. Beyond the pointwise analysis developed in the studies mentioned above, here we tackle the problem from a functional angle, considering the entire ROC curve or parts of it. This viewpoint indeed appears as particularly relevant in scoring applications. Although the asymptotic results established in this paper are of a theoretical nature, they are considerably meaningful from a computational perspective. It turns out indeed that smoothing is the key ingredient for the boot-strap confidence band to be accurate, whereas a naive bootstrap approach would yield bands of low coverage probability in this case and should be consequently avoided by practicioners for analyzing ROC curves.

The rest of the paper is organized as follows. In Section 2, notations are first set out and certain key notions of ROC analysis are briefly recalled. The choice of an adequate (pseudo-) metric on the ROC space, a crucial point of the analysis, is also considered. The smoothed bootstrap algorithm is presented in Section 3, together with the theoretical results establishing its asymptotic accuracy. In particular, the gain in convergence rate acquired by the smoothing step is thoroughly discussed. This fact is also supported by strong empirical evidence: Section 4 is devoted to displaying preliminary numerical results, illustrating the impact of smoothing on the bootstrap performance. Eventually, technical details are postponed to the Appendix.

Background

Here we briefly recall basic concepts of the bipartite ranking problem as well as key results related to the statistical estimation of ROC curves and set out the notations that shall be needed throughout the paper. Although the results contained in this paper can be formulated without referring to the bipartite ranking framework, in the purpose of motivating the present analysis we intentionally connected them to this major statistical learning problem, which has recently revitalized the interest for the problem of assessing the accuracy of empirical ROC curves [START_REF] Cohen | Learning to order things[END_REF]Clémençon et al., 2007).

Assumptions and notation

In the bipartite ranking problem, the matter is to order all the elements X of a set X by degree of relevance, when relevancy may be observed through some binary indicator variable Y. Precisely, one has a system consisting of a binary random output Y, taking its values in {-1, 1} say, and a random input X, taking its values in a (generally high-dimensional) feature space X , which models some observation for predicting Y. The probabilistic setup is the same as for standard binary classification, except one is concerned, borrowing the jargon of information retrieval applications for instance, by ordering all the documents x of the list X by degree of relevance for a particular request, rather than simply classifying them as relevant or not. This amounts to assign to each document x in X a score s(x) indicating its degree of relevance for this specific query. The challenge is thus to build a scoring func-tion s : X → R from sampling data, so as to rank the observations x by increasing order of their score s(x) as accurately as possible: the higher the score s(X) is, the more likely one should observe Y = +1.

True ROC curves. A standard way of measuring the ranking performance consists of plotting the ROC curve, namely the graph of the mapping

ROC s : α ∈ (0, 1) → 1 -(G s • H -1 s )(1 -α), (1)
where G s (respectively H s ) denotes s(X)'s cdf conditioned on Y = +1 (resp. conditioned on Y = -1) and

F -1 (α) = inf{x ∈ R/ F(x
(α) = 1 -G s • H -1 s (1 -α) (2)
from empirical cdf estimates:

G s (x) = 1 n+ n i=1 I {Y i =+1} K(x -s(X i )) H s (x) = 1 n- n i=1 I {Y i =-1} K(x -s(X i ))
,

where n + = n i=1 I{Y i = +1} = n -n -is the (ran- dom)
number of positive instances among the sample (distributed as the binomial Bin(n, p)) and K(u) denotes the step function I {u≥0} . In order to obtain smoothed versions G s (x) and F s (x) of the latter cdf's, a typical choice consists of picking instead a function

K(u) of the form v≥0 K h (u -v)dv, with K h (u) = h -1 K(h -1 • u)
where K ≥ 0 is a regularizing Parzen-Rosenblatt kernel (i.e. a bounded square integrable function such that K(v)dv = 1) and h > 0 is the smoothing bandwidth. Here and throughout, I {A} denotes the indicator function of any event A.

Metrics on the ROC space. When it comes to measure closeness between curves in the ROC space, various metrics may be used [START_REF] Flach | The geometry of roc space: understanding machine learning metrics through roc isometrics[END_REF]. Viewing the ROC space as a subset of the Skorohod's space D([0, 1]) 1 , essential cases include the standard metric induced by the sup norm ||.|| ∞ and the Hausdorff distance 2 d H . If the former choice appears natural when considering continuous curves, the latter is more pertinent for analyzing stepwise graphs, such as empirical ROC curves. Indeed, equipped with this metric, two piecewise constant ROC curves may be close to each other, even if their jumps do not exactly match. In particular, the topology induced by d H (which is finer than the one related to the sup norm since

d H (f 1 , f 2 ) ≤ ||f1 -f 2 || ∞ )
is much more appropriate for describing the fluctuations of the empirical ROC curve (and the deviation between the latter and its bootstrap counterpart as well). As shall be seen below, asymptotic arguments for grounding the bootstrapping of the empirical ROC curve fluctuations, when measured in terms of the sup norm ||.|| ∞ , are rather straightforward. However, given the geometry of empirical ROC curves, this metric is not convenient for our purpose and produce very wide, and thus noninformative, bands. Besides, providing a theoretical basis in the case of the Hausdorff distance is very challenging and will not be addressed in this paper, owing to space limitations. Instead, we shall consider the closely related pseudometric defined as follows:

∀(f 1 , f 2 ) ∈ D([0, 1]) 2 , d B (f 1 , f 2 ) sup t∈[0,1] d B (f 1 , f 2 ; t),
where d B (f 1 , f 2 ; t) denotes the minimum among these three quantities

|f 1 (t) -f 2 (t)|, |f -1 2 • f 1 (t) -t|, and |f -1 1 • f 2 (t) -t|.
We clearly have

d H (f 1 , f 2 ) ≤ d B (f 1 , f 2 ) ≤ ||f 1 -f 2 || ∞ .
The major advantage of considering this pseudo-metric is that it provides a control on vertical and horizontal jumps of ROC curves both at the same time, treating both types of error in a symmetric fashion. This way, d B permits to contruct builds bands of reasonable size, well adapted to the stepwise shape of empirical ROC curves, with better coverage probabilities.

1 Recall that D([0, 1]) is the set of càd-làg functions f : [0, 1] → R, i.e. functions f such that limu→x+ f(u) = f(x) and lim sup u→ x-f(u) < ∞ for all x ∈ (0, 1)
2 The Hausdorff distance between two elements f 1 and f 2 in D([0, 1]) is simply the distance between the completed graphs, obtained by connecting (x,

f i (x-)) to (x, f i (x)) with a line segment at any jump point x ∈ [0, 1]. If f i is a stepwise constant function with jump points {α (i) k } for i ∈ {1, 2}, we have d H (f 1 , f 2 ) = max i min j m i,j , where m i,j = min{|α (1) i -α (2) j |, |f 1 (α (1) i ) -f 2 (α (2) j )|}.
As the goal pursued in the present paper is to build, in the ROC space viewed as a subspace of the Skorohod's space D([0, 1]) equipped with a proper (pseudo-) metric, a confidence band for the ROC curve of a given diagnosis test statistic s(X), we shall omit to index by s the quantities considered and denote by Z the r.v. s(X) (and by Z i , 1 ≤ i ≤ n, the s(X i )'s) for notational simplicity. Throughout the paper, we assume that H(dx) and G(dx) are continuous probability distributions, with densities h(x) and g(x) respectively. Eventually, denote by P the joint distribution of (Z, Y) on R × {-1, +1} and by P n its empirical version based on the sample

D n = {(Z i , Y i )} 1≤i≤n .
Equipped with the notations above, one may write P(dz, y) = pI {y=+1} H(dz) + (1p)I {y=-1} G(dz).

Asymptotic law -Gaussian approximation

in the situation described above, the next theorem establishes the strong consistency of the empirical ROC curve (2) in sup norm and provides a strong approximation at the rate 1/ √ n, up to logarithmic factors, for the fluctuation process:

r n (α) = √ n( ROC n (α) -ROC(α)), α ∈ [0, 1]. (3)
This (gaussian) approximation plays a crucial role in understanding the asymptotic behavior of the empirical ROC curve and of its bootstrap counterpart. The following technical assumptions are also required.

H 1 The slope of the ROC curve is bounded on [0, 1]:

sup α∈[0,1] g(H -1 (α)) h(H -1 (α)) < ∞. (4)
H 2 H is twice differentiable on [0, 1]. Furthermore, ∀α ∈ [0, 1], h(α) > 0, (5)
and there exists γ > 0 such that sup

α∈[0,1] α(1 -α) d log(h • H -1 (α)) dα ≤ γ < ∞. (6)
Theorem 1 (Functional limit theorem) Suppose that assumptions H 1 -H 2 are fulfilled. Then, (i) the empirical ROC curve is strongly consistent

sup α∈[0,1] | ROC n (α) -ROC(α)| → 0 a.s. as n → ∞,
(ii) there exist a sequence of two independent brownian bridges {(B

(n) 1 (α), B (n) 2 (α))} α∈[0,1] such that
we almost surely have, uniformly over [0, 1],

r n (α) = z (n) (α) + o( (log log n) ρ 1 (γ) log ρ 2 (γ) n √ n ), (7)
where

z (n) (α) = (1 -p) -1/2 g(H -1 (1 -α)) h(H -1 (1 -α)) B (n) 1 (α) + p -1/2 B (n) 2 (ROC(α)). and    ρ 1 (γ) = 0, ρ 2 (γ) = 1, if γ < 1 ρ 1 (γ) = 0, ρ 2 (γ) = 2, if γ = 1 ρ 1 (γ) = γ, ρ 2 (γ) = γ -1 + ε, ε > 0, if γ > 1
.

These results may be immediately derived from classical strong approximations for the empirical and quantile processes (see [START_REF] Csorgo | Strong approximations in probability and statistics[END_REF], and also [START_REF] Hsieh | Nonparametric and semi-parametric statistical estimation of the ROC curve[END_REF]). Incidentally, we mention that the approximation rate is not always log 2 (n)/ √ n, contrarily to what is claimed in the paper by [START_REF] Hsieh | Nonparametric and semi-parametric statistical estimation of the ROC curve[END_REF].

Remark 1 (On assumptions H 1 -H 2 ) Hypothesis H 2 is standard for obtaining strong approximation of the quantile process (see [START_REF] Csorgo | Strong approximations of the quantile process[END_REF] or [START_REF] Csorgo | Strong approximations in probability and statistics[END_REF]). It is possible to weaken conditions (4)-( 6) in a way that they are required to hold on any compact subinterval of ]0, 1[ solely. However, in this case, one can only guarantee (ii) to be valid uniformly over compact subintervals of ]0, 1[.

Bootstrapping empirical ROC curves

Beyond consistency of the empirical curve in sup norm and the asymptotic normality of the fluctuation process, we now tackle the question of constructing confidence bands for the true ROC curve via the bootstrap approach introduced by (Efron, 1979) (extending pointwise results established in [START_REF] Hall | Nonparametric confidence intervals for receiver operating characteristic curves[END_REF]). The latter suggests to consider, as an estimate of the law of the fluctuation process r n = {r n (α)} α∈[0,1] , the conditional law given D n of the bootstrapped fluctuation process

r * n = { √ n(ROC * (α) -ROC(α))} α∈[0,1] , (8)
where ROC * is the ROC curve corresponding to a sample

D * n = {(Z * i , Y * i )}

The Algorithm

Here we describe the algorithm for building a confidence band at level 1ǫ in the ROC space from sampling data

D n = {(Z i , Y i ); 1 ≤ i ≤ n}. Set n + = 1≤i≤n I {Y i =1} = n -n -.
It is performed in four steps as follows.

Before turning to the theoretical properties of this algorithm and related numerical experiments, a few remarks are in order.

Remark 2 (Monte-Carlo approximation) From a computational angle, the true smoothed bootstrap distribution must be approximated in its turn, using a Monte-Carlo approximation scheme. A convenient way of doing this in practice, while reproducing theoretical advantages of smoothing, consists of drawing B bootstrap samples, of size n, with replacement in the original data and then perturbating each drawn data by independent centered gaussian random variables of variance h 2 (this procedure is equivalent to drawing bootstrap data from a smooth estimate P n (dz, dy) computed using a gaussian kernel K h (u) = (2πh 2 ) -1/2 exp(-u 2 /(2h 2 ))), see [START_REF] Silverman | The bootstrap: to smooth or not to smooth?[END_REF]. Regarding the choice of the number of bootstrap replications, picking B = n does not modify the rate of convergence. However, choosing B of magnitude comparable to n so that (1 + B)ǫ is an Algorithm 1 smoothed ROC bootstrap 1. Based on D n , compute the empirical class cdf estimates G and H, as well as their smoothed versions G and H. Plot the ROC curve estimate:

ROC(α) = 1 -G • H -1 (1 -α), α ∈ [0, 1].
2. From the smooth distribution estimate

P n (dz, y) = n - n I {y=+1} H(dz)+ n + n I {y=-1} G(dz), draw a bootstrap sample D * n = {(Z * i , Y * i )} 1≤i≤n conditioned on D n .

Based on D *

n , compute the bootstrap versions of the empirical class cdf estimates G * and H * . Plot the bootstrap ROC curve

ROC * (α) = 1 -G * • H * -1 (1 -α), α ∈ [0, 1].
4. Eventually, get the bootstrap confidence bands at level 1ǫ defined by the ball of center ROC and radius

δ ǫ / √ n in D([0, 1])
, where δ ǫ is defined by integer may be more appropriate: the ǫ-quantile of the approximate bootstrap distribution is the uniquely defined and this will not modify the rate of convergence neither [START_REF] Hall | On the number of bootstrap simulations required to construct a confidence interval[END_REF].

P * (||r * n || ∞ ≤ δ ǫ ) = 1 -ǫ in
Remark 3 (On tuning parameters) The primary parameters of Algorithm 1 are those related to the smoothing stage. When using a gaussian regularizing kernel, one should typically choose a bandwidth h n of order n -1/5 in order to minimize the mean square error.

Remark 4 (On recentering) From the asymptotic analysis viewpoint, it would be fairly equivalent to recenter by a smoothed version of the original empirical curve ROC(.) = 1 -G • H -1 (1 -.) in the computation of the bootstrap fluctuation process. However, numerically speaking, computing the sup norm of the estimate ( 8) is much more tractable, insofar as it solely requires to evaluate the distance between piecewise constant curves over the pooled set of jump points. It should also be noticed that smoothing the original curve, as proposed in [START_REF] Hall | Nonparametric confidence intervals for receiver operating characteristic curves[END_REF], should be also avoided in practice, since it hides the jump locations, which constitute the essential part of the information.

Remark 5 (Subsampling) A possible alternative for improving the accuracy of the bootstrap approximation could consist of using subsampling, namely drawing m bootstrap data without replacement in the original data, with m much smaller than the original sample size n [START_REF] Arcones | On the asymptotic accuracy of the bootstrap under arbitrary resampling size[END_REF]. Due to space limitations, we leave aside this question for further research.

Asymptotic analysis

We now investigate the accuracy of the bootstrap estimate output by Algorithm 1. The result stated in the next theorem extend those established in [START_REF] Hall | Nonparametric confidence intervals for receiver operating characteristic curves[END_REF] in the pointwise framework. The functional nature of the approximation result below is essential, since it should be enhanced that, in most ranking applications, assessing the uncertainty about the whole estimated ROC curve, or some part of it at least, is what really matters. In the sequel, we assume that the kernel K used in the smoothing step is "pyramidal"3 (e.g. gaussian or of the form I {u∈[-1,+1]} ).

Theorem 2 (Asymptotic accuracy) Suppose that the hypotheses of Theorem 1 are fulfilled. Assume further that smoothed versions of the cdf 's G and H are computed at step 1 using a scaled kernel K hn (u) with h n ↓ 0 as n → ∞ in a way that nh 3 n → ∞ and nh 5 n log 2 n → 0. Then, the bootstrap distribution estimates output by Algorithm 1 are such that

sup t∈R |P * (||r * n || ∞ ≤ t) -P(||r n || ∞ ≤ t)| = o P ( log(h -1 n ) √ nh n ), sup t∈R |P * (d * n ≤ t) -P(d n ≤ t)| = o P ( log(h -1 n ) √ nh n ).
Hence, up to logarithmic factors, choosing h n ∼ 1/(log 2+η n 1/5 ) with η > 0 yields an approximation error of order n -2/5 for the bootstrap estimate. Although its rate is slower than the one of the gaussian approximation (7), the smoothed bootstrap method remains very appealing from a computational perspective, the construction of confidence bands from simulated brownian bridges being very difficult to implement in practice. As shall be seen below, the rate reached by the smoothed bootstrap distribution is nevertheless a great improvement, compared to the naive bootstrap approach (see the discussion below).

Remark 6 (Bootstrapping summary statistics) From Theorem 1 above, asymptotic validity of the smooth bootstrap method for estimating the distribution of the fluctuations of a functional Φ( ROC) of the empirical ROC curve may be deduced, as soon as the function Φ defined on D([0, 1]) is sufficiently smooth (namely continuously Hadamard differentiable). For instance, it could be applied to summary statistics involving a specific piece of the ROC curve only in order to focus on the "best instances" (Clémençon & Vayatis, 2007), or more classically to the area under the ROC curve (AUC). However, in the latter case, due to the fact that this particular summary statistic is of the form of a U-statistic (Clémençon et al., 2007), the naive bootstrap rate is faster than the one we obtained here (of order n -1 ).

Discussion

Let us now give an insight into the reason why the smoothed bootstrap procedure outperforms the nonsmoothed one. In most statistical problems where the nonparametric bootstrap is useful, there is no particular reason for implementing it from a smoothed version of the empirical df rather from the raw empirical distribution itself [START_REF] Silverman | The bootstrap: to smooth or not to smooth?[END_REF]. However, in the present case, smoothing affects the rate of convergence.

Theorem 3 (Pointwise convergence rate) Suppose that the bootstrap process ( 8) is built by drawing from the raw cdf 's G and H instead of their smoothed versions at step 2 of Algorithm 1. Then, for any α ∈]0, 1[,

sup t∈R |P * (r * n (α) ≤ t) -P(r n (α) ≤ t)| = O P (n -1/4 ). (9)
Hence, the naive bootstrap rate induces an error of order O(n -1/4 ) cannot be improved, whereas it may be shown that the rate n -2/5 is attained by the smoothed bootstrap (in a similar fashion to the functional setup), provided that the amount of smoothing is properly chosen. Heuristically, this is a consequence of the oscillation behavior of the deviation between the bootstrap quantile H * -1 (1α) and its expected value H -1 (1α) given the data D n , due to the fact that the step cdf H is not regular around H -1 (1α): this corresponds to a jump with probability one.

Higher-order accuracy. A classical way of improving the pointwise approximation rate (9) consists of bootstrapping a standardized version of the r.v. r n (α).

It is natural to consider, as standardization factor, the square root of an estimate of the asymptotic variance:

σ 2 (α) = var(z (n) (α)) = α(1 -α) 1 -p g(H -1 (1 -α)) 2 h(H -1 (1 -α)) 2 + ROC(α)(1 -ROC(α)) p . ( 10
)
An estimate σ 2 n of plug-in type could be considered, obtained by plugging n + /n, ROC and smoothed density estimators h = H ′ and g = G ′ into (10) instead of their (unknown) theoretical counterparts. More interestingly, from a computational viewpoint, a bootstrap estimator of the variance could also be used. Following the argument used in [START_REF] Hall | Nonparametric confidence intervals for receiver operating characteristic curves[END_REF] for a smoothed original estimate of the ROC curve, one may show that a smoothed bootstrap of the studentized statistic r n (α)/σ n (α) yields a better pointwise rate of convergence than 1/ √ n, the one of the gaussian approximation in the Central Limit Theorem. Precisely, for a given α ∈]0, 1[, if the bandwidth used in the computation of σ 2 n (α) is chosen of order n -1/3 , we have: n (α). Notice that the bandwidth used in the standardization step (i.e. for estimating the variance) is not the same as the one used at the resampling stage of the procedure. This is a key point for achieving second-order accuracy. This time, the smoothed (studentized) bootstrap method widely outperforms the gaussian approach, when the matter is to build confidence intervals for the ordinate ROC(α) of a point of abciss α on the empirical ROC curve. However, it is not clear yet, whether this result remains true for confidence bands, when considering the whole ROC curve (this would actually require to establish an Edgeworth expansion for the supremum ||r n / σ n || ∞ ). This will be the scope of further research.

Simulation results

The striking advantage of the smoothed bootstrap is the improved rate of convergence of the resulting estimator. Furthermore, choosing d B for measuring the magnitude order of curve fluctuations has an even larger impact on the accuracy of the empirical bands.

  the case of the sup norm or by P * (d * n ≤ δ ǫ ) = 1-ǫ, when considering the d B distance, denoting by P * (.) the conditional probability given the original data D n .

  's bootstrap counterpart (i.e. the variance estimator computed from bootstrap data) by σ * 2

  The difficulty is twofold. Firstly, the target of the bootstrap procedure is here a distribution on a path space, the ROC space being viewed as a subspace of D n ([0, 1]), equipped with either ||.|| ∞ ) or else d B (., .). Secondly, both r n and d n are functionals of the quantile process { H -1 (α)} α∈[0,1] . It is well-known that the naive bootstrap (i.e. resampling from the raw empirical distribution) generally provides bad approximations of the distribution of empirical quantiles in practice: the rate of convergence for a given quantile is indeed of order O P (n -1/4 )[START_REF] Falk | Weak convergence of smoothed and nonsmoothed bootstrap quantile estimates[END_REF], whereas the rate of the gaussian approximation is n -1/2 . As shall be seen below (see §3.3), the same phenomenon may be naturally observed for ROC curves. In a similar fashion to what is generally recommended for empirical quantiles, we suggest to implement a smoothed version of the bootstrap algorithm in order to improve the approximation rate of ||r n || ∞ 's distribution, respectively of d n 's distribution . In short, this boils down to resampling the data from a smoothed version of the empirical distribution P n .

1≤i≤n of i.i.d. random pairs with a common distribution P n close to P n . We shall also consider d * n = √ nd B (ROC * , ROC), whose random fluctuations (given D n ) are expected to mimic those of d n = √ nd B ( ROC, ROC).

According to the hypothesis introduced in[START_REF] Giné | Weighted uniform consistency of kernel density estimators[END_REF], a kernel K is pyramidal if in the linear span of functions k ≥ 0 with a subgraph that can be represented as a finite number of Boolean operations among sets of the form {(s, u) : P(s, u) ≥ φ(u)}, where P is a polynomial on R 2 and φ is an arbitrary real function.

As an illustration of this theoretical result, we now display simulation results, emphasizing the gain acquired by smoothing and considering the pseudo-metric d B .

We present confidence bands for a single trajectory and the estimation of the coverage probability of the bands for a simple binormal model :

where ε and X are independent standard normal r.v.'s.

In this example, the scoring function s(x) is the maximum likelihood estimator of the probit model on the training set. We choose here β 0 = β 1 = 1, n = 1000, B = 999 and γ = 0.95 for the targeted coverage probability. Coverage probabilities are obtained over 2000 replications of the procedure, using the package ROCR of statistical software R. As mentioned before, choosing ||.|| ∞ yields very large bands with coverage probability close to 1! Though still large, bands based on the pseudo-metric d B are clearly much more informative (see Fig. 1). It should be noticed that the coverage improvement obtained by smoothing is clearer in the pontwise estimation setup (here α = 0.2) but more difficult to evidence for confidence bands.

Appendix -Technical details

Theorem 1 is a slight variant of Theorems 2.1 and 2.2 in [START_REF] Hsieh | Nonparametric and semi-parametric statistical estimation of the ROC curve[END_REF]. The proof mainly relies on r n (α)'s Taylor expansion (holding true almost surely, uniformly over

), combined with standard results in strong approximation theory [START_REF] Csorgo | Strong approximations in probability and statistics[END_REF].

Turning to Theorem 2, the proof follows from the strong approximation ( 7), combined with a standard coupling argument. Here we give a brief sketch of the latter. Note firstly that, conditioned on D n , by applying Theorem 1's second part with ROC as a target curve (instead of ROC), one gets that r * n (α) is, uniformly over [0, 1], almost surely equivalent to

with a remainder of order o(

under the stipulated conditions [START_REF] Giné | Rates of strong uniform consistency for multivariate kernel density estimators[END_REF]. Furthermore, from standard result on the modulus of continuity of the brownian bridge [START_REF] Shorack | Empirical processes with applications to statistics[END_REF], we almost surely have

Applying then the LIL to (n -/n) and (n -/n), it follows that almost surely, uniformly in α,

Since this results holds uniformly in α, by the continuous mapping theorem applied to the function sup α∈[0,1] (.), we also get the result in term of distribution up to the given almost sure order.

Eq. ( 11) is a direct consequence of the results established in Falk & Reiss (1994) applied to the second term in the expansion. The studentized version may be obtained by following the arguments of [START_REF] Hall | Nonparametric confidence intervals for receiver operating characteristic curves[END_REF] with the minor modification that the ROC curve itself is not smoothed.