Etienne Tanré

Pierre Vallois

Range of Brownian motion with drift

Keywords: Range Process, Enlargement of filtration, Brownian motion with drift MSC 2000 : 60E10, 60F05, 60G17, 60G40, 60G44, 60J10, 60J60, 60J65

Range of Brownian motion with drift

Etienne Tanré, Pierre P.

Introduction

1) The range concerning one-dimensional Markov chains or random walks has been already investigated by [START_REF] Glynn | On the range of a regenerative sequence[END_REF] and [START_REF] Feller | The asymptotic distribution of the range of sums of independent random variables[END_REF]. See for instance [START_REF] Vallois | The range of a simple random walk on Z[END_REF] for a short survey. The aim of this paper is to study the range of a Brownian motion with drift ; this process being the prototype of transient diffusions. The range (R (X) (t) ; t ≥ 0) associated with a continuous process (X t ; t ≥ 0) is the process :

R (X) (t) = sup 0≤u,v≤t (X v -X u ) = sup 0≤u≤t X u -inf 0≤u≤t X u .
(1.1)

When (X t ) is a one-dimensional Brownian motion started at 0, Feller [START_REF] Feller | The asymptotic distribution of the range of sums of independent random variables[END_REF] has computed the density function of R (X) (t), using the fact that the joint distribution of sup 0≤u≤t X u and inf 0≤u≤t X u is explicitly known. Unfortunately the result is expressed as the sum of a series, and the result cannot be generalized to diffusions since the joint distribution of the maximum and the minimum is in general unknown.

To go further we observe that t → R (X) (t) is a non-decreasing, continuous function starting at 0. Therefore we can define its right continuous inverse : θ (X) (a) = inf t ≥ 0; R (X) (t) > a .

(1.2)

Formally it is equivalent to deal with (θ (X) (a) ; a ≥ 0) either (R (X) (t) ; t ≥ 0), since we have R (X) (t) < a = θ (X) (a) > t .

(1.

3)

It is actually more convenient to work with (θ (X) (a) ; a ≥ 0). For instance, if (X t ) is a Brownian motion, then for any a > 0, the Laplace transform of the r.v. θ (X) (a) can be computed, see [START_REF] Imhof | On the range of Brownain motion and its inverse process[END_REF], [START_REF] Vallois | Diffusion arrêtée au premier instant où l'amplitude atteint un niveau donné[END_REF]. Moreover (θ (X) (a) ; a ≥ 0) has independent increments.

Assume that (X t ) is a diffusion process and for simplicity X 0 = 0. It is proved (Theorem 4, [START_REF] Vallois | Diffusion arrêtée au premier instant où l'amplitude atteint un niveau donné[END_REF]) that the process (X t ; 0 ≤ t ≤ θ (X) (U a)) is distributed as (X t ; 0 ≤ t ≤ T (X) (aU )∧T (X) (a(U -1)), where a > 0, T (X) (c) is the first hitting time of level c, and U denotes a r.v. uniformly distributed on [0, 1], independent of the underlying process (X t ). This property has been generalized by [START_REF] Salminen | On first range times of linear diffusions[END_REF].

The Laplace transform of T (X) (c) ∧ T (X) (d) can be expressed besides through eigenfunctions associated with the generator of (X t ). Consequently previous result gives an analytic expression for the Laplace transform of θ (X) (a).

2) Let (B(t) ; t ≥ 0) be a one-dimensional Brownian motion started at 0, and (S t ; t ≥ 0) its unilateral maximum : S t = sup 0≤u≤t B(u). The process (S t ; t ≥ 0) gives rise to an excursion theory. More precisely it can be proved (see for instance [START_REF] Revuz | Continuous martingales and Brownian motion, volume 293 of Grundlehren der Mathematischen Wissenschaften [Fundamental Principles of Mathematical Sciences[END_REF], Chap. XII, section 2.) that the family of processes

x -B(T (B) (x -) + t) ; 0 ≤ t ≤ T (B) (x) -T (B) (x -) ; x > 0 is Poisson point process whose characteristic measure is explicit.

It seems natural to ask for a similar question replacing the one-sided maximum by the range. Results in this direction has been obtained in [START_REF] Hsu | Brownian excursions from extremes[END_REF], [START_REF] Imhof | A construction of the Brownian path from BES 3 pieces. Stochastic Process[END_REF] and [START_REF] Vallois | Decomposing the Brownian path via the range process[END_REF]. In these decompositions, the fact that the Brownian motion is recurrent plays a crucial role. We shall be a little bit more precise below. This leads naturally to investigate the case of transient diffusions. To obtain explicit results, we restrict our attention to the case of Brownian motion with drift.

Let (B δ (t) ; t ≥ 0) be the Brownian motion with drift δ : B δ (t) = B(t) + δ t. For simplicity, we note (R δ (t) ; t ≥ 0) the range of (B δ (t)) instead of (R (B δ ) (t) ; t ≥ 0) :

R δ (t) = sup 0≤s,u≤t (B δ (s) -B δ (u)) = sup 0≤s≤t B δ (s) -inf 0≤s≤t B δ (s).
(1.4)

The right continuous inverse of (R δ (t) ; t ≥ 0) is noted (θ δ (a) ; a ≥ 0) :

θ δ (a) = inf {t ≥ 0 ; R δ (t) > a} .
(1.5)

3) Let us explain why the transient case is drastically different from the recurrent case. Let Σ δ be the random set :

Σ δ = {a ≥ 0 ; B δ (θ δ (a))B δ (θ δ (a-)) < 0} (1.6)
where θ δ (a-) is the left limit of θ δ at a.

When δ = 0, for any a > 0, the two sets Σ 0 ∩ (0, a] and Σ 0 ∩ [a, +∞) are infinite. This generates a difficulty to enumerate the points in Σ 0 . However a decomposition of the Brownian path (B 0 (t) ; t ≥ 0) has been given [START_REF] Imhof | A construction of the Brownian path from BES 3 pieces. Stochastic Process[END_REF][START_REF] Vallois | Decomposing the Brownian path via the range process[END_REF] through Σ 0 , and three dimensional Bessel processes.

Suppose δ > 0. Since B δ (t) drifts to infinity, as t → +∞, for any a > 0 the random set Σ δ ∩ [a, +∞) is finite and enables an explicit description of Σ δ . Let us first introduce the minimum -S 1 of (B δ (t) ; t ≥ 0), and ρ 1 the last time such that B δ (ρ 1 ) = -S 1 . Secondly let S 2 be the maximum of B δ (t), t running over [0, ρ 1 ], and ρ 2 be the last time in [0, ρ 1 ] such that B δ (ρ 2 ) = S 2 , and so on, see (2.1) for details. Then

Σ δ = {A n ; n ≥ 1} (1.7)
where

A n = S n + S n+1 .
We compute the density function of (S 1 , • • • , S n ) in Proposition 2.1. We also obtain an enlighted description of the law of (A 1 , • • • , A n ) (cf Theorem 2.2). Let ψ : R → (-∞, 1) be the function :

ψ(x) = e x -1 -x e x -1 if x = 0, and ψ(0) = 0. (1.8)
Then ψ is one-to-one from R to (-∞, 1) and :

ψ(2δA 1 ), ψ(-2δA 2 ) ψ(-2δA 1 ) , . . . , ψ((-1) n+1 2δA n ) ψ((-1) n+1 2δA n-1 ) d = (U 1 , U 2 , . . . , U n ) (1.9)
where U We are also able to determine the law of the sequence of processes (|B

δ (t + ρ n ) -B δ (ρ n )| ; 0 ≤ t ≤ ρ n-1 -ρ n ) n≥2 .
The result being technical, we do not state it in the Introduction. A complete formulation can be found in Theorem 2.4. We would like to note that three dimensional Bessel processes coming from the decomposition of the Brownian motion have to be replaced by non negative diffusions distributed as the process Z (δ) where

Z (δ) (t) = B t + δ t 0 coth(δZ (δ) (s)) ds.
(1.10)

Heuristically, (Z (δ) (t) ; t ≥ 0) is the process (B δ (t) ; t ≥ 0) conditioned to be positive, see [START_REF] Williams | Path decomposition and continuity of local time for onedimensional diffusions[END_REF]. Taking formally the limit δ → 0 in (1.10), we recover the three-dimensional Bessel process started at 0.

4)

In Section 4, we focus on the law of R δ (t) (resp. θ δ (a)) where t > 0 (resp. a > 0) is fixed. We compute the two distribution functions of R δ (t) and θ δ (a). We partially recover the result of [START_REF] Chong | The ruin problem and cover times of asymmetric random walks and Brownian motions[END_REF]. A path decomposition of (B δ (t) ; 0 ≤ t ≤ θ δ (a)) (Proposition 4.1) allows us to determine the Laplace transform of θ δ (a). Then, by a straightforward approach, it is easy to obtain the asymptotic behaviour of θ δ (a), a → +∞. We prove two limit results : the first one looks like a Law of Large Numbers and the second one is similar to the Central Limit Theorem. More precisely :

θ δ (a) a a.s. -→ a→∞ 1 δ , (1.11)
√ a θ δ (a) a - 1 δ d -→ a→∞ N 0, 1/δ 3 , (1.12)
where N 0, 1/δ 3 denote the Gaussian distribution with 0-mean and variance 1/δ 3 .

Notations and main results

We keep the notations introduced in the Introduction. In particular (B ±δ (t) ; t ≥ 0) is a Brownian motion with constant drift ±δ starting at 0. Throughout this paper, we suppose that δ > 0.

Let -S 1 denote the absolute (random) minimum of (B δ (t) ; t ≥ 0) and ρ 1 the last (random) time when it is reached. The decomposition of the Brownian path t → B δ (t) may be continued by induction on k ≥ 2, as follows :

   S k = sup t∈[0,ρ k-1 ] (-1) k B δ (t), ρ k = sup{t ∈ [0, ρ k-1 ] , B δ (t) = (-1) k S k }. (2.1)
Note that if we set ρ 0 = ∞, then relations (2.1) are still valid for k = 1.

Since B δ (t) goes to +∞, as t → ∞, and t → B δ (t) is a continuous function, then the random times (ρ k ) k≥1 are well defined. This does not mean that {t ∈ [0, ρ k-1 ]; B δ (t) = (-1)

k S k } is reduced to the singleton {ρ k }. We actually prove that this property holds.

We start with the law of (S 1 ,

• • • , S k ) ; k ≥ 1.
Proposition 2.1 Suppose δ > 0. Let (S k ; k ≥ 1) be defined by (2.1). Then 1. The law of S 1 is exponential of parameter 2δ (i.e. with density function 2δe -2δx 1 {x≥0} ).

2. Conditionally to {S 1 = x 1 }, S 2 has density function :

δe δx1 sh(δx 1 ) 1 sh 2 (δ(x + x 1 )) 1 {x≥0} . (2.2) 3. For every k ≥ 3, conditionally to {S 1 = x 1 , . . . , S k-1 = x k-1 }, S k has density function : δ sh δx k-1 sh δ(x k-1 + x k-2 ) sh δx k-2 1 sh 2 δ(x k-1 + x) 1 {0≤x≤x k-2 } .
Proof The proof of Proposition 2.1 is postponed in Section 3.

Note that x k + x k-1 appears in the conditional density function of S k . This leads us to introduce

A k = S k + S k+1 ; k ≥ 1.
(2.3) (A k ) k≥1 is the sequence of maximal ranges associated with (B δ (t) ; t ≥ 0).

Theorem 2.2 Suppose δ > 0. Then 1. (A 1 , • • • , A n ) has a density function given by δ n 2 n-1 k=1 e (-1) k δa k sh(δa k ) e (-1) n 2δan -1 + (-1) n+1 2δa n sh 2 δa n 1 {0≤an≤•••≤a1} . (2.4)
2. Let : ψ : R → (-∞, 1)

ψ(x) = e x -1 -x e x -1 ; x = 0 and ψ(0) = 0.
ψ is one-to-one from R to (-∞, 1) and we have

ψ(2δA 1 ), ψ(-2δA 2 ) ψ(-2δA 1 ) , . . . , ψ((-1) n+1 2δA n ) ψ((-1) n+1 2δA n-1 ) d = (U 1 , U 2 , . . . , U n ) (2.5) where U 1 , • • • , U n are i.i.d. r.v.'s, uniformly distributed on [0, 1].
Proof The proof is postponed in Section 3 Remark 2.3 1) Recall that Σ δ is defined by (1.6). In [START_REF] Vallois | Decomposing the Brownian path via the range process[END_REF], it is proved that if δ = 0, then Σ 0 is a one-dimensional Poisson point process (P.p.p.) with characteristic measure ν(da) = 1 a 1 {a>0} da. Concerning P.p.p. we refer to [START_REF] Revuz | Continuous martingales and Brownian motion, volume 293 of Grundlehren der Mathematischen Wissenschaften [Fundamental Principles of Mathematical Sciences[END_REF] (chapter XII).

If δ > 0, we claim that Σ δ is not a P.p.p.. Using (2.5), by tedious calculations we prove

E   n≥1 f (A n )   = δ 2 ∞ 0 f (x) x sh 2 (δx) dx, (2.6)
for any positive Borel function f . Suppose that Σ δ is a P.p.p. with characteristic measure ν δ , (2.6) implies that :

ν δ (dx) = δ 2 x sh 2 (δx) 1 {x>0} dx.
(2.7)

A straightforward calculation gives :

ν δ ([b, +∞)) = δb coth(δb) -log(sh(δb)) -log 2 , b > 0.
(2.8)

Consequently, N = n≥1 1 {An≥b} is a Poisson random variable with parameter ν δ ([b, +∞)).
But n → A n is a decreasing sequence, then

P(N = 0) = P(A 1 < b) = Ψ(2δb).
This generates a contradiction since

P(N = 0) = exp(-ν δ ([b, +∞))) = Ψ(2δb).
2) Due to (2.4), it is easy to check that (A 2n-1 , A 2n ; n ≥ 1) is a Markov chain, which takes its values in (x, y) ∈ R 2 ; 0 < y < x , with transition probability density function :

