Eva Löcherbach
email: eva.loecherbach@u-cergy.fr

Victor Rabiet
email: victor.rabiet@ens-cachan.fr

Ergodicity for multidimensional jump diffusions with position dependent jump rate

Keywords: Diffusions with jumps, Harris recurrence, Nummelin splitting, continuous time Markov processes, Additive functionals

We consider a jump type diffusion X = (X t ) t with infinitesimal generator given by

where µ is of infinite total mass. We prove Harris recurrence of X using a regeneration scheme which is entirely based on the jumps of the process. Moreover we state explicit conditions in terms of the coefficients of the process allowing to control the speed of convergence to equilibrium in terms of deviation inequalities for integrable additive functionals.

Introduction

Let N (ds, dz, du) be a Poisson random measure on R + × R d × R + , defined on a probability space (Ω, A, P ) with intensity measure dsµ(dz) du, where µ is a σ-finite measure on R d having infinite total mass. We consider a process X = (X t ) t≥0 , X t ∈ R d , solution of

X t = x+ t 0 g(X s ) ds+ t 0 σ(X s ) dW s + [0,t] R d ×R +
c(z, X s-)1 u≤γ(z,X s-) N (ds, dz, du), (1.1) x ∈ R d , where W is an m-dimensional Brownian motion. The associated infinitesimal generator is given for smooth test functions by

Lψ(x) = 1 2
1≤i,j≤d

a ij (x) ∂ 2 ∂x i ∂x j ψ(x) + g(x)∇ψ(x) + R d (ψ(x + c(z, x)) -ψ(x))γ(z, x)µ(dz) (1.2)
where a = σσ * . Notice that the jump rate at time t of process depends on the position of the process X t itself, i.e. the intensity measure in the infinitesimal operator L is γ(z, x)µ(dz). Moreover, since µ has infinite total mass, jumps occur with infinite activity, i.e. the process possesses infinitely many small jumps during any finite time interval [0, T ].

The principal aim of the present paper is to give easily verifiable conditions on the coefficients g, σ, c and γ under which the process is recurrent in the sense of Harris and satisfies the ergodic theorem, without imposing any non-degeneracy condition on the diffusive part. Recall that a process X is called recurrent in the sense of Harris if it possesses an invariant measure m such that any set A with m(A) > 0 is visited infinitely often almost surely (see Azéma, Duflo and Revuz [START_REF] Azéma | Mesure invariante des processus de Markov récurrents[END_REF] (1969)): For all x ∈ R d ,

P x ∞ 0 1 A (X s ) ds = ∞ = 1.
For classical jump diffusions there starts to be a huge literature on this subject, see e.g. Masuda [START_REF] Masuda | Ergodicity and exponential β-mixing bounds for multidimensional diffusions with jumps[END_REF] (2007) who works in a simpler situation where the "censure" (or "intensity") term 1 u≤γ(z,X s-) is not present. In order to prove recurrence, Masuda follows the Meyn and Tweedie approach developed in [START_REF] Meyn | Markov Chains and Stochastic Stability[END_REF] or [START_REF] Meyn | Stability of Markovian processes III : Foster-Lyapunov criteria for continuous-time processes[END_REF]. Kulik [START_REF] Alexey | Exponential ergodicity of the solutions to SDE's with a jump noise[END_REF] (2009) uses the stratification method in order to prove exponential ergodicity of jump diffusions, but the models he considers do not include the censured situation neither. Finally, let us mention Duan and Qiao [START_REF] Duan | Stationary measures for stochastic differential equations with jumps[END_REF] (2014) who are interested in solutions driven by non-Lipschitz coefficients.

On the contrary to the above mentioned papers, in our model, jumps occur at a given intensity depending both on the current position of the process and on an additional noise z. The presence of this intensity term γ(z, X s-) in (1.1) is in fact quite natural, but it implies that the study of X is technically more involved than the non-censured situation when γ is lower-bounded and strictly positive.

The aim of the present paper is to show that the jumps themselves can be used in order to generate a splitting scheme which implies the recurrence of the process. The method we use is the so-called regeneration method which we apply to the big jumps. More precisely, for some suitable set E such that µ(E) < ∞, we cut the trajectory of X into excursions between successive jumps appearing due to "noise" z belonging to E. In spirit of the splitting technique introduced by Nummelin [START_REF] Nummelin | A splitting technique for Harris recurrent Markov chains[END_REF] (1978) and Athreya and Ney [START_REF] Athreya | A new approach to the limit theory of recurrent Markov chains[END_REF] (1978), we state a non-degeneracy condition which guarantees that the jump operator associated to the big jumps possesses a Lebesgue absolutely continuous component. This amounts to imposing that the partial derivatives of the jump term c with respect to the noise z are sufficiently non-degenerate, see (2.7) and (2.8) below. We stress that we do not need any non-degeneracy condition for the diffusion coefficient σ.

In this way we are able to construct a sequence of regeneration times R n such that the trajectories (X (Rn+t),t<R n+1 -Rn ) n≥1 are i.i.d. In particular, "regeneration generates independence immediately", i.e. at each regeneration time R n , the "future" X Rn+t , t ≥ 0 is independent of the past σ{X s , s < R n }, without imposing any time lag as usual in the study of processes in continuous time.

Notice that we do not apply the splitting technique to an extracted sampled chain nor to the resolvent chain as in Meyn and Tweedie [START_REF] Meyn | Stability of Markovian processes III : Foster-Lyapunov criteria for continuous-time processes[END_REF] (1993); the loss of memory needed for regeneration is produced only by big jumps. This approach is very natural in this context, but does not seem to be used so far in the literature, except for Xu [START_REF] Xu | Exponential mixing for SDE's forced by degenerate Lévy noises[END_REF] (2011), who works in a very specific frame where the jumps do not depend on the position of the process.

Our paper is organized as follows. In Section 2 we state our main assumptions, prove a lower bound which is of local Doeblin type and state our main results on Harris recurrence and speed of convergence to equilibrium of the process. Section 3 introduces the regeneration technique based on big jumps and proves the existence of certain (polynomial) moments of the associated regeneration times. Section 4 is devoted to an informal discussion on explicit and easily verifiable conditions stated in terms of the coefficients g, σ, c and γ which imply the Harris recurrence. Finally, we give in Section 5 a proof of the local Doeblin condition.

Notations

Consider a Poisson random measure N (ds, dz, du) on R + × R d × R + , defined on a probability space (Ω, A, P ), with intensity measure dsµ(dz) du, where µ is a σ-finite measure on

(R d , B(R d )) of infinite total mass. Let X = (X t ) t≥0 , X t ∈ R d , be a solution of X t = x+ t 0 g(X s ) ds+ t 0 σ(X s ) dW s + [0,t] R d ×R + c(z, X s-)1 u≤γ(z,X s-) N (ds, dz, du), (2.3)
x ∈ R d , where W is an m-dimensional Brownian motion. Write IF = (F t ) t≥0 for the canonical filtration of the process given by

F t = σ{W s , N ([0, s] × A × B), s ≤ t, A ∈ B(R d ), B ∈ B(R + )}.
Throughout this paper, for any x ∈ R d , |x| will denote the Euclidean norm on R d . Moreover, for d × -d matrices A, A denotes the associated operator norm.

Assumptions

In order to grant existence and uniqueness of the above equation, throughout this article, we impose the following conditions on the coefficients g, σ, c and γ.

Assumption 2.1

1. g and σ are globally Lipschitz continuous; σ is bounded.

2. c and γ are Lipschitz continuous with respect to x, i.e.

|c(z, x) -c(z, x )| ≤ L c (z)|x -x | and |γ(z, x) -γ(z, x )| ≤ L γ (z)|x -x |,
where L c , L γ : R d → R + .

sup

x R d (L c (z)γ(z, x) + L γ (z)|c(z, x)|)µ(dz) < ∞. 4. sup x R d γ(z, x)|c(z, x)|µ(dz) < ∞. 5. lim Γ→∞ sup x R d γ(z, x)|c(z, x)|1 {Γ≤γ(z,x)} µ(dz) = 0.
Under these assumptions, Theorem 3 of Bally and Rabiet [START_REF] Bally | Asymptotic behavior for multi scale PDMP[END_REF] (2015) implies that (2.3) admits a unique non-explosive adapted solution which is Markov, having càdlàg trajectories.

Notice that our assumptions do not imply that there exists a finite total jump rate

R d γ(z, x)µ(dz)
for any x ∈ R d . In other words, the above integral might be equal to +∞, and jumps occur with infinite activity. We also stress that due to the presence of the censure term 1 u≤γ(z,X s-)

in equation ( 2.3) we are not in the classical frame of jump diffusions where the jump term depends in a smooth manner on z and x.

In the present article we are seeking for conditions ensuring that the process X is recurrent in the sense of Harris without using additional regularity of the coefficients, based on some minimal non-degeneracy of the jumps and without imposing any non-degeneracy condition on σ. Recall that a process X is called recurrent in the sense of Harris if it possesses an invariant measure m such that any set A of positive m-measure m(A) > 0 is visited infinitely often by the process almost surely (see Azéma, Duflo and Revuz [START_REF] Azéma | Mesure invariante des processus de Markov récurrents[END_REF]

(1969)): For all x ∈ R d , P x ∞ 0 1 A (X s ) ds = ∞ = 1.
We will prove Harris recurrence by introducing a splitting scheme that is entirely based on the "big" jumps of X. In order to do so, we introduce the following additional assumption. Assumption 2.2 1. Writing the Lebesgue decomposition µ = µ ac + µ s , with µ ac (dz) = h(z) dz, for some measurable function h ≥ 0 ∈ L 1 loc (λ), λ the Lebesgue measure on R d , we suppose that there exists z 0 ∈ R d and R > 0 such that inf

z∈R d :|z-z 0 |≤R h(z) > 0.
2.There exists a non-decreasing sequence (E n ) n of subsets of R d and an increasing sequence of positive numbers γn with γn ↑ +∞ as n → ∞, such that

E n = R d , En γ(z, x)µ(dz) =: γn (x) ≤ γn < ∞ (2.4)
for all n.

A useful lower bound

We fix some n. Thanks to (2.4), we can couple the process X t with a rate γn -Poisson process

N [n] = (N [n]
t ) t≥0 such that jumps of X t produced by noise z ∈ E n ,

∆X t = En ∞ 0 c(z, X t-)1 u≤γ(z,X t-) N (dt, dz, du),
can only occur at the jump times

T [n] k , k ≥ 1, of N [n] . Let Π(x, dy) = L(X T [n] k |X T [n]
k -= x)(dy) be the associated transition kernel. Our aim is to obtain a local Doeblin condition of the type

Π(x, dy) ≥ 1 C (x)βν(dy), (2.5)
for a suitable measurable set C, some β ∈]0, 1[ and a suitable probability measure ν.

First notice that γn in (2.4) is only an upper bound on the total jump rate produced by noise belonging to E n . As a consequence, for any k ≥ 1 and on the event that X T [n] k -= x, jumps are only accepted with probability γn(x) γn . Moreover, it is easy to see that the following lower bound holds. Write K = B(z 0 , R) with z 0 and R chosen according to Assumption 2.2 item 1. Then

Π(x, V ) ≥ 1 γn En∩K γ(z, x)1 V (x + c(z, x))µ(dz) ≥ 1 γn En∩K γ(z, x)1 V (x + c(z, x))h(z) dz, (2.6)
where h is the Lebesgue density of the absolute continuous part of µ. It is natural to use a change of variables in the r.h.s. of the above lower bound, i.e. to replace, for fixed initial position x, the argument x + c(x, z)) by y = y(z), on suitable subsets of R d where z → x + c(x, z) is a diffeomorphism. The difficulty is to control the dependence on the starting point x, since we are seeking for uniform lower bounds (2.5), uniform in x ∈ C. This uniform control is achieved in the following proposition.

