Maximum likelihood drift estimation
for a threshold diffusion

arXiv:1803.05408v1 [math.PR] 14 Mar 2018

Antoine Lejay* and Paolo Pigato†
March 15, 2018

We study the maximum likelihood estimator of the drift parameters
of a stochastic differential equation, with both drift and diffusion coefficients constant on the positive and negative axis, yet discontinuous
at zero. This threshold diffusion is called the drifted Oscillating Brownian motion. The asymptotic behaviors of the positive and negative
occupation times rule the ones of the estimators. Differently from most
known results in the literature, we do not restrict ourselves to the ergodic
framework: indeed, depending on the signs of the drift, the process may
be ergodic, transient or null recurrent. For each regime, we establish
whether or not the estimators are consistent; if they are, we prove the
convergence in long time of the properly rescaled difference of the estimators towards a normal or mixed normal distribution. These theoretical
results are backed by numerical simulations.
Keywords: Threshold diffusion, Oscillating Brownian motion, maximum likelihood estimator, null recurrent process, ergodic process, transient process, mixed
normal distribution
AMS 2010: primary: 62M05; secondary: 62F12; 60J60.

Acknowledgement: P. Pigato gratefully acknowledges financial support from
ERC via Grant CoG-683164.
*

Université
de
Lorraine,
CNRS,
Inria,
IECL,
F-54000
Nancy,
France;
Antoine.Lejay@univ-lorraine.fr
†
Weierstrass Institute for Applied Analysis and Stochastics, Mohrenstrasse 39, Berlin, 10117,
Germany; E-mail: paolo.pigato@wias-berlin.de

1

1 Introduction
We consider the process, called a drifted Oscillating Brownian motion (DOBM),
which is the solution to the Stochastic Differential equation (SDE)
𝜉𝑡 = 𝜉0 +
with

⎧
⎨𝜎

+ > 0
𝜎(𝑥) =
⎩𝜎− > 0

∫︁ 𝑡
0

𝜎(𝜉𝑠 ) d𝑊𝑠 +

if 𝑥 ≥ 0,
if 𝑥 < 0

∫︁ 𝑡
0

𝑏(𝜉𝑠 ) d𝑠,

and 𝑏(𝑥) =

⎧
⎨𝑏

∈R
⎩𝑏− ∈ R
+

(1)

if 𝑥 ≥ 0,
if 𝑥 < 0.

(2)

The strong existence to (1) follows for example from the results of [28]. Separately
on R+ and R− , the dynamics of such process is the one of a Brownian motion with
drift, with threshold and regime-switch at 0, consequence of the discontinuity of the
coefficients.
This model can be seen as an alternative to the model studied in [36], which
is a continuous time version of the Self-Exciting Threshold Autoregressive models
(SETAR), a subclass of the TAR models [45, 46].
The practical interest of such processes are numerous. In finance, we show in [31]
that an exponential form of this process generalizes the Black & Scholes model in
a way to model leverage effects. Moreover, the introduction of a piecewise constant
drift such as the one in (2) is a straightforward way to produce a mean-reverting
process, if 𝑏+ < 0 and 𝑏− > 0. In [31], we find some evidence on empirical financial
data that this may be the case. This corroborates other studies with different
models [35, 38, 44].
Still in finance, the solution to (1) models other quantities than stocks. In [13],
Eq. (1) with constant volatility serves as a model for the surplus of a company after
the payment of dividends, which are payed only if the profits of the company are
higher than a certain threshold. Similar threshold dividend pay-out strategies are
considered in [3]. In these works, the behavior of the process at the discontinuity
is referred to as “refraction”. SETAR models have also applications to deal with
transaction costs or regulator interventions [49], to interest and exchange rates [6,9],
...
More general discontinuous drifts and volatilities arise in presence of Atlas models
and other ranks based models [17]. SDE with discontinuous coefficients have also
numerous applications in physics [39,42], meteorology [15] and many other domains.
In [43, 44], F. Su and K.-S. Chan study the asymptotic behavior of the quasilikelihood estimator of a diffusion with piecewise regular diffusivity and piecewise
affine drift with an unknown threshold. The quasi-likelihood they use is based on the
Girsanov density where the diffusivity is replaced by 1. In particular, they construct
some hypothesis test to decide whether or not the drift is affine or piecewise affine
in the ergodic situation.
In [26], Y. Kutoyants consider the estimation of a threshold 𝑟 of a diffusion with a
known or unknown drift switching at 𝑟. His results are then specialized to OrnsteinUhlenbeck type processes. Also this framework assumes that the diffusion is ergodic.
In the present paper, we derive some maximum likelihood estimators for the drift
parameters 𝑏− and 𝑏+ from continuous observations. We study their asymptotic

2

𝛽𝑇− − 𝑏−

𝛽𝑇+ − 𝑏+
(E)
(N0)
(N1)

𝑏+ < 0, 𝑏− > 0

ergodic

≈

√1
𝑇

√︂

1−

𝑏+
𝜎 𝒩
𝑏− +

√︂

≈

𝑏+ = 0, 𝑏− > 0 null recurrent

√1 𝜎+ 𝒩 +
𝑇

(T0) 𝑏+ > 0, 𝑏− ≥ 0

transient

≈

√1 𝜎+ 𝒩
𝑇

(T1)

transient

≈

√1 𝜎+ 𝒩
𝑇

1−

𝑏−
𝜎 𝒩
𝑏+ −

√1 𝛽 −
𝑇 1

√1 𝛽 +
𝑇 1

𝑏+ = 0, 𝑏− = 0 null recurrent

𝑏+ > 0, 𝑏− < 0

≈

√1
𝑇

√

≈

1
𝑇 1/4

𝑏−
𝜎− √𝜎+ √𝒩

−

|𝒩 |

ℛT0 as 𝑇 → ∞
ℛ+
T1 as 𝑇 → ∞

Table 1: Asymptotic behavior of estimators, where 𝒩 , 𝒩 + and 𝒩 − are independent,
unit Gaussian variables. The law of (𝛽1− , 𝛽1+ ) in case (N0) is given in (25).
The r.v.s ℛT0 and ℛ+
T1 follow the law in (22). Results of both sides in
(T1) are wrt to P+ (cf. Proposition 6), which intuitively can be thought
as conditioning to the process diverging towards postive infinity.
behavior as the time tends to infinity in order to derive some confidence intervals
when available. This article completes [30], where we estimate (𝜎− , 𝜎+ ) for highfrequency data. We use our estimators on financial historical data in [31].
Our estimators of 𝑏±
𝑇 are
𝛽𝑇± = ±

(±𝜉𝑇 ) ∨ 0 − (±𝜉0 ) ∨ 0 − 𝐿𝑇 (𝜉)/2
,
𝑄±
𝑇

−
where 𝑄+
𝑇 (resp. 𝑄𝑇 ) is the occupation time of the positive (resp. negative) side
of the real axis, and 𝐿𝑇 (𝜉) is the symmetric local time of 𝜉 at 0. Estimators for
𝑄± and 𝐿𝑇 (𝜉) are quite straightforward to implement from discrete observations of
a trajectory of 𝜉, and so are estimators for 𝑏± . As for the estimators of (𝜎− , 𝜎+ )
in [30], the local time and the occupation times play a central role in the study of
the estimators of (𝑏− , 𝑏+ ).
The long time asymptotic regime of the process depends on the respective signs of
the coefficients (𝑏− , 𝑏+ ). Using symmetries, this leads 5 different cases in which the
process may be ergodic, null recurrent or transient and the estimators have different
asymptotic behaviors. In some situations, the estimators are not convergent. In
others, we establish consistency as well as Central Limit Theorems, with speed 𝑇 1/2
or 𝑇 1/4 , depending again on the signs of 𝑏± . We summarize in Table 1 the various
asymptotic behaviors. We are in a situation close to the one encountered by M. Ben
Alaya and A. Kebaier in [1] for estimating square-root diffusions, where several
situations shall be treated. The work [26, 43, 44] mentioned above only consider
ergodic situations. Non-parametric estimation of the drift in the recurrent case is
considered in [4].
Finally, we develop in Section 7.1 a hypothesis test for the value of the drift which
is based on the Wilk’s theorem, which relates asymptotically the log-likelihood to
a 𝜒2 distribution with 2 degrees of freedom. Besides, we show in Section 7.2 the Local
Asymptotic Normality (LAN [27, 29]) and the Local Asymptotic Mixed Normality
(LAMN [22]) in the ergodic case and the null recurrent case with non vanishing drift.

3

These LAN/LAMN properties are related to the efficiency of the operators. The
Wilk’s as well as the LAN/LAMN properties are proved by combining the quadratic
nature of the log-likelihood with our martingale central limit theorems.
Outline. In Section 2, we present the maximum likelihood estimator, which is
based on the Girsanov transform. In Section 3, we characterize the different regimes
of the process accordingly to the signs of the drifts. Our main results are presented
in Section 4. The limit theorems that we use are presented in Section 5. The proofs
for each cases are detailed in Section 6. We present the Wilk theorem and the
LAN/LAMN property in Section 7. Finally, in Section 8, we conclude this article
with numerical experiments.

2 The maximum likelihood estimator
In this section, we propose and discuss an estimator for the parameters (𝑏− , 𝑏+ ) of
the drift coefficient of 𝜉 from continuous time observations
Data 1. We observe of a path (𝜉𝑡 )𝑡∈[0,𝑇 ] on the time interval [0, 𝑇 ] of the solution
to (1), together with its negative and positive occupation times
𝑄−
𝑇

=

∫︁ 𝑇
0

1𝜉𝑠 <0 d𝑠 and

𝑄+
𝑇

=

∫︁ 𝑇
0

1𝜉𝑠 >0 d𝑠,

as well as its symmetric local time
1 ∫︁ 𝑇
1−𝜖≤𝜉𝑠 ≤𝜖 d𝑠.
𝜖→0 2𝜖 0

𝐿𝑇 (𝜉) = lim
The coefficients (𝜎− , 𝜎+ ) are known.

+
Remark 1. Approximations of (𝑄−
𝑇 , 𝑄𝑇 , 𝐿𝑇 (𝜉)) are easy to construct from (𝜉𝑡 )𝑡∈[0,𝑇 ]
so that observing (𝜉𝑡 )𝑡∈[0,𝑇 ] is sufficient to build approximations of our estimator.
This is detailed in Section 8.1.

The Girsanov weight of the distribution of (1), with respect to the one of the
̃︁ , for a Brownian motion 𝑊
̃︁ , is
solution to d𝜉𝑡 = 𝜎(𝜉𝑡 ) d𝑊
𝐺(𝑏− , 𝑏+ ) = exp

(︃∫︁
𝑇
0

)︃

𝑏(𝜉𝑠 ) ̃︁
1 ∫︁ 𝑇 𝑏2 (𝜉𝑠 )
d𝑊𝑠 −
d𝑠 .
𝜎(𝜉𝑠 )
2 0 𝜎 2 (𝜉𝑠 )

(3)

A reasonable way to set up an estimator of (𝑏− , 𝑏+ ) is to consider 𝐺(𝑏− , 𝑏+ ) as
a likelihood and to optimize this quantity. This is how estimators for the drift are
classically constructed [25, 33].
Notation 1. To avoid confusion with the + and − used as indices, we write L𝑥M+ :=
max{𝑥, 0} for the positive part of 𝑥 and L𝑥M− := max{−𝑥, 0} ≥ 0 for the negative
part.
Proposition 1. The likelihood 𝐺(𝑏− , 𝑏+ ) is maximal at (𝛽𝑇− , 𝛽𝑇+ ) given by
𝛽𝑇± = ±