P ((u, v); (x, y)) = 4δ 2 Ψ(2δv) 1 e 2δx -1 e 2δy e 2δy -1 Ψ(2δy)1 {0<y<x<v} . (2.9)
Note that this quantity does not depend on u.

3) If we set

X n = ψ (-1) n+1 A n ; n ≥ 1, then (X n ; n ≥ 1
) is a Markov chain on (-∞, 1), with initial distribution the uniform distribution on [0, 1], with transition probability kernel :

K(x, f ) = E f (U ψ(-ψ -1 (x))) (2.10)
where U denotes a r.v. uniformly distributed in [0, 1].

We now give the law of the sequence of processes ((

B δ (t + ρ k ) -B δ (ρ k ) ; t ∈ [0, ρ k-1 -ρ k ])) k≥1 .
The result makes use of the process Z (δ) defined by (1.10).

Theorem 2.4 Let δ > 0 and k ≥ 2. Conditionally to S 1 = x 1 , . . . , S k = x k , 1) (B δ (t) ; 0 ≤ t ≤ ρ k ), (B δ (t + ρ 1 ) -B δ (ρ 1 ) ; t ≥ 0), (B δ (t + ρ 2 ) -B δ (ρ 2 ) ; 0 ≤ t ≤ ρ 1 -ρ 2 ), ..., (B δ (t + ρ k ) -B δ (ρ k ) ; 0 ≤ t ≤ ρ k-1 -ρ k ) are independent ; 2) (B δ (t) ; 0 ≤ t ≤ ρ k )
has the law of Brownian motion with drift (-1) k δ stopped at the first hitting time of level (-1) k x k , and conditioned not to hit

(-1) k+1 x k-1 ; 3) For any 2 ≤ l ≤ k, (|B δ (t + ρ l ) -B δ (ρ l )| ; 0 ≤ t ≤ ρ l-1 -ρ l
) is a process with the same law as Z (δ) stopped at the first moment it reaches the level

x l + x l-1 .
3 Proofs of Theorems 2.2, 2.4 and Proposition 2.1

We keep the notations introduced in Sections 1 and 2.

3.1 Proof of Theorem 2.2

1) Formula (2.4) is a direct consequence of Proposition 2.1. Indeed, let F be a test function. We have :

E [F (A 1 , . . . , A n )] = E [F (S 1 + S 2 , . . . , S n + S n+1 )] = R n+1 F (x 1 + x 2 , . . . , x n + x n+1 )2δ n+1 e -δx1 sh δx n n-1 k=1 sh δ(x k+1 + x k ) × 1 sh 2 δ(x n+1 +x n ) 1 {x1≥x3≥•••≥0} 1 {x2≥x4≥•••≥0} dx 1 . . . dx n+1 .
We use the following change of variables :

       x 1 = x 1 a 1 = x 1 + x 2 . . . a n = x n + x n+1 . (3.1)
Then,

x k -x k+2 = a k -a k+1 ; 1 ≤ k ≤ n -1 x n = α + (-1) n+1 x 1 x n+1 = a n -(α + (-1) n+1 x 1 ) where α = a n-1 -a n-2 + • • • + (-1) n a 1 .
Consequently,

E [F (A 1 , . . . , A n )] = R n F (a 1 , . . . , a n ) 2δ n n-1 k=1 sh(δa k ) sh 2 δa n ×1 {0 ≤ a n ≤ a n-1 ≤ • • • ≤ a 1 } I(α)da 1 . . . da n ,
where

I(α) = ∞ 0 δe -δx1 sh δ(α + (-1) n+1 x 1 )1 {0≤α+(-1) n+1 x1≤an} dx 1 .
Suppose for instance that n is even. Setting y = α -x 1 in order to obtain Then (2.4) follows immediately.

I(α) = e -
2) The density function

h 1 of A 1 is h 1 (x) = δ 2 sh 2 δx e -2δx -1 + 2δx . (3.2) Let H 1 denote the distribution function of A 1 H 1 (x) = x 0 h 1 (y) dy. (3.3)
We have

x

0 e 2y -1 -2y sh 2 y dy = e x sh x e -2x -1 + 2x . (3.4) Consequently H 1 (x) = ψ(2δx). (3.5) So A 1 d = 1 2δ ψ -1 (U 1 ), (3.6)
where U 1 is a uniform random variable on [0, 1].

We compute now the density h

(x) 2 of A 2 conditionally to A 1 = x. h (x) 2 (y) = δ 2 2 e 2δy -1 -2δy sh 2 δy e -δx sh δx 1 {0≤y≤x} 2 sh 2 δx δ 1 e -2δx -1 + 2δx = δ e -δx sh δx e -2δx -1 + 2δx e 2δy -1 -2δy sh 2 δy 1 {0≤y≤x} . (3.7) Let H (x) 2
denote the associated distribution function

H (x) 2 (y) = e -δx sh δx e -2δx -1 + 2δx δ y 0 e 2δt -1 -2δt sh 2 δt dt. (3.8)
Relation (3.4) implies that :

H (x) 2 (y) = e -δx sh δx e -2δx -1 + 2δx e δy sh δy e -2δy -1 + 2δy = ψ(-2δy) ψ(-2δx) .
So, we have proved that ψ(2δA 1 ),

ψ(-2δA 2 ) ψ(-2δA 1 ) is distributed as (U 1 , U 2 ),
where

U 1 and U 2 are two independent r.v.'s, uniformly distributed on [0, 1].
Reasoning by induction, identity (2.5) can be proved by the same way via (2.4).

Proofs of Theorem 2.4 and Proposition 2.1

In our approach, Theorem 2.4 and Proposition 2.1 are a direct consequence of Proposition 3.2 stated below. Therefore, we first focus on this key result, and we prove Theorem 2.4 and Proposition 2.1 at the end of this subsection.

The description of the laws of (B δ (t) ; 0 ≤ t ≤ ρ 1 ) and (B δ (t) -B δ (ρ 1 ) ; t ≥ ρ 1 ) is given by the well-known theorem of Williams [START_REF] Williams | Path decomposition and continuity of local time for onedimensional diffusions[END_REF].

Proposition 3.1 ([17])
1) The law of S 1 is exponential with parameter 2δ (i.e. its density function is

2δe -2δx 1 {x≥0} ). 2) Conditionally to S 1 = x 1 , a. (B δ (t) ; t ≤ ρ 1 ) and (B δ (t+ρ 1 )-B δ (ρ 1 ) ; t ≥ 0) are independent processes b. (B δ (t) ; t ≤ ρ 1
) is a process with the same law as a B -δ stopped at its first hitting time of -x

1 c. (B δ (t + ρ 1 ) -B δ (ρ 1 ) ; t ≥ 0) is a (positive) process distributed as Z (δ) .
To obtain the decomposition of the Brownian motion with drift given in Theorem 2.4, Proposition 3.1 leads us to study (B -δ (t) ; 0

≤ t ≤ T -δ -a ) conditionally to max 0≤t≤T -δ -a B -δ (t) where : T -δ -a = inf {t ≥ 0 , B -δ (t) < -a} , a > 0. (3.9)
Since (-B -δ (t) ; 0 ≤ t ≤ T -δ -a ) and (B δ (t) ; 0 ≤ t ≤ T δ a ) have the same law, it is equivalent to determine the distribution of (B δ (t) ; 0 ≤ t ≤ T δ a ) conditionally to

S = -inf t∈[0,T δ a ] B δ (t), (3.10)
T δ a being the first hitting time of level a :

T δ a = inf {t ≥ 0 , B δ (t) > a} , a > 0. Let ρ be the random time :

ρ = sup{t ∈ [0, T δ a ], B δ (t) = -S}. (3.11)
Time ρ plays a central role in our approach, as shows the following proposition. (3.12)

2) Conditionally to S = b, a. (B δ (t) ; 0 ≤ t ≤ ρ) and (B δ (t + ρ) -B δ (ρ) ; 0 ≤ t ≤ T δ a -ρ) are inde- pendent, b. (B δ (t) ; 0 ≤ t ≤ ρ
) is a process with the same law as a Brownian motion with drift -δ, stopped at its first hitting time of level -b and conditioned to stay less than a, c. (B δ (t + ρ) -B δ (ρ) ; 0 ≤ t ≤ T δ a -ρ) is a (positive) process distributed as the process Z (δ) , stopped at its first hitting time of level a + b Before starting the proof of Proposition 3.2, we would like to explain our approach.

Let (F t ) be the natural filtration generated by (B(t) ; t ≥ 0). Unfortunately the random time ρ is not stopping time. However ρ is a last exit time :

ρ = sup{t ∈ [0, T δ a ]; B δ (t) = sup 0≤u≤t B δ (u)}.
This leads us to apply the theory of enlargement of filtrations (See Protter [START_REF] Ph | Stochastic integration and differential equations[END_REF] Ch. VI).

Let (F ρ t ) be the smallest right-continuous filtration containing (F t ) (i.e. F t ⊂ F ρ t , for any t ≥ 0) such that ρ is a (F ρ t )-stopping time. Let :

Y ρ t = P (ρ ≤ t|F t ) , (3.13)
(Y ρ t ; t ≥ 0) is the optional projection of 1 {ρ≤t} ; t ≥ 0 (cf [9] p 371).
It is easy to check that (Y ρ t ; t ≥ 0) is a (F t )-submartingale. Let :

Y ρ t = M ρ t + A ρ t , (3.14)
be its Doob-Meyer decomposition, where (M ρ t ; t ≥ 0) denotes the martingale part, (A ρ t ; t ≥ 0) is a non-decreasing and adapted process such that A ρ 0 = 0. The processes (Y ρ t ; t ≥ 0) and (A ρ t ; t ≥ 0) will be given in Lemma 3.3. To determine the law of B δ (s); s ≤ ρ , we need the following result. Let (U t ; t ≥ 0) be a non-negative and (F t )-adapted process then ( [START_REF] Ph | Stochastic integration and differential equations[END_REF] p 371) :

E[U ρ ] = E ∞ 0 U t dA ρ t . (3.15)
Note that S = -

min 0≤s≤ρ B δ (s). Taking U t = F (B δ (s); s ≤ t) f (-min 0≤s≤t B δ (s)),
where F and f are measurable and non negative, we get :

E [F (B δ (s); s ≤ ρ) f (S)] = E ∞ 0 F (B δ (s); s ≤ t) f (-min 0≤s≤t B δ (s))dA ρ t .
Since (A ρ t ; t ≥ 0) is explicitly known, previous identity allows us to determine the law of B δ (s); s ≤ ρ conditionally to S, see Lemma 3.5.

To obtain the law of (B δ (t + ρ) -B δ (ρ) ; 0 ≤ t ≤ T δ a -ρ), we use ([9] Theorem 18, p 375) the following property :

B(t) = B(t) -χ(t) is a (F ρ t ) -Brownian motion (3.16)
where :

χ(t) = - t∧ρ 0 1 1 -Y ρ s d B, M ρ s + t 0 1 {ρ<s} 1 Y ρ s d B, M ρ s . (3.17)
This allows us to prove that (B δ (t + ρ) -B δ (ρ) ; 0 ≤ t ≤ T δ a -ρ) solves a stochastic differential equation of type (1.10).

Lemma 3.3 Let (Y ρ t ) t≥0 be the process defined by (3.13). 1. We have :

Y ρ t = e 2δB δ (t∧T δ a ) -e -2δB δ (t∧T δ a )
e 2δB δ (t∧T δ a ) -e -2δa

; t ≥ 0, (3.18)

where B δ (t) = -inf s≤t B δ (s).
2. (Y ρ t ) t≥0 is a (F t ) sub-martingale with Doob-Meyer decomposition (3.14) and

A ρ t = ln e 2δB δ (t∧T δ a ) -e -2δa 1 -e -2δa
(3.19)

M ρ t = 2δ t∧T δ a 0 1 e 2δB δ (s) -e -2δa
e -2δB δ (s) dB(s). b) Let t > 0 fixed. We have :