Proposition 2.3 Grant Assumption 2.2. Suppose moreover that there exist x 0 ∈ R d and r > 0 such that for all x ∈ B(x 0 , r), i) there exists A > 0 with

|∇ z c(z 0 , x)h| ≥ A|h|, ∀h ∈ R d , (2.7)
ii) there exists K > 0 such that for all z ∈ B(z 0 , R),

∇ z c(z 0 , x) -1 i,j ∂ 2 c ∂z i ∂z j (z, x) ≤ K d , (2.8)
iii)

inf z:|z-z 0 |≤R,x:|x-x 0 |≤r γ(z, x)h(z) = ε > 0, S = sup z:|z-z 0 |≤R sup x:|x-x 0 |≤r sup i |∂ z i c(z, x)| < ∞, (2.9)
where h(z) is the Lebesgue density of the absolutely continuous part of µ.

Fix n 0 with B(z 0 , R) ⊂ E n 0 . Then there exist η > 0 and some ball B ⊂ R d such that for all n ≥ n 0 ,

inf x∈B(x 0 ,η) P [X T [n] k ∈ V |X T [n] k -= x] ≥ 1 γn S d ελ(V ∩ B).
(2.10)

As a consequence of Proposition 2.3, the local Doeblin condition (2.5) holds with C = B, β = λ(B)ε γnS d ∧ 1 and ν(dy) = 1 λ(B) 1 B (y)dy. Notice that the set C is not a "petite" set in the sense of Meyn and Tweedie (1993) [START_REF] Meyn | Stability of Markovian processes III : Foster-Lyapunov criteria for continuous-time processes[END_REF].

The main ingredient of the proof of Proposition 2.3 is the following result.

Lemma 2.4 Let Ψ x (z) = x + c(z, x), K = B(z 0 , R), Ψ x (K) = {Ψ x (z), z ∈ K} and a x = x + c(z 0 , x) = Ψ x (z 0 ). Put ρ = A 2 R ∧ 1 2K . ( 2

.11)

Then there exists η > 0 such that

B(a x 0 , ρ 2 ) ⊂ x∈B(x 0 ,η) Ψ x (K).
(2.12)

Moreover, for all x ∈ B(x 0 , r), B(a x , ρ) ⊂ Ψ x (K), and there exists

K x ⊂ K such that z → Ψ x (z) is a C 1 -diffeomorphism from K x to B(a x , ρ).
The proof of this Lemma and of Proposition 2.3 is given in Section 5 below.

Remark 2.5

The ball B appearing in (2.10) can be chosen as B = B(a x 0 , ρ/2) with ρ as in (2.11) and a x 0 = x 0 + c(z 0 , x 0 ). If moreover

L c = sup z:|z-z 0 |≤R L c (z) < ∞, then we can choose η = ρ 2(1 + L c ) ∧ r = A R ∧ 1 2K 4(1 + L c ) ∧ r.
We close this section with two examples where the ball B and the radius η are explicitly given.

Example 2.6 We consider the one-dimensional case, with µ(dz) = dz. Throughout this example, f will be a bounded

1-Lipschitz function such that |f (x)| ≥ f > 0 for all x ∈ B(x 0 , r). 1. Suppose that c(z, x) = e -|z| f (x) for all z ∈ B(z 0 , R) and that |z 0 | ≥ R + a, a > 0. Then for all x ∈ B(x 0 , r), |∇ z c(z 0 , x)h| = f (x)e -|z 0 | h ≥ A|h| with A = f e -|z 0 | . Moreover ∇ z c(z 0 , x) -1 ∂ 2 c ∂z 2 (z, x) = e |z 0 | |f (x)| |f (x)|e -|z| ≤ e |z 0 |-a =: K, ∀z ∈ B(z 0 , R).
Recall that B = B(a x 0 , ρ/2) where a

x 0 = x 0 + c(z 0 , x 0 ) = x 0 + e -|z 0 | f (x 0 ). We have ρ 2 = A 4 R ∧ 1 2K = f e -|z 0 | 4 R ∧ e a-|z 0 | 2 .
Finally, since

L c = sup z:|z-z 0 |≤R L c (z) ≤ e -a , η = ρ 2(1 + L c ) ∧ r ≥ f e -|z 0 | 4(1 + e -a ) R ∧ e a-|z 0 | 2 ∧ r. 2. Suppose now that c(z, x) = f (x) 1 + z 2 and that |z 0 | ≥ R + a, a > 0. Then for all x ∈ B(x 0 , r), |∇ z c(z 0 , x)h| = f (x) 2z 0 (z 2 0 +1) 2 h ≥ A|h| with A = 2f |z 0 | (z 2 0 +1) 2 . Moreover ∇ z c(z 0 , x) -1 ∂ 2 c ∂z 2 (z, x) = (z 2 0 + 1) 2 2|f (x)||z 0 | × 2|f (x)| |3z 2 -1| (z 2 + 1) 3 ≤ (z 2 0 + 1) 2 |3(|z 0 | + R) 2 + 1| |z 0 |(a 2 + 1) 3 = K, ∀z ∈ B(z 0 , R).
In this case,

ρ 2 = f |z 0 | 2(z 2 0 + 1) 2 R ∧ |z 0 |(a 2 + 1) 3 2(z 2 0 + 1) 2 |3(|z 0 | + R) 2 + 1| . Since L c = sup z:|z-z 0 |≤R L c (z) ≤ 1 1+a 2 , we have η = ρ 2(1 + L c ) ∧ r ≥ f |z 0 | 2(1 + a 2 )(z 2 0 + 1) 2 R ∧ |z 0 |(a 2 + 1) 3 2(z 2 0 + 1) 2 |3(|z 0 | + R) 2 + 1| ∧ r.

Drift criteria

The set C = B(x 0 , η) appearing in the local Doeblin condition (2.5) and (2.10) will play the role of a small set in the sense of Nummelin [START_REF] Nummelin | A splitting technique for Harris recurrent Markov chains[END_REF] (1978) and Meyn-Tweedie [START_REF] Meyn | Markov Chains and Stochastic Stability[END_REF] (1993). In order to be able to profit from the lower bound (2.10), we have to show that (X T [n] k -) k comes back to the set C i.o. For that sake, we introduce a drift condition in terms of the continuous time process, inspired by Douc, Fort and Guillin [START_REF] Douc | Subgeometric rates of convergence of ergodic strong Markov processes[END_REF] (2009).

Assumption 2.7 There exists a continuous function