L𝜉𝑇 M± − L𝜉0 M± − 𝐿𝑇 (𝜉)/2
.
𝑄±
𝑇

4

(4)

Proof. Let us denote by 𝜉 the canonical process. Let us consider the measure P
̃︁ , 𝑡 ∈ [0, 𝑇 ] for a Brownian
such that 𝜉 is the unique solution to d𝜉𝑡 = 𝜎(𝜉𝑡 ) d𝑊
𝑡
̃︁ .
motion 𝑊
Under the distribution Q, with density 𝐺(𝑏− , 𝑏+ ) given by (3), with respect to P,
the process 𝜉 is solution to
d𝜉𝑡 = 𝜎(𝜉𝑡 ) d𝑊𝑡 + 𝑏(𝜉𝑡 ) d𝑡
̃︁ − 𝑡 𝑏 (𝜉 ) d𝑠. Under Q, the process 𝑊 is a Brownian
for 𝑊 defined by 𝑊𝑡 = 𝑊
𝑡
0 𝜎 𝑠
motion.
We define 𝛾(𝑥) = 𝑏(𝑥)/𝜎 2 (𝑥) and 𝐹 (𝑥) = 𝛾(𝑥)𝑥 for 𝑥 ∈ R. The function 𝐹 is
piecewise linear with
∫︀

𝐹 ′ (𝑥) = 𝛾(𝑥) for 𝑥 ̸= 0 and 𝛽 = 𝛾(0+) − 𝛾(0−) =

𝑏−
𝑏+
− 2.
2
𝜎+ 𝜎−

From the Itô-Tanaka formula1 [23, Theorem 7.1, p. 218],
𝐹 (𝜉𝑡 ) − 𝐹 (𝜉0 ) =

∫︁ 𝑡
0

𝐹−′ (𝜉𝑠 ) d𝜉𝑠

∫︁ 𝑡
𝛽
𝛽
+ 𝐿𝑡 (𝜉) =
𝛾(𝜉𝑠 ) d𝜉𝑠 + 𝐿𝑡 (𝜉), 𝑡 ∈ [0, 𝑇 ].
2
2
0

̃︁ ,
Since 𝛾(𝜉𝑠 ) d𝜉𝑠 = 𝑏(𝜉𝑠 )/𝜎(𝜉𝑠 ) d𝑊
𝑠
∫︁ 𝑇
0

𝑏(𝜉𝑠 ) ̃︁
𝛽
d𝑊𝑠 = 𝐹 (𝜉𝑇 ) − 𝐹 (𝜉0 ) − 𝐿𝑇 (𝜉).
𝜎(𝜉𝑠 )
2

Injecting this in the formula (3),
𝛽
1 ∫︁ 𝑇 𝑏2 (𝜉𝑠 )
log 𝐺(𝑏+ , 𝑏− ) = 𝐹 (𝜉𝑇 ) − 𝐹 (𝜉0 ) − 𝐿𝑇 (𝜉) −
d𝑠
2
2 0 𝜎 2 (𝜉𝑠 )
𝑏2− −
𝑏2
𝛽
𝑄𝑇 . (5)
= 𝐹 (𝜉𝑇 ) − 𝐹 (𝜉0 ) − 𝐿𝑇 (𝜉) − +2 𝑄+
𝑇 −
2
2
2𝜎+
2𝜎−
Maximizing (5) over 𝑏+ and 𝑏− leads to (4).
The proof of the next lemma is a direct consequence of the Itô-Tanaka formula.
It is the key to study the asymptotic behavior of 𝛽± .
Lemma 1. For any 𝑇 ≥ 0,
𝛽𝑇± = 𝑏± +

𝑀𝑇±
,
𝑄±
𝑇

2 ±
where 𝑀 ± :=± 0· 𝜎± 1±𝜉𝑠 ≥0 d𝐵𝑠 are continuous time martingales with ⟨𝑀 ± ⟩ = 𝜎±
𝑄
+
−
and ⟨𝑀 , 𝑀 ⟩ = 0.

∫︀

The occupation time is non decreasing. For the sake of simplicity, let us write
±
𝑄±
∞ = lim 𝑄𝑇 ∈ R+ ∪ {∞}.
𝑇 →∞

1

Our local time 𝐿(𝜉) is twice the one that appear in Theorem 7.1 of [23].

5

3 Analytic characterization of the regime of the
process
3.1 Scale function and speed measure
A well known fact [19, 23, 41] states that the infinitesimal generator (ℒ, Dom(ℒ)) of
the process 𝜉 solution to (1) may be written as
(︃

)︃

∫︁ 𝑥
2𝑏(𝑦)
1 2
−ℎ(𝑥) d
ℎ(𝑥) d𝑓 (𝑥)
ℒ𝑓 = 𝜎 (𝑥)𝑒
𝑒
with ℎ(𝑥) =
d𝑦
2
d𝑥
d𝑥
0 𝜎 2 (𝑦)
for all 𝑓 ∈ Dom(ℒ) = {𝑓 ∈ 𝒞0 (R) | ℒ𝑓 ∈ 𝒞0 (R)}.

The process 𝑋 is fully characterized by its speed measure 𝑀 with a density 𝑚 and
its scale function 𝑆 with
𝑚(𝑥) :=

∫︁ 𝑥
2
exp(ℎ(𝑥))
and
𝑆(𝑥)
:=
exp(−ℎ(𝑦)) d𝑦.
𝜎(𝑥)2
0

(6)

3.2 The regimes of the process
The diffusion 𝑋 is either recurrent or transient. If lim𝑥→+∞ 𝑆(𝑥) = +∞ and
lim𝑥→−∞ 𝑆(𝑥) = −∞, then the process is (positively or null) recurrent. Otherwise,
it is transient [19, 23].
When 𝑏(𝑥) = 𝑏+ for 𝑥 ≥ 0,
⎧
⎪
⎪
𝑥 (︃
⎪
⎪
⎪
2
⎪
⎪
⎨ 𝜎+

if 𝑏+ = 0,
(︃

)︃)︃

2𝑏+ 𝑥
1 − exp − 2
𝑆(𝑥) = 2𝑏+
𝜎)︃+
⎪
(︃
(︃
)︃
⎪
2
⎪
⎪
𝜎
2|𝑏
|𝑥
+
⎪
+
⎪
−1
⎪
⎩ 2|𝑏 | exp
2
𝜎+
+

if 𝑏+ > 0,
if 𝑏+ < 0.

Similar formulas hold for 𝑏− . Hence, the process 𝜉 is transient if only if 𝑏+ > 0
or 𝑏− < 0.
A recurrent process is either null recurrent
or positive recurrent. The process
∫︀
is positive recurrent if and only if 𝑀 (R) := R 𝑚(𝑥) d𝑥 < +∞, in which case it is
actually ergodic. Therefore, the process 𝜉 is ergodic if and only if 𝑏+ < 0 and 𝑏− > 0.
Otherwise, the process 𝜉 is only null recurrent.
When the process is ergodic (𝑏+ < 0, 𝑏− > 0), its invariant measure is
−2𝑥|𝑏

⎧
⎪
⎪
⎪
⎪
⎨

+
𝑏−
𝜎2
+
𝑒
𝑚(𝑥)
d𝑥 = 𝑏− + |𝑏+ | 2𝑥𝑏−
⎪
𝑀 (R)
|𝑏+ |
⎪
2
⎪
⎪
𝑒 𝜎−
⎩
𝑏− + |𝑏+ |

|

if 𝑥 ≥ 0,
if 𝑥 < 0.

Therefore, the regimes of 𝜉 depends only on the respective signs of 𝑏+ and 𝑏− .
Nine combinations are possible. As some cases are symmetric, we actually consider
five cases exhibiting different asymptotic behaviors of 𝑄±
𝑇 , hence of the estimators.
This is summarized in Table 2.
These cases are:

6

𝑏− > 0
𝑏− = 0
𝑏− < 0

𝑏+ > 0

𝑏+ = 0

𝑏+ < 0

transient T0
transient T0
transient T1

null recurrent N1
null recurrent N0
transient T0

ergodic E
null recurrent N1
transient T0

Table 2: Recurrence and transience properties of 𝜉.

E) Ergodic case 𝑏+ < 0, 𝑏− > 0.
N0) Null recurrent case 𝑏+ = 0, 𝑏− = 0.
N1) Null recurrent case 𝑏+ = 0, 𝑏− > 0.
T0) Transient case 𝑏+ > 0, 𝑏− ≥ 0.
T1) Transient case 𝑏+ > 0, 𝑏− < 0.
Case T0 corresponds to two entries of table 2. The case 𝑏+ < 0, 𝑏− = 0 is symmetric
to N1. Case 𝑏+ ≤ 0, 𝑏− < 0 is symmetric to T0.

4 Asymptotic behavior of the estimators
In this section, we state our main results on the asymptotic behavior of the occupation times of the process and the corresponding ones of the estimators, for each of
the 5 cases.
Proposition 2 (Ergodic case E). If 𝑏+ < 0, 𝑏− > 0, then
(︃

−
𝑄+
𝑇 𝑄𝑇
,
𝑇 𝑇

)︃

a.s.

−−−→
𝑇 →∞

(︃

)︃

|𝑏+ |
|𝑏− |
,
.
|𝑏− | + |𝑏+ | |𝑏− | + |𝑏+ |

(7)

In addition,
a.s.

(𝛽𝑇+ , 𝛽𝑇− ) −−−→ (𝑏+ , 𝑏− )
√

𝑇 →∞

⎛

⎞

𝜎+
𝜎−
(𝛽𝑇+ − 𝑏+ , 𝛽𝑇− − 𝑏− ) −−−→ ⎝ √︁
𝒩 + , √︁
𝒩 −⎠ ,
and √︁
𝑇 →∞
|𝑏− | + |𝑏+ |
|𝑏− |
|𝑏+ |
𝑇

law

where 𝒩 + and 𝒩 − are two independent, unit Gaussian random variables.
Proposition 3 (Null recurrent case with vanishing drift N0). Assume 𝑏+ = 𝑏− = 0.
Assume 𝜉0 = 0. Then
(︃

−
𝑄+
𝑇 𝑄𝑇
,
𝑇 𝑇

)︃

law

= (Λ, 1 − Λ) for all 𝑇 > 0,

where Λ follows a law of arcsine type with density
𝑝Λ (𝑢) :=
Besides,

1
1
𝜎+ /𝜎−
√︁
𝜋 𝑢(1 − 𝑢) 1 − (1 − (𝜎+ /𝜎− )2 )𝑢
√
law
𝑇 (𝛽𝑇+ , 𝛽𝑇− ) = (𝛽1+ , 𝛽1− )

7

for 0 < 𝑢 < 1.

(8)