1 {ρ≤t} = 1 {ρ≤t≤T δ a } + 1 {T δ a ≤t} . (3.22)
If t ≤ T δ a , ρ ≤ t means that after the time t, B δ hits a before -B δ (t). Consequently :

P(ρ ≤ t ≤ T δ a |F t ) = e 2δB δ (t) -e -2δB δ (t) e 2δB δ (t) -e -2δa 1 {t≤T δ a } . (3.23)
Y ρ denotes the process :

Y ρ t = e 2δB δ (t) -e -2δB δ (t) e 2δB δ (t) -e -2δa for t ≥ 0. (3.24) Note that Y ρ T δ a = 1, therefore : Y ρ t = Y ρ t∧T δ a .
c) We know that e -2δB δ (t) = e -2δB(t)-2δ 2 t is a martingale. So using the classical rules of stochastic calculus we get :

d Y ρ t = 1 e 2δB δ (t) -e -2δa
e -2δB δ (t) 2δdB(t) + e -2δB δ (t) -e -2δa e 2δB δ (t) -e -2δa 2 2δe 2δB δ (t) dB δ (t).

(3.25)

Since the support of the random measure dB δ is included in {t ≥ 0 ; B δ (t) = -B δ (t)}, then (3.19) and (3.20) follow immediately.

The following result will be useful in the proof of Lemma 3.5 below.

Proposition 3.4 Let x > 0 and y > 0.

Then (B δ (t) ; 0 ≤ t ≤ T δ x ) condi- tioned by {T δ x < T δ -y } is distributed as (B -δ (t) ; 0 ≤ t ≤ T -δ x ) conditioned by {T -δ x < T -δ -y }.
Proof Let F be a test function. We have :

A = E F B(t) + δt ; 0 ≤ t ≤ T δ x |T δ x < T δ -y = E F B(t) + δt ; , 0 ≤ t ≤ T δ x 1 {T δ x <T δ -y } 1 P(T δ x < T δ -y ) . The r.v. F B(t) + δt ; 0 ≤ t ≤ T δ x 1 {T δ x <T δ -y } being F T δ
x measurable, and T δ x < ∞ a.s., Girsanov's theorem and (3.21) imply :

A = E F B(t), 0 ≤ t ≤ T 0 x 1 {T 0 x <T 0 -y } exp(δB(T 0 x ) - δ 2 2 T 0 x ) e -2δx -e 2δy 1 -e 2δy , = E F B(t), 0 ≤ t ≤ T 0 x 1 {T 0 x <T 0 -y } exp(- δ 2 2 T 0 x ) e -δx -e 2δy+δx
1 -e 2δy .

Replacing δ by -δ, we obtain similarly

à = E F B(t) -δt, 0 ≤ t ≤ T -δ x |T -δ x < T -δ -y , = E F B(t), 0 ≤ t ≤ T 0 x 1 {T 0 x <T 0 -y } exp(-δB(T 0 x ) - δ 2 2 T 0 x ) × e 2δx -e -2δy 1 -e -2δy , = E F B(t), 0 ≤ t ≤ T 0 x 1 {T 0 x <T 0 -y } exp(- δ 2 2 T 0 x )
e δx -e -2δy-δx 1 -e -2δy .

We are now able to describe the law of (B δ (t); 0 ≤ t ≤ ρ) conditionally to S.

Lemma 3.5 1) The r.v. S has a density function ϕ given by (3.12).

2) Conditionally to S = b, (B δ (t) ; 0 ≤ t ≤ ρ) is a process with the same law as a Brownian motion with drift -δ, stopped at its first hitting time of level -b and conditioned to stay less than a.

Proof of Lemma 3.5 1) Let f and F be two non-negative test functions and

∆ = E [f (S)F (B δ (s), 0 ≤ s ≤ ρ)].
We have :

∆ = E f (B δ (ρ))F (B δ (s), 0 ≤ s ≤ ρ) . Setting U t = f (B δ (t))F (B δ (s), 0 ≤ s ≤ t), we have : ∆ = E [U ρ ]. Since (U t , t ≥ 0
) is a predictable process, property (3.15) implies that :

∆ = E ∞ 0 U t dA ρ t .
Using (3.19), we obtain

∆ = 2δE T δ a 0 f (B δ (t))F (B δ (s), 0 ≤ s ≤ t) e 2δB δ (t) e 2δB δ (t) -e -2δa
dB δ (t) .

In particular, if F = 1, the change of variable : "x = B δ (t)" yields to :

E [f (S)] = 2δE S 0 f (x) e 2δx e 2δx -e -2δa dx , = 2δE ∞ 0 f (x)
e 2δx e 2δx -e -2δa P(S ≥ x) dx .

Φ(x) = ∞ x ϕ(s) ds = 1 -e -2δa e 2δx -e -2δa .
(3.28)

This implies (3.12).

2) We now study the law of (B δ (t) ; 0 ≤ t ≤ ρ). (T δ -x ; x ≥ 0) being the right continuous inverse of (B δ (t) ; t ≥ 0), then

∆ = 2δ ∞ 0 f (x)E F (B δ (s), 0 ≤ s ≤ T δ -x )1 {T δ -x <T δ a }
e 2δx e 2δx -e -2δa dx,

= 2δ ∞ 0 f (x)E F (B δ (s), 0 ≤ s ≤ T δ -x )|T δ -x < T δ a × e 2δx e 2δx -e -2δa P(T δ -x < T δ a ) dx, = ∞ 0 f (x)E F (B δ (s), 0 ≤ s ≤ T δ -x )|T δ -x < T δ a ϕ(x) dx.
We now use the result of Proposition 3.4 to obtain the law of (B δ (t) ; 0 ≤ t ≤ ρ) conditionally to S.

Proof of Proposition 3.2

It remains to determine the law of (B δ (t + ρ) -B δ (ρ), 0 ≤ t ≤ T δ a -ρ) and to prove that this process is independent of (B δ (t); 0 ≤ t ≤ ρ).

Let (F ρ t ) t≥0 denote the smallest filtration including (F t ) t≥0 for which ρ is a stopping time. ρ is an honest time since, for any t > 0, on {ρ < t}, ρ coincides with a F t -measurable r.v. (see [START_REF] Ph | Stochastic integration and differential equations[END_REF] p 370) That allows to use (3.16). Combining (3.17), (3.20) with (3.18), we obtain χ(t) = -t∧ρ 0 (e 2δB δ (s) -e -2δa )2δe -2δB δ (s) (e -2δB δ (s) -e -2δa )(e 2δB δ (s) -e -2δa ) ds

+ t∧T δ a 0 1 {ρ<s}
(e 2δB δ (s) -e -2δa )2δe -2δB δ (s) (e 2δB δ (s) -e -2δB δ (s) )(e 2δB δ (s) -e -2δa ) ds,

χ(t) = - t∧ρ 0 2δe -2δB δ (s) e -2δB δ (s) -e -2δa ds + t∧T δ a 0 1 {ρ<s} 2δe -2δB δ (s) e 2δB δ (s) -e -2δB δ (s)
ds.

Since B δ (s) = S, then for any ρ ≤ s ≤ T δ a , we have

χ(t + ρ) -χ(ρ) = 2δ (t+ρ)∧T δ a ρ e -2δB δ (s) e 2δS -e -2δB δ (s) ds; 0 ≤ t ≤ T δ a -ρ.
Let (W (t)) and Z(t) be the processes

W (t) = B(t+ρ)-B(ρ) ; t ≥ 0, Z(t) = B δ (ρ+t)-B δ (ρ) = B δ (ρ+t)+S ; 0 ≤ t ≤ T δ a -ρ.
According to (3.16), for any t ∈ [0, T δ a -ρ] we have :

Z(t) = W (t) + δt + 2δ t 0 e -2δB δ (s+ρ) e 2δS -e -2δB δ (s+ρ) ds. = W (t) + δ t 0 coth(δ Z(s))ds. We know that B(t) is a F ρ ρ+t -Brownian motion. Therefore (W (t)) is in- dependent of F ρ ρ . Since S = -B δ (ρ) and (B δ (t); 0 ≤ t ≤ ρ) are F ρ ρ -measurable, then (W (t)) is independent of S and (B δ (t); 0 ≤ t ≤ ρ). Furthermore, T δ a -ρ = inf{t ≥ 0 ; Z(t) = a + b}, with b = S.
To conclude the proof of Proposition 3.2 (part 1 and 3), we use the fact that stochastic differential equation (1.10) admits a unique non negative strong solution. We can prove this property with a similar method as for a Bessel process (we apply Yamada Watanabe Theorem to the square of Z (δ) ).

Remark 3.6 Conditionally to {S ≤ z}, the density of S is :

δ sh δa sh δ(z + a) sh δz 1 sh 2 δ(x + a) 1 {0≤x≤z} .
(3.29) Formula (3.29) is really a direct consequence of (3.28).

Proofs of Theorem 2.4 and Proposition 2.1

The processes (-B δ (t) ; t ≥ 0) and (B -δ (t) ; t ≥ 0) have the same law. Consequently Proposition 3.2 admits a version where (B δ (t) ; 0

≤ t ≤ T δ a ) (resp. S) is replaced by (B -δ (t) ; 0 ≤ t ≤ T -δ -a ) (resp. sup 0≤t≤T -δ -a
B -δ (t)). Then Theorem 2.4

and Proposition 2.1 can be proved by induction on k.

4 First range time, range at a fixed time, asymptotic results

In this section, we focus on B δ stopped at its first range time θ δ (a) defined by (1.5). The first result concerns the law of (B δ (t) ; 0 ≤ t ≤ θ δ (a)). We determine the joint law of (B δ (t) ; 0 ≤ t ≤ θδ (a)) and (B δ (t + θδ (a))

-B δ ( θδ (a)) ; 0 ≤ θ δ (a) -θδ (a)) where θδ (a) =        sup{t ≤ θ δ (a) ; B δ (t) = inf 0≤u≤θ δ (a) B δ (u)} if B δ (θ δ (a)) > 0 sup{t ≤ θ δ (a) ; B δ (t) = sup 0≤u≤θ δ (a)
B δ (u)} otherwise . A similar result was obtained by Vallois [START_REF] Vallois | Diffusion arrêtée au premier instant où l'amplitude atteint un niveau donné[END_REF] when δ = 0.

f (u) = 2δ (1 -e -2δa ) 2 (e 2δu -e -2δa )1 {-a≤u<0} + 2δ (e 2δa
As for the proof of Proposition 3.2, a straightforward approach may be developed using enlargement of filtrations, since θδ (a) is an honest time. More precisely, the optional projection of 1 { θδ (a)≤t} ; t ≥ 0 can be decomposed as :

P θδ (a) ≤ t|F t = M t + A t ,
where (M t ; t ≥ 0) is a martingale and :

A t = 2δe -δa e δa -e -δa B δ (t ∧ θ(a)) - 2δe δa e δa -e -δa B δ (t ∧ θ(a)). (4.3)
The details are left to the reader. Proposition 4.1 allows us to determine the Laplace transform of θ δ (a). Let σ δ (a) be the stopping time

σ δ (a) = inf s ≥ 0 , Z (δ) (s) > a . (4.4)
Proposition 4.2 The Laplace transforms of θδ (a), σ δ (a) and θ δ (a) are given by :

E e -λ θδ (a) = δ λ 2λ + δ 2 coth(a 2λ + δ 2 ) coth(δa) -δ - √ 2λ + δ 2 sh(a √ 2λ + δ 2 ) sh(δa) , (4.5
)

E e -λσ δ (a) = √ 2λ + δ 2 δ sh(δa) sh(a √ 2λ + δ 2 ) , (4.6
)

E e -λθ(a) = √ 2λ + δ 2 λ √ 2λ + δ 2 ch(a √ 2λ + δ 2 ) ch(δa) sh 2 (a √ 2λ + δ 2 ) - δ sh(δa) sh(a √ 2λ + δ 2 ) - √ 2λ + δ 2 sh 2 (a √ 2λ + δ 2 ) (4.7)
where λ ≥ 0.

Remark 4.3 1) In [START_REF] Vallois | Range reliability in random walks[END_REF], using approximation by random walks, it is proved :

E [θ δ (a)] = a 2 2 f (δa), (4.8)
where

f (x) = 1 x 2 (x -x coth x + 1)(x + x coth x -1) ; x > 0. (4.9)
Moreover f is decreasing, in particular :

E [θ δ (a)] ≤ E [θ 0 (a)] = a 2 2 . (4.10)
The variance of θ δ (a) is computed :

Var [θ δ (a)] = a 4 12 g 1 (δa)g 2 (δa), (4.11)
with

g 1 (x) = 3 sh 2 x -x 2 x 2 sh 2 x ; g 2 (x) = x 2 coth 2 x + 4x coth x -5 -x 2 x 2 ; x > 0.
Moreover,

Var (θ δ (a)) ≤ Var (θ 0 (a)) = a 4 12 .
2) It is possible to prove Propositions 4.1 and 4.2 using Theorems 2.2, 2.4 and the identity : Ω = ∪ n≥0 ρ n+1 ≤ θδ (a) < ρ n , where the sequence (ρ n ) n≥0 is defined at the beginning of Section 2. It is however easier and shorter to prove directly Propositions 4.1 and 4.2.