V : R d → [1, ∞[, an increasing concave positive function Φ : [1, ∞[→ (0, ∞) and a constant b < ∞ such that for any s ≥ 0 and any x ∈ R d , E x (V (X s )) + E x s 0 Φ • V (X u ) du ≤ V (x) + bE x s 0 1 C (X u ) du , (2.13)
where

C = B(x 0 , η 2 ), η as in Proposition 2.3.
This drift condition ensures that the process comes back to the set C = B(x 0 , η 2 ) infinitely often. The choice of η/2 is on purpose and will be explained by Proposition 3.3 below.

If V ∈ D(L) belongs to the domain of the extended generator L of the process X, then Theorem 3.4 of Douc, Fort and Guillin [START_REF] Douc | Subgeometric rates of convergence of ergodic strong Markov processes[END_REF] (2009) shows that the following condition

LV (x) ≤ -Φ • V (x) + b1 B (x) (2.14)
implies the above Assumption 2.13. We discuss in Section 4 examples where (2.13) or (2.14) are verified. Under Assumption 2.7, Douc, Fort and Guillin [START_REF] Douc | Subgeometric rates of convergence of ergodic strong Markov processes[END_REF] (2009) give estimates on modulated moments of hitting times. Modulated moments are expressions of the type

E x τ 0 r(s)f (X s ) ds,
where τ is a certain hitting time, r a rate function and f any positive measurable function. Knowledge of the modulated moments permits to interpolate between the maximal rate of convergence (taking f ≡ 1) and the maximal shape of functions f that can be taken in the ergodic theorem (taking r ≡ 1). In the present chapter we are interested in the maximal rate of convergence and hence we shall always take f ≡ 1.

For the function Φ of (2.13) put

H Φ (u) = u 1 ds Φ(s) , u ≥ 1, (2.15) and r Φ (s) = r(s) = Φ • H -1 Φ (s). (2.16)
If for instance Φ(v) = cv α with 0 ≤ α < 1, this gives rise to polynomial rate functions

r(s) ∼ Cs α 1-α ; α = 1 yields r(s) = ce cs .
In most of the cases, we will deal with the case Φ(v) = cv α , 0 ≤ α < 1 and thus work in the context of polynomial rates of convergence. In this situation, the most important technical feature about the rate function is the following sub-additivity property r(t + s) ≤ c(r(t) + r(s)), (2.17) for t, s ≥ 0 and c a positive constant. We shall also use that

r(t + s) ≤ r(t)r(s),
for all t, s ≥ 0.

Main results

Theorem 2.8 Grant the assumptions of Proposition 2.3, Assumptions 2.1, 2.2 and 2.7. 1. The process X is recurrent in the sense of Harris having a unique invariant probability measure m such that Φ • V ∈ L 1 (m). The invariant probability measure m is the unique solution of

R d   1 2 1≤i,j≤d a ij (x) ∂ 2 ψ(x) ∂x i ∂x j + g(x)∇ψ(x)   m(dx) = R d R d (ψ(x) -ψ(x + c(z, x)))γ(z, x)µ(dz) m(dx), (2.18) for all ψ ∈ C 2 (R d ) being of compact support. 2.
Moreover, for any measurable function f ∈ L 1 (m), we have

1 t t 0 f (X s ) ds → m(f ) as t → ∞, P x -almost surely for any x ∈ R d .
The above ergodic theorem is an important tool e.g. for statistical inference based on observations of the process X in continuous time. In this direction, the following deviation inequality is of particular interest. Recall that ν is the measure given in the local Doeblin condition (2.5).

Theorem 2.9 Grant the assumptions of Proposition 2.3, Assumptions 2.1, 2.2 and 2.7 with

Φ(v) = cv α , 0 ≤ α < 1. Put p = 1/(1 -α). Let f ∈ L 1 (m) with f ∞ < ∞, x be any initial point and 0 < ε < f ∞ .
Then for all t ≥ 1 the following inequality holds:

P x 1 t t 0 f (X s ) ds -m(f ) > ε ≤ K(p, ν, X) V (x) t -(p-1) × × 1 ε 2(p-1) f 2(p-1) ∞ if p ≥ 2 1 ε p f p ∞ if 1 < p < 2 . ( 2

.19)

Here K(p, ν, X) is a positive constant, different in the two cases, which depends on p, ν and on the process X, but which does not depend on f , t, ε.

Finally, we obtain the following quantitative control of the convergence of ergodic averages.

Proposition 2.10 Grant the assumptions of Proposition 2.3, Assumptions 2.1, 2.2 and 2.7 with Φ(v) = cv α , 0 ≤ α < 1. Then for any x, y ∈ R d ,

1 t t 0 P s (x, •) ds - 1 t t 0 P s (y, •)ds T V ≤ C 1 t (V (x) (1-α) + V (y) (1-α) ), (2.20)
where C > 0 is a constant. In particular, if α ≥ 1 2 , then

1 t t 0 P s (x, •) ds -m T V ≤ C 1 t V (x) (1-α) . (2.21)
The proof of Theorems 2.8 and 2.9 and of Proposition 2.10 relies on the regeneration method that we are going to introduce now.

3 Regeneration for the chain of big jumps

Regeneration times

We show how the lower bound (2.10) on the jump kernel (2.5) allows us to introduce regeneration times for the process X.

We start "splitting" the jump transition kernel Π(x, dy) of (2.5) in the following way. Since Π(x, dy) ≥ β1 C (x)ν(dy), we may introduce a split kernel Q((x, u), dy), which is a transition kernel from

R d × [0, 1] to R d , defined by Q((x, u), dy) =    ν(dy) if (x, u) ∈ C × [0, β] 1 1-β (Π(x, dy) -βν(dy)) if (x, u) ∈ C×]β, 1] Π(x, dy) if x / ∈ C. (3.22)
Notice that

1 0 Q((x, u), dy)du = Π(x, dy);
it is in this sense that Q((x, u), dy) can be considered as "splitting" the original transition kernel Π by means of the additional "color" u.

We now show how to construct a version of the process X recursively over time intervals [T

[n] k , T [n]
k+1 [, k ≥ 0. We start at time t = 0 with X 0 = x and introduce the process Z t defined by

Z t = x + t 0 g(Z s ) ds + t 0 σ(Z s ) dW s + t 0 E c n ∞ 0 c(z, Z s-)1 u≤γ(z,Z s-) N (ds, dz, du).
For t < T

[n]

1 , we clearly have Z t = X t . Notice also that T

[n] 1

is independent of the rhs of the above equation and exponentially distributed with parameter γn . We put X T

[n] 1 -:= Z T [n] 1 -(notice that Z T [n] 1 = Z T [n] 1 -, since Z almost surely does not jump at time T [n] 1 ). On X T [n]
1 -= x , we do the following. 1. We choose a uniform random variable U 1 ∼ U (0, 1), independently of anything else.

2. On U 1 = u, we choose a random variable V 1 ∼ Q((x , u), dy) and we put

X T [n] 1 := V 1 . (3.23)
We then restart the above procedure with the new starting point V 1 instead of x.

We will write X t for the process with additional color U k , defined by

X t = k≥0 1 [T [n] k ,T [n] k+1 [ (t)(X t , U k ).
Remark 3.1 Notice that the above splitting procedure does not even use the strong Markov property of the underlying process. It only uses the independence properties of the driving Poisson random measure.

This new process is clearly Markov with respect to its filtration, and by abuse of notation we will not distinguish between the original filtration IF introduced in Section 2 and the canonical filtration of X t . In this richer structure, where we have added the component U k to the process, we obtain regeneration times for the process X. More precisely, write

A := C × [0, β] and put R 0 := 0, R k+1 := inf{T [n] m > R k : X T [n] m -∈ A}.
(3.24)

Then we clearly have

Proposition 3.2 a) X R k ∼ ν(dx)U (du) on R k < ∞, for all k ≥ 1. b) X R k +• is independent of F R k -on R k < ∞, for all k ≥ 1. c) If R k < ∞ for all k, then the sequence (X R k ) k≥1 is i.i.d.
It is clear that in this way the speed of convergence to equilibrium of the process is determined by the moments of the extended stopping times R k . In the next section we show that the drift condition of Assumption 2.7 ensures in particular that R k < ∞ P x -almost surely for any x.

Existence of moments of the regeneration times

Recall the local Doeblin condition (2.5) , the definition of the set C and of C = B(x 0 , η/2). Let τ C = inf{t ≥ 0 : X t ∈ C } be the first hitting time of C . It is known (Douc, Fort and Guillin [START_REF] Douc | Subgeometric rates of convergence of ergodic strong Markov processes[END_REF] (2009)) that the condition (2.13) implies that

E x τ C 0 r(s) ds ≤ V (x), (3.25)
where r is given as in (2.16).

Return times to C

In particular, equation (3.25) implies that τ C < ∞ P x -surely for all x. We show that this implies that the regeneration times R k introduced in (3.24) above are finite almost surely. Recall that T

[n]

k are the successive jump times of the dominating Poisson point process N [n] having rate γn . The regeneration times R k are expressed in terms of the jump chain

X T [n]
k -, k ≥ 0. We have to ensure that the control of return times to C for the continuous time process imply analogous moments for the jump chain.

Before stating the first result going into this direction, we have to introduce the following objects. Let σ ∞ be the sup-norm of the diffusion coefficient σ and let B be such that |g(x)| ≤ B(1 + |x|), ∀x ∈ R d . Since g is supposed to be globally Lipschitz continuous, such a constant B clearly exists. Finally, we choose n sufficiently large such that γn > B (3.26) (recall that γn → ∞ as n → ∞) and such that

σ ∞ √ π 2 γn (γ n -B) 3 2 + B η γn (γ n -B) 2 < η 4 , (3.27)
where Proof Recall the process Z t defined by

B η = sup x R d c(z, x)γ(z, x) dµ(z) + B 1 + |x 0 | + η 2 . Proposition 3.3
Z t = x + t 0 g(Z s ) ds + t 0 σ(Z s ) dW s + t 0 E c n ∞ 0 c(z, Z s-)1 u≤γ(z,Z s-) N (ds, dz, du)
and recall that for any t < T

[n]

1 , Z t = X t . Recall also that T

[n] 1

is independent of the rhs of the above equation, exponentially distributed with parameter γn . Now let x ∈ C and upper-bound

P x [X T [n] 1 -/ ∈ C] = P x [Z T [n] 1 -/ ∈ C]. Clearly, P x |Z t -x| ≥ η 2 ≤ 2 η E x [|Z t -x|]. Let T > 0. Then, with Z * T = sup t∈[0,T ] |Z t | and m ∈ N, E x [|Z t -x|1 Z * T <m ] ≤ E x | t 0 σ(Z s ) dW s |1 Z * T <m + E x t 0 |g(Z s )| ds1 Z * T <m + E x t 0 E c n ∞ 0 |c(z, Z s-)|1 u≤γ(z,Z s-) N (ds, dz, du)1 Z * T <m ≤ E x | t 0 σ(Z s ) dW s | + E x t 0 |g(Z s )| ds1 Z * T <m + E x t 0 E c n ∞ 0 |c(z, Z s-)|1 u≤γ(z,Z s-) N (ds, dz, du)
with (using the Itô isometry and the fact that σ is bounded)

E x | t 0 σ(Z s ) dW s | ≤ E x | t 0 σ(Z s ) dW s | 2 ≤ σ ∞ √ t. Moreover, for x ∈ B(x 0 , η 2 ) = C , E x t 0 |g(Z s )| ds1 Z * T <m ≤ E x t 0 B(1 + |Z s |) ds1 Z * T <m = Bt + B t 0 E x |Z s |1 Z * T <m ds ≤ Bt(1 + |x|) + B t 0 E x |Z s -x|1 Z * T <m ds ≤ Bt 1 + |x 0 | + η 2 + B t 0 E x |Z s -x|1 Z * T <m ds.
We upper bound

E x t 0 E c n ∞ 0 |c(z, Z s-)|1 u≤γ(z,Z s-) N (ds, dz, du) = E x t 0 E c n ∞ 0 |c(z, Z s-)|1 u≤γ(z,Z s-) dz dµ(z) du ≤ t sup x E c n c(z, x)γ(z, x) dµ(z) ≤ t sup x R d c(z, x)γ(z, x) dµ(z)
and put

B η = sup x R d c(z, x)γ(z, x) dµ(z) + B 1 + |x 0 | + η 2 . Then E x [|Z t -x|1 Z * T <m ] ≤ σ ∞ √ t + B η t + B t 0 E x |Z s -x|1 Z * T <m ds.
Then Gronwall's lemma (see Proposition 5.5 in the appendix) implies that

E x [|Z t -x|1 Z * T <m ] ≤ ( σ ∞ √ t + B η t)
E x [|Z t -x|] ≤ ( σ ∞ √ t + B η t)e Bt , (3.29)
for all t ≤ T. In the above rhs, the constants do not depend on T, hence (3.29) is actually true for any t ≥ 0. Furthermore

E x [|Z t--x|] ≤ ( σ ∞ √ t + B η t)e Bt , (3.30)
which can be seen as follows. Using (3.29), we have, for s < t,

E x [|Z s -x|1 Z * t <m ] ≤ E x [|Z s -x|] ≤ ( σ ∞ √ s + B η s)e Bs ≤ ( σ ∞ √ t + B η t)e Bt
so, using dominated convergence when s tends to t from inferior values,

E x [|Z t--x|1 Z * t <m ] ≤ ( σ ∞ √ t + B η t)e Bt ,
and letting m → +∞, monotone convergence gives (3.30). Now, T

[n] 1

is independent from (Z t ) t , exponentially distributed with parameter γn . By choice of n, γn > B. Then

E x [|Z T [n] 1 --x|] ≤ +∞ 0 ( σ ∞ √ t + B η t)e Bt γn e -γnt dt = γn σ ∞ +∞ 0 √ te -(γn-B)t dt + B η +∞ 0 te -(γn-B)t dt = γn σ ∞ Γ( 3 2 ) (γ n -B) 3 2 + B η Γ(2) (γ n -B) 2 = γn σ ∞ √ π 2 1 (γ n -B) 3 2 + B η (γ n -B) 2
for every x ∈ B(x 0 , η 2 ). Since n was chosen to have

σ ∞ √ π 2 γn (γ n -B) 3 2 + B η γn (γ n -B) 2 < η 4 , we obtain sup x∈B(x 0 , η 2 ) P x [X T [n] 1 -/ ∈ C] ≤ sup x∈B(x 0 , η 2 ) P x [|Z T [n] 1 --x| ≥ η 2 ] ≤ sup x∈B(x 0 , η 2 ) 2 η E x [|Z T [n] 1 --x|] < 1 2 .

•

The above arguments imply the following statement.

Corollary 3.5 Let S 1 = inf{T [n] k , k ≥ 1 : X T [n] k -∈ C}. Then P x (S 1 < ∞) = 1 for all x.
Proof We introduce the following sequence of stopping times.

t 1 = τ C , s 1 = inf{T [n] k > t 1 }, . . . , t l = inf{s ≥ s l-1 : X s ∈ C }, s l = inf{T [n] k > t l }.
The above stopping times are all finite almost surely. We put τ * = inf{l : X s l -∈ C}.

Then, using (3.28), for any x ∈ R d ,

P x (τ * ≥ n 0 ) ≤ 1 2 n 0
, which shows that τ * < ∞ P x -almost surely for all x. In particular, S 1 ≤ s τ * < ∞ P x -almost surely for all x.

•

The above proof shows in particular that the polynomial control obtained for the first entrance time in C , obtained in (3.25) remains true for S 1 . Moreover we have the following control on polynomial moments of the regeneration times.

Proposition 3.6 Grant Assumption 2.7 with Φ(v) = cv α , 0 ≤ α < 1. Let p = 1 1-α .
Then there exists a constant c such that

E x (S p 1 ) ≤ cV (x). ( 3
E x s 1 0 r(s) ds = E x τ C 0 r(s) ds + E x s 1 τ C r(s) ds ≤ V (x) + E x s 1 τ C r(s) ds,
where we have used (3.25). Now, using that s 1 -τ C is independent of F τ C , exponentially distributed with parameter λ, we upper-bound

E x s 1 τ C r(s) ds = E x E Xτ C s 1 -τ C 0 r(τ C + s) ds ≤ E x (r(τ C ))E(r(s 1 -τ C )) = cE x (r(τ C )), since E Xτ C (r(s 1 -τ C )) = ∞ 0 λe -λt r(t) dt < ∞ does not depend on X τ C . Using that r(t) ≤ c + t 0 r(s) ds, (3.32)
we obtain

E x (r(τ C )) ≤ c + E x τ C 0 r(s) ds ≤ c + V (x).
Therefore,

E x s 1 0 r(s) ds ≤ c + cV (x) ≤ cV (x), (3.33)
where we have used that V (x) ≥ 1.

2. We now use r(t + s) ≤ r(t)r(s) in order to obtain a control of E x tτ * 0 r(s) ds. We certainly have

E x tτ * 0 r(s) ds = E x t 1 0 r(s) ds + n≥1 E x t n+1 tn r(s) ds1 {n<τ * } ≤ V (x) + n≥1 E x 1 {n-1<τ * } r(t n ) t n+1 -tn 0 r(s) ds = V (x) + n≥1 E x 1 {n-1<τ * } r(t n )E Xt n t 1 0 r(s) ds ≤ V (x) + n≥1 E x 1 {n-1<τ * } r(t n )V (X tn ) , (3.34)
where we have used (3.25) and the fact that 1 {n-1<τ * } is F s n-1 -measurable. Now, X tn belonging to C , we can upper-bound V (X tn ) ≤ V C = c, and obtain

E x tτ * 0 r(s) ds ≤ V (x) + c n≥1 E x 1 {n-1<τ * } r(t n ) . (3.35)
We use r(t + s) ≤ r(t)r(s) and the Markov property with respect to t 1 to obtain

E x 1 {n-1<τ * } r(t n ) ≤ E x r(t 1 ) sup y∈C E y (r(t n-1 )1 {n-2<τ * } ).
Using (3.32), the first factor can be treated as follows

E x r(t 1 ) ≤ c + E x t 1 0 r(s) ds ≤ c + V (x) ≤ cV (x), since V (x) ≥ 1.
Further, let p ∈] 1-α α ∨ 1, 1 α [, and q ≥ 1 with 1 p + 1 q = 1. Then, using that

P (n -2 < τ * ) ≤ ( 1 2 ) n-2 , E y 1 {n-2<τ * } r(t n-1 ) ≤ E y (r p (t n-1 )) 1/p 1 2 (n-2)/q .
We have, by definition of r that r p (t) ≤ ct α 1-α p , where α 1-α p > 1 by choice of p. Using Jensen's inequality we obtain

r p (t n-2 ) ≤ c(n -2) α 1-α p-1 (t α 1-α p 1 + . . . + (t n 2 -t n-3 ) α 1-α p ).
We now use, by choice of p, that α

1-α p -1 ≤ α 1-α 1 α -1 = α 1-α
, and therefore

t α 1-α p = 1 α 1-α p -1 t 0 s α 1-α p-1 ds ≤ c t 0 r(s) ds.
This allows to rewrite

r p (t n-2 ) ≤ c(n -2) α 1-α p-1 t 1 0 r(s) ds + . . . + tn 2 -t n-3 0 r(s) ds .
Using successively the Markov property at times t 1 , t 2 , . . . , t n-3 , we obtain

E y r p (t n-2 ) ≤ c(n -2) α 1-α p-1 (n -2) sup z∈C E z t 1 0 r(s) ds.
Finally, by (3.25), sup z∈C E z t 1 0 r(s) ds ≤ sup z∈C V (z) = c, and therefore

(E y r p (t n-2 )) 1/p ≤ c(n -2) α 1-α .
Coming back to (3.35) we conclude that

E x t τ * 0 r(s) ds ≤ V (x) + cV (x) n≥1 1 2 n-2 q (n -2) α 1-α ≤ cV (x).
3. We now argue as follows.

E x s τ * 0 r(s) ds = E x t τ * 0 r(s) ds + E x s τ * t τ * r(s) ds ≤ cV (x) + n≥1 E x 1 {τ * =n} sn tn r(s) ds ≤ cV (x) + n≥1 E x 1 {τ * >n-1} r(t n ) sn-tn 0 r(s) ds ≤ cV (x) + n≥1 E x [1 {τ * >n-1} r(t n ) sup y∈C s 1 0 r(s) ds] ≤ cV (x) + c n≥1 E x 1 {τ * >n-1} r(t n ),
where we have used the Markov property with respect to t n and (3.33). The last sum is treated as (3.35), which concludes our proof, since r(s) ≥ cs α 1-α .

•

The above result implies an analogous control for moments of the regeneration times R k of (3.24). More precisely, we can now define

S l = inf{T [n] k > S l-1 : X T [n] k -∈ C}, l ≥ 2, and let R 1 = inf{S l : U l ≤ β}, R k+1 = inf{S l > R k : U l ≤ β}. (3.36)
An analogous argument as the one used in the proof of Proposition 3.6 then implies Theorem 3.7 Grant Assumption 2.7 with

Φ(v) = cv α , 0 ≤ α < 1 and let p = 1/(1 -α). Then E x R p 1 ≤ cV (x). (3.37)
We are now ready to prove Theorems 2.8 and 2.9.

3.3 Proof of Theorems 2.8 and 2.9

Proof of Theorem 2.8. Let

m(O) := E R 2 R 1 1 O (X s ) ds,
for any measurable set O. By the strong law of large numbers, any set O with m(O) > 0 is visited i.o. P x -almost surely by the process X, for any starting point (

x, u) ∈ R d × [0, 1].
Hence, the process is recurrent in the sense of Harris, and by the Kac occupation time formula, m is the unique invariant measure of the process (unique up to multiplication with a constant). Now, recall that ν is of compact support, hence V ∈ L 1 (ν). Using (3.37) in the case α = 0, we obtain

m(R d × [0, 1]) = E(R 2 -R 1 ) = E ν R 1 ≤ cν(V ) < ∞. This implies that X is positive recurrent.
The invariant measure m of the original process X is the projection onto the first coordinate of m. In particular, X is also positive Harris recurrent, and m can be represented as

m(f ) = E R 2 R 1 f (X s ) ds.
The ergodic theorem is then simply a consequence of the positive Harris property of X. Finally, the fact that Φ • V ∈ L 1 (m) is an almost immediate consequence of (2.13), based on Dynkin's formula.

• Proof of Theorem 2.9. Theorem 2.9 follows from Theorem 5.2 of Löcherbach and Loukianova (2013) in [START_REF] Löcherbach | Polynomial deviation bounds for recurrent Harris processes having general state space[END_REF] together with Proposition 3.6.

•

We finally proceed to the proof of Proposition 2.10.

Proof of Proposition 2.10. Let X and Y be a copies of the process, issued from x (from y respectively) at time 0. Let R 1 and R 1 be the respective regeneration times. Using the same realization V k for X and for Y (recall (3.23)), it is clear that R 1 and R 1 are shift-coupling epochs for X and for Y, i.e. X R 1 +• = Y R 1 +• . If follows then from Thorisson [START_REF] Thorisson | Analysis, algebra, and computers in mathematical research[END_REF] (1994), see also Roberts and Rosenthal [START_REF] Gareth | Quantitative bounds for convergence rates of continuous time Markov processes[END_REF] (1996), Proposition 5, that

1 t t 0 P s (x, •) ds - 1 t t 0 P s (y, •)ds T V ≤ C 1 t (E x (R 1 ) + E y (R 1 )). (3.38)
Recall that p = 1/(1 -α). Then

E x (R 1 ) ≤ (E x R p 1 ) 1/p ≤ c(V (x)) (1-α) . Now, if α ≥ 1 2 , then 1 -α ≤ α and therefore, E x (R 1 ) ≤ c Φ • V (x) ∈ L 1 (m).
In this case, we can integrate (3.38) against m(dy) and obtain the second part of the assertion. •

Discussing the drift condition

In this section, we discuss in an informal way several easily verifiable sufficient conditions implying Assumption 2.7 with Φ(v) = cv α , 0 < α ≤ 1. These conditions will involve different coefficients of the process. Recall that the infinitesimal generator L of the process X is given for every C 2 -function ψ with compact support on R d by

Lψ(x) = 1 2 i,j a ij (x) ∂ 2 ∂x i ∂x j ψ(x) + b(x)∇ψ(x) + R d [ψ(x + c(z, x)) -ψ(x)] K(x, dz),
where a = σσ * and K(x, dz) = γ(z, x)µ(dz). In order to grant Assumption 2.7, we are seeking for conditions implying that

LV ≤ -cV α (x) + b1 C (x), (4.39)
for some 0 ≤ α ≤ 1, with C = B(x 0 , η 2 ) and b, c > 0.

Example 4.1 If we choose for instance V (x) = |x -x 0 | 2 and α =1 2 it suffices to impose that for all x ∈ R d \ C ,

Tr(σσ * ) + 2 g(x), x -x 0 + R d 2(x -x 0 ) + c(z, x), c(z, x) γ(z, x)µ(dz) ≤ -c|x -x 0 |. (4.40)
We now discuss several concrete sufficient conditions implying (4.39). In this context, it is interesting to notice that the influence of the different coefficients can be quite different. Some coefficients can work in a favorable way in order to ensure (4.39). In that case we will say that they are "pushing" the diffusion into the set C . Other coefficients might play a neutral role or even work against (4.39). Since we have three natural parts of coefficients (diffusion part, drift and the jump part), we will discuss here the following cases: "pushing" with the jumps only, "pushing" with jumps and drift together 1 and "pushing" with the drift only.

Pushing with the jumps

Consider first a pure jump process, i.e. the case when a = g = 0. We choose V (x) = |x -x 0 |2 and propose the following conditions. 1. Global condition with respect to z. ∀z ∈ R d , ∀x ∈ B x 0 , η 2 c , c(z, x) + 2(x -x 0 ), c(z, x) ≤ 0.

(4.41)

2. Local conditions with respect to z on some set K. There exists a set K such that the following holds.

1. there exists ξ > 0 such that for all x ∈ B x 0 , η On the one hand, this implies that (x -x 0 ), c(z, x) ≤ -ζ 2 |x -x 0 ||c(z, x)|, which means that c(z, x) belongs to the convex cone of direction (x -x 0 ) and angle arccos -ζ 2 . On the other hand, using (4.44), the following condition

2 (x -x 0 ), c(z, x) + |c(z, x)||x -x 0 | ≤ -ζ|x -x 0 ||c(z, x)|
is a sufficient (but not necessary!) condition which leads to (x -x 0 ), c(z, x) ≤ -(1+ζ) 2 |xx 0 ||c(z, x)|. In other words, it suffices that c(z, x) belongs to the convex cone of direction (x -x 0 ) and angle arccos -(1+ζ) 2 .

The above conditions on the jump mechanism are naturally quite restrictive since they ensure that from everywhere in R d \ C , the jumps force the process into the set C . Nevertheless, this example is useful, and we will come back to these arguments later when discussing the influence of the drift coefficient.

In a next step, let us suppose that σ = 0. Then under the above conditions, for all

x ∈ B x 0 , η 2 c , LV (x) = Tr(σ(x)σ * (x)) + E c(z, x) + 2(x -x 0 ), c(z, x) γ(z, x)µ(dz) ≤ Tr(σ(x)σ * (x)) -ζ|x -x 0 |ξ. Let Σ = sup x∈B x 0 , η 2 c |Tr(σ(x)σ * (x))| |x-x 0 |
and suppose that K is such3 that ζξ > Σ. Then LV (x) ≤ -c(V (x))

1 2
with c = ζξ -Σ. Finally, if g = 0, the minimal additional condition x-x 0 , g(x) < 0, for all x ∈ B x 0 , η 2 c , ensures that the above result will remain true.

Pushing with both jumps and drift part

The conditions we made on the jump mechanism in the above paragraph are of course very strong. In this paragraph, we will therefore consider that these conditions hold only for x belonging to some set E 1 . Moreover, we will suppose that the drift coefficient contributes to force the diffusion into C when x belongs to another set E 2 . More precisely, we suppose that

E 1 ⊂ B x 0 , η 2 c and put E 2 = B x 0 , η 2
(σσ * ) + 2 g(x), x -x 0 ≤ -c|x -x 0 |.
Such a condition is true for example if

g(x) = - 1 2 (c + Σ) x -x 0 |x -x 0 | ,
where we recall that Σ = sup

x∈B x 0 , η 2 c |Tr(σ(x)σ * (x))| |x-x 0 | . Example 4.3
We continue Example 2.6 item 1. and consider the one-dimensional case with µ(dz) = dz and c(z, x) = e -|z| f (x). We suppose that x 0 = 0 and let

E 1 = [-M -η 2 ] ∪ [ η 2 + M ].
Moreover we choose K = [a, a + 2R] in such a way that K e -|z| dz = 1 2 . Finally we will suppose that for all (z, x)

∈ K × E 1 γ(z, x) ≥ γ > 0 and f (x) ≥ f > 0 (4.45) with f • γ > 2 Σ ζ . ( 4

Pushing partially with both

In the last paragraph we supposed that on the subset E 2 where the drift is driving the process towards C , the jumps do not act in a contradictory way -this is actually ensured by the condition (4.41). Notice that it is a priori not possible to weaken this assumption on E 2 . Indeed, without condition (4.41) we have the following structural problem: we cannot even be sure that

E |c(z, x)| 2 γ(z, x)µ(dz) < +∞. (4.47)
Moreover if we do not suppose a global condition as (4.41), it will be necessary to compensate the possible non-negative part E\K c(z, x) + 2(x -x 0 ), c(z, x) γ(z, x)µ(dz) due to jumps in order to obtain a suitable control for LV (x).

Pushing only with the drift

If we decide to ensure the Lyapunov condition by means of the drift coefficient g only, in the same spirit as above, we could take E 1 = ∅, but would have to keep global conditions, like the condition (4.41), if we use the same Lyapunov function. However, if we choose another Lyapunov function, the situation might be more favorable as we are going to explain now. Let for example V

(x) = |x| (= x 2 1 + • • • + x 2 d ) for x ∈ B x 0 , η 2 c . Then ∇V (x) = x |x| , ∂ 2 ∂ i ∂ j V (x) = δ ij |x| - x i x j |x| 3 .
Let D be such that for some 0 < α ≤ 1. Then it is immediate to see that

LV (x) ≤ D 2 d |x| + 1≤i,j≤d |x i ||x j | |x| 3 + x, g(x) |x| + E |c(z, x)|γ(z, x)µ(dz) ≤D = D 2 d |x| + |x| 2 1 |x| 3 + x, g(x) |x| + D ≤ D |x| + x, g(x) |x| + D ≤ D |x| - 1 |x| (c|x| 1+α + D|x| + D) + D ≤ -c|x| α .
5 Proofs

Proof of Proposition 2.3

Proof We first admit Lemma 2.4 and we put K = B(z 0 , R). As a consequence, there exists a ball B(x 0 , η) such that for all x ∈ B(x 0 , η), B(a x 0 , ρ 2 ) ⊂ Ψ x (K). Choose K ⊂ K such that Ψ x : K → B(a x 0 , ρ 2 ) is a C 1 -diffeomorphism for all x ∈ B(x 0 , η). 4 Since for all (z, x) ∈ K × B(x 0 , η), γ(z, x)h(z) ≥ ε, we now have This last result highlights two facts. First, h is an 1 2 -contraction from the complete space B(0, R) into itself, so the fixed-point theorem gives us the existence of z ∈ B(0, R) such that h(z) = z, and, secondly, the range of h defined on B(0, R) is B(0, R), so we have in fact the existence of z ∈ B(0, R) such that h(z) = z, or equivalently, g(z) = y, which ends the proof.

• Remark 5.2 1. g is in fact a C 1 -diffeomorphism from V = B(0, R) ∩ g -1 B 0, R 2 to B 0,

•

We are now able to prove Lemma 2.4. Proof [of Lemma 2.4] 1) Let x ∈ B(x 0 , r). We can apply Proposition 5.4 with a = z 0 , f = Ψ x which gives ρ = A 2 R ∧ Since lim n→+∞ G (Bt) n n! = 0, the assertion follows.

•

Remark 3 . 4

 34 For any n verifying (3.26) and (3.27), inf x∈C P x (X T [n]1 -∈ C) The choice 1 2 in the above lower bound is arbitrary, by choosing larger values of n, we could achieve any bound 1 -on the right hand side of (3.28).

. 31 ) 1 0

 311 Proof We adopt the notation of the proof of Corollary 3.5.1. In what follows, c will denote a constant that might change from line to line. We start by studying E x s r(s) ds, where r is as in (2.16) ands 1 = inf{T n k > τ C }. Let λ = γnbe the rate of the Poisson process associated to T[n] k , k ≥ 1.Then by definition of s 1 ,

2 ( 4 .42) 2 . 2 with c = ζξ > 0 .Remark 4 .2 1 . 2 c 3 .

 242204123 x)|γ(z, x)µ(dz) > ξ. There exists ζ ∈ (0, 1] such that for all z ∈ K and for all x ∈ B x 0 , η2 c . c(z, x) + 2(x -x 0 ), c(z, x) ≤ -ζ|c(z, x) + 2(x -x 0 )||c(z, x)|,(4.43)3. For all z ∈ K and for allx ∈ B x 0 , η 2 c , |c(z, x)| ≤ |x -x 0 |.(4.44)Notice that this last condition implies in particular that |c(z, x)+ 2(x -x 0 )| ≥ |x -x 0 |. Then under the above conditions, for all x ∈ B x 0 , x + c(z, x)) -V (x))γ(z, x)µ(dz) = E c(z, x) + 2(x -x 0 ), c(z, x) γ(z, x)µ(dz) ≤ -ζ K |c(z, x) + 2(x -x 0 )||c(z, x)|γ(z, x)µ(dz) ≤ -ζ|x -x 0 | K |c(z, x)|γ(z, x)µ(dz) ≤ -ζ|x -x 0 |ξ = -c(V (x))1 Using the Cauchy-Schwarz inequality, (4.41) implies that for all z ∈ R d and for allx ∈ B x 0 , η , |c(z, x)| ≤ 2|x -x 0 |. In particular for all x ∈ B x 0 , η 2 c , sup z∈R d |c(z, x)| < +∞.2. The condition (4.44) is a natural condition to force the process to enter into the set B x 0 , η 2 . There is a simple geometric interpretation of the conditions (4.43) and (4.44). Indeed, they lead to the (effective) condition c(z, x) + 2(x -x 0 ), c(z, x) ≤ -ζ|x -x 0 ||c(z, x)| or 2 (x -x 0 ), c(z, x) + |c(z, x)| 2 ≤ -ζ|x -x 0 ||c(z, x)|.