where the explicit joint density of (𝛽1+ , 𝛽1− ) is given by (25) below. In particular,
(𝛽𝑇+ , 𝛽𝑇− ) converges almost surely to (𝑏+ , 𝑏− ) = (0, 0).
Proposition 4 (Null recurrent case with non-vanishing drift N1). Assume 𝑏+ = 0,
𝑏− > 0. Then
𝑄+
a.s.
a.s.
𝑇
−−−→ 1 and (𝛽𝑇+ , 𝛽𝑇− ) −−−→ (𝑏+ , 𝑏− ).
𝑇
→∞
𝑇
→∞
𝑇
In addition, there exists three independent unit Gaussian random variables 𝒩 − , 𝒩 +
and 𝒩 such that
⎛
⎞
(︃ −
)︃
√
−
𝑄𝑇 √
𝜎
𝑏
𝒩
law
+
−
⎠.
√ , 𝑇 (𝛽𝑇+ − 𝑏+ ), 𝑇 1/4 (𝛽𝑇− − 𝑏− ) −−−→ ⎝ |𝒩 |, 𝜎+ 𝒩 + , 𝜎− √︁ · √︁
𝑇 →∞
𝑏−
𝑇
𝜎+
|𝒩 |
(9)
Proposition 5 (Transient case for upward drift T0). Assume 𝑏+ > 0, 𝑏− ≥ 0 so
that the process 𝜉 is transient and lim𝑇 →∞ 𝜉𝑇 = +∞. Then 𝑄+
𝑇 /𝑇 converges almost
surely to 1 as 𝑇 → ∞ and
√
law
a.s.
(10)
𝛽𝑇+ −−−→ 𝑏+ and 𝑇 (𝛽𝑇+ − 𝑏+ ) −−−→ 𝜎+ 𝒩 +
𝑇 →∞

𝑇 →∞

for a unit Gaussian random variable 𝒩 + . Let ℓ0 be the last passage time to 0, which
is almost surely finite. Assume 𝜉0 = 0. We have
𝛽𝑇− 1𝑇 >ℓ0 = ℛT0 1𝑇 >ℓ0 and lim 𝛽𝑇− = ℛT0 a.s. with ℛT0 :=
𝑇 →∞

𝐿∞ (𝜉)
𝐿∞ (𝜉)
. (11)
− =
2𝑄ℓ0
2𝑄−
∞

The density of ℛT0 is given by (22) below. The case 𝑏+ ≤ 0, 𝑏− < 0 is treated by
symmetry.
Proposition 6 (Transient case for diverging drift T1). Assume 𝑏+ > 0, 𝑏− < 0 so
that the process 𝜉 is transient. Assume that 𝜉0 = 0. Then there exists a Bernoulli
random variable ℬ ∈ {0, 1} such that
𝜎− 𝑏+
,
𝜎+ 𝑏− + 𝜎− 𝑏+
)︃
(︃ −
)︃
(︃ +
𝑄𝑇
𝑄𝑇
P+
−−−→ 1 = 1 and P−
−−−→ 1 = 1
𝑇 𝑇 →∞
𝑇 𝑇 →∞
with P+ (·) = P(· | ℬ = 1) and P− (·) = P(· | ℬ = 0).
P(ℬ = 1) = 1 − P(ℬ = 0) =

On the event {ℬ = 1} (resp. {ℬ = 0}), 𝛽𝑇+ (resp. 𝛽𝑇− ) converges almost surely
to 𝑏+ (resp. 𝑏− ) while 𝛽𝑇− − 𝑏− (resp. 𝛽𝑇+ − 𝑏+ ) is the ratio of two a.s. finite random
variables.
In addition, for unit Gaussian random variables 𝒩 + and 𝒩 − ,
√
law
𝑇 (𝛽𝑇+ − 𝑏+ ) −−−→ 𝜎+ 𝒩 + under P+ ,
(12)
𝑇 →∞
√
law
𝑇 (𝛽𝑇− − 𝑏− ) −−−→ 𝜎− 𝒩 − under P− .
(13)
𝑇 →∞

8

In addition,
𝐿∞ (𝜉)
,
2𝑄−
∞
𝐿∞ (𝜉)
=
.
2𝑄+
∞

−
lim 𝛽𝑇− = ℛ−
T1 a.s. under P+ with ℛT1 =

𝑇 →∞

+
lim 𝛽𝑇+ = −ℛ+
T1 a.s. under P− with ℛT1

𝑇 →∞

(14)
(15)

+
The distribution of ℛ−
T1 is that of of ℛ given (22) below. That of ℛT1 is found by
symmetry.

5 Auxiliary tools
In this section, we give first some results on a martingale central limit theorem that
will be used constantly. To deal with the transient or null recurrent cases, we make
use of some analytic properties of one-dimensional diffusions.

5.1 Limit theorems on martingales
The following result follows immediately from [32, Proposition 1, p. 148; Theorem 1,
p. 150].
Proposition 7 (A criterion for convergence). Under the true probability P,
(i) as 𝑇 → ∞, 𝛽𝑇+ (resp. 𝛽𝑇− ) converges a.s. to 𝑏+ (resp. 𝑏− ) on the event {𝑄+
∞ =
−
+∞} (resp. {𝑄∞ = +∞}).
(ii) as 𝑇 → ∞, 𝑀𝑇+ (resp. 𝑀𝑇− ) converges a.s. to a finite value on the event {𝑄+
∞ <
±
+∞} (resp. {𝑄−
<
+∞}).
In
other
words,
𝛽
is
not
a
consistent
estimator
on
𝑇
∞
±
{𝑄∞ < +∞}.
We now state an instance of a Central Limit theorem for martingales which follows
from [7]. This theorem will be used to deal with the cases E, N1 and T1. Let us
start by recalling the notion of stable convergence introduced by A. Rényi [21, 40].
Definition 1 (Stable convergence). A sequence (𝑋𝑛 )𝑛∈N on a probability space
(Ω, ℱ, P) is said to converge stably with respect to a 𝜎-algebra 𝒢 ⊂ ℱ if for any
bounded, continuous function 𝑓 and any bounded, 𝒢-measurable random variable 𝑌 ,
E(𝑓 (𝑋𝑛 )𝑌 ) −−−→ E(𝑓 (𝑋)𝑌 ).
𝑛→∞

Proposition 8 (A central limit theorem for martingales). Let (Ω, ℱ, P) be the underlying probability space of the process 𝜉 with a filtration (ℱ𝑡 )𝑡≥0 . If for some
constants 𝑐+ , 𝑐− > 0,
𝑄−
𝑄+
P
P
𝑇
−−−−→ 𝑐+ and 𝑇 −−−−→ 𝑐− ,
𝑇
→+∞
𝑇
→+∞
𝑇
𝑇
then for the martingales 𝑀 ± defined in Lemma 1,
(︃

𝑀𝑇+ 𝑀𝑇−
√ ,√
𝑇
𝑇

)︃

√
√
ℱ∞ −stably
−−−−−−→ (𝜎+ 𝑐+ 𝒩 + , 𝜎− 𝑐− 𝒩 − ),
𝑇 →∞

9

(16)

on a probability space (Ω′ , ℱ ′ , P′ ) extending (Ω, ℱ, P) and containing two independent unit Gaussian random variables 𝒩 + , 𝒩 − , themselves independent from 𝜉. In
addition,
(︃ +
)︃
(︃
)︃
√
𝑀𝑇 𝑀𝑇− ℱ∞ −stably
𝜎+ + 𝜎− −
𝑇
,
(17)
, − −−−−−−→ √ 𝒩 , √ 𝒩
𝑇 →∞
𝑄𝑇
𝑐+
𝑐−
𝑄+
𝑇
Proof. Set
]︃
[︃
]︃
√
2 +
1/ 𝑇
0
𝜎
𝑄
0
𝑇
√
𝑎𝑇 :=
and 𝑞𝑇 = ⟨𝑀, 𝑀 ⟩𝑇 = +
2 − .
0
𝜎−
𝑄𝑇
0
1/ 𝑇
[︃

Thus,
⎡

𝑎𝑇 𝑞𝑇 𝑎′𝑇

=

2
⎣𝜎+

𝑄+
𝑇
𝑇

0

⎤

[︃

]︃

𝑐
0 ⎦ P
−−−→ + .
−
𝑄
𝜎 2 𝑇 𝑇 →∞ 𝑐−

(18)

− 𝑇

Theorem 2.2 in [7] yields (16). Besides,
√ 𝑀𝑇±
𝑇
𝑀𝑇±
𝑇 ± = ±×√ .
𝑄𝑇
𝑄𝑇
𝑇

(19)

If a sequence (𝑋𝑛 )𝑛 converges ℱ∞ -stably and a sequence (𝑌𝑛 )𝑛 of ℱ∞ -measurable
random variables converges in probability, then (𝑋𝑛 , 𝑌𝑛 )𝑛 converges ℱ∞ -stably. Using the property in (19) and (17) yields (17).

5.2 The fundamental system
Along with the characterization through the scale function and the speed measure,
much information on the process can be read from the so-called fundamental system
[11, 19, 41]: For any 𝜆 > 0, there exists some functions 𝜑𝜆 and 𝜓𝜆 such that
∙ 𝜓𝜆 and 𝜑𝜆 are continuous, positive from R to R with 𝜑𝜆 (0) = 𝜓𝜆 (0) = 1.
∙ 𝜓𝜆 is increasing with lim𝑥→−∞ 𝜓𝜆 (𝑥) = 0, lim𝑥→∞ 𝜓𝜆 (𝑥) = +∞.
∙ 𝜑𝜆 is decreasing with lim𝑥→−∞ 𝜑𝜆 (𝑥) = +∞, lim𝑥→∞ 𝜑𝜆 (𝑥) = 0.
∙ 𝜑𝜆 and 𝜓𝜆 are solutions to ℒ𝑓 = 𝜆𝑓 .
In the case of piecewise constant coefficients with one discontinuity at 0, these
solutions may be computed as linear combinations of the minimal functions for
constant coefficients. Using the fact that 𝜑, 𝜑′ , 𝜓 and 𝜓 ′ are continuous at 0,
√
⎧
)︂
(︂
2𝜆
−𝑏− + 𝑏2− +2𝜎−
⎪
⎪
if 𝑥 < 0
⎨exp 𝑥
2
𝜎−
√2
√2
)︂
(︂
)︂
(︂
(20)
𝜓𝜆 (𝑥) =
2
2
−𝑏+ + 𝑏+ +2𝜎+ 𝜆
−𝑏+ − 𝑏+ +2𝜎+ 𝜆
⎪
⎪
⎩𝜅+ exp 𝑥
+
𝛿
exp
𝑥
if
𝑥
≥
0,
2
2
+
𝜎+
𝜎+
√2
√2
⎧
(︂
)︂
(︂
)︂
2
2𝜆
−𝑏− − 𝑏− +2𝜎− 𝜆
−𝑏− + 𝑏− +2𝜎−
⎪
⎪
⎨𝜅− exp 𝑥
+ 𝛿− exp 𝑥
if 𝑥 < 0,
2
2
𝜎−
√ 2 𝜎− 2 )︂
(︂
𝜑𝜆 (𝑥) =
(21)
−𝑏+ − 𝑏+ +2𝜎+ 𝜆
⎪
⎪
⎩exp 𝑥
if
𝑥
≥
0
𝜎2
+

10

with
√︁

√︁

𝜅+ :=

2
2
2
2
2
2
−𝑏− 𝜎+
+ 𝑏+ 𝜎 −
+ 𝜎−
𝑏2+ + 2𝜆𝜎+
+ 𝜎+
𝑏2− + 2𝜆𝜎−

√︁

2
2
2𝜎−
𝑏2+ + 2𝜆𝜎+

√︁

√︁

𝛿+ :=

2
2
2
2
2
2
+ 𝜎−
− 𝑏+ 𝜎−
𝑏− 𝜎 +
𝑏2+ + 2𝜆𝜎+
− 𝜎+
𝑏2− + 2𝜆𝜎−

√︁

2
2
2𝜎−
𝑏2+ + 2𝜆𝜎+

2
2
2
2
2
2
𝑏2+ + 2𝜆𝜎+
𝑏2− + 2𝜆𝜎−
+ 𝜎+
− 𝜎−
+ 𝑏+ 𝜎 −
−𝑏− 𝜎+

√︁

2
2
2𝜎+
𝑏2− + 2𝜆𝜎−

,

√︁

√︁

𝛿− :=

,

√︁

√︁

𝜅− :=

,

2
2
2
2
2
2
𝑏− 𝜎 +
− 𝑏+ 𝜎−
+ 𝜎−
𝑏2+ + 2𝜆𝜎+
+ 𝜎+
𝑏2− + 2𝜆𝜎−

√︁

2
2
2𝜎+
𝑏2− + 2𝜆𝜎−

.

We also define the quantities [37]
√︁

2
𝜆
−𝑏− + 𝑏2− + 2𝜎−
1 𝜓𝜆′ (0)
̂︀
𝜓(𝜆)
:=
=
≥0
2
2 𝜓𝜆 (0)
2𝜎−

and

̂︀
𝜑(𝜆)
:= −

1 𝜑′𝜆 (0)
2 𝜑𝜆 (0)

=

𝑏+ +

√︁

2
𝑏2+ + 2𝜎+
𝜆

≥ 0.

2
2𝜎+

In particular,
𝑏+
̂︀
̂︀
𝜓(0)
= 0 and 𝜑(0)
= 2 when 𝑏− ≥ 0 and 𝑏+ ≥ 0.
𝜎+

5.3 Last passage time and occupation time for the transient
process
When 𝜉 is a transient process, the last passage time ℓ0 = sup{𝑡 ≥ 0 | 𝜉𝑡 = 0} of 𝜉
at 0 is almost surely finite. Its Laplace transform is (See (53) in [37]):
̂︀
̂︀
𝜓(0)
+ 𝜑(0)
E0 [exp(−𝜆ℓ0 )] = ̂︀
.
̂︀
𝜓(𝜆) + 𝜑(𝜆)

Let us now assume that 𝑏+ > 0 and 𝑏− ≥ 0. This is the transient case T0 where
−
the process ends up almost surely in the positive semi-axis. Thus, 𝑄−
ℓ0 = 𝑄∞ and
𝐿ℓ0 (𝜉) = 𝐿∞ (𝜉).
2
Let us write ̂︀𝑏± := 𝑏± /𝜎±
. From Corollary 5 in [37],
E0 [exp(−𝛼𝐿∞ (𝜉) − 𝜆𝑄−
∞ )] =

̂︀
̂︀
̂︀
𝜓(0)
+ 𝜑(0)
𝑏+
=
.
√︁
̂︀
̂︀
̂︀
2
𝛼 + 𝜑(0) + 𝜓(𝜆)
𝛼 + ̂︀𝑏+ − 𝑏2− + 2𝜎12 𝑏2− + 2𝜎−
𝜆
−

Let 𝑝𝐿∞ (𝜉) (𝑡) be the density of 𝐿∞ (𝜉). Setting 𝜆 = 0, we see that 𝐿∞ (𝜉) is distributed
according to an exponential distribution of rate ̂︀𝑏+ . Thus, 𝑝𝐿∞ (𝜉) (𝑡) = ̂︀𝑏+ exp(−̂︀𝑏+ 𝑡).
A conditioning shows that
E0 [exp(−𝛼𝐿∞ (𝜉) − 𝜆𝑄−
∞ )] =

∫︁ +∞
0

(︁

⃒

)︁

⃒
exp(−𝛼𝑡)E0 exp(−𝜆𝑄−
∞ ) ⃒ 𝐿∞ (𝜉) = 𝑡 𝑝𝐿∞ (𝑡) d𝑡.

11

By inverting the Laplace transform with respect to 𝛼, since 𝑝𝐿∞ (𝑡) = ̂︀𝑏+ exp(−̂︀𝑏+ 𝑡),
⎛

⃒
⃒
E0 exp(−𝜆𝑄−
)
∞ ⃒ 𝐿∞ (𝜉)
(︁

⎛

⎯
⎸

⎞⎞

2
1 ⎸
−
⎷ 𝑏− + 𝜆⎠⎠ .
= 𝑡 = exp ⎝−𝑡 ⎝− + √
2
2
2𝜎− 2𝜎−

)︁

̂︀
𝑏

Inverting the latter Laplace transform with respect to 𝜆, the density 𝑝𝑄−∞ (𝑠|𝑡) of 𝑄−
∞
given {𝐿∞ (𝜉) = 𝑡} is
𝑏2− 𝑠
𝑡𝑏−
𝑡2
𝑡
√
−
.
𝑝𝑄∞ (𝑠|𝑡) =
exp
− 2 − 2
2
2𝜎−
2𝜎− 8𝜎− 𝑠
𝜎− 2 2𝜋𝑠3/2
)︃