Proof of Proposition 4.2

We only give the main ideas of the proof, the details being left to the reader.

Let f λ be the function :

f λ (x) = δ √ 2λ + δ 2 sh(x √ 2λ + δ 2 ) sh(δx) ; x > 0. (4.12)
It is easy to check that f λ is an eigenfunction of the infinitesimal generator associated with Z (δ) :

1 2 f λ (x) + δ coth(δx)f λ (x) = λf λ (x) ; x > 0. (4.13)
The function f λ being locally bounded on [0, +∞), then

E 0 [exp(-λσ δ (a))] = f λ (0) f λ (a) . (4.14)
This proves (4.6).

Recall (cf formulas 2.1.4 and 2.2.4 p. 295 of [START_REF] Borodin | Handbook of Brownian motion-facts and formulae[END_REF])

E e -λT δ u |T δ u < T δ u+a = sh((a + u) √ 2λ + δ 2 ) sh(a √ 2λ + δ 2 ) sh(δa) sh(δ(a + u)) , (4.15
)

E e -λT δ u |T δ u < T δ u-a = sh((a -u) √ 2λ + δ 2 ) sh(a √ 2λ + δ 2 )
sh(δa) sh(δ(a -u))

. By inversion of the Laplace transform of θ δ (a), the authors in [START_REF] Chong | The ruin problem and cover times of asymmetric random walks and Brownian motions[END_REF] have computed the probability density function of this r.v.

We develop an alternative approach based on the knowledge (c.f. [START_REF] Borodin | Handbook of Brownian motion-facts and formulae[END_REF], formula 1.15.8 (1) p. 271) of the joint distribution of (B δ (t), R δ (t)). By tedious calculations (see Section 5), we determine the probability distribution function of R δ (t).

Relation (1.3) allows us to obtain the probability distribution function of θ δ (a) and the rate of decay of P(θ δ (a) > t), as t → ∞. Proposition 4.4 Let a > 0, t > 0 and C k = k 2 π 2 + a 2 δ 2 , k ∈ N. Then :

P (R δ (t) < a) = ∞ k=1 4k 2 π 2 C 2 k exp - C k t 2a 2 1 -(-1) k ch(δa) × 1 + k 2 π 2 t a 2 - 4a 2 δ 2 C k -(-1
) k aδ sh(δa) .

(4.17)

P (θ δ (a) > t) ∼ t→∞ 4π 4 (1 + ch(δa)) a 2 (π 2 + a 2 δ 2 ) t exp - π 2 + a 2 δ 2 2a 2 t . (4.18)
Remark 4.5 1) Formula (4.17) has been obtained in [START_REF] Siebenaler | Hitting times of the range of a birth and death random walk[END_REF], the approach being different.

2) Taking the a-derivative in (4.17) gives the density function of R δ (t).

3) Relation (1.3) implies that :

P (θ δ (a) > t) = ∞ k=1 4k 2 π 2 C 2 k exp - C k t 2a 2 1 -(-1) k ch(δa) × 1 + k 2 π 2 t a 2 - 4a 2 δ 2 C k -(-1) k aδ sh(δa) . ( 4
θ δ (a) a a.s. -→ a→∞ 1 δ , (4.20)
R δ (t) t a.s. -→ t→∞ δ, (4.21)
δ 3/2 √ a θ δ (a) a - 1 δ d -→ a→∞ N (0, 1) , (4.22)
√ t R δ (t) t -δ d -→ t→∞ N (0, 1). ( 4

.23)

Proof 1) We first examine (4.20). By (4.7), we have :

E e -λ a θ δ (a) = 2λ a + δ 2 λ a 2λ a + δ 2 ch(a 2λ a + δ 2 ) ch(δa) sh 2 (a 2λ a + δ 2 ) - δ sh(δa) sh(a 2λ a + δ 2 ) - 2λ a + δ 2 sh 2 (a 2λ a + δ 2 ) . E e -λθ δ (a) a = a δ + λ δa + o 1 a λ δ + λ δa + o 1 a ch(δa + λ δ + o(1)) ch(δa) sh 2 (δa + λ δ + o(1)) - δ sh(δa) sh(δa + λ δ + o(1)) + δ + λ δa + o 1 a sh 2 (δa + λ δ + o(1))
.

The limit of the two first terms in the bracket is easy to determine : sh(δa) sh(δa + As for the third term, it may be neglected (being equivalent to 4δ exp(-2δa)).

Hence θ δ (a)/a converges in distribution to the constant 1/δ. This implies that θ δ (a)/a converges in probability to 1/δ. Since (R δ (t) ; t ≥ 0) is the right continuous inverse of (θ δ (a) ; a ≥ 0), R δ (t)/t converges in probability to δ.

(R δ (t) ; t ≥ 0) is a subadditive process (cf [START_REF] Durrett | Probability : theory and examples[END_REF], example 6.2 p.320, in the discrete case). The subadditive ergodic theorem implies that R δ (t)/t converges a.s., as t → ∞. As a result, R δ (t)/t converges a.s. towards δ when t → ∞.

2) To prove (4.22), we use the characteristic function of θ δ (a). This function can be explicitly determined through (4.7) and analytic continuation argument. In the discrete case (cf proof of Theorem 20 of [START_REF] Vallois | The range of a simple random walk on Z[END_REF]) a detailed approach is developed.

3

) We claim that R δ (t) t -δ √ t converges in distribution to N (0, 1) when t → +∞.
Let a > 0 and t > 0 such that a √ t + δt > 0. We have :

p = P R δ (t) t -δ √ t < a = P (R δ (t) < s) ,
where s = a √ t + δt.

Property (1.3) implies that

p = P (θ δ (s) > t) = P θ δ (s) s - 1 δ √ s > u(t) ,
where

u(t) = t s - 1 δ √ s = - a δ δ + a √ t .
Since t → +∞ implies s → +∞, and u(t) ∼ t→∞ -a δ 3/2 , (4.23) follows immediately.

of [START_REF] Borodin | Handbook of Brownian motion-facts and formulae[END_REF] : We decompose C i (δ) + C i (-δ), 1 ≤ i ≤ 3, as follows :

P (R δ (t) < a) = 1 √ 2πt e - δ 2 t 2 k∈Z µ k , with µ k = a -a 2k
C 1 (δ) + C 1 (-δ) = -aA 1 + (-1) k A 2 , C 2 (δ) + C 2 (-δ) = a 2 A 1 + a 2 A 3 -(-1) k a 2 A 4 , C 3 (δ) + C 3 (-δ) = -a 3 A 3 -(-1) k a 3 A 4 -2a 3 A 5 + (-1) k 2a 3 A 6 ,
where we have set :

A 1 = 1 δa + ikπ + 1 -δa + ikπ = -2ikπ C k

22

Coming back to (5.4), we easily obtain successively : In particular :

I k exp k 2 π 2 t 2a 2 = ikπt a 2 - k 2 π
(I k + I -k ) exp k 2 π 2 t 2a 2 = k 2 π 2 C 2 k (1 -(-1) k ch(δa)) × 4t(k 2 π 2 + δ 2 a 2 ) a 2 -4 (a 2 δ 2 -k 2 π 2 )t a 2 + 8k 2 π 2 -24a 2 δ 2 C k -(-1) k 8 k 2 π 2 C 2 k aδ sh(δa).
This achieves the proof of Proposition 4.4.

δα δ an 0 e

 0 δy sh δy dy = e -δα 4 e 2δan -1 -2δa n .

Proposition 3 . 2

 32 Let S and ρ be r.v.'s defined by (3.10) and (3.11). Then, 1) S has a density function ϕ ϕ(x) = δe δa sh δa 1 sh 2 δ(x + a) 1 {x≥0} .

(3. 20 )

 20 Proof a) Recall a classical result concerning hitting times of Brownian motion with drift (see for instance Borodin and Salminen[START_REF] Borodin | Handbook of Brownian motion-facts and formulae[END_REF] formula 2.1.2 p. 295)P T δ b < T δ a B δ (0) = x) = e -2δx -e -2δae -2δb -e -2δa , for x between a and b.(3.21)

  P(S ∈ dx) = 2δe 2δx e 2δx -e -2δa P(S ≥ x) dx. Hence, S has a density function ϕ and :ϕ(x) = 2δe 2δx e 2δx -e -2δa ∞ x ϕ(s) ds.(3.27)(3.27) may be written as a linear differential equation with respect to Φ(x) = ∞ x ϕ(s) ds. It is easy to solve explicitly since Φ(0) = 1 :

(4. 1 ) 4 . 1 1.

 141 Proposition The r.v. B δ ( θδ (a)) has density function f :

( 4 . 2 ) 2 .

 422 -1) 2 (e 2δa -e 2δu )1 {0≤u≤a} . Conditionally to {B δ ( θδ (a)) = u}, a. the processes (B δ (s) ; 0 ≤ s ≤ θδ (a)) and (B δ (s + θδ (a)) -u ; 0 ≤ s ≤ θ(a) -θδ (a)) are independent. b. (B δ (s) ; 0 ≤ s ≤ θδ (a)) has the same law as a Brownian motion with drift δ stopped at its first hitting time of level u, conditioned by hitting u before u -sgn(u)a. c. (|B δ (s + θδ (a)) -u| ; 0 ≤ s ≤ θ(a) -θδ (a)) has the same law as the process Z (δ) defined by (1.10), stopped at its first hitting time of level a.

  (4.16) (4.5) is a direct consequence of (3.21) and Proposition 4.1.

  λ δ + o(1)) = e δa -e -δae δa+ λ δ +o(1) -e -δa-

2 + 2a 3 (- 1 ) k e δa - 1 (aδ + ikπ) 3 ., I k exp k 2 π 2 t 2a 2 = - ikπt a 3 (

 231133 ConsequentlyC 1 (δ) + C 1 (-δ)) -k 2 π 2 t a 4 (C 2 (δ) + C 2 (-δ)) + ikπ a 3 (C 3 (δ) + C 3 (-δ)).(5.6)

2 k A 4 2 k A 5 3 k A 6 = (2δ 3 a 3 -

 2425363 + ikπ) 2 = 2(a 2 δ 2 -k 2 π 2 ) C = e δa (δa + ikπ) 2 + e -δa (-δa + ikπ) 2 = 2 ch(δa)(a 2 δ 2 -k 2 π 2 ) -4iaδkπ sh(δa) C = 1 (δa + ikπ) 3 + 1 (-δa + ikπ) 3 = -6δ 2 a 2 ikπ + 2ik 3 π 3 C = e δa (δa + ikπ) 3 + e -δa (-δa + ikπ) 3 , 6k 2 π 2 δa) sh(δa) + i(-6δ 2 a 2 kπ + 2k 3 π 3 ) ch(δa) C 3 k .

  .19)Again, taking the t-derivative in (4.19), we obtain the density function of θ δ (a). However that series expension is more complicated than (4.19) (c.f. also Theorem 9 of[START_REF] Chong | The ruin problem and cover times of asymmetric random walks and Brownian motions[END_REF]) : it is more convenient to use probability distribution function instead of probability density function.The law of θ δ (a) or R δ (t) being complicated, it seems natural to consider the asymptotic behaviour of θ δ (a) (resp. R δ (t)) when a goes to +∞ (resp. t → +∞). The asymptotic comportments of θ δ (a) and R δ (t) are given by the following convergence results :

	Proposition 4.6

  2 t a 2 A 1 -(-1) k ikπt a 2 A 2 -k 2 π 2 t a 2 + ikπ A 3 + (-1) k k 2 π 2 t a 2 -ikπ A 4 -2ikπA 5 + (-1) k 2ikπA 6 = Re I k exp k 2 π 2 t 2a 2 = 2k 2 π 2 t C k a 2 -(-1) k 2k 2 π 2 t ch(δa) C k a 2 -2 (a 2 δ 2 -k 2 π 2 )k 2 π 2 t C 2 k a 2 +(-1) k 2 ch(δa)(a 2 δ 2 -k 2 π 2 )k 2 π 2 t -4k 2 π 2 a 3 δ sh(δa) C 2 k a 2 + 4k 4 π 4 -12k 2 π 2 a 2 δ 2 C 3 k +(-1) k (12δ 2 a 2 k 2 π 2 -4k 4 π 4 ) ch(δa) C 3

k

.

INRIA, Projet OMEGA,

route des Lucioles, BP93, 06902 Sophia-Antipolis, France, Etienne.Tanre@sophia.inria.fr 2 Institut Élie Cartan, BP 239, 54506 Vandoeuvre-lès-Nancy Cedex, France, Pierre.Vallois@iecn.u-nancy.fr

Proof ofProposition 4.4 In order to calculate the probability distribution function of R δ (a), we perform the following calculations. Beginning with the formula 1.15.8 (1) p. 271


Alban Quadrat
email: alban.quadrat@inria.fr

Towards an effective study of the algebraic parameter estimation problem

Keywords: Parameter estimation problem, algebraic systems theory, symbolic computation, annihilators, orthogonal polynomials, expansion into a basis, ring of differential operators