E

  |c(z, x)|γ(z, x)µ(dz) ≤ D and |a ij | < D, where a = σσ * . With γ > 0 such that γ|x| 1 ≤ |x| (where |x| 1 = |x 1 |+• • •+|x d |) and D def = D 2 d + 1 γ 2 we assume that g verifies, for every x ∈ B x 0 , η 2 c , x, g(x) ≤ -c|x| 1+α -D|x| -D, (4.48)

En 1 VK 1 V 2 ) 1 V 2 ,

 11212 (ψ x (z))γ(z, x) dµ(z) ≥ ε (ψ x (z)) dz = ε B(ax 0 , ρ (y)|J ψ -1x (y)| dy.Put z = ψ -1 x (y), then|J ψ -1 x (y)| = 1 |J ψx (z)| = 1 |det(∇ z c(z, x))|and, using Hadamard's Inequality,|det(∇ z c(z, x))| ≤ d i=1 |∂ z i c(z, x)|.As a consequence, we obtainEn 1 V (ψ x (z))γ(z, x) dµ(z) ≥ ε S d λ V ∩ B a x 0 ,with (2.6), ends the proof. • It remains to give a proof of Lemma 2.4. This proof goes through several intermediate steps which are given now. Lemma 5.1 Let g : R d→ R d be a C 2 -function such that 1. g(0) = 0, 2. dg 0 = Id, 3. there exist R, K > 0 such that for all z ∈ B(0, R), i,j,k ∂ 2 g k ∂z i ∂z j (z) ≤ K . Put R = R ∧ 1 2K . Then B 0, R 2 ⊂ g B(0, R) .Proof The third condition allows to apply the Mean Value Inequality to z → dg z since|||d(dg) z ||| ≤ K, ∀z ∈ B(0, R). Therefore, with R = R ∧ 1 2K , dg z -Id = dg z -dg 0 ≤ K|z| ≤ 1 2 , ∀z ∈ B(0, R). Let now y ∈ B 0, R 2 and set h : B(0, R) → R d , z → h(z) := y + z -g(z). We havedh z = Id -dg z ≤ 1 ∀z ∈ B(0, R).Using again the Mean Value Inequality, we obtain for all z, z ∈ B(0, R),|h(z) -h(z )| ≤ 1 2 |z -z |.In particular |h(z)| ≤ 1 2 |z -z | + |h(z )|, so |h(z)| ≤ 1 2 |z| + |h(0)| = 1 2 |z| + |y| < R, for all z ∈ B(0, R).

R 2 . 2 .Lemma 5 . 3 1 . 5 . 4 Proof 1 ) 1 a

 225315411 We could have taken, of course, R = R ∧ 1-ε K for any ε ∈]0, 1[. Let A be a d × d matrix such that ∀h ∈ R d , |Ah| ≥ K|h|. Then B(Au, K R) ⊂ A(B(u, R)). Proof Notice first that A is clearly invertible. Let now y ∈ B(Au, K R). Then for v ∈ R d , |v| = |A(A -1 v)| ≥ K|A -1 v|, so, with v = y -Au, K R ≥ |y -Au| ≥ K|A -1 (y -Au)| = K|A -1 y -u|,or, equivalently, R ≥ |A -1 y -u| implying that A -1 y ∈ B(u, R) and y ∈ A(B(u, R)).•We now have the following extension of Lemma 5.Proposition Let f : R d → R d a C 2 -function and a ∈ R d such that 1. |df a h| ≥ A|h| for all h ∈ R d ,2. there exist R, K > 0 such that for all y ∈ B(a, R), We use Lemma 5.1 withg(z) = df -1 a (f (a + z) -f (a)).All hypotheses needed in Lemma 5.1 are satisfied since ∂ 2 g ∂z i ∂z j (z) = df -

2 ) 2 ⊂ 2 ⊂

 222 Since f (y) = df a g(y -a) + f (a), using Lemma 5.3, df a g B(0, R) , where we have used the preceding step in order to obtain the last inclusion. Therefore, B f (a), A R f B(a, R) .

  .46) It is clear that (4.42) is verified for all (z, x) ∈ K × E 1 , and, moreover, that the jumps are strong enough to ensure the drift condition even in presence of the Brownian part.If we impose moreover that for all x ∈ B 0, η

	2 Adding finally the condition that for all x ∈ B 0, η 2	c , |f (x)| ≤ |x|, then (4.44) is satisfied. c , sgn(f (x)) = -sgn(x), (4.41) is true as
	well and (4.43) follows with ζ = 1.

  Proof It is easy to obtain by induction that, for every n ∈ N * ,g(t) ≤ A 1 + ) dt n • • • dt 1 dt ,

	n-1 k=1 g(t n which implies (Bt) k k! t t 1 t n-1 + B n 0 0 ... 0
	g(t) ≤ A 1 +	n-1 k=1	(Bt) k k!	+ G	(Bt) n n!	.

1 2K such that B(a x , ρ) ⊂ Ψ x B(z 0 , 2ρ A ) ⊂ Ψ x (K),

this will be the most interesting case

This condition has to be seen in relation with Condition (2.9).

Actually it is always possible to multiply γ(z, x) with a sufficiently large constant ensuring that ξ > Σ/ζ.

Indeed, from Lemma 2.4, there exists K ⊂ K such that Ψx : K → B(ax, ρ) is a C 1 -diffeomorphism, and since B(ax 0 , ρ 2 ) ⊂ B(ax, ρ), there exists K ⊂ K ⊂ K such that Ψx :K → B(ax 0 , ρ 2 ) is a C 1 -diffeomorphism.

Acknowledgments

We thank Vlad Bally for many nice discussions on the subject.

This research has been conducted as part of the project Labex MME-DII (ANR11-LBX-0023-01).

where we recall that K = B(z 0 , R). Since our conditions are uniform in x, the radius ρ will be the same for all x ∈ B(x 0 , r).

2) The previous point implies in particular that