(︃

Hence, the distribution 𝑝(𝑄−∞ ,𝐿∞ (𝜉)) (𝑠, 𝑡) of the pair (𝑄−
∞ , 𝐿∞ (𝜉)) is
𝑡𝑏+
√
𝑝(𝑄−∞ ,𝐿∞ (𝜉)) (𝑠, 𝑡) = 2
exp
𝜎+ 𝜎− 2 2𝜋𝑠3/2

(︃(︃

𝑏2− 𝑠
𝑏+
𝑡2
𝑏−
.
−
𝑡
−
−
2
2
2
2
𝑠
2𝜎−
𝜎+
2𝜎−
8𝜎−
)︃

)︃

The density 𝑝ℛT0 (𝑟) of the random variable ℛT0 := 𝐿∞ (𝜉)/2𝑄−
∞ is then
𝑝ℛT0 (𝑟) = 2

∫︁ +∞
0

𝑠 · 𝑝(𝑄−∞ ,𝐿∞ (𝜉)) (𝑠, 2𝑟𝑠) d𝑠
(︃

𝑟𝑏+
2𝑟𝑏+ (𝑟 − 𝑏− )2
= 2 √
+
2
2
2𝜎−
𝜎+ 𝜎− 2 𝜎+

)︃−3/2

, 𝑟 > 0.

The distribution of 𝐿∞ (𝜉)/2𝑄+
∞ when 𝑏− < 0, 𝑏+ ≤ 0 is found by symmetric
arguments.
Let us now assume that 𝑏+ > 0 and 𝑏− < 0. This is the transient case T1 where
the process can end up in both semi-axis. The Laplace transform is
E0 [exp(−𝛼𝐿∞ (𝜉) − 𝜆𝑄−
∞ )] =

̂︀
̂︀
̂︀
𝜓(0)
+ 𝜑(0)
𝑏+ − ̂︀𝑏−
.
=
√︁
̂︀𝑏−
̂︀
̂︀
1
2
2
̂︀
𝛼 + 𝜑(0)
+ 𝜓(𝜆)
𝛼 + 𝑏+ − 2 + 2𝜎2 𝑏− + 2𝜎− 𝜆
−

With analogous computations as before we get that the density 𝑝ℛ− (𝑟) of ℛ−
T1 is
T1

𝑝ℛ− (𝑟) =
T1

𝑟
√

𝜎−

(︃

𝑏+
𝑏−
− 2
2
2 𝜎+ 𝜎−

)︃(︃

2𝑟𝑏+ (𝑟 − 𝑏− )2
+
2
2
𝜎+
2𝜎−

)︃−3/2

, 𝑟 > 0.

Considering also the previous case, we can write the following formula, holding for
the density of ℛ = ℛT0 or ℛ = ℛ−
T1 in both cases T0 and T1:
𝑝ℛ (𝑟) =

𝑟
√

𝜎−

𝑏+
L𝑏− M−
+
2
2
𝜎−
2 𝜎+
(︃

)︃(︃

2𝑟𝑏+ (𝑟 − 𝑏− )2
+
2
2
𝜎+
2𝜎−

)︃−3/2

, 𝑟 > 0.

(22)

Notice that this is the density of a positive random variable which is not integrable.
This gives the limit behavior of the estimator 𝛽𝑇− of 𝑏− . The behavior of 𝛽𝑇+ in the
corresponding cases can be found by symmetric arguments.

12

6 Proofs of the asymptotic behavior of the estimators
6.1 Asymptotic behavior for the ergodic case (E)
The ergodic case is the most favorable one.
The process 𝜉 is ergodic, so that
∫︀
for any bounded, measurable function 𝑓 , 𝑇1 0𝑇 𝑓 (𝜉𝑠 ) d𝑠 converges almost surely to
∫︀
𝑚(𝑥)
d𝑥.
𝑓 (𝑥) 𝑀
(R)
With the explicit expression of 𝑀 that follows from (6),
𝑀 (R+ ) =

−1
1
|𝑏+ 𝑏− |
, 𝑀 (R− ) =
and 𝑀 (R) =
.
𝑏+
𝑏−
|𝑏− | + |𝑏+ |

From the ergodic theorem, since 𝑄±
𝑡 =
𝑄±
𝑇
𝑇

∫︀ 𝑡
0

1±𝜉𝑠 ≥0 d𝑠,
𝑀 (R± )
𝑀 (R)

a.s.

−−−→
𝑇 →∞

so that

𝑄+
|𝑏− |
𝑄− a.s.
|𝑏+ |
a.s.
𝑇
and 𝑇 −−−→
.
−−−→
𝑇 𝑇 →∞ |𝑏− | + |𝑏+ |
𝑇 𝑇 →∞ |𝑏− | + |𝑏+ |
The rate of convergence follows from Proposition 8 and (23).

(23)

6.2 Asymptotic behavior for the null recurrent case with
vanishing drift (N0)
When 𝑏− = 𝑏+ = 0, the process 𝜉 is an Oscillating Brownian motion (OBM, introduced first in [24], see also [30]). Supposing 𝜉0 = 0, using the scaling relation [30, Remark 3.7], for any 𝑇 > 0,
(︃

L𝜉𝑇 M+ L𝜉𝑇 M− 𝐿𝑇 (𝜉) 𝑄+
√ , √ , √ , 𝑇
𝑇
𝑇
𝑇
𝑇

)︃

law

= (L𝜉1 M+ , L𝜉1 M− , 𝐿1 (𝜉), 𝑄+
1 ).

Therefore,

√
√ ⎞
⎛
⎞
+
+
1
1
L𝜉
M
/
𝑇
−
𝐿
(𝜉)/
𝑇
L𝜉
M
−
𝐿
(𝜉)
𝑇
𝑇
1
1
2
2
(︃ )︃
(︃ )︃
⎟
⎜
⎜
⎟
+
√
+
⎟ law 𝛽1+
⎜
⎜
⎟
𝛽𝑇+
𝑄
/𝑇
𝑄
𝑇
⎜
⎟
⎜
⎟
1
√
√
=
𝑇
−⎟.
− = ⎜1
− = ⎜1
−
⎟
𝛽1
𝛽𝑇
⎝ 2 𝐿1 (𝜉) − L𝜉1 M ⎠
⎝ 2 𝐿𝑇 (𝜉)/ 𝑇 − L𝜉𝑇 M / 𝑇 ⎠
𝑄−
𝑄−
1
𝑇 /𝑇
⎛

We recall now that 𝑋 = Φ(𝜉) := 𝜉/𝜎(𝜉) is a Skew Brownian motion [10, 30]. An
explicit formula for the density for the position a Skew Brownian motion, its local
time and its occupation time is known [2, 12]. Since the transform Φ is piecewise
linear, one easily recover the one of an OBM, its local and occupation times. Hence,
the density of (𝜉1 , 𝐿1 (𝜉), 𝑄+
1 (𝜉)) is
𝑝(𝜉1 ,𝐿1 (𝜉),𝑄+ ) (𝜌, 𝜆, 𝜏 )
1

⎧
(𝜆/2 + 𝜌)𝜆/2
⎪
⎪
⎪
⎪
3
⎪
⎨ 2𝜋𝜎− 𝜎+
(1 − 𝜏 )3/2 𝜏 3/2

=⎪
⎪
⎪
⎪
⎪
⎩

(︃

(𝜆/2)2
(𝜆/2 + 𝜌)2
exp − 2
−
2
2𝜎− (1 − 𝜏 )
2𝜎+
𝜏
(︃

(𝜆/2 − 𝜌)𝜆/2
(𝜆/2)2
(𝜆/2 − 𝜌)2
exp
−
−
3
2
2
2𝜋𝜎+ 𝜎−
(1 − 𝜏 )3/2 𝜏 3/2
2𝜎+
𝜏
2𝜎−
(1 − 𝜏 )

13

)︃

for 𝜌 ≥ 0,
(24)

)︃

for 𝜌 < 0.

The change of variable in the density suggested by
(︃

(𝛽1+ , 𝛽1− , 𝑄+
1)

=

L𝜉1 M+ − 𝐿1 (𝜉)/2 𝐿1 (𝜉)/2 − L𝜉1 M− +
,
, 𝑄1
𝑄+
1 − 𝑄+
1
1

)︃

gives
𝑝(𝛽 + ,𝛽 − ,𝑄+ ) (𝑎, 𝑏, 𝛿)
1

1

1

= 2𝛿(1 − 𝛿)𝑝(𝜉1 ,𝐿1 (𝜉),𝑄+ ) (𝑎𝛿 + 𝑏(1 − 𝛿), |𝑎𝛿 + 𝑏(1 − 𝛿)| − 𝑎𝛿 + 𝑏(1 − 𝛿), 𝛿)
1

and then, since 𝑄+
1 ∈ [0, 1],
𝑝(𝛽 + ,𝛽 − ) (𝑎, 𝑏)
1

=

1

∫︁ 1
0

2𝛿(1 − 𝛿)𝑝(𝜉1 ,𝐿1 (𝜉),𝑄+ ) (𝑎𝛿 + 𝑏(1 − 𝛿), |𝑎𝛿 + 𝑏(1 − 𝛿)| − 𝑎𝛿 + 𝑏(1 − 𝛿), 𝛿) d𝛿.
1

(25)

6.3 Asymptotic behavior for the null recurrent case with
non-vanishing drift (N1)
We consider 𝑏+ = 0, 𝑏− > 0. The particle is then pushed upward when its position
is negative. Yet the process is only null recurrent. The measure 𝑀 satisfies
𝑀 (R− ) =

1
and 𝑀 (R+ ) = +∞.
𝑏−

(26)

Using 9) in [19, Section 6.8, p. 228] or [34, 47], with (26),
𝑄−
𝑀 (R− )
a.s.
𝑇
= 0.
−−−→
𝑇
→∞
𝑇
𝑀 (R)
−
+
Since 𝑄+
𝑇 + 𝑄𝑇 = 𝑇 , it holds that 𝑄𝑇 /𝑇 converges almost surely to 1.
Using Proposition 8 on 𝑀 + only, we obtain that