The paper aims at developing the first steps toward a symbolic computation approach to the algebraic parameter estimation problem defined by Fliess and Sira-Ramirez and their coauthors. In this paper, within the algebraic analysis approach, we first give a general formulation of the algebraic parameter estimation for signals which are defined by ordinary differential equations with polynomial coefficients such as the standard orthogonal polynomials (Chebyshev, Jacobi, Legendre, Laguerre, Hermite, . . . polynomials). Based on a result on holonomic functions, we show that the algebraic parameter estimation problem for a truncated expansion of a function into an orthogonal basis of L 2 defined by orthogonal polynomials can be studied similarly. Then, using symbolic computation methods such as Gröbner basis techniques for (noncommutative) polynomial rings, we first show how to compute ordinary differential operators which annihilate a given polynomial and which contain only certain parameters in their coefficients. Then, we explain how to compute the intersection of the annihilator ideals of two polynomials and characterize the ordinary differential operators which annihilate a first polynomial but not a second one. These results, which are at the core of the algebraic parameter estimation, are implemented in the Non-A package built upon the OreModules software.

INTRODUCTION

In this paper, we consider the problem of estimating constant parameters θ of a signal x when the signal y(t) = x(θ, t) + γ(t) + (t)

(1) is observed, and where γ is a perturbation and a zero-mean noise. Many different approaches have been developed for this problem in the literature of signal processing and control theory. In this paper, we follow the algebraic approach developed by Fliess, Sira-Ramirez, Mboup, and their coauthors. See [START_REF] Fliess | An algebraic framework for linear identification[END_REF]; [START_REF] Sira-Ramírez | Algebraic Identification and Estimation Methods in Feedback Control Systems[END_REF]; [START_REF] Mboup | Numerical differentiation with annihilators in noisy environment[END_REF]; [START_REF] Ushirobira | Algebraic parameter estimation of a biased sinusoidal waveform signal from noisy data[END_REF]Ushirobira et al. ( , 2013)), and the references therein. In the absence of noise, i.e., = 0, the algebraic approach provides explicit and exact formulas for the constant parameters θ in terms of integrals of the observed signal y. The use of integrals helps to filter the noise . The explicit expressions of the parameters θ can then be used to do real-time estimation.

In the algebraic approach to parameter estimation, the dynamics of the signal x is defined by a time-invariant ordinary differential (OD) system [START_REF] Fliess | An algebraic framework for linear identification[END_REF]; [START_REF] Sira-Ramírez | Algebraic Identification and Estimation Methods in Feedback Control Systems[END_REF]; [START_REF] Mboup | Numerical differentiation with annihilators in noisy environment[END_REF]; [START_REF] Ushirobira | Algebraic parameter estimation of a biased sinusoidal waveform signal from noisy data[END_REF]Ushirobira et al. ( , 2013))). In [START_REF] Fliess | Algebraic change-point detection[END_REF], there is a remark about the possibility to extend this framework to ODEs with polynomial coefficients (see Remark 1 in page 133). Being not aware of this remark but inspired by works developed in the symbolic computation community [START_REF] Chyzak | Effective algorithms for parametrizing linear control systems over Ore algebras[END_REF]), in [START_REF] Ushirobira | Algebraic estimation of a biased and noisy continuous signal via orthogonal polynomials[END_REF], we show that the algebraic approach can be extended to signal x which is defined by an OD system with polynomial coefficients. Indeed, within the algebraic approach to parameter estimation, the signal is studied in the frequency domain by means of the operational calculus. The Laplace transform can still be used for OD systems with polynomial coefficients since the time variable t is then transformed into -∂ s , where ∂ s = d ds is the derivation with respect to the Laplace variable s. Within the algebraic analysis approach [START_REF] Kashiwara | Foundations of Algebraic Analysis[END_REF]; [START_REF] Chyzak | Effective algorithms for parametrizing linear control systems over Ore algebras[END_REF]), the Laplace transform is an automorphism of the Weyl algebra A 1 (k) = k t, ∂ t | ∂ t t = t ∂ t + 1 of OD operators with polynomial coefficients in t over a field k of characteristic zero (e.g., k = Q, R), i.e., the k-linear map L : A 1 (k) -→ A 1 (k) defined by L (t) = -∂ s and L (∂ t ) = s is an automorphism of k-algebras, where the target Weyl algebra is in the complex variable s, i.e., A 1 (k) = k s, ∂ s | ∂ s s = s ∂ s + 1 . The algebraic approach to parameter estimation can then be developed similarly to the linear time-invariant OD systems. In [START_REF] Ushirobira | Algebraic estimation of a biased and noisy continuous signal via orthogonal polynomials[END_REF], we initiated the study of signals x which are defined by dynamical systems with polynomial coefficients or signals which can be expanded into orthogonal bases of L 2 defined by orthogonal polynomials such as Hermite, Jacobi, Laguerre, . . . , polynomials [START_REF] Abramowitz | Handbook of mathematical functions with formulas, graphs, and mathematical tables[END_REF]; [START_REF] Chyzak | Effective algorithms for parametrizing linear control systems over Ore algebras[END_REF]). Note that the algebraic parameter estimation problem was studied for Taylor expansion series [START_REF] Mboup | Numerical differentiation with annihilators in noisy environment[END_REF]) and initiated for Fourier expansions, i.e., expansions into the orthogonal basis {e i n θ } n∈Z of L 2 (T), where T = {z ∈ C | |z| = 1} is the unit torus of C [START_REF] Ushirobira | Algebraic parameter estimation of a biased sinusoidal waveform signal from noisy data[END_REF](Ushirobira et al. ( , 2013))).

The first goal of this paper is to give a formulation of the algebraic parameter estimation for signals x which satisfy an OD equation with polynomial coefficients. Then, we show that the case of a signal x which is a truncated expansion of a function into a basis defined by a family of orthogonal polynomials can similarly be studied. Finally, using computer algebra techniques (e.g., Gröbner basis techniques), we initiate an effective study of the computations of annihilators of polynomials, a problem which plays a key role in the algebraic parameter estimation.

ALGEBRAIC PARAMETER ESTIMATION

Classes of signals x under study

Let a signal x satisfy the following OD equations (ODE)

n i=0 a i (t) x (i) (t) = 0, ( 2
)
where the a i 's are polynomials in t and with coefficients in a commutative ring K of constants (namely, ċ = 0 for c ∈ K), i.e., a i ∈ K[t] for i = 0, . . . , n. Using the following standard results on the Laplace transform L

• L f (n) (s) = s n L (f ) (s) - n-1 i=0 s n-i-1 f (i) (0), • L (t n f ) (s) = (-1) n ∂ n s (L (f ) (s)
), with the notation x = L (x), we obtain:

n i=0 a i (-∂ s )   s i x(s) - i-1 j=0 s i-j-1 x (j) (0)   = 0.
The above identity can be rewritten as:

n i=0 a i (-∂ s ) s i x(s)- n i=0 i-1 j=0 a i (-∂ s ) s i-j-1 x (j) (0) = 0.
Expending the second term of the above identity, we get

n i=0 i-1 j=0 a i (-∂ s ) s i-j-1 x (j) (0) = n j=1   n i=j a i (-∂ s ) s i-j   x (j-1) (0) = n-1 k=0 n i=k+1 a i (-∂ s ) s i-k-1 x (k) (0),
and thus, we obtain:

n i=0 a i (-∂ s ) s i x(s) - n-1 k=0 n i=k+1 a i (-∂ s ) s i-k-1 x (k) (0) = 0. If we note R = n i=0 a i (-∂ s ) s i , S k = - n i=k+1 a i (-∂ s ) s i-k-1 ,
for k = 0, . . . , n -1, then the above identity becomes:

R(s, ∂ s ) x(s) + n-1 k=0 S k (s) x (k) (0) = 0.
(3)

Let D = K[s] ∂ s be the noncommutative polynomial ring of OD operators with coefficients in the commutative ring A = K[s], namely, the ring of noncommutative polynomials of the form n i=0 a i ∂ i s , where a i ∈ A for i = 0, . . . , n, which satisfy the following relation:

∀ a ∈ A, ∂ s a = a ∂ s + da ds .
Now, using the operator identity ∂ s s = s ∂ s + 1 in D, which corresponds to the standard Leibniz identity, i.e.,

(∂ s s)(z(s)) = d ds (s z(s)) = s d ds z(s)+z(s) = (s ∂ s +1) z(s),
the terms a i (-∂ s ) s k can be rewritten as an OD operator of the form m l=0 b l (s) ∂ l s , i.e., as an element of D. For more details, see, e.g., [START_REF] Kashiwara | Foundations of Algebraic Analysis[END_REF][START_REF] Mcconnell | Noncommutative Noetherian Rings[END_REF]. Remark 1. Within the framework developed in [START_REF] Fliess | An algebraic framework for linear identification[END_REF], the algebraic estimation problem is stated for a general time-invariant linear control system defined by:

n i=0 a i x (i) (t) = n-1 j=0 b j u (j) (t).
The presence of an input u does not bring substantial differences from the input free case (u = 0). Indeed, an extra term coming from the contribution of u has to be added. To simplify, in this paper, we shall only consider the case u = 0 as well as the single-output case. For the general case, we refer the reader to Quadrat (2017).

Let us now consider the signal z = x + γ, where γ is a perturbation which admits a Laplace transform γ = L (γ). Then, combining z = x + γ with (3), we obtain:

R(s, ∂ s ) z(s) + n-1 k=0 S k (s) x (k) (0) -R(s, ∂ s ) γ(s) = 0.
To simplify the notation, let ϑ k = x (k) (0), k = 0, . . . , n-1, and

S = n-1 k=0 S k (s) ϑ k ∈ K[s]. Hence, we obtain: R(s, ∂ s ) z(s) + S(s) -R(s, ∂ s ) γ(s) = 0. ( 4
)
When γ satisfies an ODE, then the perturbation γ is called structured (see [START_REF] Fliess | An algebraic framework for linear identification[END_REF]). If γ satisfies an ODE of order m with polynomial coefficients, i.e., m i=0 a i (t) γ (i) (t) = 0, where a i ∈ K[t] for i = 0, . . . , m, then, similarly as above, γ satisfies the following equation R (s, ∂ s ) γ(s) + S (s) = 0, (5) where R = m i=0 a i (-∂ s ) s i and:

S k = - m i=k+1 a i (-∂ s ) s i-k-1 , S = m-1 k=0 S k (s) γ (k) (0). Example 1. If γ is a constant biais, i.e., γ(t) = γ H(t),
where H is the Heaviside distribution (i.e., H(t) = 1 for t > 0 and 0 for t < 0) and γ is a real constant, then, within algebraic analysis, it is well-known that this distribution is the solution of the ODE with polynomial coefficients t ∂ t (γ(t)) = 0. This can be checked again by using the Laplace transform which yields ∂ s (s γ(s)) = 0 and its solution γ(s) = γ/s is the Laplace transform of γ(t). Similarly, if γ(t) is an impulse, i.e., γ(t) = γ δ(t), where δ is the Dirac distribution and γ is a real constant, then, within algebraic analysis, it is known that this distribution satisfies the equation t γ(t) = 0. This can be checked again since, using the Laplace transform, we get ∂ s γ(s) = 0, i.e., γ(s) = γ, which is the Laplace transform of γ(t) = γ δ(t).

To simplify, in this paper, we shall only consider the perturbation γ to be a constant biais (see Example 1). For the general case, see Quadrat (2017). Using the Laplace transform of γ, i.e., γ(s) = γ/s, and (4), we then get:

R(s, ∂ s ) z(s) + S(s) -R(s, ∂ s ) γ s = 0. ( 6
)
We can clean the denominator of ( 6) to obtain an OD operator with polynomial coefficients in s, i.e., an element of the ring D. We first have:

R(s, ∂ s ) γ s = γ n i=1 a i (-∂ s ) s i-1 + a 0 (-∂ s ) 1 s .
If the degree of a 0 (-∂ s ) in ∂ s is n 0 , then using the identity

d n0 ds n0 s -1 = -n 0 ! s -(n0+1) , we get that a 0 (-∂ s ) s -1 = d 0 (s)/s n0+1 , where d 0 ∈ k[s].
Therefore, multiplying (6) by s n0+1 , we obtain:

s n0+1 R(s, ∂ s ) z(s) + s n0+1 S(s) -γ s n0+1 n i=0 a i (-∂ s ) s i-1 = 0. If we note                  P = s n0+1 R(s, ∂ s ) ∈ D = K[s] ∂ s , Q = s n0+1 S(s) = s n0+1 n-1 k=0 S k (s) ϑ k ∈ K[ϑ 1 , . . . , ϑ n-1 , s], Q = -s n0+1 n i=0 a i (-∂ s ) s i-1 γ ∈ K[γ, s],
(7) then we finally obtain that the signal z satisfies:

P (s, ∂ s ) z(s) + Q(s) + Q(s) = 0. ( 8
) Example 2. If x = m k=0 a k
k! t k is a polynomial of degree m, then we have x (m+1) (t) = 0, and thus n = m + 1, a n = 1, a i = 0 for i = 0, . . . , n-1, ϑ k = x (k) (0) for k = 0, . . . , n-1, n 0 = 0, and:

P = s m+2 , Q = -s m k=0 s m-k ϑ k , Q = -γ s m+1 . Example 3. If x = e -i 2 π
T n t is an element of the Fourier orthogonal basis {e -i 2 π T n t } n∈Z of L 2 (R), then x satisfies the first order ODE ẋ + i ν n x = 0, where ν = 2 π/T , i.e., a 1 = ∂ s , a 0 = i ν n and n 0 = 0. Thus, we get (8), where:

       ϑ 0 = x(0), P = s (s + i ν n) , Q = -s ϑ 0 , Q = -γ (s + i ν n) . Example 4. Let x be a signal satisfying a second ODE a 2 (t) ẍ(t) + a 1 (t) ẋ(t) + a 0 (t) x(t) = 0, (9) where a i ∈ K[t], i = 0, 1, 2. If deg t a 0 = n 0 , then (7) gives:                ϑ 0 = x(0), ϑ 1 = ẋ(0), P = s n0+1 a 0 (-∂ s ) 1 + a 1 (-∂ s ) s + a 2 (-∂ s ) s 2 , Q = -s n0+1 (Q 0 (s) ϑ 0 + Q 1 (s) ϑ 1 ) , = -s n0+1 ((a 1 (-∂) 1 + a 2 (-∂ s ) s) ϑ 0 + a 2 (-∂ s ) ϑ 1 ) , Q = -γ s n0+1 (a 0 (-∂ s ) s -1 + a 1 (-∂ s ) 1 + a 2 (-∂ s ) s).
Now, let us consider the following case    a 2 (t) = a 22 t 2 + a 21 t + a 20 , a 1 (t) = a 11 t + a 10 , a 0 (t) = a 00 , where a 22 , a 21 , a 20 , a 11 , a 10 , a 00 ∈ K. Then, we have n 0 = deg t a 0 = 0 and: 2021(1) If x = A sin(ω t+φ), then a 2 = 1, a 1 = 0 and a 0 = ω 2 , and thus we obtain ( 8), where:

         P = s (a 22 s 2 ∂ 2 s + s (-a 21 s + 4 a 22 -a 11 ) ∂ s + a 20 s 2 + (a 10 -2 a 21 ) s + (a 00 -a 11 + 2 a 22 ), Q = -s ((a 20 s + a 10 -a 21 ) ϑ 0 + a 20 ϑ 1 ), Q = -(
         ϑ 0 = x(0), ϑ 1 = ẋ(0), P = s (s 2 + ω 2 ), Q = -s (s ϑ 0 + ϑ 1 ) , Q = -(s 2 + ω 2 ) γ. (2) If x = T n is the n th Chebyshev polynomial of the first
kind, then we have a 2 = -t 2 +1, a 1 = -t and a 0 = n 2 , and thus we get ( 8), where:

         ϑ 0 = x(0), ϑ 1 = ẋ(0), P = -s s 2 ∂ 2 s + 3 s ∂ s -s 2 -n 2 + 1 , Q = -s (s ϑ 0 + ϑ 1 ) , Q = -(s 2 + n 2 ) γ. If x = U n
is the n th Chebyshev polynomial of the second kind, then we have a 2 = -t 2 + 1, a 1 = -3 t and a 0 = n (n + 2), and thus we get (8), where:

         ϑ 0 = x(0), ϑ 1 = ẋ(0), P = -s s 2 ∂ 2 s + s ∂ s -s 2 -(n + 1) 2 , Q = -s (s ϑ 0 + ϑ 1 ) , Q = -(s 2 + n (n + 2)) γ. (3) If x = P α,β n
is the n th Jacobi polynomial which depends on the parameters α and β, then we have

a 2 = -t 2 + 1, a 1 = -(α + β + 2) t + β -α and a 0 = n (n + α + β + 1)
, and thus we get (8), where:

             ϑ 0 = x(0), ϑ 1 = ẋ(0), P = s (-s 2 ∂ 2 s + (β + α -2) s ∂ s + s 2 + (β -α) s + (n + 1) (n + α + β)), Q = -s (ϑ 0 (s + β -α) + ϑ 1 ), Q = -s 2 + (β -α) s + n (n + α + β + 1) γ.
The n th Legendre polynomial is a particular case of the n th Jacobi polynomial for which α = β = 0, i.e.:

         ϑ 0 = x(0), ϑ 1 = ẋ(0), P = s -s 2 ∂ 2 s -2 s ∂ s + s 2 + (n + 1) n , Q = -s (ϑ 0 s + ϑ 1 ) , Q = -s 2 + n (n + 1) γ. (4) If x = L (α) n
is the n th Laguerre polynomial which depends on α, then we have a 2 = t, a 1 = -t + α + 1 and a 0 = n, and thus we get (8), where:

       ϑ 0 = x(0), ϑ 1 = ẋ(0), P = s (s (1 -s) ∂ s + (α -1) s + n + 1), Q = -s α ϑ 0 , Q = -(α s + n) γ.
(5) If x = H n is the n th Hermite polynomial, then a 2 = 1, a 1 = -2 t, a 0 = 2 n, and we get (8), where:

         ϑ 0 = x(0), ϑ 1 = ẋ(0), P = s (2 s ∂ s + s 2 + 2 (n + 1)) Q = -s (s ϑ 0 + ϑ 1 ) Q = -(s 2 + 2 n) γ. Remark 2.
The signal x in (1) could also be a wavelet. This case will be studied in a future publication. Wavelets satisfy functional equations such as dilation equations. Classes of functional equations can be studied by means of the so-called Ore algebras. For more details, see [START_REF] Chyzak | Effective algorithms for parametrizing linear control systems over Ore algebras[END_REF] and the references therein.

Estimation of an expansion of x into a basis

Till now, we have consider x to be the solution of an ODE with polynomial coefficients. More generally, x can be the output of a linear OD system with polynomial coefficients. An intermediate case is a signal x = N k=0 λ k x k which is a linear combination of signals x k which satisfy ODEs of the form (2) and the λ k 's are constants. This case corresponds to the case of uncoupled linear OD system since the x k 's are independent from one another. For instance, {x k } k=0,...,N can be the first N terms of the Taylor basis, i.e., x k = t k /k!. They can also be the first generators of a basis of L 2 such as the Fourier orthogonal basis or an orthogonal basis defined by a family of orthogonal polynomials. In this case, x is an approximation of the expansion of a function into an orthogonal basis. An important remark is that x then satisfies an ODE of the form of (2). This result a consequence of the concept of holonomic functions developed in algebraic analysis [START_REF] Kashiwara | Foundations of Algebraic Analysis[END_REF]). See the next paragraph for a simple explanation in the particular of ODEs with rational functions coefficients. As a consequence for the algebraic parameter estimation problem, we can always assume that x is defined by (2), i.e., we do not have to make a distinction between the case of a single signal x and of a finite linear combination of signals x k , when each signal x k satisfies an ODE with polynomial coefficients.

Let us suppose that K is a field and let D = K(t) ∂ be the noncommutative ring of OD operators in ∂ t with coefficients in the field K(t) of rational functions in t with coefficients in K. An element of D is of the form [START_REF] Chyzak | Effective algorithms for parametrizing linear control systems over Ore algebras[END_REF]). If x is a D-finite function, then the set ∂ i t x = x (i) i=0,...,n of cardinality n+1 admits at least one relation over K(t), i.e., there exist a i ∈ K(t) for i = 0, . . . , n such that n i=0 a i (t) ∂ i t x(t) = 0, i.e., x satisfies an ODE with rational function coefficients, and thus an ODE of the form (2). Conversely, if x satisfies an ODE of the form (2), then we get x (n) = n-1 i=0 ai an x (i) , which shows that D x is a K(t)-vector space of dimension n with the basis x (i) i=0,...,n-1 . If x 1 and x 2 are two D-finite functions and λ 1 and λ 2 two constants, then we clearly have

n i=0 a i (t) ∂ i t , where a i ∈ K(t) for i = 0, . . . , n. A function x is said to be D-finite if the left D-module defined by x, namely, D x = {d x | d ∈ D}, is a finite-dimensional K(t)- vector space
I = D (λ 1 x 1 + λ 2 x 2 ) ⊆ D x 1 + D x 2 ,
which shows that I is a finite-dimensional K(t)-vector space since D x 1 and D x 2 are both finite-dimensional over K(t), which implies that λ 1 x 1 + λ 2 x 2 satisfies an ODE with rational function coefficients, and thus an ODE of the form of (2). For more details, see [START_REF] Chyzak | Effective algorithms for parametrizing linear control systems over Ore algebras[END_REF]).

Note that the order n of the ODE for λ 1 x 1 + λ 2 x 2 , which is equal to the dimension of I, is at most n 1 + n 2 , where n i is the order of the ODE defining x i for i = 1, 2.

A natural problem is then to compute the ODEs that x satisfies when the ODEs for the x i 's are known. This can be easily done by considering the successive derivatives of x = N i=0 λ i x i and to search for K(t)-linear combinations among them [START_REF] Chyzak | Effective algorithms for parametrizing linear control systems over Ore algebras[END_REF]).

Example 5. Let x j = e i (ωj t+φj ) for j = 1, 2, 3. Clearly, x j satisfies ∂ t x j -i ω j x j = 0, i.e., x j is annihilated by ∂ t -i ω j . Let us compute the annihilator of the linear combination x = λ 1 x 1 + λ 2 x 2 + λ 3 x 3 , where the λ i 's are constants. Differentiating x with respect of t, we get

        x ẋ ẍ x (3) . . .         =         λ 1 λ 2 λ 3 i ω 1 λ 1 i ω 2 λ 2 i ω 3 λ 3 -ω 2 1 λ 1 -ω 2 2 λ 2 -ω 2 3 λ 3 -i ω 3 1 λ 1 -i ω 3 2 λ 2 -i ω 3 3 λ 3 . . . . . . . . .           x 1 x 2 x 3   .
Let L be the second matrix appearing in the above equation. Linear relations between the rows of L yield OD operators annihilating x. Indeed, if J = (1 ∂ t ∂ 2 t . . .) T and r is a row vector with entries in K satisfying r L = 0, then the above equation can be rewritten as J x = L x, and thus we then get (r J) x = (r L) x = 0, i.e., the OD operator r J annihilates x. We can check that there are no relations between the first three rows of L but r = (i ω 1 ω 2 ω 3 -(ω 1 ω 2 +ω 1 ω 3 +ω 2 ω 3 ) -i (ω 1 +ω 2 +ω 3 ) 1) is in the left kernel of L, which yields that the OD operator

∂ 3 t -i (ω 1 +ω 2 +ω 3 ) ∂ 2 t -(ω 1 ω 2 +ω 1 ω 3 +ω 2 ω 3 ) ∂ t +i ω 1 ω 2 ω 3 annihilates x = λ 1 x 1 + λ 2 x 2 + λ 3 x 3 .
Example 6. Let x 1 (resp., x 2 ) be the n th (resp., m th ) Hermite polynomial, λ 1 and λ 2 two constants, and the signal x = λ 1 x 1 + λ 2 x 2 . The functions x 1 and x 2 are respectively defined by the following ODE: ẍ1 (t) -2 t ẋ1 (t) + 2 m x 1 (t) = 0, ẍ2 (t) -2 t ẋ2 (t) + 2 n x 2 (t) = 0. Let us compute the ODE satisfies by x. To do that, we differentiate 2 + 2 = 4 times x and we get

(z ż z z (3) z (4) ) T = L (x 1 x 2 ẋ1 ẋ2 ) T ,
where the matrix L is defined by:

L =     λ 1 λ 2 0 0 0 0 λ 1 λ 2 -2 m λ 1 -2 n λ 1 2 t λ 1 2 t λ 2 -4 m (2 t 2 -m) λ 1 -4 n (2 t 2 -n) λ 2 8 t (t 2 -m) λ 1 8 t (t 2 -m) λ 2     .
The left kernel of L is generated by the row vector:

r = 4 m n -4 (m + n) t 4 t 2 + 2 (m + n) -4 t 1 .
Therefore, we get d(t, ∂ t ) z(t) = 0, where d is the OD operator defined by

d = r 1 ∂ t ∂ 2 t ∂ 3 t ∂ 4 t T , i.e.: ∂ 4 t -4 t ∂ 3 t + (4 t 2 + 2 (m + n)) ∂ 2 t -4 (m + n) t ∂ t + 4 m n.

Statement of the parameter estimation problem

To do study the algebraic parameter estimation problem, we start with (8), i.e., P (s, ∂ s ) z(s) + Q(s) + Q(s) = 0. Let θ = {θ 1 , . . . , θ r } be a set of independent variables and let us suppose that K = k[θ 1 , . . . , θ r ], simply denoted by