(5.50) Therefore,

so it is sufficient to prove that

which can be seen as follows. Let y ∈ B(a x 0 , ρ 2 ) and x ∈ B(x 0 , η), then

so y ∈ B(a x , ρ), for every x ∈ B(x 0 , η) and the statement is proved.

•

Proof [of Remark 2.5] Recall that we have imposed the additional hypothesis

in order to grant (5.50). •

Appendix

In this paper, we have used the following version of Gronwall's lemma.


Wuwei He

Daniel Sidobre

Improving Human-Robot Object Exchange by Online Force Classification

Keywords: Robotic, Physical interaction, Object exchange, Relevance Vector Machine

Handing an object over to a human is a challenging task for a robot to perform, especially when the human partner has no experience interacting with robots. This paper presents our work to enable a robot to learn how to achieve this task with wrist force/torque sensing. Firstly, we present a device to record the data, then we discuss the techniques used for the teaching. We choose to focus on the classification problem defined to enable the robot to find the finger opening movement. The main challenge is that the classification should be run online, at a comparable rate to the controller. To achieve a computationally efficient classifier, the Wavelet Packet Transformation is used for feature extraction, and then we used the Fisher criterion to reduce the dimension of features. A Relevance Vector Machine is used for the continuous classification procedure mainly for its sparsity. Some recorded data and the results from dimension reduction are shown, then we discuss the accuracy and sparsity of the classification by Relevance Vector Machine in this application. The software of continuous classification on forces is then tested on the robot for interactive object exchange between human and robot, which gives promising results.