𝑀𝑇+ law
√ −−−→ 𝜎 + 𝒩 +
𝑇 𝑇 →∞
for a Gaussian random variable 𝒩 + ∼ 𝒩 (0, 1) and then that

(27)
√

law

𝑇 (𝛽𝑇+ − 𝑏+ ) −−−→
𝑇 →∞

𝜎+𝒩 +.
The asymptotic behavior of 𝑄−
𝑇 is more delicate to deal with as the process 𝜉 is
only null recurrent. For this, we use the results of [14] which extends the one of
D.A. Darling and M. Kac [8] on additive and martingale additive functionals.
The Green kernel with respect to the invariant measure 𝑀 of ℒ is given by [11,
19, 41]
⎧

1 ⎨𝜓𝜆 (𝑥)𝜑𝜆 (𝑦)
𝑔𝜆 (𝑥, 𝑦) :=
𝑊𝜆 ⎩𝜑𝜆 (𝑥)𝜓𝜆 (𝑦)

if 𝑥 < 𝑦,
𝜓 ′ (0)𝜑𝜆 (0) − 𝜓𝜆 (0)𝜑′𝜆 (0)
with 𝑊𝜆 = 𝜆
.
𝑆 ′ (0)
if 𝑥 ≥ 𝑦

14

The Wronskian 𝑊𝜆 is then equal to
√︁
√
2
𝑏2− + 2𝜆𝜎−
2𝜆
𝑏−
𝑊𝜆 =
+
− 2.
2
𝜎+
𝜎−
𝜎−
√
√
In particular, 𝑊𝜆 / 𝜆 converges to 2/𝜎+ as 𝜆 converges to 0.
On the other hand, it follows from (20) and (21) that

𝜓𝜆 (𝑥) −−→ 𝜓0 (𝑥) := 1 and 𝜑𝜆 (𝑥) −−→ 𝜑0 (𝑥) := 1 when 𝑥 > 0
𝜆→0

𝜆→0

while
√
𝜓𝜆 (𝑥) = exp(𝑥 2𝜆) −−→ 𝜓0 (𝑥) := 1 and 𝜑𝜆 (𝑥) = 𝜑0 (𝑥) := 1 when 𝑥 ≥ 0.
𝜆→0

For a measurable function 𝑓 : R → R+ such that
convergence results imply that

∫︀
R

𝑓 d𝑀 < +∞, the above

√ ∫︁
𝜎+ ∫︁
√
𝜆 𝑔𝜆 (𝑥, 𝑦)𝑓 (𝑥)𝑚(𝑦) d𝑦 −−→
𝑓 (𝑦)𝑚(𝑦) d𝑦, ∀𝑥 ∈ R.
𝜆→0
R
2 R
√
We then define ℓ(𝜆):= 2/𝜎+ which is a constant function, and 𝛼:=1/2, the exponent
of 𝜆.
Let (ℳ𝑡 )𝑡≥0 be a Mittag-Leffler process of index 1/2 (it is the inverse of an
increasing stable process of index 1/2). The process 2−1/2 ℳ is equal in distribution
to the running maximum of a Brownian motion, or equivalently, to the local time
of a Brownian motion [14, Remark 2.9, p. 21].
2 −
𝑄𝑡 , 𝑡 ≥ 0
From Theorem 3.1 and Corollary 3.2 in [14, p. 26], since ⟨𝑀 − ⟩𝑡 = 𝜎−
and 𝑀 is continuous,
(︃√︁

−
𝑀𝑛𝑡
𝑄−
ℓ(𝑛) 1/4 , ℓ(𝑛) √𝑛𝑡
𝑛
𝑛

)︃
𝑡∈[0,1]

)︁
(︁ √
law
−
−
−
−
→
𝜎
− 𝜈𝐵 (ℳ𝑡 ), 𝜈ℳ𝑡
𝑛→∞

𝑡∈[0,1]

(28)

with respect to the uniform topology, where 𝐵 is a Brownian motion independent
from ℳ and
(︂∫︁ 1
)︂
∫︁
1
𝜈 = 𝑚(𝑦)E𝑦
1𝜉𝑠 ≤0 d𝑠 d𝑦 = 𝑀 (R− ) = .
𝑏−
R
0
From the reflection principle,√the distribution of ℳ1 is the same as the one of a
truncated normal distribution 2𝒯 where 𝒯 := |𝒢| with 𝒢 ∼ 𝒩 (0, 1). Setting 𝑡 = 1
in (28) and using the scaling property of the Brownian motion 𝐵 − ,
(︃

𝑀𝑇− 𝑄−
𝑇
, 1/2
1/4
𝑇
𝑇

)︃

)︃
√ √
𝜎+
𝜎
+
−−−→ 𝜎− √
𝒯 · 𝒩 −,
𝒯
𝑇 →∞
𝑏−
𝑏−
law