K = k[θ],
where k is a subfield of K. In other words, we suppose that the polynomials a i 's which define the dynamics of the signal x, i.e., (2), depend on certain constant parameters θ. Similarly, let ϑ = {ϑ 1 , . . . , ϑ s } be another set of independent variables and Θ = θ ∪ ϑ. We simply denote k[θ 1 , . . . , θ r , ϑ 1 , . . . , ϑ s ] by k[Θ]. Now, using (7), we obtain that:

P ∈ k[θ] ∂ , Q ∈ k[θ, ϑ, s] = k[Θ, s], Q ∈ k[θ, γ, s].
Note that the ϑ i 's (resp., γ) appear linearly in Q (resp., Q). Since P = s n0+1 n i=0 a i (-∂ s ) s i (see ( 7)), the coefficients of the polynomials a i 's appear linearly in P . Hence, if the parameters θ i 's appear linearly in the a j 's, then they also appear linearly in P and Q. In what follows, as usually done in the literature, we assume that this condition holds.

The parameter estimation problem aims at "estimating" the parameters Θ (or a given subset of Θ) from the observation of y = z + . In the algebraic parameter estimation problem, "estimating" means, in the noiseless case, i.e., when = 0, exactly and explicitly determine the parameters Θ in terms of y and of its integrals. For more details, see [START_REF] Fliess | An algebraic framework for linear identification[END_REF]; Sira-Ramírez et al. ( 2014); [START_REF] Mboup | Numerical differentiation with annihilators in noisy environment[END_REF]; [START_REF] Ushirobira | Algebraic parameter estimation of a biased sinusoidal waveform signal from noisy data[END_REF]Ushirobira et al. ( , 2013) ) and the references therein. The use of integrals of y is made to filter the effects of the noise . Example 7. If we consider again Example 2, then we can take ϑ = {ϑ k } k=0,...,n-1 and Θ = ϑ. Hence, the algebraic estimation problem aims at estimating the ϑ k 's or equivalently the a k 's from the measure of y = z + . Example 8. If we consider again 1 (resp., 2, resp. 3) of Example 4, then, for instance, we can take θ = {ω}, ϑ = {ϑ 0 , ϑ 1 }, and thus Θ = {ω, ϑ 0 , ϑ 1 } (resp., Θ = {n, ϑ 0 , ϑ 1 }, Θ = {α, β, ϑ 0 , ϑ 1 }). The algebraic estimation problem aims at estimating Θ from y = z + .

To get rid off γ, we have to find annihilators of Q, namely OD operators Π(s, ∂ s ) satisfying Π(s, ∂ s ) Q = 0, so that:

Π(s, ∂ s ) P (s, ∂ s ) z(s) + Π(s, ∂ s ) Q(s) = 0. ( 10
)
In particular, we want to find annihilators Π of Q such that the parameters Θ can be obtained explicitly in terms of z and its derivatives. In another words, we want to find annihilators Π j of Q, j = 1, . . . , r, such that the OD system Π j (s, ∂ s )

P (s, ∂ s ) z(s) + Π j (s, ∂ s ) Q(s) = 0, j = 1, . . . , r, (11) yields R 1 Θ = R 2 , where R 1 = R 11 (s, ∂ s ) z(s) + R 10 (s), R 2 = R 21 (s, ∂ s ) z(s) + R 20 (s)
, and R 1 is generically a non-singular matrix so that we get: Θ = R -1 1 R 2 . We can divise the both R 1 and R 2 by s N for a large enough integer N so that R 1 = s -N R 1 and R 2 = s -N R 2 are polynomials only in s -1 . Doing that (i.e., adding integrators in the corresponding expressions) helps to filter the noise while coming back to the time domain by means of the inverse Laplace transform. Thus, we can consider Θ = R 1 -1 R 2 . Note also that the coefficients of the annihilators Π j 's cannot depend on the parameters to be estimated till they are unknown. For more details, see [START_REF] Fliess | An algebraic framework for linear identification[END_REF]; [START_REF] Sira-Ramírez | Algebraic Identification and Estimation Methods in Feedback Control Systems[END_REF]; [START_REF] Mboup | Numerical differentiation with annihilators in noisy environment[END_REF]; [START_REF] Ushirobira | Algebraic parameter estimation of a biased sinusoidal waveform signal from noisy data[END_REF]Ushirobira et al. ( , 2013)). Remark 3. If the parameters θ i 's are not assumed to appear linearly in the a j 's, i.e., in P and Q, then, from (11), algebraic equations have to be solved to get the θ i 's.

To initiate an effective study of the algebraic estimation parameter problem, we propose to investigate the following three problems by means of computer algebra methods:

(1) Let Q ∈ k[θ, γ, s] be the polynomial defined in (7) and let us note , s] is the polynomial defined in ( 7), then compute ann E (.q) ∩ ann E (.Q) and: ann E (.q)/ (ann E (.q) ∩ ann E (.Q)) .

q = Q/γ ∈ k[θ, s]. Compute a set of generators of the following left ideal of D = k[θ, s] ∂ s ann D (.q) = {d ∈ D | d q(s) = 0}. (2) If θ is a subset of θ (e.g., if θ = ∅ or θ = θ) and E = k[θ , s] ∂ s , from Point 1, deduce the left ideal ann D (.q) ∩ E of E. (3) If Q ∈ k[θ, ϑ
The first problem aims at computing the annihilators of the polynomial Q to get (11). Note that the annihilators of Q that we can be used do only have to contain parameters that are already estimated. Thus, the second problem solves this point by controlling the coefficients that can appear in the annihilators. Finally, the last problem solves the problem of recognizing whether or not an annihilator of Q also annihilates the polynomial Q which contains the parameters ϑ, i.e., the initial conditions of (2). For instance, if all annihilators of Q are also annihilators of Q, then no parameters ϑ can algebraically be estimated.

COMPUTATION OF ANNIHILATORS

Let K be a noetherian ring, namely, every ideal of K can be generated by a finite set of generators (see, e.g., Rotman ( 2009)), and p ∈ A = K[s] a polynomial of degree q -1 (q ≥ 1) in s with coefficients in K. Let D = A ∂ s be the ring of OD operators in ∂ s = d ds with coefficients in A. Let us also consider the following column vector:

J q = (1 ∂ s . . . ∂ q s ) T ∈ D q×1 .
(12) Applying J q to p, we get the following polynomial vector: R = J q p ∈ A q×1 . Since A is a noetherian ring, the following A-module ker A (.R) = {µ ∈ A 1×q | µ R = 0} formed by all the A-linear relations among of the rows of R is noetherian, and thus is finitely generated, i.e., can be generated by a finite generating set (see Rotman (2009)). Let {S i• } i=1,...,r be a set of generators of ker A (.R), where S i• ∈ A 1×q for i = 1, . . . , r, and S = (S T 1• . . . S T r• ) T ∈ A r×q the matrix whose rows are the S i• 's. Let us consider the following A-homomorphisms (i.e., A-linear maps):

.R :

A 1×q -→ A µ -→ µ R, .S : A 1×r -→ A 1×q ν -→ ν S. Then, we obtain ker A (.R) = im A (.S). If M = A/(A 1×q R)
is the A-module finitely presented by R, i.e., the factor module of A by the ideal of A generated by the entries of the column vector R, and π the A-homomorphism which maps 1 to its residue class π(1), then we have the following exact sequence of A-modules (see, e.g., Rotman ( 2009))

A 1×r .S / / A 1×q .R / / A π / / M / / 0, namely, π is surjective and ker A (.R) = im A (.S). Remark 4. If K = k[x 1 , . . . ,
x n ] is a polynomial ring over a computable field k (e.g., k = Q, Q(y 1 , . . . , y m )), then Gröbner basis techniques can be used to explicitly compute the matrix S (see, e.g., [START_REF] Becker | Gröbner Bases. A Computational Approach to Commutative Algebra[END_REF]. For instance, it can be computed by the SyzygyModule command of the OreModules package [START_REF] Chyzak | OreModules: A symbolic package for the study of multidimensional linear systems[END_REF]).

We have (S J q ) p = S R = 0, i.e., Q = S J q ∈ D r×1 is a column of OD operators which annihilates p, i.e., Q p = 0.

Since the degree of p is q -1, we get ∂ q s p = 0, which shows that the last entry of R is 0, and thus that the vector (0 . . . 0 1) ∈ ker A (.R) = im A (.S). Thus, there exists ω ∈ A 1×r such that (0 . . . 0 1) = ω S, and thus we obtain that ∂ q s = (0 . . . 0 1) J q = ω S J q = ω Q. Let us introduce the left ideal of D defined by all the OD operators which annihilate p, i.e., its annihilator over D:

ann D (.p) = {d ∈ D | d p = 0}.
If we write Q = (q 1 . . . q r ) T , where q i ∈ D, then we obtain that q i ∈ ann D (p.), i.e., r i=1 D q i ⊆ ann D (p.). Let us prove that we have ann

D (.p) = r i=1 D q i . If d ∈ ann D (.p), then d = t i=0 a i ∂ i s , where a i ∈ A. We can write d = c ∂ q s + d , where d ∈ D is of degree strictly less than q and c ∈ D. Let us write d = q-1 i=0 b i ∂ i
s and let β = (b 0 . . . b q-1 ) ∈ A 1×q . Since the degree of the polynomial p is q -1, then ∂ q s p = 0, and thus we get

q-1 i=0 b i ∂ i s p = c ∂ q s p + q-1 i=0 b i ∂ i s p = d p = 0,
which shows that β R = (β J q ) p = 0, i.e., β ∈ ker A (.R) = im A (.S). Hence, there exists γ ∈ A 1×r such that β = γ S, which yields β J q = γ S J q = γ Q, and thus we obtain

d = c ∂ q s + γ Q = (c ω + γ) Q, where c ω + γ ∈ D 1×r , which shows that d ∈ r i=0 D q i . Lemma 1.
With the above notations, we have:

ann D (.p) = {d ∈ D | d p = 0} = r i=1 D q i .
In other words, the annihilator ann D (.p) of p over D is generated by the entries of the column vector Q. Remark 5. If k is a field of characteristic 0 (e.g., k = Q, R), then a classical result of J. T. Stafford asserts that any left/right ideal of A 1 (k) can be generated by two elements (see [START_REF] Stafford | Module structure of the Weyl algebra[END_REF]). For an implementation of this result in the Stafford package, see Quadrat (2007). This result extends to a polynomial extension of the Weyl algebra, i.e., to A 1 (k) [y], where y is a commuting indeterminate. If we consider B 1 (k), then B 1 (k) is a principal left/right ideal domain, namely, every left/right ideal of B 1 (k) can be generated by one element. This generator can simply be computed by means of gcd computations.

Example 9. Let p = s 2 + ω 2 = -Q/γ ∈ A = Q(ω)[s],
where Q is defined in 1 of Example 4. The degree of p being 2, q = 3, and we get ker A (.R) = im A (.S), where:

R = J 3 p = p dp ds d 2 p ds 2 0 T = (s 2 + ω 2 2 s 2 0) T , S =   -2 s ω 2 0 0 -1 s 0 0 0 0 1   , Q = S J 3 =   ω 2 ∂ 2 s + s ∂ s -2 s ∂ 2 s -∂ s ∂ 3 s   . If D = A 1 (Q(ω))
, then, we obtain:

ann D (.p) = D (ω 2 ∂ 2 s +s ∂ s -2)+D (s ∂ 2 s -∂ s )+D ∂ 3 s . (13)
The generators of ann D (.p) do not form of Gröbner basis [START_REF] Becker | Gröbner Bases. A Computational Approach to Commutative Algebra[END_REF]). Computing a Gröbner basis, we get:

ann D (.p) = D (ω 2 ∂ 2 s + s ∂ s -2) + D ((s 2 + ω 2 ) ∂ s -2 s). If E = B 1 (Q(ω)), then we can easily check that ω 2 ∂ 2 s + s ∂ s -2 = 1 (s 2 + ω 2 ) (ω 2 ∂ s + 2) ((s 2 + ω 2 ) ∂ s -2 s), which shows that ann E (.p) = E ((s 2 + ω 2 ) ∂ s -2 s). Now, if we consider A = Q[ω, s],
then we obtain

S =      -2 s ω 2 0 -2 s s 2 + ω 2 0 0 0 -1 s 0 0 0 0 1      , (14)
and thus ann D (.p) is defined by the three generators of ( 13) and (s 2 + ω 2 ) ∂ s -2 s, and a Gröbner basis then yields: [START_REF] Ushirobira | Algebraic estimation of a biased and noisy continuous signal via orthogonal polynomials[END_REF]. Since the degree of p is 3, we have q = 4 and we get ker A (.R) = im A (.S), where the matrix R is defined by

ann D (.p) = D ((s 2 + ω 2 ) ∂ s -2 s) + D (s ∂ 2 s -∂ s ) + D ∂ 3 s . Finally, we can check that s ∂ 2 s -∂ s = 1 2 ∂ 2 s (s 2 + ω 2 ) ∂ s -2 s) - 1 2 (s 2 + ω 2 ) ∂ 3 s , which shows that ann D (.p) = D ((s 2 + ω 2 ) ∂ s -2 s) + D ∂ 3 s . Example 10. Let us consider K = Q[λ 0 , λ 1 , λ 2 ] and: p = λ 0 s 3 + 2 λ 1 s 2 -2 λ 2 s (s 2 -4) ∈ A = K[s]. See
R = J 4 p = p dp ds d 2 p ds 2 d 3 p ds 3 0 T ,
and S is given in Figure 1. Then, ann D (.p) is generated by the entries of Q = S J 4 ∈ D 7×1 . If Q = (q 1 . . . q 7 ) T , then, e.g., q 2 = s 3 ∂ 3 s -3 s 2 ∂ 2 s + 6 s ∂ s -6 and q 7 = ∂ 4 s .