Introduction

The exchange of objects between human and robot is one of the elementary tasks that an interactive mobile manipulator has to perform. Either in the daily life or in the industrial contexts, natural and intuitive object exchange is preferable for the human partner. To achieve this goal, a device is proposed to teach robots the ability to exchange objects naturally as the humans do.

To achieve safe object hand-over, one of the important challenges for the robot is to determine if it is safe to open the fingers, and the right moment to open the fingers.

In this paper, we show how to learn this skill from force signals recorded by an instrumented device during object exchange between people. This device, named Bidule, uses a 6D force/torque sensor to record interaction forces during the exchanges. We model this force signal in classes that we identify by classification techniques based on a Relevance Vector Machine applied to the information extracted from the force signal. The underlying objective is to identify information that is useful to people during human to human hand-over tasks.

In order to detect grasping, we use wavelet analysis and classification techniques to study object exchange between people. We also propose solutions to reduce the dimension of the classification model. This model, learned from human experience, is then used to teach the robot the capability to know if the human has grasped the manipulated object and to decide when to react.

We present first experimental results obtained with our robots at LAAS-CNRS implementing a full exchange task that integrates the results from human aware motion planning, trajectory generation and control, and safety issues such as online collision checking and external force monitoring on robot joints [START_REF] Sidobre | Humanrobot interaction[END_REF].

The objective of this paper is to present the object designed and the set of techniques used to extract efficiently patterns from the force signals during object exchange. Even if it is still a primarily study, we show the possibility of learning a model from force signals during human interaction aiming to improve physical Human Robot Interaction. We paid a particular attention to choose efficient tools to build our demonstrator and show the capacities of the approach to extract the relevant information from the transmitted force, and so enhance object exchanges.

Related works are presented in Section 2. The Bidule object, designed to study the interaction forces during exchange, is presented in Section 3. The main part of the paper is Sections 4 and 5, which explain the techniques for feature extraction and classification. Then we give some results obtained with a robot in Section 6, followed by the conclusion.

Related Work

A solution to exchange objects between a human and a four-fingered robot hand was presented in [START_REF] Nagata | Delivery by hand between human and robot based on fingertipforce-torque information[END_REF]. Authors used a 6-axis force/torque (f/t) sensor mounted on each fingertip to control the grasping force and evaluate modifications in human grasp condition. Similarly, in the domain of cooperative manipulation with humans [START_REF] Aggarwal | Cooperation between a 4 DOF robotic hand and a human for carrying an object together[END_REF][START_REF] Takubo | Human-robot cooperative manipulation using a virtual nonholonomic constraint[END_REF], researchers try to detect the different stages of cooperation like contact and slip.

The study of forces during object exchange can be done with a sensorized object [START_REF] Kim | Control of a robot hand emulating human's hand-over motion[END_REF][START_REF] Nakazawa | Force control of a robot gripper based on human grasping schemes[END_REF][START_REF] Nakasawa | Force control of a robot hand Emulating human's grasping Motion[END_REF]. Their object was designed to be grasped by two fingers, which are the thumb and the forefinger so the forces are only measured in two directions. When they designed the object, they proposed to answer a set of questions. The same questions can define the objective of our learning based approach too:

• "How does the giver know the receiver's contact?"

• "Which one starts to act first for the smooth hand-over?"

• "How do both of them control the grasp forces during the hand-over?" And they established that "the giver may feel slight vibrations on the fingertip and a change in the weight of the object as the receiver contacts the object". They also showed that the grasp forces are adapted to the weight of the object.

A recent work using a sensorized baton to teach a Willow Garage PR2 robot demonstrates the validity of the approach [START_REF] Chan | A humaninspired object handover controller[END_REF]. The authors built four controllers from the results of a study on human-human hand-overs using the baton. Then they asked participants to compare the controllers. The results show that the controller directly derived from learned elements is the one preferred by participants. But, the conclusion of the paper states that the tests were limited to one object exchanged in a vertical orientation. As the manual construction of a general controller will be very tedious, we propose here to use analysis and classification tools. Likewise the measure of the grasping force is difficult in the general case, in [START_REF] Chan | A humaninspired object handover controller[END_REF] a particular grasp corresponding to the position of a Force Sensitive Resistor sensor is used. In our case, the robot does not control the grasp force and simply open the fingers of the gripper, so the users can use more natural and varied grasps.

Another important consideration is related to the localization of the information in both time and frequency domains. Vibration condition, defined as a threshold on the high-frequency hand acceleration signal, can be used to detect the contact with environment when the robot places an object [START_REF] Romano | Human-inspired robotic grasp control with tactile sensing[END_REF]. Here again, the threshold is determined by manual trial and error, which would be especially difficult when multiple types of contact should be classified.

The methods used in this paper are wavelet analysis based. Readers can refer to [START_REF] Mallat | A wavelet tour of signal processing: The sparse way[END_REF] for a detailed discussion on wavelets and the comparison to Fourier analysis and its variants. The Linear Discriminant Analysis is used to select features, a nice explanation of which can be found in the book [START_REF] Bishop | Pattern recognition and machine learning[END_REF]. Relevance Vector Machine, used as a classifier in this work, was firstly introduced in [START_REF] Tipping | Sparse bayesian learning and the relevance vector machine[END_REF]. Frequency information based classification, although with different feature formation, is also reported with successful application on touch classification [START_REF] Sato | Touche: enhancing touch interaction on humans, screens, liquids, and everyday objects[END_REF], and especially on image classification, such as [START_REF] Rehman | Image classification based on complex wavelet structural similarity[END_REF].

Bidule: a Device for Manipulation Learning

The first objective of this research was to design an autonomous object being able to record contact forces and motion. We designed an object, called Bidule, to learn useful information during object exchange like force and motion signals. The goal is that the learned information can be directly used on robots equipped with a wrist f/t sensor. The f/t sensor is common for mobile manipulators and for robots without a f/t sensor, the wrist f/t can be often estimated. The device was also designed to help to learn the trajectories and the force control policies, however the details of these last two parts will not be covered in this paper.

Figure 1: Left: Bidule design with the 6D force/torque sensor (orange), the handle (blue) are fixed on rings screwed to the sensor, the computer is on the left side and the battery with the power supply board on the right side. Right: Bidule, the device used during this research.

Bidule is built up around an embedded Gumstix microcomputer, a WiFi link, and a 6D force/torque sensor. An Inertial Measurement Unit (IMU) is used to compensate inertial forces and to record the motions of Bidule. The shape and the size were designed to fit better the human hand (Figure 1). Two tubes can be grasped by humans and transmit the external forces to the sensor.

The force signals are amplified, filtered and re-scaled before an Analog to Digital Converter (ADC) digitizes these signals. The link between the computer board and the acquisition board is done via a SPI link running at 3M bits/s available on a 60 pins connector. Two additional bronze rings, each of mass 0.227 kg, can be added around the Bidule object of mass 0.529 kg. This enables us to study different exchange forces when the mass of exchanged objects varies.

During manipulation, the f/t sensor records not only the contact forces of Bidule with the environment, but also the inertia of Bidule. Dynamic compensation needs to estimate firstly the inertial parameters (the mass, center of mass, and all elements in the inertia matrix) of the upper and bottom parts of Bidule. This problem can be illustrated by Figure 2. The whole set of inertial parameters of the upper part of Bidule is named

G 1 : G1 = (m cx cy cz Ixx Ixy Ixz Iyy Iyz Izz) T (1)
Where m is the mass of this part, vector (c x c y c z ) T represents the position of the center of mass in the frame F S of the force sensor and (I xx I xy I xz I yy I yz I zz ) T includes all the terms in the inertial matrix I, which can be written as:

I =   Ixx Ixy Ixz Ixy Iyy Iyz Ixz Iyz Izz   (2)
Inertia matrix I is represented in the frame F 1 located at the center of mass (see Figure 2). Details of the inertial model and how Kalman Filters are used to achieve the online estimation can be found in [START_REF] He | Reactive control and sensor fusion for mobile manipulators in human robot interaction[END_REF]. The same method is used for the

O G 2 G 3 F 1 F s G 1 Figure 2:
The frame F 1 is located at the center of mass of the upper part of the Bidule object. The frame F s is located at the center of the force/torque sensor. G 1 represents the inertial model of the upper part and G 2 that of the bottom part.

bottom part G 2 . It should be noted that G 3 and different coordination transformation need to be handled carefully. We will omit the details in this paper because of limited space. Kubus et al. proposed a recursive approach in [START_REF] Kubus | Improving force control performance by computational elimination of non-contact forces/torques[END_REF] to solve a similar problem. Although we use an Unscented Kalman Filter instead of Recursive Least Square, the readers can refer to the paper for more discussion on the topic. Once the whole set of inertial parameters are estimated, the inertial forces (f inertia and τ inertia ) are simulated and removed from the sensored signals (f sensor and τ sensor ), and we obtain contact forces (f contact and τ contact ).

fcontact = fsensor -finertia (3) τcontact = τsensor -τinertia (4)
The structure of this calculation can be summarized as: Bidule is designed to learn three aspects of object manipulation for a service robot:

• How to classify events from the rich information coming from the f/t sensor during exchange or when the robot end effector contacts environment, either for human grasp or collision with different obstacles.

• The trajectories for hand-over motions adapted to the mass and types of the manipulated object.

• The force control policies for robot during object exchange with human partner.

In this paper, we focus on the first part, while the learning of the trajectory and of the force control policies will be the subject for future research. The control of our robot is centered around a trajectory controller based on Online Trajectory Generation (OTG) [START_REF] Broquère | Soft motion trajectory planner for service manipulator robot[END_REF] [START_REF] Kroger | Online trajectory generation: basic concepts for instantaneous reactions to unforeseen events[END_REF] [START_REF] Broquère | Planification de trajectoire pour la manipulation d'objets et l'interaction homme-robot[END_REF]. During the hand-over, several situations concerning the contact forces will occur:

• To give an object to a person, the robot approaches and tracks the exchange point, which is computed in a frame attached to the human partner. When the person grasps the object, the robot detects the event and opens the gripper. If it is still moving, the robot stops its movement and then it retracts the arm along the same path. This enables a dynamic and reactive object exchange.

• During the hand-over motion, the manipulated object can collide with the human hand or with the environment. The human can also fail to grasp the object. For safety issue, the robot should react differently, depending on the type of collision or failed grasps.

Two such events are shown in Figure 4. Our first objective for learning is to use continuous classification to monitor the force signal history, the previous second for example, to identify what situation it is. We want the robot to be able to use classification to distinguish different contact types and, more importantly, to open the fingers or to stop the motion at the right moment. We define firstly three classes: collision, grasp, and all the rest cases: pure noise, or grasp/collision but not yet the moment to react. In this paper, we limit environment to rigid bodies, such as a table, chair, and other rigid objects in the working space of the robot. The signals would be different collision happened between compliant bodies, and the study for more environment types will be part of our future work.