(︃

for a Gaussian random variable 𝒩 − ∼ 𝒩 (0, 1) independent from 𝒯 .
It remains to show the independence of 𝒩 + , 𝒩 − and 𝒯 .
2 ±
Since ⟨𝑀 ± ⟩𝑡 = 𝜎±
𝑄𝑡 and ⟨𝑀 + , 𝑀 − ⟩𝑡 = 0 for any 𝑡 ≥ 0, the Knight theorem [18, Theorem 7.3’, p. 92] implies that there exists on an extension of (Ω, ℱ, P)
2 ±
a 2-dimensional Brownian motion (𝐵 + , 𝐵 − ) such that 𝑀𝑡± = 𝐵 ± (𝜎±
𝑄𝑡 ) for any

15

√
𝑡 ≥ 0. Let us set 𝐵𝑛+ (𝑡) = 𝑛−1/2 𝐵 + (𝑛𝑡) and 𝐵𝑛− (𝑡) = 𝑛−1/4 𝐵 − ( 𝑛𝑡) for any integer
𝑛 and any 𝑡 ≥ 0. From the scaling property, (𝐵𝑛+ , 𝐵𝑛− ) is still a 2-dimensional Brownian motion which converges in distribution to a 2-dimensional Brownian motion
+
−
(𝐵∞
, 𝐵∞
) in the space 𝒞([0, 1], R2 ) of continuous functions.
+
+
For any 0 ≤ 𝑠 ≤ 𝑡 ≤ 1, 𝑄+
𝑛𝑡 − 𝑄𝑛𝑠 ≤ 𝑛(𝑡 − 𝑠) so that (𝑄𝑛𝑡 /𝑛)𝑡∈[0,1] is tight in the
space of continuous functions. Hence, (𝑄+
𝑛𝑡 /𝑛)𝑡∈[0,1] converges in probability to the
identity map 𝑡 ↦→ 𝑡 in the space of continuous function
𝒞([0, 1], R).
(︁
)︁
+
−1/2 −
Combining this result with (28), it holds that 𝐵𝑛 (𝑡), 𝐵𝑛− (𝑡), 𝑛−1 𝑄+
𝑄𝑛𝑡
𝑛𝑡 , 𝑛

𝑡∈[0,1]

is tight in 𝒞([0, 1], R4 ) and then necessarily
(︃

𝐵𝑛+ (𝑡), 𝐵𝑛− (𝑡),

−
𝑄+
𝑛𝑡 𝑄𝑛𝑡
,√
𝑛
𝑛

)︃

law

𝑡∈[0,1]

−
+
(𝑡), 𝑡, 𝜈ℳ𝑡 )𝑡∈[0,1]
(𝑡), 𝐵∞
−
−−→ (𝐵∞
𝑛→∞

in the space of continuous functions 𝒞([0, 1], R4 ). Being the inverse of a 1/2-stable
+
−
process, hence a pure jump process, ℳ is independent from (𝐵∞
, 𝐵∞
) for the arguments presented in [14, p. 38] or [18, Theorem 6.3, p. 77].
For any 𝑡 ∈ [0, 1], it holds that
+
+
𝑀𝑛𝑡
1
2 +
2 𝑄𝑛𝑡
√ = √ 𝐵 + (𝜎+
𝑄𝑛𝑡 ) = 𝐵𝑛+ 𝜎+
𝑛
𝑛
𝑛

(︃

)︃

−
+
𝑀𝑛𝑡
1
2 −
2 𝑄𝑛𝑡
√
and 1/4
= √ 𝐵 − (𝜎−
𝑄𝑛𝑡 ) = 𝐵𝑛+ 𝜎−
.
𝑛
𝑛
𝑛

(︃

)︃

−1/2 −
Using (𝑛−1 𝑄+
𝑄𝑛𝑡 )𝑡∈[0,1] as random time changes, we deduce from
𝑛𝑡 )𝑡∈[0,1] and (𝑛
the results in [5, p. 144] that

(︃

−
−
+
𝑀𝑛𝑡
𝑄+
𝑀𝑛𝑡
𝑛𝑡 𝑄𝑛𝑡
√
,
,
,
𝑛1/2 𝑛1/4 𝑛
𝑛

)︃

law

2
2
−−−→ (𝐵 + (𝜎+
𝑡), 𝐵 − (𝜎−
𝜈ℳ𝑡 ), 𝑡, 𝜈ℳ𝑡 )𝑡∈[0,1] .
𝑡∈[0,1]

𝑛→∞

√
√
2
2
)/𝜎+ in (27), 𝒯 := ℳ1 / 2 and 𝒩 − := 𝐵 − (𝜎−
𝜈)/𝜎− 𝜈, this
By setting 𝒩 + := 𝐵 + (𝜎+
proves (9) using 𝑡 = 1 in the above limit.

6.4 Asymptotic behavior for the transient case (T0)
We recall that we consider only 𝑏+ > 0 and 𝑏− ≥ 0.
It is known that the last passage time ℓ0 to 0 is almost surely finite so that
−
+
−
𝑄∞ < ∞ almost surely. Since 𝑄+
𝑇 = 𝑇 − 𝑄𝑇 , we obtain that 𝑄𝑇 /𝑇 converges
+
almost surely to 1. The convergence results regarding 𝛽𝑇 follows from Proposition 8
applied only to one component.
The asymptotic behavior of 𝛽𝑇+ follows from Proposition 7(ii).
When 𝑇 > ℓ0 and 𝜉0 = 0, then 𝜉𝑇 > 0 and thus 𝛽𝑇− = 𝐿𝑇 (𝜉)/2𝑄−
𝑇 . Yet the local
−
time 𝐿𝑇 (𝜉) and the occupation times 𝑄𝑇 are constant when 𝑇 > ℓ0 . The result
follows by the computations of Section 5.3.

6.5 Asymptotic behavior for the transient case generated by
diverging drift (T1)
When 𝑏− < 0 and 𝑏+ > 0, the process is also transient as 𝑆(+∞) < +∞ and
2
𝑆(−∞) > −∞. The scale function 𝑆 map R to (−𝛾− , 𝛾+ ) with 𝛾± = |𝜎±
|2𝑏± . From

16

the Feller test [23, Theorem 5.29], the process does not explode. Thus, as 𝜉0 = 0, it
follows from [48] or [23, Proposition 5.5.22, p. 354] that
(︂

)︂

(︂

)︂

𝑝 := P lim 𝜉𝑇 = +∞ = 1 − P lim 𝜉𝑇 = −∞ =
𝑇 →∞

𝑇 →∞

𝜎 2 𝑏+
𝛾−
= 2 − 2 .
𝛾− + 𝛾+
𝜎+ 𝑏− + 𝜎− 𝑏+

Then event that {lim𝑇 →∞ |𝜉𝑇 | = +∞} arise when the process starts an excursion
with infinite lifetime, thus after the last passage time to 0. We denote by 𝐴± the
+
−
+
−
event {lim𝑇 →∞ 𝑄±
𝑇 /𝑇 = 1}, so that 𝐴 ∩ 𝐴 = ∅. Hence, P(𝐴 ) = 1 − P(𝐴 ) = 𝑝.
Using the same arguments as in the case T0, given 𝐴± , 𝛽𝑇± converges almost
surely to 𝑏± while 𝛽𝑇∓ = ±𝐿∞ (𝜉)/𝑄∓
ℓ0 .
√
For the Central Limit Theorem, we apply Corollary 2.3 in [7] on 𝑀𝑇+ / 𝑇 . As
𝑄+
𝑇 > 0 a.s. as soon as 𝑇 > 0 since 𝜉0 = 0 and
{ℬ = 𝑝} = 𝐴+ , a.s. for 𝑝 = 0, 1,
it follows that for a normal distribution 𝒩 ,
√︃

𝑇
2 +
𝜎+
𝑄𝑇

𝑀𝑇 ℱ∞ −stably
· √ −−−−−−→ 𝒩 under P(· | ℬ = 1).
𝑇 𝑁 →∞

It follows that
√ 𝑀𝑇
𝑇
− 𝑏+ = 𝑇 + = 𝜎 +
𝑄𝑇
𝑄+
𝑇
√︃

𝛽𝑇+

√︃

𝑇
2 +
𝜎+
𝑄𝑇

𝑀 ℱ∞ −stably
√ 𝑇 −−
−−−−→ 𝜎+ 𝒩 under P(· | ℬ = 1).
𝑇 𝑁 →∞

Hence the result.

7 Wilk’s theorem and LAN property
Owing to the quadratic nature of the log-likelihood, we easily deduce both a Wilk
theorem, on which a hypothesis test may be developed, as well as the Local Asymptotic Normality (LAN) property, which proves that our estimators are asymptotically
efficient.

7.1 Wilk’s theorem and a hypothesis test
The log-likelihood log 𝐺(𝑏+ , 𝑏− ) can be computed from the data using (5). Moreover,
this function is quadratic in 𝑏+ and 𝑏− . We have then
⎤

⎡

𝐷(𝑏+ , 𝑏− ) := ∇ log 𝐺(𝑏+ , 𝑏− ) =

L𝜉0 M+
L𝜉𝑇 M+
−
− 2𝜎12 𝐿𝑇 (𝜉) − 𝜎𝑏+2 𝑄+
2
2
𝑇 ⎥
𝜎+
⎢ 𝜎+
+
+
⎦
⎣ L𝜉 M−
−
L𝜉
M
𝑏
− 𝜎𝑇2 + 𝜎02 + 2𝜎12 𝐿𝑇 (𝜉) − 𝜎−2 𝑄−
𝑇
−
−
−
−

⎡

𝑄+

and 𝐻(𝑏+ , 𝑏− ) := Hess log 𝐺(𝑏+ , 𝑏− ) = ⎣

+

𝑇
⎢− 𝜎2

0

⎤

0 ⎥
𝑄−

− 𝜎2𝑇

−

17

⎦.

Therefore, around any point (𝑏+ , 𝑏− ),
log 𝐺(𝑏+ + Δ𝑏+ , 𝑏− + Δ𝑏− ) = log 𝐺(𝑏+ , 𝑏− )
𝑄+
𝑄−
Δ𝑏+
+ 𝐷(𝑏+ , 𝑏− ) ·
− 𝑇2 Δ𝑏2+ − 𝑇2 Δ𝑏2− . (29)
Δ𝑏−
2𝜎+
2𝜎−
[︃

]︃

In particular, we prove a result in the asymptotic behavior of the log-likelihood for
the ergodic case E or the null recurrent case N1. A similar result can be given for
the null recurrent cases N0 with a different limit distribution that can be identified
with the density given in Section 6.2.
Proposition 9 (Wilk’s theorem; Ergodic case E or null recurrent case N1). Denote
true
by (𝑏true
+ , 𝑏− ) be the real parameters. Then
law

true
+ 2
− 2
true
− 2 log 𝐺(𝑏true
+ , 𝑏− ) −−−→ 𝜒2 := (𝒩 ) + (𝒩 ) under P(𝑏true
+ ,𝑏− )
𝑇 →∞

(30)

for two independent, unit Gaussian random variables 𝒩 + and 𝒩 − . Besides, when
true
(𝑏+ , 𝑏− ) ̸= (𝑏true
+ , 𝑏− ), then
a.s.

true .
− log 𝐺(𝑏+ , 𝑏− ) −−−→ + ∞ under P(𝑏true
+ ,𝑏− )

𝑇 →∞

(31)

Proof. Considering (29) at (𝑏+ , 𝑏− ) = (𝛽𝑇+ , 𝛽𝑇− ) since 𝐷(𝛽𝑇+ , 𝛽𝑇− ) = (0, 0), for any
parameter (𝑏+ , 𝑏− ) and any 𝛼+ , 𝛼− > 0,
log 𝐺(𝑏+ , 𝑏− ) = −

(︁ 𝛼+
)︁
)︁
(︁ 𝛼−
𝑄−
𝑄+
𝑇
𝑇
+ 2
− 2
2 (𝑏
2 (𝑏
−
.
𝑇
−
𝛽
)
−
𝛽
)
𝑇
+
−
𝑇
𝑇
2
2
2𝑇 𝛼+ 𝜎+
2𝑇 𝛼− 𝜎−

When the process is ergodic (Case E), we set (𝛼+ , 𝛼− ) = (1, 1). It follows from
±
true
Proposition 2 that (30) holds. When (𝑏+ , 𝑏− ) ̸= (𝑏true
+ , 𝑏− ), then 𝛽𝑇 − 𝑏± does not
+
converge to 0 while 𝑄𝑇 converges a.s. to infinity. This proves (31). The result is
similar in the case N1 with (𝛼+ , 𝛼− ) = (1, 1/2).
A hypothesis test can be developed from Proposition 9. The null hypothesis
is (𝑏+ , 𝑏− ) = (𝑏0+ , 𝑏0− ) for a given drift (𝑏0+ , 𝑏0− ) while the alternative hypothesis is
(𝑏+ , 𝑏− ) ̸= (𝑏0+ , 𝑏0− ).
Using (5), we compute log 𝐺(𝑏0+ , 𝑏0− ). The null hypothesis is rejected with a confidence level 𝛼 if − log 𝐺(𝑏0+ , 𝑏0− ) > 𝑞𝛼 where 𝑞𝛼 is the 𝛼-quantile P[𝜒2 ≤ 𝑞𝛼 ] = 𝛼
while 𝜒2 follows a 𝜒2 distribution with 2 degrees of freedom.

7.2 The LAN property
The LAN (local asymptotic normal) property, introduced by L. Le Cam in [29],
characterizes the efficiency of the estimator (See also [16, 27] among many other
references). It was extended as the LAMN (local asymptotic mixed normal) to deal
with a mixed normal limits by P. Jeganathan in [22].
The quadratic nature of the log-likelihood as well as our limit theorems implies
that the LAN (resp. LAMN) property is verified in the ergodic case E (resp. null
recurrent case N1).

18

Proposition 10 (LAN property; Ergodic case E). In the ergodic case E,√the LAN√
2
true
2
property holds for the likelihood at (𝑏true
+ , 𝑏− ) with rate of convergence (𝜎+ / 𝑇 , 𝜎− / 𝑇 )
and asymptotic Fisher information
]︃