TOWERS OF ANNIHILATORS

Let θ be a subset of θ and

K = k[θ ] the subring of K = k[θ].
S =              -9 λ 0 + 18 λ 2 3 λ 0 s -6 λ 2 s -2 λ 1 2 λ 1 s + 4 λ 2 4 λ 2 s 0 -6 6 s -3 s 2 s 3 0 6 λ 1 -4 λ 1 s s (λ 1 s -4 λ 2 ) 4 λ 2 s 2 0 0 -6 λ 0 + 12 λ 2 3 s (λ 0 -2 λ 2 ) 2 λ 1 s + 8 λ 2 0 0 2 λ 1 -2 λ 1 s -4 λ 2 s (λ 1 s + 4 λ 2 ) 0 0 0 -3 λ 0 + 6 λ 2 3 λ 0 s -6 λ 2 s + 2 λ 1 0 0 0 0 0 1              ∈ A 7×5 .
Fig. 1. Matrix S for p = λ 0 s 3 + 2 λ 1 s 2 -2 λ 2 s (s 2 -4) monomial order which eliminates the elements of θ \ θ .

For more details, see, e.g., [START_REF] Chyzak | Effective algorithms for parametrizing linear control systems over Ore algebras[END_REF].

Example 11. We consider again Example 9. Let K = Q be the subring of

K = Q[ω] (i.e., θ = ∅, θ = {ω}), A = K[s], B = K [s] = Q[s], D = A ∂ s and E = B ∂ s . We obtain: ann E (.p) = ann D (.p) ∩ E = E (s ∂ 2 s -∂ s ) + E ∂ 3 s . Example 12. We consider again Example 10. Let K = Q be the subring of K = Q[λ 0 , λ 1 , λ 2 ] (i.e., θ = {λ 0 , λ 1 , λ 2 }, θ = ∅), A = K[s], B = K [s], D = A ∂ s and E = B ∂ s .
We get that ann E (.p) = E e 1 + E e 2 , where e 1 = ∂ 4

s and

e 2 = s 3 ∂ 3 s -3 s 2 ∂ 2 s + 6 s ∂ s -6. Now, if θ = {λ 0 , λ 2 }, K = Q[λ 0 , λ 2 ], B = K [
s] and E = B ∂ s , then we obtain that ann E (.p) = E e 1 + E e 2 + E e 3 , where:

e 3 = -8 λ 2 s ∂ 3 s + (λ 0 -2 λ 2 ) (3 s 2 ∂ 2 s -12 s ∂ s + 18).

INTERSECTION OF ANNIHILATORS

With the notations of Section 4, let p 1 , p 2 ∈ B and ann E (.p i ) for i = 1, 2. Let us show how to effectively compute the left ideal I = ann E (.p 1 ) ∩ ann E (.p 2 ) of E.

As explained in Section 3 with the ring B instead of A, let deg s p i = q i -1, J i = (1 . . . ∂ qi s ) T , R i = J i p i and S i ∈ B ri×qi be such that ker B (.R i ) = im B (.S i ). If Q i = S i J i = (q 1,i . . . q ri,i ) T , then, by Lemma 1, we have:

ann E (.p i ) = ri j=1 E q j,i , i = 1, 2.
An element p ∈ I is an element of E which satisfies p = ri j=1 e j,i q j,i for i = 1, 2 and for certain q j,i 's in E.

If we note Q = (Q T 1 Q T 2 ) T ∈ E (r1+r2)×1 , then we have: (e 1,1 . . . e 1,r1 -e 1,2 . . . -e r2,2 ) Q = 0. Hence, if T ∈ E s×(r1+r2) is such that ker E (.Q) = im E (.T ), and T = (T 1 -T 2 ), where T 1 ∈ E s×r1 and T 2 ∈ E s×r2 , then we have T 1 Q 1 = T 2 Q 2 ,
and thus we obtain: ann

E (.p 1 ) ∩ ann E (.p 2 ) = E 1×s (T 1 Q 1 ) = E 1×s (T 2 Q 2 ).
The computation of the matrix T can be obtained by the SyzygyModule command of OreModules.

Finally, let us explicit characterize the left E-module: P = ann E (.p 1 )/ (ann E (.p 1 ) ∩ ann E (.p 2 )) . Using Noether's second isomorphism theorem (see, e.g., Rotman ( 2009)), we have:

P ∼ = P = (ann E (.p 1 ) + ann E (.p 2 )) /ann E (.p 2 ) =   r1 j=1 E q j,1 + r2 j=1 E q j,2   /   r2 j=1 E q j,2   = E 1×(r1+r2) Q / E 1×r2 Q 2 .
We note that e ∈ ann E (.p 1 ) is of the form of e = η Q 1 , where η ∈ E 1×r1 . If τ : ann E (.p 1 ) -→ P (resp., κ : E 1×(r1+r2) Q -→ P ) is the left E-homomorphism defining the canonical projection onto P (resp., P ), then the above isomorphism ψ : P -→ P is defined by: ψ (τ (η Q 1 )) = κ(η Q 1 ). Now, using Lemma 3.1 on pages 349-350 of [START_REF] Cluzeau | Factoring and decomposing a class of linear functional systems[END_REF], we obtain the following isomorphism:

P ∼ = E 1×(r1+r2) / E 1×(r1+s) 0 I r2 T 1 -T 2 ∼ = P = E 1×r1 / E 1×s T 1 .
More precisely, if σ : E 1×r1 -→ P is the isomorphism defining the canonical projection onto P , then we have:

φ : P -→ P κ((η ζ) Q) -→ σ(η), φ -1 : P -→ P σ(η) -→ κ(η Q 1 ).
We obtain the isomorphism φ • ψ : P -→ P defined by: ∀ η ∈ E 1×r1 , (φ • ψ)(τ (η Q 1 )) = σ(η). Using the isomorphism φ • ψ, testing if e = η Q 1 ∈ ann E (.p 1 ) belongs to ann E (.p 2 ), i.e., if τ (e) = 0, is equivalent to testing if σ(η) = 0, i.e., if η = λ T 1 for a certain λ ∈ E 1×r1 . This condition can be checked by means of the Factorize command of the OreModulespackage. We also note that P ∼ = P = 0 iff the matrix T 1 admits a left inverse U 1 ∈ E r1×s , i.e., U 1 T 1 = I r1 . This condition can be checked by means of the LeftInverse command of the OreModules package. If P = 0, then we have ann E (.p 1 ) ⊆ ann E (.p 2 ), i.e., we cannot find an annihilator of p 1 which is not an annihilator of p 2 . Example 13. Let us consider again Example 9 Applying d 1 = ∂ 3 s ∈ ann E (.p) ∩ ann E (.Q) to (8) to get rid off of both Q and Q, we obtain d 1 (∂ s ) P (s) z(s) = 0, i.e., c 1 θ 1 = c 2 , where θ 1 = ω 2 and:    c 1 = s z (3) (s) + 3 z (2) (s), c 2 = -s 3 z (3) (s) + 9 s 2 z (2) (s) + 18 s z (1) (s) + 6 z(s) .

. If p = Q/γ, D = Q[ω, ϑ 0 , ϑ 1 ] ∂ s and E = A 1 (Q), then we have: ann D (.p) = D ((s 2 + ω 2 ) ∂ s -2 s) + D (s ∂ 2 s -∂ s ) + D ∂ 3 s , ann E (.p) = E (s ∂ 2 s -∂ s ) + E ∂ 3 s , ann D (.Q) = D (s 2 ∂ 2 s -2 s ∂ s + 2) + D (ϑ 1 s ∂ 2 s + 2 ϑ 0 s ∂ s -4 ϑ 0 ) + D ((2 ϑ 0 s + ϑ 1 ) ∂ 2 s -2 ϑ 0 ∂ s ) + D ∂ 3 s , ann E (.Q) = E (s 2 ∂ 2 s -2 s ∂ s + 2) + E ∂
Let Q 1 = (s ∂ 2 s -∂ s ∂ 3 s ) T be the column vector formed by the two generators of ann E (.p). Since the residue class of η = (1 0) in P is not zero, i.e., η / ∈ D 1×2 T 1 , we get d 2 = η Q 1 = s ∂ 2 s -∂ s ∈ ann E (.p) \ ann E (.Q). Hence, applying d 2 to (8), we get c 3 θ 1 + ϑ 1 = c 4 , where:

   c 3 = s 2 z (2) (s) + s z (1) (s) -z(s), c 4 = -s 4 z (2) (s) + 5 s 3 z (1) (s) + 3 s 2 z(s) .

We can easily solve these two linear equations to get:

       θ 1 = c 2 c 1 = s -4 c 2 s -4 c 1 , ϑ 1 = c 1 c 4 -c 2 c 3 c 1 = s -5 (c 1 c 4 -c 2 c 3 ) s -5 c 1 .
We are now left with the identification of the last parameter ϑ 0 appearing in Q as the coefficient of s 2 . Using the results of Section 4, we obtain ann D (.p) ∩ ann D (.s 2 ) = D (s ∂ 2 s -∂ s ) + D ∂ 3 s . Hence, considering the first generator d 3 = (s 2 + ω 2 ) ∂ s -2 s of ann D (.p) and apply it to (8), we get 2 θ 1 ϑ 0 = s ϑ 1 + c 5 θ 2 1 + c 6 θ 1 + c 7 , with the notations:      c 5 = s z (2) (s) + 2 z (1) (s), c 6 = s 4 z (2) (s) + 7 s 2 z (1) (s) + 5 s z(s), c 7 = s 4 z (1) (s) + s 3 z(s). Therefore, we obtain:

ϑ 0 =
s ϑ 1 + c 5 θ 2 1 + c 6 θ 1 + c 7 2 θ 1 .

Rewriting the above fraction as a quotient of polynomials in s -1 and coming to the time domain by means of the inverse Laplace transform, we obtain explicit formulas for the parameters θ 1 , ϑ 0 and ϑ 1 (see [START_REF] Ushirobira | Algebraic parameter estimation of a biased sinusoidal waveform signal from noisy data[END_REF]).

The algorithms for solving the three problems are implemented in the Non-A package dedicated to the algebraic estimation problem built upon the OreModules package. 4) and (5). The algebraic estimation problem can then be directly studied by means of (15). The main advantage of this approach is that it does not require the integration in closed-form solutions of the dynamics of the perturbation γ, which allows us to consider general type of structured perturbation γ. For more details, See Quadrat (2017).

  From the computation of the annihilators of p ∈ K[s] over the ring A = K[s] = k[θ, s] (see Section 3), in this section, we show how to obtain the annihilators of p with coefficients in the subringB = K [s] = k[θ , s] of A. In other words, if E = B ∂ s (resp., D = A ∂ s ),is the ring of OD operators in ∂ s with coefficients in B (resp., A), then, from the knowledge ann D (.p), we explain how to compute the left E-module ann E (.p) = {e ∈ E | e p = 0}. By Lemma 1, the annihilators of p over D are completely determined by the left D-module D 1×r Q = r i=1 D q i . In other words, the entries of Q generate ann D (.p). Since B ⊆ A, we get E ⊆ D, and thus we have: ann E (.p) = ann D (.p) ∩ E. Computing ann D (.p) ∩ E can be obtained by means of a noncommutative Gröbner basis computation with a

  3 s . Using the OreModules package, we obtain:ann E (.p) ∩ ann E (.Q) = E ∂ 3 s , ann E (.p)/ (ann E (.p) ∩ ann E (.Q)) ∼ = P = E 1×2 / D 1×2 ∂ s 0 0 1 .Let us show how to use the above results to solve the parameter estimation problem for x = A sin(ω t + φ) defined in 1 of Example 4[START_REF] Ushirobira | Algebraic parameter estimation of a biased sinusoidal waveform signal from noisy data[END_REF]).

  It is well-known that the ring D = A ∂ s of OD operators is a left Ore domain, namely, for all d 1 , d 2 ∈ D \ {0}, there exist e 1 , e 2 ∈ D \ {0} such that e 1 d 1 = e 2 d 2 (see, e.g., McConnell et al. (2000)). With the notations of Section 2.1, if d 1 = R and d 2 = R , then there exist T, T ∈ D \ {0} such that T R = T R . A set of generators of the OD operators T and T satisfying the above equality can be computed by means of the SyzygyModule command of the OreModules package. Using (4) and (5), we get T (s, ∂ s ) R(s, ∂ s ) z(s) + T (s, ∂ s ) S(s) + T (s, ∂ s ) S (s) = 0, (15) i.e., γ has been eliminated from (