The examination of the time series of forces/torques shows they depend on the human partner and the manipulated object. For example, a larger object will be handled more slowly and grasped more firmly, so the contact forces should be bigger in magnitude and last longer. Although some simple statistics on the f/t signals could also gain some useful results, we prefer to machine learning techniques to take advantage of the rich time-frequency information in the signals and obtain robust control policies for the robot.

In Figure 5, we see several examples of the forces/torques produced when the object is grasped by a person or during typical collision between Bidule and a table. They show different forms of vibration. We can see that grasp exhibits more variation and collision has a higher-frequencies pattern. For either of the two classes, neither the magnitude nor the duration of the vibration could be easily used to identify the signal. However some characteristics can be used to achieve the classification, which is what we try to extract by machine learning. The choices of techniques are oriented to obtain a computationally efficient and hence fast classifier so it can be run at a high rate, together with other algorithms such as Online Trajectory Generation and online collision checking.

5 Feature Extraction and Learning

Wavelet Packet Transformation

The contact forces are transient signals, for which wavelet based methods have been established to be a reasonable choice to extract time-frequency information. Comparing to Fourier analysis, wavelets capture not only the frequency information, but also the time localization. The Wavelet Packet Transformation (WPT) used for feature extraction is reported in various work relating to myoelectric signal [START_REF] Englehart | Classification of the myoelectric signal using time-frequency based representations[END_REF] or underwater acoustic signal [START_REF] Learned | A wavelet packet approach to transient signal classification[END_REF]. Wavelet Packet Transformation is an algorithm to decompose discrete signals into subbands with different time-frequency resolutions [START_REF] Mallat | A wavelet tour of signal processing: The sparse way[END_REF]. Each step of decomposition is implemented as being filtered by a two-channel filter bank, together with downsampling operations. Figure 6 shows an example of a signal of size sixteen decomposed by a Haar filter bank (A. Jensen & Cour-Harbo, 2001). The original signal is represented at level zero of the tree. At level j = 1, the blue part is obtained by lowpass-downsampling operations of the signal, while the pink part is obtained by highpass-downsampling operations. The same computation is then executed on these two signals, decomposing the signal into four subbands (level j = 2). The signal is decomposed until j = 4, which is the maximum levels J for a signal of size 16. If the Notice that in Discrete Wavelet Transformation (DWT), only high frequency part is recursively decomposed. Then it can be seen as one subset of the WPT. The problem with DWT is it may lose the possibility to highlight some features. As reported in [START_REF] Coifman | Entropy-based algorithms for best basis selection[END_REF], only some signals have the information focused in this fixed representation. And [START_REF] Learned | A wavelet packet approach to transient signal classification[END_REF] shows a clear example of a transient signal where some patterns are captured by a WPT, but they will not appear in the results of a DWT.

In our system, as the transformation should be computed in real-time, we choose to implement each step of the wavelet transformation in lifting scheme. The lifting scheme reduces the computation complexity, requires less memory, and maintains always perfect reconstruction. For the implementation details, readers may refer to [START_REF] Sweldens | The lifting scheme: A new philosophy in biorthogonal wavelet constructions[END_REF].

To illustrate the time-frequency pattern of a signal, we computed WPT of a linear chirp signal y(t) = sin(kπt 2 ) by Daubechies wavelets (db4). Then the fifth level is chosen from the tree and ordered by frequency. The results are presented in Figure 7, in which the coordinates of time and frequency is without units because they lost the physical meaning, and magnitude plots the function log(1 + x 2 (i)), in which x(i) are the wavelet coefficients in level 5. We can see clearly the pattern of an increasing frequency. The ridges that appear perpendicular to the main diagonal line of increasing frequency are artifacts, which will not influence the results of pattern recognition.

For this signal with a clear increasing frequency, it is easy to illustrate its pattern by constructing a time-frequency plot, but to classify complex and variant contact forces we need tools to help us finding the best pattern representation from the WPT.

Feature Formation and Selection

From the WPT tree, several kinds of feature representations exist: part of the decomposition coefficients, statistical information of the coefficients, or the energy map [START_REF] Fatourechi | A wavelet-based approach for the extraction of event related potentials from eeg[END_REF] [START_REF] Birbaumer | The thought translation device (ttd) for completely paralyzed patients[END_REF] [START_REF] Xue | Wavelet packet transform for feature extraction of eeg during mental tasks[END_REF]. The whole coefficients set captures all the in-class variations and noise, which would need a large data set to achieve convergence of the optimization process of learning.

Figure 7: Linear chirp (top) and its time-frequency pattern (bottom), constructed from level 5 of a WPT tree structure. We see that the WPT tree captures the information in time and in frequency.

Feature selection is a critical step to minimize the classification error [START_REF] Peng | Feature selection based on mutual information criteria of max-dependency, max-relevance, and min-redundancy[END_REF]. In many applications, increasing the data set would also minimize the classification error, but in our case, building a large data set from manipulations with different candidates and different situations is a tedious task. On the other hand, the classifier will run online at a comparable frequency as the controller, so reduce the dimension of the learned model is also important. Given these reasons, we decide to differentiate in time and in frequency separately by combining a part of the coefficients and the energy map, instead of keeping all the coefficients. As shown in the following, feature selection can reduce significantly the feature size.

The Relative Energy Map is used to capture the frequency information. Each level of the WPT decomposes the signal into a time-frequency representation with different time-frequency precision. For each subband in the tree, an energy value is computed from all the coefficients in the subband. If we name one subband b and its coefficients x(k), then its energy e b is computed as:

e b = k (x(k)) 2 (5)
By this definition, e 0 , which is computed on the subband S 0 in Figure 8, is the total en-ergy of the signal, and e L 1 the energy in subband S L 1 . The energy is calculated through the whole tree. Then they are all normalized to the total energy of the signal in the time window, which means:

E b = e b e 0 (6)
For a sample to learn, we write this feature set as E = {E b }, which includes the relative energy in all the selected subbands. We notice that a large part of the feature produced by this method are zero, which is no surprise since wavelets produce a sparse representation of the original signal (see [START_REF] Mallat | A wavelet tour of signal processing: The sparse way[END_REF] for more details about the sparsity of WPT). We should then choose the most discriminative features from E to further reduce the size of the feature vector. The energy map is independent of the magnitude and the duration of vibration in the original force signal, and the time information is lost. It can be noted that the relative energy map is already scaled to between 0 and 1.

S LL 2 S LH 2 S HL 2 S HH 2 S L 1 S H 1 S 0
Figure 8: Subbands for the WPT tree of total level 3. L means that the subband has been filtered by a low-pass filter and then downsampled, H means filtered by a highpass filter and downsampled. The subbands in blue are the low-frequency part at each level.

Coefficients of One Subband are selected to capture the global form of the signal. The objective is twofold: to keep the information of the changing in weight when the receiver starts to hold the exchanged object, and to deduce the right time to open fingers or to react for the robot. The robot should differentiate two signals as shown in Figure 9 where the same signal is shifted in time. We choose the coefficients from the first subband of a level in the WPT tree. It contains information of low frequency part of the signal and with a reduced size. We write this feature set as S. In Figure 8, the candidates are in the blue color. Since each level contains such a low frequency subband, the choice of the level is decided by compromising the classification error and the feature size produced. The higher level we choose for this subband, the less it is precise in time. The coefficients are then normalized to between 0 and 1 as the relative energy map.

The choice of the first subband of a level is based on prior knowledge that for signals in all classes of the problem at hand, this subband does indeed capture well the peaks of the original force signals, as shown by Figure 10. For some signals that contain only high frequency components, the choice of subband needs to be made by feature selection algorithms.

Fisher Linear Discriminant Analysis is then used to select features and to reduce further the dimension of the classification model. Fisher-LDA is a tool used for classi-

F x F y F z
Figure 9: The classifier should be able to make difference between the two signals above to find the right time to react to an event detected by the f/t sensor.

Figure 10: Part of the WPT tree plotted for the signal F (F = F x + F y + F z ) in Figure 9. The first subbands (highlighted by red rectangle) capture the peak in signal (the peaks are in white color). The WPT is computed with Symlets 5 wavelets (sym5). The zeros in the tree show also the sparsity of the WPT. The Unit for horizontal axis is omitted since time measurement unit varies for different levels in the WPT. And level 0 is the original time signal. fication and dimension reduction of features [START_REF] Welling | Fisher linear discriminant analysis[END_REF]. Detailed discussion for this tool can be found in [START_REF] Bishop | Pattern recognition and machine learning[END_REF]. Fisher-LDA considers finding the best weight w by maximizing the objective:

J(w) = w T S B w w T S W w (7)
where S B is the "between classes scatter matrix", and S W is the "within classes scatter matrix" [START_REF] Welling | Fisher linear discriminant analysis[END_REF]. The value of J means the separability of classes for a set or a subset of features. Here, we use the criterion to select the best features in our feature set of the relative energy map E. The criterion is reformulated into definition defined as (8). It is calculated for every feature in E to evaluate the separability measure [START_REF] Gu | Linear discriminant dimensionality reduction[END_REF] of that feature. As WPT produces a sparse representation of the signal, which means many terms in the tree are nearly zero, it is reasonable not to include those features, which are invariant between classes.

J = tr(S -1 W S B ) (8)
During the object exchange manipulation, the orientations of f/t changes in different situations: the posture of the robot and of people, the exchange positions, and more. Not including orientations for the event detection simplifies the problem. In fact, the f/t are expressed in the local frame of the sensor, which makes the orientations more difficult to use for the classifier. Learning a model with orientations will also be much more complex and produce a model of much higher dimension, resulting in a slower classifier. While the orientations of f/t are excluded from the feature set, the feature set building process can be summarized as:

• The WPT tree representation is built from force F and from torque T , F is computed as F = ||Fx + Fy + Fz||, and T is computed as T = ||Tx + Ty + Tz||.

• The relative energy map is computed for the WPT.

• It is then combined with a subband from the tree. The subband is also scaled to range from 0 and 1, matching the same as relative energy map.

• Fisher LDA criterion is used to reduce the size of the feature set of energy map E and the level from the subband S selected as a trade-off between the classification accuracy and feature size.

RVM and Classification

Once the features are well defined, the problem is formulated as classification of multipleclasses. In the previous section, the f/t signals are transformed into the feature space of the energy map and a subband from the WPT tree. If we write the whole feature vector as x, then this problem is: given (x i , d i ) i=1,2,...N , a set of N training data where d i is the class label, how to determine a classifier y(x) that correctly classifies the force signals. The Support Vector Machine (SVM) established itself as the state-of-the-art algorithm with strong theoretical foundations for classification and regression [START_REF] Vapnik | The nature of statistical learning theory[END_REF]) [START_REF] Burges | A tutorial on support vector machines for pattern recognition[END_REF]. However, the new Relevance Vector Machines (RVM) approach is really promising, and we choose it for this work. Compared to SVM, RVM has some advantages, as discussed in [START_REF] Tipping | Sparse bayesian learning and the relevance vector machine[END_REF], which will not be repeated in this paper.

We show how a sparse model is achieved by RVM. For clarity, the simple case of a binary classification problem is discussed, which means d i ∈ {0, 1}. Given a new input vector x * , the probability of its class label is given as:

p(d|x * ) = 1 1 + exp(-y(x * )) (9)
the RVM classifier function y is given by:

y(x) = N i=1 w i K(x * , x i ) (10)
In which, K(•, •) is a kernel function, x i (i = 1, 2, 3..., N ) are the N training data. As discussed in [START_REF] Tipping | Sparse bayesian learning and the relevance vector machine[END_REF], the parameters w i in (10) are determined by Bayesian estimation, and a sparse prior is introduced on these parameters, forcing them to be highly concentrated around 0. Then they are computed by maximizing the posterior distribution of the class labels given in the inputs. A very few nonzero terms of w i means that a very few samples in the training data are used in the classifier function given by ( 10), achieving a sparse model. Like for SVM, the kernel trick is used to expand the basis functions for y(x). In this work, the RBF kernel function is used, which is:

K(x * , x) = exp(- ||x * -x|| 2 2σ 2 ) (11)
In which, the optimized σ > 0 defines the kernel width. The choice of RBF kernel is made by trying different kernels and RBF gave satisfying results. Readers can refer to [START_REF] Bishop | Pattern recognition and machine learning[END_REF] and [START_REF] Tipping | Sparse bayesian learning and the relevance vector machine[END_REF] to see how RVM is generalized for multiclass problems, and for information for other kernels. RVM has been successfully applied for fault diagnosis [START_REF] Widodo | Fault diagnosis of low speed bearing based on relevance vector machine and support vector machine[END_REF], supervised hyperspectral classification [START_REF] Mianji | Robust hyperspectral classification using relevance vector machine[END_REF], and for recovering 3D human body pose by regression [START_REF] Agarwal | 3d human pose from silhouettes by relevance vector regression[END_REF].

6 Data Acquisition and Results

The Experimental Protocol

The experiments are carried out with nine pairs of volunteers, all adults but including both genders. We selected people who have worked around robots and also four pairs who have not had the experience. Three qualitative velocities were chosen: slow, normal and fast and explained to each participant. First, Bidule was presented to volunteers as a precious and fragile object that needs to be handle with care. These exchanges define the slow velocity exchanges. Then the real solidity of Bidule was revealed and they were asked to carry out natural exchanges, namely normal velocity exchanges.

Then the volunteers were asked to do fast exchanges. We defined four different orientations for Bidule, the two vertical ones (antenna up and upside down), one horizontal and one tilted. We did one third of these experiments with an additional mass. Figure 11 presents the seventeen different cases of exchange, while Figure 12 shows one such exchange experiment. Each couple of volunteers performed two times each case of exchange, one as giver and one as receiver. Each record lasts three to five seconds, depending on the speed of object exchange.

The Data Set

Firstly the time window is fixed to approximately 1 second because we observed that the slowest manipulation lasts about 1 second, which gives 2048 data points as the acquisition rate is set to 1kHz (1024 for force data points and 1024 data points for torque). By discarding failed data records because of various reasons such as communication failure, 213 experiments for each type of events ("collision" and "grasp") are selected. For class "collision", we recorded data of class collision between Bidule and different environment: human hand, table, bottles, and video tapes. We limit the environment to the working space of the robot, which includes only rigid bodies. Including soft bodies in the "collision" database would change how the classifier identifies events, and how that works on the robot needs further investigation. As we need to model the correct time to release the object for the giver, the samples for class "noise" are simply cut from the data, which could be pure noise, or just signal at different instant from the reaction time, shifted by 0.05, 0.1, and 0.15 second, etc. We used a random process to segment data samples of "noise" from recorded noise signals, with the smallest time difference of 0.05 second. The total data set for training is of size 639. Then 240 data of the three classes are also selected in the same way to test the model.

Feature Selection Results

In this section, we will discuss the choice of features by some results. Table 1 shows the result of classification with different feature set: the full energy map (E), the subband from level 1 S, and the two combined. The choice to combine two different feature sets is justified by the results in the classification errors. To reduce the feature size of the relative energy map E and the subband S, we did two studies: Firstly, we evaluated E by the Fisher LDA criterion, then we fixed E to evaluate the different choices on S by directly comparing the classification error. Firstly, we calculated the energy map from the WPT until level 6, which produces a total of 126 features for force and 126 for torque. The relative energy of level 0 (original signal) is omitted since it is 1 for all instances. We examine the separability J of the relative energy map feature set E. The results can be seen in Figure 13, here for the forces F . Several terms in the feature set dominate the discriminative measure, which enables us to reduce largely the size of the feature set without losing too much the discriminative information of the feature set. The final set of features is chosen from the most discriminating terms from the combination the WPT results of force and of torque.

Figure 13: The separability measurement of features in the relative energy map E of forces. J is calculated by Equation ( 8) through all features in E. The zeros in this measurement are largely due to the sparsity of the WPT tree. Similar results are obtained for torques too. The features are sequenced by the separability measurement J in a decreasing order.

For the choice of the level from which to extract the low frequency subband, we use RVM with RBF kernel and compare directly the classification error rate. During this study the energy map is reduced to the first 13 most discriminative features. Firstly we choose S at the first level, then evaluate the classification error. Then we evaluate at level 2 and until that the error rate of classification increased significantly at level 6.

Figure 14: The trade-off between feature size (blue curve) and classification error (red curve). The small feature size means a faster model for online prediction. Level 5 is chosen as it gives a satisfying balance between number of features and error rate.

It can be explained by the losing of time precision in the higher level of the WPT tree. The subband is fixed to level 5 to keep error rate low with low number of features. By this study, we finally chose the first 13 features out of 252 of the relative energy map E. The subband of level 5 gives 64 features. After the dimension reduction, a classification accuracy rate of 97.5% is obtained, with a compact model having dimension of 9 (the nonzero weights in RVM), while feature size is 77. When the whole WPT tree is used directly to train a RVM, the learning does not converge because of the in-class variations and due to the small size of the training data set, and the whole WPT tree decomposed until level 6 is of size 12288 (each level has 1024 elements for forces and 1024 elements for torques).

Comparison to SVM

When the same feature set is used for a SVM, it gives a similar classification accuracy. As shown in Table 2, RVM produces a much sparser model (hence faster classification). In the table, the dimensions are the number of Support Vectors or of Relevance Vectors. Since we want to run the classifier online at a comparable speed of the controller with other resource consumming computations, reduced complexity in model is important. Notice that only the used features but not the whole WPT tree need to be computed for online classification. The results for SVM is obtained by running the libsvm (Chang & Lin, 2011) (RBF kernel).

The 2.5% classification error happens all when the classifier gives label "collision" when the object is grasped by a person. In this case the robot does not open the gripper when the object is grasped. Although we believe that 2.5% error rate is acceptable, it will be interesting to investigate how we can improve this in the future.

Figure 15: Time delay between the moment when the object comes in contact (shown by "event", and classified as type "collison") and the moment when the robot reacts (shown as "detect", and the reaction is "stop"). The event happens at time a, and the robot reacts at time b. At time b, the system monitors the past 1 second of f/t. This time window that the system monitors is illustrated by the black line.

Firsts Results with a Robot

The classifier is implemented on the robot for the wrist f/t sensor, which is capable of using the learned model from human experience. The classifier runs at 50Hz, together with other software of perception and planning on the robot. Figure 15 demonstrates the interaction forces and the time delay before the robot react. For this experiment, the object was held by the robot, and a person grasped the object, but not firmly. The event given by the classifier is not "grasp" but "collision", showing the ability to tell if the object is properly grasped. The figure shows the time delay between the grasp event and the time for the robot to react (here, is not to open the gripper but stop the motion, since "collision" is detected, instead of "grasp"), with detection by force signals.

The time delay is normally less than half a second between the event happening and the reaction, with an average of 0.37 seconds. So the choice of a one second observation window seems adequate. One important aspect is that the result does not show the fastest speed the classifier can detect the event, but the right moment to react. This delay is learned by the model, which is decided by the instances acquired by human to human manipulations. It is possible that volunteers tend to exchange Bidule slower than how they do in daily life with familiar objects, but we do not have strong results about this hypothesis. The classifier, on the other hand, can achieve faster detection if only fast training data is selected. We achieved a model which gives 0.2 second reaction by a smaller set of selected fast data, but we decide that the reaction speed should be determined by user study instead of hand selection of training data.

We would like to mention that external force monitoring is implemented on the joint space when collision detection is needed for safety issues, so this reaction time does not introduce safety issue for the robot.

First manipulations between robot and people have been carried out with promising results. Figure 16 shows how the software, from signal processing to classification, are used for manipulations. During this manipulation, the robot decides to open the gripper. A user study of this approach to see how it can improve the quality of object exchange for human users in interactive manipulation is still to finish. The user study will include many aspects, such as human preference for the speed of reaction of the robot. Another example would be to compare the false negative and false positive for grasp detection and see which one has the greater influence on user's experience.

Conclusion and Future Work

The device proposed to teach the robot how to exchange object like humans is very promising as rich information can be extracted. Force signal classification is studied in this paper. Wavelet Packet Transformation is chosen, from which the energy map and a subband are used to extract the features, and Fisher Linear Discriminant Criterion is used to reduce the dimension of the feature set. Furthermore, the Relevance Vector Machine is employed as classifier because it produces a sparse model. Some results for feature selection and for classification are shown. It can be noted that by constructing a classifier, robust control policies can be achieved with no thresholds to choose by trial and error, which becomes especially difficult when more classes should be considered.

Our experimental device called Bidule has been designed to learn more than what we discussed in this paper, like force control policies during exchange and the movements. Those topics deserve further research and another paper to discuss. For the classification problem, more classes can be added to make a difference between the collisions (soft surface against hard surface), or separate more on the different grasp types. On the other hand, the user study about the event detection is the next work to finish. We need also to investigate the possibilities of using the same set of classification techniques for other manipulations, like putting objects on a table or grasping a handle. At the end, we want to mention that the collision detection by f/t sensor is also to be used for other tasks, such as for the robot to put properly an object on the table. For example, touching of an object with a table by its bottom surface of by an edge would surely produce different patterns on the force signals. And in this case, it would be reasonable to include the torques into the feature set too. With these new methods discussed in this paper to build compact and efficient models, online classification can improve many aspects of the robotic manipulations.

Figure 3 :

 3 Figure 3: Contact Forces/Torques Computation. The gravitational vector g, acceleration a and angular velocity ω can be computed from an inertial measurement unit (IMU). f s and τ s are measured forces/torques vectors. The same strategy can be used for a wrist force sensor on the robot arm to achieve online estimation of the contact forces/torques.

Figure 4 :

 4 Figure 4: Left: the manipulated object collides with the green box. Right: human takes the object given by the robot.

Figure 5 :

 5 Figure 5: Forces/torques signals during several exchange and collision events. On the left, we see three different object exchanges between two people, when the object mass changes too. On the right, f/t of three different collisions between object and rigid environment happened. The scale for time, force and torque are set to be the same for all the plots.

Figure 6 :

 6 Figure 6: Wavelet Packet Transformation decomposes the original signal recursively into a tree structure, here through the Haar wavelet filter bank. Each subband is decomposed into a high-frequency part (pink) and a low-frequency part (blue).

Figure 11 :

 11 Figure 11: The configurations of Bidule used for the experiments

Figure 16 :

 16 Figure 16: Software tested during manipulation. (a): Robot Jido, the human partner and objects to manipulate. (b): The robot decides to take the object on the table. (c): The person asks for the object. (d): Robot gives the object to him, and he takes the object. (e): Detecting that the object is taken by the person, the robot opens the gripper. (f): Object exchange finished with success.

Table 1 :

 1 Error rate for different feature set. RVM with a RBF kernel is used as the classifier

	Features	E	S	{E S}
	Error(%)	15.42	21.67	2.5

Table 2 :

 2 For this application, RVM produces a much sparser model than SVM, and a comparable classification accuracy.

		Dimension Accuracy(%)
	SVM	151	97.92
	RVM	9	97.5

Acknowledgements

This work has been supported by the European Community's Seventh Framework Program FP7/2007-2013 SAPHARI under grant agreement no. 287513 and by the French National Research Agency project ANR-10-CORD-0025 ICARO.

Authors' names and contact information: Wuwei He and Daniel Sidobre CNRS, LAAS, 7 avenue du colonel Roche, F-31400 Toulouse, France; and Univ de Toulouse; UPS, LAAS; F-31400 Toulouse, France. Email: wuwei.he@laas.fr and daniel.sidobre@laas.fr.