[︃

2
1
|𝑏−|
0
𝜎+
.
Γ :=
2
0
𝜎− |𝑏+|
|𝑏− | + |𝑏+ |
true
Proof. At (𝑏true
+ , 𝑏− ), the gradient may be written
𝑄+ √
√ 𝑇
𝑇 (𝛽𝑇+
2
⎢ 𝑇 𝜎+
⎣ 𝑄− √
√ 𝑇
𝑇 (𝛽𝑇−
𝑇 𝜎2

⎡
true
𝐷(𝑏true
+ , 𝑏− ) =

⎤

− 𝑏true
+ )⎥

−

⎦.

− 𝑏true
− )

Using (29),
(︁

𝑅(𝑇 ) := log

2 Δ𝑏
2 Δ𝑏
√ + true
√−
𝐺 𝑏true
+ + 𝜎+ 𝑇 , 𝑏− + 𝜎− 𝑇

)︁

true
𝐺(𝑏true
+ , 𝑏− )
]︃ [︃
]︃
2 −
2 +
+
true
𝜎−
𝑄𝑇
𝑄𝑇
𝜎+
Δ𝑏+
𝑄+
2
𝑇 (𝛽𝑇 − 𝑏+ )
·
−
=
Δ𝑏
Δ𝑏2−
−
−
−
true
+
Δ𝑏−
𝑄𝑇 (𝛽𝑇 − 𝑏− ).
2𝑇
2𝑇
]︃
]︃
[︃
]︃
[︃
[︃
]︃ [︃
1 Δ𝑏+
1 𝑀𝑇+
Δ𝑏
Δ𝑏+
+
(32)
· ⟨𝑀 + , 𝑀 − ⟩𝑇
−
=√
− ·
Δ𝑏−
Δ𝑏−
2𝑇 Δ𝑏−
𝑇 𝑀𝑇

[︃

With 𝑐 = |𝑏+ | + |𝑏− |, Proposition 2 implies that (𝑇 −1/2 𝑀 + ,−1/2 𝑀 − ) converges in
distribution to 𝒢 ∼ 𝒩 (0, Γ) and 𝑇 −1 ⟨𝑀 + , 𝑀 − ⟩𝑇 converges to the diagonal, definite
positive matrix Γ. This proves the LAN property.
Proposition 11 (LAMN property; Null recurrent case N1). In the null recurrent
true
case N1, the LAMN property holds for the likelihood at (𝑏true
+ , 𝑏− ) with rate of
2
1/2
2
1/4
convergence (𝜎+ /𝑇 , 𝜎− /𝑇 ) and asymptotic (random) Fisher information
𝜎2
0
Γ := + 2 𝜎+
0 𝜎− 𝑏− |𝒩 |
[︃

]︃

for a normal random variable 𝒩 ∼ 𝒩 (0, 1).

8 Simulation study
8.1 From continuous to discrete data
In this section, we apply the estimator (4) to simulated processes. We test whether
the results are good or not depending on sign and magnitude of involved quantities
(cf. [36, Section 4.1]).
+
Estimator (4) involves (𝜉0 , 𝜉𝑇 ) which are observable, as well as (𝐿𝑇 (𝜉), 𝑄−
𝑇 , 𝑄𝑇 )
which are not.
Imagine that we observe 𝜉𝑡 on a discrete time grid {𝑘𝑇 /𝑛; 𝑘 = 0, . . . , 𝑛}. The
time step between two observations is Δ𝑡 = 𝑇 /𝑛. In [30] we discuss, in the case of

19

the Oscillating Brownian Motion, the estimator for the occupation time 𝑄+
𝑇 given
the Riemann sums:
𝑛
̂︀ + = Δ𝑡
𝑄
𝑇,𝑛

∑︁

1𝜉𝑘𝑇 /𝑛 ≥0 .

𝑘=1

√
We proved that the speed of convergence is strictly better than 𝑛, meaning that
√ ̂︀ +
P
𝑛(𝑄𝑇,𝑛 − 𝑄+
−−→ 0. This result can be extended to the drifted process 𝜉 via
𝑇) −
𝑛→∞
the Girsanov theorem.
Again in [30], we propose the following approximation of the local time of 𝜉 at 0
̂︀
𝐿
𝑇,𝑛

−3
=
2

√︂

𝑛 (︁
)︁ (︁
)︁
𝜋 𝜎+ + 𝜎− ∑︁
L𝜉𝑘𝑇 /𝑛 M+ − L𝜉(𝑘−1)𝑇 /𝑛 M+ · L𝜉𝑘𝑇 /𝑛 M− − L𝜉(𝑘−1)𝑇 /𝑛 M− .
2Δ𝑡 𝜎+ 𝜎− 𝑘=1
(33)

Up to some constants, this is essentially an approximation of the covariation
between L𝜉𝑡 M+ and L𝜉𝑡 M− . In [30] we also proved the consistence of such estimator.
For the Brownian motion, this estimator converges at speed 𝑛1/4 . The proof has not
been adapted to our case because of the technical difficulties due to the discontinuity
of the coefficients in 0. Anyway we conjecture a rate of 1/4. An alternative estimator
can be constructed counting the number of crossings, as in [20]. We use (33) for
empirical reasons, since it performs better on simulated trajectory.
At this point, we also notice that estimator (33) involves 𝜎± , which are quantities not directly observable on the trajectory; let us mention that the same is true
for the classic estimator of local time for SDEs with differentiable coefficients [20].
Therefore, if 𝜎± are known a priori, (33) can be implemented using such quantities. Otherwise, 𝜎± must be estimated on discrete time observations of 𝜉 as well, a
problem which has been thoroughly investigated in [31].
We do not push further in the present paper the theoretical discussion on the quality of these discrete time approximations. Some more insights, based on numerical
results, are given in the following section.

8.2 Implementation and simulation
The aim of the following section is to show on figures the numerical evidence of the
central limit theorems stated in Section 4. We also mean to say something more on
the choice of the step of the time grid in relation to the quality of the estimation of the
local time. The code used for the following simulations have been implemented using
the software R. We will consider time grids of the type {0, 𝑇 /𝑁, 2𝑇 /𝑁, . . . , 𝑇 } ⊂ N,
so with 𝑁 ∈ N, and use as approximation of local and occupation time the following:
̂︀
̂︀ +
̂︀ +
𝐿
𝑇,𝑁 and 𝑄𝑇 := 𝑄𝑇,𝑁 .

Summing up, as estimators of 𝑏± we use the approximation
̂︀
𝑏±

𝑇,𝑁

̂︀
L𝜉𝑇 M± − L𝜉0 M± − 𝐿
𝑇,𝑁 /2
.
=±
±
̂︀
𝑄𝑇,𝑁

20

(34)

These are our discrete times approximations of the continuous time estimators for
which the convergence results have been proved. In practice, in some cases the
estimator does not really depend on the local time, but is essentially determined by
the final value of the process and the occupation times. In those cases, the quality
of the discrete time approximation of the local time does not really matter, and
therefore we can take 𝑁 = 1. In most cases, anyway, a good approximation of the
local time is needed in order to observe on simulations the theoretical central limit
behavior expected from Section 4. In these cases 𝑁 ∈ N must be taken large enough.
In what follows, the choice of the parameters is detailed for every figure. The
diffusion parameter is taken constant 𝜎+ = 𝜎− = 0.01, the same one for all the
different simulations, and supposed known a priori. In such manner the CLTs that
we want to test can be better observed. To use estimator (34) on real data, with 𝜎±
not known, one must first use the estimators for 𝜎± in [31] and then use such
estimations in the estimator (33) for the local time. We indicate with [+] and [−]
estimation on positive and negative semiaxis, e.g. (N1)[+] stands for “estimation
of 𝑏+ ” in case N1.
(E)[−]

15
0

0

5

5

10

T (β−T − b−)

15
10

T (β+T − b+)

20

20

25

25

30

(E)[+]

−0.04

0.00

0.04

−0.06

−0.02

0.02

0.06

Figure 1: (E). SDE parameters: 𝜎± = 0.01, 𝑏− = 0.004, 𝑏+ = −0.003. Simulation
parameters: 𝑇 = 1000, 𝑁 = 100 000. We show both √
sides (positive and
negative) of the estimation, displaying the density of 𝑇 (𝛽𝑇± − 𝑏± ). The
CLT in Proposition (7) is accurate for large 𝑇 and time step 𝑇 /𝑁 small,
since the quality of the estimation of the local time is key in this case. The
limit behavior is Gaussian.

References
[1] Mohamed Alaya Ben and Ahmed Kebaier. Parameter estimation for the squareroot diffusions: ergodic and nonergodic cases. Stoch. Models, 28(4):609–634,
2012.

21

−0.04

−0.02

0.00

0.02

0.04

−0.02

0.00

0.02

0.04

40
30
10
0

10
0
−0.04

20

T (β−T − b−)

30

40

(T1)[−] wrt P−

20

T (β+T − b+)

30
0

10

20

T (β+T − b+)

30
20
0

10

T (βT+ − b+)

(T1)[+] wrt P+

40

(T0)[+]

40

(N1)[+]

−0.04

−0.02

0.00

0.02

0.04

−0.04

−0.02

0.00

0.02

0.04

Figure 2: (N1)[+], (T0)[+], (T1)[+] w.r.t P+ , (T1)[−] w.r.t P− . SDE parameters:
𝜎± = 0.01; in case N1: 𝑏− = 0.004, 𝑏+ = 0; in case T0: 𝑏− = 0.004, 𝑏+ =
0.006; in case T1: 𝑏− = −0.004, 𝑏+ = 0.003.√Simulation parameters:
𝑇 = 1000, 𝑁 = 1000. We
the density of 𝑇 (𝛽𝑇+ − 𝑏+ ) in cases N1
√ display
and T0. We also show 𝑇 (𝛽𝑇± −𝑏± ) in case T1, but the density is w.r.t P±
(cf. (6)). This is approximated computing the estimator on trajectories
such that 𝜉𝑇 is larger or respectively smaller than 0. In all these cases the
CLT is Gaussian and we do not need to have a fine discretisation/time
grid, since the local time is asymptotically negligible and the quantities
which matter in the estimator are 𝜉𝑇 and the occupation times. This
accounts of (9)-positive part, (10), (12) and (13).

[2] T. Appuhamillage, V. Bokil, E. Thomann, E. Waymire, and B. Wood. Occupation and local times for skew Brownian motion with applications to dispersion
across an interface. Ann. Appl. Probab., 21(1):183–214, 2011.
[3] Soren Asmussen and Michael Taksar. Controlled diffusion models for optimal
dividend pay-out. Insurance: Mathematics and Economics, 20(1):1–15, 1997.
[4] Federico M. Bandi and Peter C. B. Phillips. Fully nonparametric estimation of
scalar diffusion models. Econometrica, 71(1):241–283, 2003.
[5] Patrick Billingsley. Convergence of probability measures. Wiley Series in Probability and Statistics: Probability and Statistics. John Wiley & Sons, Inc., New
York, second edition, 1999. A Wiley-Interscience Publication.
[6] Michael P Clements and Jeremy Smith. A Monte Carlo study of the forecasting
performance of empirical SETAR models. Journal of Applied Econometrics,
14(2):123–141, March 1999.
[7] Irene Crimaldi and Luca Pratelli. Convergence results for multivariate martingales. Stochastic Process. Appl., 115(4):571–577, 2005.
[8] D. A. Darling and M. Kac. On occupation times for Markoff processes. Transactions of the American Mathematical Society, 84(2):444–458, 1957.

22

30
20
0

10

T(1 4)(β−T − b−)

40

50

(N1)[−]

−0.04

−0.02

0.00

0.02

0.04

Figure 3: (N1)[−]. SDE parameters: 𝜎± = 0.01; in case N1: 𝑏− = 0.004; 𝑏+ = 0.
Simulation parameters: 𝑇 = 1000, 𝑁 = 100 000. The CLT in (9)-negative
part, is accurate for large 𝑇 and time step 𝑇 /𝑁 small, since the quality of
the estimation of the local time is key in this case. Remark that in this
case (null recurrent), the CLT has speed of convergence 𝑇 1/4 and the limit
law is not Gaussian. This accounts of (9)-negative part.

[9] Marc Decamps, Marc Goovaerts, and Wim Schoutens. Self exciting threshold
interest rates models. Int. J. Theor. Appl. Finance, 9(7):1093–1122, 2006.
[10] Pierre Étoré. On random walk simulation of one-dimensional diffusion processes
with discontinuous coefficients. Electron. J. Probab., 11:no. 9, 249–275, 2006.
[11] William Feller. Generalized second order differential operators and their lateral
conditions. Illinois J. Math., 1:459–504, 1957.
[12] A. Gairat and V. Shcherbakov. Density of skew Brownian motion and its functionals with application in finance. Mathematical Finance, 26(4):1069–1088,
2016.
[13] Hans U. Gerber Gerber and Elias S.W. Shiu. On optimal dividends: From
reflection to refraction. Journal of Computational and Applied Mathematics,
186(1):4 – 22, 2006. Special Issue: Jef Teugels.
[14] R. Höpfner and E. Löcherbach. Limit theorems for null recurrent Markov processes. Mem. Amer. Math. Soc., 161(768), 2003.
[15] Scott Hottovy and Samuel N. Stechmann. Threshold models for rainfall and
convection: Deterministic versus stochastic triggers. SIAM Journal on Applied
Mathematics, 75(2):861–884, 2015.

23

− β−T

0.00

0.01

0.02

0.03

0.04

0.05

10
0

0

0

5

20

20

10

30

15

β+T

40

β−T

40

20

50

60

25

60

70

(T1)[−], wrt P+

30

(T1)[+], wrt P−

80

(T0)[−]

0.00

0.01

0.02

0.03

0.04

0.05

0.00

0.01

0.02

0.03

0.04

0.05

Figure 4: (T0)[−], (T1)[+] w.r.t to P− and (T1)[−] w.r.t to P+ . SDE parameters:
𝜎± = 0.01; in case T0: 𝑏− = 0.004, 𝑏+ = 0.003; in case T1: 𝑏− =
−0.003, 𝑏+ = 0.01. In case T0: simulation parameters: 𝑇 = 20, 𝑁 =
1000. We display the density of 𝛽𝑇− . In case T1: simulation parameters:
𝑇 = 20, 𝑁 = 3 × 107 . We display the density of 𝛽𝑇± w.r.t P∓ (cf. (6)).
This is approximated computing the estimator on trajectories such that
𝜉𝑇 is smaller or respectively larger than 0. In these cases the estimator is
not consistent, so what we show is not actually a CLT but the convergence
of the estimators towards the law (22) (cf. results (11), (14), (15)). This
convergence is accurate for large 𝑇 but also depends on the time step
𝑇 /𝑁 . Moreover, we see that the theoretical distribution of 𝛽𝑇− in case T1
is almost singular at the origin, and therefore the exact behavior near the
origin is hard to catch on simulated trajectories. This can be improved
using different kernels (instead of the Gaussian one) in the estimation of
the density. This can be easily done with the function “density” in R.
The limit behavior is better approximated when 𝑏+ and 𝑏− have similar
magnitude. We chose here to display the case 𝑏− = −0.003; 𝑏+ = 0.01
to mention this critical behavior. Anyway, this feature does not really
matter in statistical application, because in this case the estimator not
only is not consistent, but does not even guess the correct sign of the
parameter. Indeed, we are here in the very critical case of a transient
process generated by diverging drift T1.

[16] I. A. Ibragimov and R. Z. Has’minskiı̆. Statistical estimation, volume 16 of
Applications of Mathematics. Springer-Verlag, New York-Berlin, 1981.
[17] Tomoyuki Ichiba, Ioannis Karatzas, and Mykhaylo Shkolnikov. Strong solutions
of stochastic equations with rank-based coefficients. Probab. Theory Related
Fields, 156(1-2):229–248, 2013.
[18] Nobuyuki Ikeda and Shinzo Watanabe. Stochastic differential equations and
diffusion processes, volume 24 of North-Holland Mathematical Library. North-

24

(N0)[−]

0

0

5

5

10

Tβ−T
10

Tβ+T

15

15

20

20

25

(N0)[+]

−0.3

−0.2

−0.1

0.0

0.1

−0.1

0.0

0.1

0.2

0.3

Figure 5: (N0). SDE parameters: 𝜎± = 0.01, 𝑏− = 0 𝑏+ = 0;. Simulation parameters: 𝑇 = 10, 𝑁 = 1000; 𝑇 = 100, 𝑁 = 10 000; 𝑇 = 1000, 𝑁 = 100 000.
Differently from before, we do not show the convergence to the scaled limit
law (25), but the scaling relation (8), for three different final times. We
show it on both positive and negative semiaxes. Because of the estimation
of the local time, this also depends on the choice of 𝑁 .

Holland Publishing Co., Amsterdam; Kodansha, Ltd., Tokyo, second edition,
1989.
[19] K. Itô and H.P.J. McKean. Diffusion Processes and their Sample Paths: Reprint
of the 1974 Edition. Classics in Mathematics. Springer Berlin Heidelberg, 1996.
[20] Jean Jacod. Rates of convergence to the local time of a diffusion. Ann. Inst.
H. Poincaré Probab. Statist., 34(4):505–544, 1998.
[21] Jean Jacod and Albert N. Shiryaev. Limit theorems for stochastic processes,
volume 288 of Grundlehren der Mathematischen Wissenschaften. SpringerVerlag, Berlin, second edition, 2003.
[22] P. Jeganathan. On the asymptotic theory of estimation when the limit of the
log-likelihood ratios is mixed normal. Sankhyā Ser. A, 44(2):173–212, 1982.
[23] Ioannis Karatzas and Steven E. Shreve. Brownian motion and stochastic calculus, volume 113 of Graduate Texts in Mathematics. Springer-Verlag, New York,
second edition, 1991.
[24] Julian Keilson and Jon A. Wellner. Oscillating Brownian motion. J. Appl.
Probability, 15(2):300–310, 1978.
[25] Yu. A. Kutoyants. Efficient density estimation for ergodic diffusion processes.
Statistical Inference for Stochastic Processes, 1(2):131–155, 1998.
[26] Yury A. Kutoyants. On identification of the threshold diffusion processes. Ann.
Inst. Statist. Math., 64(2):383–413, 2012.

25

[27] Lucien Le Cam and Grace Lo Yang. Asymptotics in statistics. Springer Series
in Statistics. Springer-Verlag, New York, second edition, 2000.
[28] J.-F. Le Gall. One-dimensional stochastic differential equations involving the
local times of the unknown process. Stochastic Analysis. Lecture Notes Math.,
1095:51–82, 1985.
[29] Lucien LeCam. On some asymptotic properties of maximum likelihood estimates and related Bayes’ estimates. Univ. California Publ. Statist., 1:277–329,
1953.
[30] A. Lejay and P. Pigato. Statistical estimation of the Oscillating Brownian
Motion. Bernoulli, 2017. To appear.
[31] A. Lejay and P. Pigato. A threshold model for local volatility: evidence of
leverage and mean reversion effects on historical data, 2018. Preprint.
[32] D. Lépingle. Sur le comportement asymptotique des martingales locales. In
Séminaire de Probabilités, XII (Univ. Strasbourg, Strasbourg, 1976/1977), volume 649 of Lecture Notes in Math., pages 148–161. Springer, Berlin, 1978.
[33] R.S. Lipster and A.N. Shiryaev. Statistics of random processes. II. Applications.
Springer, 2001.
[34] Gisirô Maruyama and Hiroshi Tanaka. Some properties of one-dimensional
diffusion processes. Mem. Fac. Sci. Kyusyu Univ. Ser. A. Math., 11:117–141,
1957.
[35] Michael Monoyios and Lucio Sarno. Mean reversion in stock index futures
markets: A nonlinear analysis. Journal of Futures Markets, 22(4):285–30, April
2002.
[36] Pedro P. Mota and Manuel L. Esquível. On a continuous time stock price model
with regime switching, delay, and threshold. Quant. Finance, 14(8):1479–1488,
2014.
[37] Jim Pitman and Marc Yor. Hitting, occupation and inverse local times of
one-dimensional diffusions: martingale and excursion approaches. Bernoulli,
9(1):1–24, 02 2003.
[38] James M. Poterba and Lawrence H. Summers. Mean reversion in stock prices.
Journal of Financial Economics, 22(1):27–59, 1988.
[39] J. M. Ramirez, E. A. Thomann, and E. C. Waymire. Advection–dispersion
across interfaces. Statist. Sci., 28(4):487–509, 2013.
[40] Alfréd Rényi. On stable sequences of events. Sankhyā Ser. A, 25:293–302, 1963.
[41] L. C. G. Rogers and David Williams. Diffusions, Markov processes, and martingales. Vol. 2. Cambridge Mathematical Library. Cambridge University Press,
Cambridge, 2000. Itô calculus, Reprint of the second (1994) edition.

26

[42] F. Sattin. Fick’s law and Fokker–Planck equation in inhomogeneous environments. Physics Letters A, 37(22):3941–3945, 2008.
[43] Fei Su and Kung-Sik Chan. Quasi-likelihood estimation of a threshold diffusion
process. J. Econometrics, 189(2):473–484, 2015.
[44] Fei Su and Kung-Sik Chan. Testing for threshold diffusion. J. Bus. Econom.
Statist., 35(2):218–227, 2017.
[45] Howell Tong. Threshold models in nonlinear time series analysis, volume 21 of
Lecture Notes in Statistics. Springer-Verlag, New York, 1983.
[46] Howell Tong. Threshold models in time series analysis—30 years on. Stat.
Interface, 4(2):107–118, 2011.
[47] Hisao Watanabe and Minoru Motoo. Ergodic property of recurrent diffusion
processes. J. Math. Soc. Japan, 10:272–286, 1958.
[48] Shinzo Watanabe. Generalized arc-sine laws for one-dimensional diffusion processes and random walks. In Stochastic analysis (Ithaca, NY, 1993), volume 57
of Proc. Sympos. Pure Math., pages 157–172. Amer. Math. Soc., Providence,
RI, 1995.
[49] Pradeep K Yadav, Peter F Pope, and Krishna Paudyal. Threshold autoregressive modeling in finance: the price differences of equivalent assets. Mathematical
Finance, 4(2):205–221, April 1994.

27

