FitzHugh-Nagumo.eps
A Contraction Theory Approach to Stochastic
Incremental Stability
Quang-Cuong Pham ∗
LPPA, Collège de France
Paris, France
cuong.pham@normalesup.org
Nicolas Tabareau
LPPA, Collège de France
Paris, France
tabareau.nicolas@gmail.com
Jean-Jacques Slotine
Nonlinear Systems Laboratory, MIT
Cambridge, MA 02139, USA
jjs@mit.edu
September 6, 2018
Abstract
We investigate the incremental stability properties of Itô stochastic
dynamical systems. Specifically, we derive a stochastic version of non-
linear contraction theory that provides a bound on the mean square
distance between any two trajectories of a stochastically contracting
system. This bound can be expressed as a function of the noise inten-
sity and the contraction rate of the noise-free system. We illustrate
these results in the contexts of stochastic nonlinear observers design
and stochastic synchronization.
1 Introduction
Nonlinear stability properties are often considered with respect to an equi-
librium point or to a nominal system trajectory (see e.g. [31]). By contrast,
incremental stability is concerned with the behaviour of system trajectories
with respect to each other. From the triangle inequality, global exponential
incremental stability (any two trajectories tend to each other exponentially)
is a stronger property than global exponential convergence to a single tra-
jectory.
Historically, work on deterministic incremental stability can be traced
back to the 1950’s [23, 7, 16] (see e.g. [26, 20] for a more extensive list and
historical discussion of related references). More recently, and largely inde-
pendently of these earlier studies, a number of works have put incremental
To whom correspondance should be addressed.
http://arxiv.org/abs/0704.0926v2
stability on a broader theoretical basis and made relations with more tradi-
tional stability approaches [14, 32, 24, 2, 6]. Furthermore, it was shown that
incremental stability is especially relevant in the study of such problems as
state detection [2], observer design or synchronization analysis.
While the above references are mostly concerned with deterministic sta-
bility notions, stability theory has also been extended to stochastic dynam-
ical systems, see for instance [22, 17]. This includes important recent de-
velopments in Lyapunov-like approaches [12, 27], as well as applications to
standard problems in systems and control [13, 34, 8]. However, stochastic
versions of incremental stability have not yet been systematically investi-
gated.
The goal of this paper is to extend some concepts and results in in-
cremental stability to stochastic dynamical systems. More specifically, we
derive a stochastic version of contraction analysis in the specialized context
of state-independent metrics.
We prove in section 2 that the mean square distance between any two
trajectories of a stochastically contracting system is upper-bounded by a
constant after exponential transients. In contrast with previous works on
incremental stochastic stability [5], we consider the case when the two tra-
jectories are subject to distinct and independent noises, as detailed in sec-
tion 2.2.1. This specificity enables our theory to have a number of new and
practically important applications. However, the fact that the noise does
not vanish as two trajectories get very close to each other will prevent us
from obtaining asymptotic almost-sure stability results (see section 2.3.2).
In section 3, we show that results on combinations of deterministic con-
tracting systems have simple analogues in the stochastic case. These combi-
nation properties allow one to build by recursion stochastically contracting
systems of arbitrary size.
Finally, as illustrations of our results, we study in section 4 several ex-
amples, including contracting observers with noisy measurements, stochas-
tic composite variables and synchronization phenomena in networks of noisy
dynamical systems.
2 Main results
2.1 Background
2.1.1 Nonlinear contraction theory
Contraction theory [24] provides a set of tools to analyze the incremental
exponential stability of nonlinear systems, and has been applied notably
to observer design [24, 25, 1, 21, 36], synchronization analysis [35, 28] and
systems neuroscience modelling [15]. Nonlinear contracting systems enjoy
desirable aggregation properties, in that contraction is preserved under many
types of system combinations given suitable simple conditions [24].
While we shall derive global properties of nonlinear systems, many of our
results can be expressed in terms of eigenvalues of symmetric matrices [19].
Given a square matrix A, the symmetric part of A is denoted by As. The
smallest and largest eigenvalues of As are denoted by λmin(A) and λmax(A).
Given these notations, the matrix A is positive definite (denoted A > 0) if
λmin(A) > 0, and it is uniformly positive definite if
∃β > 0 ∀x, t λmin(A(x, t)) ≥ β
The basic theorem of contraction analysis, derived in [24], can be stated
as follows
Theorem 1 (Contraction) Consider, in Rn, the deterministic system
ẋ = f(x, t) (2.1)
where f is a smooth nonlinear function. Denote the Jacobian matrix of f
with respect to its first variable by ∂f
. If there exists a square matrix Θ(x, t)
such that M(x, t) = Θ(x, t)TΘ(x, t) is uniformly positive definite and the
matrix
F(x, t) =
Θ(x, t) +Θ(x, t)
Θ−1(x, t)
is uniformly negative definite, then all system trajectories converge exponen-
tially to a single trajectory, with convergence rate | sup
x,t λmax(F)| = λ > 0.
The system is said to be contracting, F is called its generalized Jacobian,
M(x, t) its contraction metric and λ its contraction rate.
2.1.2 Standard stochastic stability
In this section, we present very informally the basic ideas of standard stochas-
tic stability (for a rigourous treatment, the reader is referred to e.g. [22]).
This will set the context to understand the forthcoming difficulties and dif-
ferences associated with incremental stochastic stability.
For simplicity, we consider the special case of global exponential stability.
Let x(t) be a Markov stochastic process and assume that there exists a non-
negative function V (V (x) may represent e.g. the squared distance of x
from the origin) such that
∀x ∈ Rn ÃV (x) ≤ −λV (x) (2.2)
where λ is a positive real number and Ã is the infinitesimal operator of the
process x(t). The operator Ã is the stochastic analogue of the deterministic
differentiation operator. In the case that x(t) is an Itô process, Ã corre-
sponds to the widely-used [27, 34, 8] differential generator L (for a proof of
this fact, see [22], p. 15 or [3], p. 42).
For x0 ∈ Rn, let Ex0(·) = E(·|x(0) = x0). Then by Dynkin’s formula
([22], p. 10), one has
∀t ≥ 0 Ex0V (x(t))− V (x0) = Ex0
ÃV (x(s))ds
≤ −λEx0
V (x(s))ds = −λ
Ex0V (x(s))ds
Applying the Gronwall’s lemma to the deterministic real-valued function
t → Ex0V (x(t)) yields
∀t ≥ 0 Ex0V (x(t)) ≤ V (x0)e−λt
If we assume furthermore that Ex0V (x(t)) < ∞ for all t, then the above
implies that V (x(t)) is a supermartingale (see lemma 3 in the Appendix for
details), which yields, by the supermartingale inequality
T≤t<∞
V (x(t)) ≥ A
≤ Ex0V (x(T ))
≤ V (x0)e
(2.3)
Thus, one obtains an almost-sure stability result, in the sense that
∀A > 0 lim
T≤t<∞
V (x(t)) ≥ A
= 0 (2.4)
2.2 The stochastic contraction theorem
2.2.1 Settings
Consider a noisy system described by an Itô stochastic differential equation
da = f(a, t)dt+ σ(a, t)dW d (2.5)
where f is a Rn × R+ → Rn function, σ is a Rn × R+ → Rnd matrix-valued
function and W d is a standard d-dimensional Wiener process.
To ensure existence and uniqueness of solutions to equation (2.5), we
assume, here and in the remainder of the paper, the following standard
conditions on f and σ
Lipschitz condition: There exists a constant K1 > 0 such that
∀t ≥ 0, a,b ∈ Rn ‖f(a, t)− f(b, t)|+ ‖σ(a, t) − σ(b, t)‖ ≤ K1‖a− b‖
Restriction on growth: There exists a constant K2 > 0
∀t ≥ 0, a ∈ Rn ‖f(a, t)‖2 + ‖σ(a, t)‖2 ≤ K2(1 + ‖a‖2)
Under these conditions, one can show ([3], p. 105) that equation (2.5)
has on [0,∞[ a unique Rn-valued solution a(t), continuous with probability
one, and satisfying the initial condition a(0) = a0, with a0 ∈ Rn.
In order to investigate the incremental stability properties of system (2.5),
consider now two system trajectories a(t) and b(t). Our goal will consist of
studying the trajectories a(t) and b(t) with respect to each other. For this,
we consider the augmented system x(t) = (a(t),b(t))T , which follows the
equation
f(a, t)
f(b, t)
σ(a, t) 0
0 σ(b, t)
dW d1
dW d2
f(x, t)dt+
σ(x, t)dW 2d (2.6)
Important remark As stated in the introduction, the systems a and
b are driven by distinct and independent Wiener processes W d1 and W
This makes our approach considerably different from [5], where the authors
studied two trajectories driven by the same Wiener process.
Our approach enables us to study the stability of the system with respect
to variations in initial conditions and to random perturbations: indeed, two
trajectories of any real-life system are typically affected by distinct “real-
izations” of the noise. In addition, it leads very naturally to nice results on
the comparison of noisy and noise-free trajectories (cf. section 2.4), which
are particularly useful in applications (cf. section 4).
However, because of the very fact that the two trajectories are driven
by distinct Wiener processes, we cannot expect the influence of the noise
to vanish when the two trajectories get very close to each other. This con-
strasts with [5], and more generally, with the standard stochastic stability
case, where the noise vanishes near the origin (cf. section 2.1.2). The con-
sequences of this will be discussed in detail in section 2.3.2.
2.2.2 The basic stochastic contraction theorem
We introduce two hypotheses
(H1) f(a, t) is contracting in the identity metric, with contraction rate λ,
(i.e. ∀a, t λmax
≤ −λ)
(H2) tr
σ(a, t)Tσ(a, t)
is uniformly upper-bounded by a constant C (i.e.
∀a, t tr
σ(a, t)Tσ(a, t)
In other words, (H1) says that the noise-free system is contracting, while
(H2) says that the variance of the noise is upper-bounded by a constant.
Definition 1 A system that verifies (H1) and (H2) is said to be stochas-
tically contracting in the identity metric, with rate λ and bound C.
Consider now the Lyapunov-like function V (x) = ‖a−b‖2 = (a−b)T (a−
b). Using (H1) and (H2), we derive below an inequality on ÃV (x), similar
to equation (2.2) in section 2.1.2.
Lemma 1 Under (H1) and (H2), one has the inequality
ÃV (x) ≤ −2λV (x) + 2C (2.7)
Proof Since x(t) is an Itô process, Ã is given by the differential operator
L of the process [22, 3]. Thus, by the Itô formula
ÃV (x) = L V (x) =
∂V (x)
f(x, t) +
σ(x, t)T
∂2V (x)
σ(x, t)
1≤i≤2n
f(x, t)i +
1≤i,j,k≤2n
σ(x, t)ij
∂xi∂xk
σ(x, t)kj
1≤i≤n
f(a, t)i +
1≤i≤n
f(b, t)i
1≤i,j,k≤n
σ(a, t)ij
∂ai∂ak
σ(a, t)kj
1≤i,j,k≤n
σ(b, t)ij
∂bi∂bk
σ(b, t)kj
= 2(a− b)T (f(a, t) − f(b, t))
+tr(σ(a, t)Tσ(a, t)) + tr(σ(b, t)Tσ(b, t))
Fix t ≥ 0 and, as in [10], consider the real-valued function
r(µ) = (a− b)T (f(µa+ (1− µ)b, t)− f(b, t))
Since f is C1, r is C1 over [0, 1]. By the mean value theorem, there exists
µ0 ∈]0, 1[ such that
r′(µ0) = r(1)− r(0) = (a− b)T (f(a)− f(b))
On the other hand, one obtains by differentiating r
r′(µ0) = (a− b)T
(µ0a+ (1− µ0)b, t)
(a− b)
Thus, one has
(a− b)T (f(a)− f(b)) = (a− b)T
(µ0a+ (1− µ0)b, t)
(a− b)
≤ −λ(a− b)T (a− b) = −2λV (x) (2.8)
where the inequality is obtained by using (H1).
Finally,
ÃV (x) = 2(a− b)T (f(a)− f(b)) + tr(σ(a, t)Tσ(a, t)) + tr(σ(b, t)Tσ(b, t))
≤ −2λV (x) + 2C
where the inequality is obtained by using (H2). �
We are now in a position to prove our main theorem on stochastic incre-
mental stability.
Theorem 2 (Stochastic contraction) Assume that system (2.5) verifies
(H1) and (H2). Let a(t) and b(t) be two trajectories whose initial condi-
tions are given by a probability distribution p(x(0)) = p(a(0),b(0)). Then
∀t ≥ 0 E
‖a(t)− b(t)‖2
+ e−2λt
‖a0 − b0‖2 −
dp(a0,b0)
(2.9)
where [·]+ = max(0, ·). This implies in particular
∀t ≥ 0 E
‖a(t)− b(t)‖2
‖a(0) − b(0)‖2
e−2λt (2.10)
Proof Let x0 = (a0,b0) ∈ R2n. By Dynkin’s formula ([22], p. 10)
Ex0V (x(t)) − V (x0) = Ex0
ÃV (x(s))ds
Thus one has ∀u, t 0 ≤ u ≤ t < ∞
Ex0V (x(t)) − Ex0V (x(u)) = Ex0
ÃV (x(s))ds
≤ Ex0
(−2λV (x(s)) + 2C)ds (2.11)
(−2λEx0V (x(s)) + 2C)ds
where inequality (2.11) is obtained by using lemma 1.
Denote by g(t) the deterministic quantity Ex0V (x(t)). Clearly, g(t) is a
continuous function of t since x(t) is a continuous process. The function g
then satisfies the conditions of the Gronwall-type lemma 4 in the Appendix,
and as a consequence
∀t ≥ 0 Ex0V (x(t)) ≤
V (x0)−
e−2λT
Integrating the above inequality with respect to x0 yields the desired
result (2.9). Next, inequality (2.10) follows from (2.9) by remarking that
‖a0 − b0‖2 −
dp(a0,b0) ≤
‖a0 − b0‖2dp(a0,b0)
‖a(0) − b(0)‖2
(2.12)
Remark Let ǫ > 0 and Tǫ =
E(‖a0−b0‖2)
. Then inequal-
ity (2.10) and Jensen’s inequality [30] imply
∀t ≥ Tǫ E(‖a(t)− b(t)‖) ≤
C/λ+ ǫ (2.13)
Since ‖a(t)−b(t)‖ is non-negative, (2.13) together with Markov inequal-
ity [11] allow one to obtain the following probabilistic bound on the distance
between a(t) and b(t)
∀A > 0 ∀t ≥ Tǫ P (‖a(t)− b(t)‖ ≥ A) ≤
C/λ+ ǫ
Note however that this bound is much weaker than the asymptotic
almost-sure bound (2.4).
2.2.3 Generalization to time-varying metrics
Theorem 2 can be vastly generalized by considering general time-dependent
metrics (the case of state-dependent metrics is not considered in this article
and will be the subject of a future work). Specifically, let us replace (H1)
and (H2) by the following hypotheses
(H1’) There exists a uniformly positive definite metric M(t) = Θ(t)TΘ(t),
with the lower-bound β > 0 (i.e. ∀x, t xTM(t)x ≥ β‖x‖2) and f(a, t)
is contracting in that metric, with contraction rate λ, i.e.
Θ(t) +Θ(t)
Θ−1(t)
≤ −λ uniformly
or equivalently
M(t) +
M(t) ≤ −2λM(t) uniformly
(H2’) tr
σ(a, t)TM(t)σ(a, t)
is uniformly upper-bounded by a constant C
Definition 2 A system that verifies (H1’) and (H2’) is said to be stochas-
tically contracting in the metric M(t), with rate λ and bound C.
Consider now the generalized Lyapunov-like function V1(x, t) = (a −
b)TM(t)(a − b). Lemma 1 can then be generalized as follows.
Lemma 2 Under (H1’) and (H2’), one has the inequality
ÃV1(x, t) ≤ −2λV1(x, t) + 2C (2.14)
Proof Let us compute first ÃV1
ÃV1(x, t) =
f(x, t) +
σ(x, t)T
σ(x, t)
= (a− b)T
(a− b) + 2(a− b)TM(t)(f(a, t) − f(b, t))
+tr(σ(a, t)TM(t)σ(a, t)) + tr(σ(b, t)TM(t)σ(b, t))
Fix t > 0 and consider the real-valued function
r(µ) = (a− b)TM(t)(f(µa+ (1− µ)b, t)− f(b, t))
Since f is C1, r is C1 over [0, 1]. By the mean value theorem, there exists
µ0 ∈]0, 1[ such that
r′(µ0) = r(1)− r(0) = (a− b)TM(t)(f(a) − f(b))
On the other hand, one obtains by differentiating r
r′(µ0) = (a− b)TM(t)
(µ0a+ (1− µ0)b, t)
(a− b)
Thus, letting c = µ0a+ (1− µ0)b, one has
(a− b)T
(a− b) + 2(a− b)TM(t)(f(a) − f(b))
= (a− b)T
(a− b) + 2(a− b)TM(t)
(c, t)
(a− b)
= (a− b)T
M(t) +M(t)
(c, t)
(c, t)
(a− b)
≤ −2λ(a− b)TM(t)(a− b) = −2λV1(x) (2.15)
where the inequality is obtained by using (H1’).
Finally, combining equation (2.15) with (H2’) allows to obtain the de-
sired result. �
We can now state the generalized stochastic contraction theorem
Theorem 3 (Generalized stochastic contraction) Assume that system
(2.5) verifies (H1’) and (H2’). Let a(t) and b(t) be two trajectories whose
initial conditions are given by a probability distribution p(x(0)) = p(a(0),b(0)).
∀t ≥ 0 E
(a(t)− b(t))TM(t)(a(t) − b(t))
+ e−2λt
(a0 − b0)TM(0)(a0 − b0)−
dp(a0,b0) (2.16)
In particular,
∀t ≥ 0 E
‖a(t)− b(t)‖2
‖a(0) − b(0)‖2
e−2λt
(2.17)
Proof Following the same reasoning as in the proof of theorem 2, one
obtains
∀t ≥ 0 Ex0V1(x(t)) ≤
V1(x0)−
e−2λt
which leads to (2.16) by integrating with respect to (a0,b0). Next, observing
‖a(t)− b(t)‖2 ≤ 1
(a(t)− b(t))TM(t)(a(t) − b(t)) = 1
EV1(x(t))
and using the same bounding as in (2.12) lead to (2.17). �
2.3 Strength of the stochastic contraction theorem
2.3.1 “Optimality” of the mean square bound
Consider the following linear dynamical system, known as the Ornstein-
Uhlenbeck (colored noise) process
da = −λadt+ σdW (2.18)
Clearly, the noise-free system is contracting with rate λ and the trace of
the noise matrix is upper-bounded by σ2. Let a(t) and b(t) be two system
trajectories starting respectively at a0 and b0 (deterministic initial condi-
tions). Then by theorem 2, we have
∀t ≥ 0 E
(a(t)− b(t))2
(a0 − b0)2 −
e−2λt (2.19)
Let us verify this result by solving directly equation (2.18). The solution
of equation (2.18) is ([3], p. 134)
a(t) = a0e
−λt + σ
eλ(s−t)dW (s) (2.20)
Next, let us compute the mean square distance between the two trajec-
tories a(t) and b(t)
E((a(t)− b(t))2) = (a0 − b0)2e−2λt +
((∫ t
eλ(s−t)dW1(s)
((∫ t
eλ(u−t)dW2(u)
= (a0 − b0)2e−2λt +
(1− e−2λt)
(a0 − b0)2 −
e−2λt
The last inequality is in fact an equality when (a0 − b0)2 ≥ σ
. Thus,
this calculation shows that the upper-bound (2.19) given by theorem 2 is
optimal, in the sense that it can be attained.
2.3.2 No asymptotic almost-sure stability
From the explicit form (2.20) of the solutions, one can deduce that the
distributions of a(t) and b(t) converge to the normal distribution N
([3], p. 135). Since a(t) and b(t) are independent, the distribution of the
difference a(t)−b(t) will then converge to N
. This observation shows
that, contrary to the case of standard stochastic stability (cf. section 2.1.2),
one cannot – in general – obtain asymptotic almost-sure incremental stability
results (which would imply that the distribution of the difference converges
instead to the constant 0).
Compare indeed equations (2.2) (the condition for standard stability, sec-
tion 2.1.2) and (2.7) (the condition for incremental stability, section 2.2.2).
The difference lies in the term 2C, which stems from the fact that the influ-
ence of the noise does not vanish when two trajectories get very close to each
other (cf. section 2.2.1). The presence of this extra term prevents ÃV (x(t))
from being always non-positive, and as a result, it prevents V (x(t)) from be-
ing always “non-increasing”. As a consequence, V (x(t)) is not – in general –
a supermartingale, and one cannot then use the supermartingale inequality
to obtain asymptotic almost-sure bounds, as in equation (2.3).
Remark If one is interested in finite time bounds then the supermartin-
gale inequality is still applicable, see ([22], p. 86) for details.
2.4 Noisy and noise-free trajectories
Consider the following augmented system
f(a, t)
f(b, t)
0 σ(b, t)
f(x, t)dt +
σ(x, t)dW2d
(2.21)
This equation is the same as equation (2.6) except that the a-system is
not perturbed by noise. Thus V (x) = ‖a − b‖2 will represent the distance
between a noise-free trajectory and a noisy one. All the calculations will be
the same as in the previous development, with C being replaced by C/2.
One can then derive the following corollary
Corollary 1 Assume that system (2.5) verifies (H1’) and (H2’). Let a(t)
be a noise-free trajectory starting at a0 and b(t) a noisy trajectory whose
initial condition is given by a probability distribution p(b(0)). Then
∀t ≥ 0 E
‖a(t)− b(t)‖2
‖a0 − b(0)‖2
e−2λt
(2.22)
Remarks
• One can note here that the derivation of corollary 1 is only permitted
by our initial choice of considering distinct driving Wiener process for
the a- and b-systems (cf. section 2.2.1).
• Corollary 1 provides a robustness result for contracting systems, in the
sense that any contracting system is automatically protected against
noise, as quantified by (2.22). This robustness could be related to the
exponential nature of contraction stability.
3 Combinations of contracting stochastic systems
Stochastic contraction inherits naturally from deterministic contraction [24]
its convenient combination properties. Because contraction is a state-space
concept, such properties can be expressed in more general forms than input-
output analogues such as passivity-based combinations [29]. The following
combination properties allow one to build by recursion stochastically con-
tracting systems of arbitrary size.
Parallel combination Consider two stochastic systems of the same di-
mension {
dx1 = f1(x1, t)dt+ σ1(x1, t)dW1
dx2 = f2(x2, t)dt+ σ2(x2, t)dW2
Assume that both systems are stochastically contracting in the same
constant metric M, with rates λ1 and λ2 and with bounds C1 and C2.
Consider a uniformly positive bounded superposition
α1(t)x1 + α2(t)x2
where ∀t ≥ 0, li ≤ αi(t) ≤ mi for some li,mi > 0, i = 1, 2.
Clearly, this superposition is stochastically contracting in the metric M,
with rate l1λ1 + l2λ2 and bound m1C1 +m2C2.
Negative feedback combination In this and the following paragraphs,
we describe combinations properties for contracting systems in constant met-
rics M. The case of time-varying metrics can be easily adapted from this
development but is skipped here for the sake of clarity.
Consider two coupled stochastic systems
dx1 = f1(x1,x2, t)dt+ σ1(x1, t)dW1
dx2 = f2(x1,x2, t)dt+ σ2(x2, t)dW2
Assume that system i (i = 1, 2) is stochastically contracting with respect to
Mi = Θ
i Θi, with rate λi and bound Ci.
Assume furthermore that the two systems are connected by negative
feedback [33]. More precisely, the Jacobian matrices of the couplings are of
the form Θ1J12Θ
2 = −kΘ2JT21Θ
1 , with k a positive constant. Hence,
the Jacobian matrix of the augmented system is given by
J1 −kΘ−11 Θ2JT21Θ
J21 J2
Consider a coordinate transform Θ =
associated with
the metric M = ΘTΘ > 0. After some calculations, one has
ΘJΘ−1
Θ1J1Θ
Θ2J2Θ
≤ max(−λ1,−λ2)I uniformly (3.1)
The augmented system is thus stochastically contracting in the metric
M, with rate min(λ1, λ2) and bound C1 + kC2.
Hierarchical combination We first recall a standard result in matrix
analysis [19]. Let A be a symmetric matrix in the formA =
A21 A2
Assume that A1 and A2 are definite positive. Then A is definite positive if
sing2(A21) < λmin(A1)λmin(A2) where sing(A21) denotes the largest singu-
lar value of A21. In this case, the smallest eigenvalue of A satisfies
λmin(A) ≥
λmin(A1) + λmin(A2)
λmin(A1)− λmin(A2)
+ sing2(A21)
Consider now the same set-up as in the previous paragraph, except that
the connection is now hierarchical and upper-bounded. More precisely, the
Jacobians of the couplings verify J12 = 0 and sing
2(Θ2J21Θ
1 ) ≤ K. The
Jacobian matrix of the augmented system is then given by J =
J21 J2
Consider a coordinate transform Θǫ =
0 ǫΘ2
associated with the
metric Mǫ = Θ
ǫ Θǫ > 0. After some calculations, one has
ΘJΘ−1
Θ1J1Θ
ǫ(Θ2J21Θ
ǫΘ2J21Θ
Θ2J2Θ
Set now ǫ =
2λ1λ2
. The augmented system is then stochastically con-
tracting in the metric Mǫ, with rate
(λ1 + λ2 −
λ21 + λ
2)) and bound
2C2λ1λ2
Small gains In this paragraph, we require no specific assumption on the
form of the couplings. Consider the coordinate transformΘ =
associated with the metric Mk = Θ
Θk > 0. Aftersome calculations, one
Θ1J1Θ
Θ2J2Θ
where Bk =
kΘ2J21Θ
Θ1J12Θ
Following the matrix analysis result stated at the beginning of the pre-
vious paragraph, if infk>0 sing
2(Bk) < λ1λ2 then the augmented system is
stochastically contracting in the metric Mk, with bound C1 + kC2 and rate
λ verifying
λ ≥ λ1 + λ2
λ1 − λ2
+ inf
sing2(Bk) (3.2)
4 Some examples
4.1 Effect of measurement noise on contracting observers
Consider a nonlinear dynamical system
ẋ = f(x, t) (4.1)
If a measurement y = y(x) is available, then it may be possible to choose
an output injection matrix K(t) such that the dynamics
˙̂x = f(x̂, t) +K(t)(ŷ − y) (4.2)
is contracting, with ŷ = y(x̂). Since the actual state x is a particular
solution of (4.2), any solution x̂ of (4.2) will then converge towards x expo-
nentially.
Assume now that the measurements are corrupted by additive “white
noise”. In the case of linear measurement, the measurement equation be-
comes y = H(t)x+Σ(t)ξ(t) where ξ(t) is a multidimensional “white noise”
and Σ(t) is the matrix of measurement noise intensities.
The observer equation is now given by the following Itô stochastic dif-
ferential equation (using the formal rule dW = ξdt)
dx̂ = (f(x̂, t) +K(t)(H(t)x −H(t)x̂))dt+K(t)Σ(t)dW (4.3)
Next, remark that the solution x of system (4.1) is a also a solution of
the noise-free version of system (4.3). By corollary 1, one then has, for any
solution x̂ of system (4.3)
∀t ≥ 0 E
‖x̂(t)− x(t)‖2
+ ‖x̂0 − x0‖2e−2λt (4.4)
where
λ = inf
∣∣∣∣λmax
∂f(x, t)
−K(t)H(t)
)∣∣∣∣
C = sup
Σ(t)TK(t)TK(t)Σ(t)
Remark The choice of the injection gain K(t) is governed by a trade-
off between convergence speed (λ) and noise sensitivity (C/λ) as quantified
by (4.4). More generally, the explicit computation of the bound on the
expected quadratic estimation error given by (4.4) may open the possibility
of measurement selection in a way similar to the linear case. If several
possible measurements or sets of measurements can be performed, one may
try at each instant (or at each step, in a discrete version) to select the
most relevant, i.e., the measurement or set of measurements which will best
contribute to improving the state estimate. Similarly to the Kalman filters
used in [9] for linear systems, this can be achieved by computing, along
with the state estimate itself, the corresponding bounds on the expected
quadratic estimation error, and then selecting accordingly the measurement
which will minimize it.
4.2 Estimation of velocity using composite variables
In this section, we present a very simple example that hopefully suggests the
many possibilities that could stem from the combination of our stochastic
stability analysis with the composite variables framework [31].
Let x be the position of a mobile subject to a sinusoidal forcing
ẍ = −U1ω2 sin(ωt) + 2U2
where U1 and ω are known parameters. We would like to compute good
approximations of the mobile’s velocity v and acceleration a using only mea-
surements of x and without using any filter. For this, construct the following
observer
−αv 1
−αa 0
(αa − α2v)x
−αaαvx− U1ω3 cos(ωt)
3 cos(ωt)
(4.5)
and introduce the composite variables v̂ = v + αvx and â = a + αax. By
construction, these variables follow the equation
v̂ − v
−U1ω3 cos(ωt)
(4.6)
and therefore, a particular solution of (v̂, â) is clearly (v, a). Choose now
αa = α
v = α
2 and let Mα =
α2 −α/2
−α/2 1
. One can then show that
system (4.6) is contracting with rate λα = α/2 in the metric Mα. Thus, by
the basic contraction theorem [24], (v̂, â) converges exponentially to (v, a)
with rate λα in the metric Mα. Also note that the β-bound corresponding
to the metric Mα is given by βα =
1+α2−
α4−α2+1
Next, assume that the measurements of x are corrupted by additive
“white noise”, so that xmeasured = x+ σξ. Equation (4.5) then becomes an
Itô stochastic differential equation
3 cos(ωt)
By definition of B, the variance of the noise in the metric Mα is upper-
bounded by α
. Thus, using again corollary 1, one obtains (see Figure 1
for a numerical simulation)
∀t ≥ 0 E
‖v̂(t)− v(t)‖2 + ‖â(t)− a(t)‖2
‖v̂0 − v0‖2 + ‖â0 − a0‖2
4.3 Stochastic synchronization
Consider a network of n dynamical elements coupled through diffusive con-
nections
dxi =
f(xi, t) +
j 6=i
Kij(xj − xi)
 dt+ σi(xi, t)dW di i = 1, . . . , n
(4.7)
x, t) =
f(x1, t)
f(xn, t)
 , ⌢σ(⌢x, t) =
σ1(x1, t) 0 0
. . . 0
0 0 σn(xn, t)
The global state
x then follows the equation
x, t)− L⌢x
x, t)dW nd (4.8)
 0  2  4  6  8  10
 0  2  4  6  8  10
Figure 1: Estimation of the velocity of a mobile using noisy measurements
of its position. The simulation was performed using the Euler-Maruyama
algorithm [18] with the following parameters: U1 = 10, U2 = 2, ω = 3,
σ = 10 and α = 1. Left plot: simulation for one trial. The plot shows
the measured position (red), the actual velocity (blue) and the estimate of
the velocity using the measured position (green). Right plot: the average
over 1000 trials of the squared error ‖v̂ − v‖2 + ‖â − a‖2 (green) and the
asymptotic bound
= 200
given by our approach (red).
In the sequel, we follow the reasoning of [28], which starts by defining an
appropriate orthonormal matrix V describing the synchronization subspace
(V represents the state projection on the subspace M⊥, orthogonal to the
synchronization subspace M = {(x1, . . . ,xn)T : x1 = . . . = xn}, see [28] for
details). Denote by
y the state of the projected system,
y = V
x. Since the
mapping is linear, Itô differentiation rule simply yields
y = Vd
x, t) −VL⌢x
x, t)dW nd
y, t) −VLVT⌢y
y, t)dW nd (4.9)
Assume now that ∂f
is uniformly upper-bounded. Then for strong
enough coupling strength, A = V ∂f
VT − VLVT will be uniformly neg-
ative definite. Let λ = |λmax(A)| > 0. System (4.9) then verifies condition
(H1) with rate λ. Assume furthermore that each noise intensity σi is upper-
bounded by a constant Ci (i.e. supx,t tr(σi(x, t)
Tσi(x, t)) ≤ Ci). Condition
(H2) will then be satisfied with the bound C =
i Ci.
Next, consider a noise-free trajectory
yu(t) of system (4.9). By theo-
rem 3 of [28], we know that
yu(t) converges exponentially to zero. Thus,
by corollary 1, one can conclude that, after exponential transients of rate λ,
E(‖⌢y(t)‖2) ≤ C
On the other hand, one can show that
‖⌢y(t)‖2 = 1
‖xi − xj‖2
Thus, after exponential transients of rate λ, we have
‖xi − xj‖2 ≤
Remarks
• The above development is fully compatible with the concurrent syn-
chronization framework [28]. It can also be easily generalized to the
case of time-varying metrics by combining theorem 3 of this paper and
corollary 1 of [28].
• The synchronization of Itô dynamical systems has been investigated
in [4]. However, the systems considered by the authors of that article
were dissipative. Here, we make a less restrictive assumption, namely,
we only require ∂f
to be uniformly upper-bounded. This enables us
to study the synchronization of a broader class of dynamical systems,
which can include nonlinear oscillators or even chaotic systems.
Example As illustration of the above development, we provide here a de-
tailed analysis for the synchronization of noisy FitzHugh-Nagumo oscillators
(see [35] for the references). The dynamics of two diffusively-coupled noisy
FitzHugh-Nagumo oscillators can be described by
dvi = (c(vi + wi − 13v
i + Ii) + k(v0 − vi))dt+ σdWi
dwi = −1c (vi − a+ bwi)dt
i = 1, 2
Let x = (v1, w1, v2, w2)
T and V = 1√
1 0 −1 0
0 1 0 −1
. The Jaco-
bian matrix of the projected noise-free system is then given by
c− c(v
− k c
−1/c −b/c
Thus, if the coupling strength verifies k > c then the projected system
will be stochastically contracting in the diagonal metric M = diag(1, c) with
rate min(k − c, b/c) and bound σ2. Hence, the average absolute difference
between the two “membrane potentials” |v1 − v2| will be upper-bounded by
min(1, c)min(k − c, b/c) (see Figure 2 for a numerical simulation).
Acknowledgments We are grateful to Dr S. Darses, Prof D. Bennequin
and Prof M. Yor for stimulating discussions, and to the Associate Editor
and the reviewers for their helpful comments.
 0  2  4  6  8  10
Figure 2: Synchronization of two noisy FitzHugh-Nagumo oscillators. The
simulation was performed using the Euler-Maruyama algorithm [18] with
the following parameters: a = 0.3, b = 0.2, c = 30, k = 40 and σ = 1. The
plot shows the “membrane potentials” of the two oscillators.
A Appendix
A.1 Proof of the supermartingale property
Lemma 3 Consider a Markov stochastic process x(t) and a non-negative
function V such that ∀t ≥ 0 EV (x(t)) < ∞ and
∀x ∈ Rn ÃV (x) ≤ −λV (x) (A.1)
where λ is a non-negative real number and Ã is the infinitesimal operator
of the process x(t). Then V (x(t)) is a supermartingale with respect to the
canonical filtration Ft = {x(s), s ≤ t}.
We need to show that for all s ≥ t, one has E(V (x(s))|Ft) ≤ V (x(t)).
Since x(t) is a Markov process, it suffices to show that
∀x0 ∈ Rn E(V (x(t))|x(0) = x0) ≤ V (x0)
By Dynkin’s formula, one has for all x0 ∈ Rn
Ex0V (x(t)) = V (x0) + Ex0
ÃV (x(s))ds
≤ V (x0)− λEx0
V (x(s))ds ≤ V (x0)
where Ex0(·) = E(·|x(0) = x0).
A.2 A variation of Gronwall’s lemma
Lemma 4 Let g : [0,∞[→ R be a continuous function, C a real number
and λ a strictly positive real number. Assume that
∀u, t 0 ≤ u ≤ t g(t)− g(u) ≤
−λg(s) + Cds (A.2)
∀t ≥ 0 g(t) ≤ C
g(0) − C
e−λt (A.3)
where [·]+ = max(0, ·).
Proof Case 1 : C = 0, g(0) > 0.
Define h(t) by
∀t ≥ 0 h(t) = g(0)e−λt
Remark that h is positive with h(0) = g(0), and satisfies (A.2) where the
inequality has been replaced by an equality
∀u, t 0 ≤ u ≤ t h(t)− h(u) = −
λh(s)ds
Consider now the set S = {t ≥ 0 | g(t) > h(t)}. If S = ∅ then the
lemma holds true. Assume by contradiction that S 6= ∅. In this case, let
m = inf S < ∞. By continuity of g and h and by the fact that g(0) = h(0),
one has g(m) = h(m) and there exists ǫ > 0 such that
∀t ∈]m,m+ ǫ[ g(t) > h(t) (A.4)
Consider now φ(t) = g(m)− λ
g(s)ds. Equation (A.2) implies that
∀t ≥ m g(t) ≤ φ(t)
In order to compare φ(t) and h(t) for t ∈]m,m+ ǫ[, let us differentiate the
ratio φ(t)/h(t).
φ′h− h′φ
−λgh+ λhφ
λh(φ − g)
Thus φ(t)/h(t) is increasing for t ∈]m,m + ǫ[. Since φ(m)/h(m) = 1, one
can conclude that
∀t ∈]m,m+ ǫ[ φ(t) ≥ h(t)
which implies, by definition of φ and h, that
∀t ∈]m,m+ ǫ[
g(s)ds ≤
h(s)ds (A.5)
Choose now t0 such that m < t0 < m+ ǫ, then one has by (A.4)
g(s)ds >
h(s)ds
which clearly contradicts (A.5).
Case 2 : C = 0, g(0) ≤ 0
Consider the set S = {t ≥ 0 | g(t) > 0}. If S = ∅ then the lemma holds
true. Assume by contradiction that S 6= ∅. In this case, let m = inf S < ∞.
By continuity of g and by the fact that g(0) ≤ 0, one has g(m) = 0 and
there exists ǫ such that
∀t ∈]m,m+ ǫ[ g(t) > 0 (A.6)
Let t0 ∈]m,m+ ǫ[. Equation (A.2) implies that
g(t0) ≤ −λ
g(s)ds
which clearly contradicts (A.6).
Case 3 : C 6= 0
Define ĝ = g − C/λ. One has
∀u, t 0 ≤ u ≤ t ĝ(t)−ĝ(u) = g(t)−g(u) ≤
−λg(s)+Cds = −
λĝ(s)ds
Thus ĝ satisfies the conditions of Case 1 or Case 2, and as a consequence
∀t ≥ 0 ĝ(t) ≤ [ĝ(0)]+e−λt
The conclusion of the lemma follows by replacing ĝ by g−C/λ in the above
equation. �
References
[1] N. Aghannan and P. Rouchon. An intrinsic observer for a class of la-
grangian systems. IEEE Transactions on Automatic Control, 48, 2003.
[2] D. Angeli. A lyapunov approach to incremental stability properties.
IEEE Transactions on Automatic Control, 47:410–422, 2002.
[3] L. Arnold. Stochastic Differential Equations : Theory and Applications.
Wiley, 1974.
[4] T. Caraballo and P. Kloeden. The persistence of synchronization under
environmental noise. Proceedings of the Royal Society A, 461:2257–
2267, 2005.
[5] T. Caraballo, P. Kloeden, and B. Schmalfuss. Exponentially stable
stationary solutions for stochastic evolution equations and their pertu-
bation. Applied Mathematics and Optimization, 50:183–207, 2004.
[6] L. d’Alto and M. Corless. Incremental quadratic stability. In Proceed-
ings of the IEEE Conference on Decision and Control, 2005.
[7] B. Demidovich. Dissipativity of a nonlinear system of differential equa-
tions. Ser. Mat. Mekh., 1961.
[8] H. Deng, M. Krstic, and R. Williams. Stabilization of stochastic nonlin-
ear systems driven by noise of unknown covariance. IEEE Transactions
on Automatic Control, 46, 2001.
[9] E. Dickmanns. Dynamic Vision for Intelligent Vehicles. Course Notes,
MIT EECS dept, 1998.
[10] K. El Rifai and J.-J. Slotine. Contraction and incremental stability.
Technical report, MIT NSL Report, 2006.
[11] W. Feller. An Introduction to Probability Theory and Its Applications.
Wiley, 1968.
[12] P. Florchinger. Lyapunov-like techniques for stochastic stability. SIAM
Journal of Control and Optimization, 33:1151–1169, 1995.
[13] P. Florchinger. Feedback stabilization of affine in the control stochastic
differential systems by the control lyapunov function method. SIAM
Journal of Control and Optimization, 35, 1997.
[14] V. Fromion. Some results on the behavior of lipschitz continuous sys-
tems. In Proceedings of the European Control Conference, 1997.
[15] B. Girard, N. Tabareau, Q.-C. Pham, A. Berthoz, and J.-J. Slotine.
Where neuroscience and dynamic system theory meet autonomous
robotics: a contracting basal ganglia model for action selection. Neural
Networks, 2008.
[16] P. Hartmann. Ordinary differential equations. Wiley, 1964.
[17] R. Has’minskii. Stochastic Stability of Differential Equations. Sijthoff
and Nordhoff, Rockville, 1980.
[18] D. Higham. An algorithmic introduction to numerical simulation of
stochastic differential equations. SIAM Review, 43:525–546, 2001.
[19] R. Horn and C. Johnson. Matrix Analysis. Cambridge University Press,
1985.
[20] J. Jouffroy. Some ancestors of contraction analysis. In Proceedings of
the IEEE Conference on Decision and Control, 2005.
[21] J. Jouffroy and T. Fossen. On the combination of nonlinear contracting
observers and uges controllers for output feedback. In Proceedings of
the IEEE Conference on Decision and Control, 2004.
[22] H. Kushner. Stochastic Stability and Control. Academic Press, 1967.
[23] D. Lewis. Metric properties of differential equations. American Journal
of Mathematics, 71:294–312, 1949.
[24] W. Lohmiller and J.-J. Slotine. On contraction analysis for nonlinear
systems. Automatica, 34:671–682, 1998.
[25] W. Lohmiller and J.-J. Slotine. Nonlinear process control using con-
traction theory. A.I.Ch.E. Journal, 2000.
[26] W. Lohmiller and J.-J. Slotine. Contraction analysis of nonlinear dis-
tributed systems. International Journal of Control, 78, 2005.
[27] X. Mao. Stability of Stochastic Differential Equations with Respect to
Semimartingales. Longman, White Plains, NY, 1991.
[28] Q.-C. Pham and J.-J. Slotine. Stable concurrent synchronization in
dynamic system networks. Neural Netw, 20(1):62–77, Jan. 2007.
[29] V. Popov. Hyperstability of Control Systems. Springer-Verlag, 1973.
[30] W. Rudin. Real and complex analysis. McGraw-Hill, 1987.
[31] J.-J. Slotine and W. Li. Applied Nonlinear Control. Prentice-Hall, 1991.
[32] E. Sontag and Y. Wang. Output-to-state stability and detectability of
nonlinear systems. Systems and Control Letters, 29:279–290, 1997.
[33] N. Tabareau and J.-J. Slotine. Notes on contraction theory. Technical
report, MIT NSL Report, 2005.
[34] J. Tsinias. The concept of “exponential iss” for stochastic systems and
applications to feedback stabilization. Systems and Control Letters,
36:221–229, 1999.
[35] W. Wang and J.-J. E. Slotine. On partial contraction analysis for cou-
pled nonlinear oscillators. Biol Cybern, 92(1):38–53, Jan. 2005.
[36] Y. Zhao and J.-J. Slotine. Discrete nonlinear observers for inertial
navigation. Systems and Control Letters, 54, 2005.
	Introduction
	Main results
	Background
	Nonlinear contraction theory
	Standard stochastic stability
	The stochastic contraction theorem
	Settings
	The basic stochastic contraction theorem
	Generalization to time-varying metrics
	Strength of the stochastic contraction theorem
	``Optimality'' of the mean square bound
	No asymptotic almost-sure stability
	Noisy and noise-free trajectories
	Combinations of contracting stochastic systems
	Some examples
	Effect of measurement noise on contracting observers
	Estimation of velocity using composite variables
	Stochastic synchronization
	Appendix
	Proof of the supermartingale property
	A variation of Gronwall's lemma

  We investigate the incremental stability properties of It\^o stochastic
dynamical systems. Specifically, we derive a stochastic version of nonlinear
contraction theory that provides a bound on the mean square distance between
any two trajectories of a stochastically contracting system. This bound can be
expressed as a function of the noise intensity and the contraction rate of the
noise-free system. We illustrate these results in the contexts of stochastic
nonlinear observers design and stochastic synchronization.

Introduction
Nonlinear stability properties are often considered with respect to an equi-
librium point or to a nominal system trajectory (see e.g. [31]). By contrast,
incremental stability is concerned with the behaviour of system trajectories
with respect to each other. From the triangle inequality, global exponential
incremental stability (any two trajectories tend to each other exponentially)
is a stronger property than global exponential convergence to a single tra-
jectory.
Historically, work on deterministic incremental stability can be traced
back to the 1950’s [23, 7, 16] (see e.g. [26, 20] for a more extensive list and
historical discussion of related references). More recently, and largely inde-
pendently of these earlier studies, a number of works have put incremental
To whom correspondance should be addressed.
http://arxiv.org/abs/0704.0926v2
stability on a broader theoretical basis and made relations with more tradi-
tional stability approaches [14, 32, 24, 2, 6]. Furthermore, it was shown that
incremental stability is especially relevant in the study of such problems as
state detection [2], observer design or synchronization analysis.
While the above references are mostly concerned with deterministic sta-
bility notions, stability theory has also been extended to stochastic dynam-
ical systems, see for instance [22, 17]. This includes important recent de-
velopments in Lyapunov-like approaches [12, 27], as well as applications to
standard problems in systems and control [13, 34, 8]. However, stochastic
versions of incremental stability have not yet been systematically investi-
gated.
The goal of this paper is to extend some concepts and results in in-
cremental stability to stochastic dynamical systems. More specifically, we
derive a stochastic version of contraction analysis in the specialized context
of state-independent metrics.
We prove in section 2 that the mean square distance between any two
trajectories of a stochastically contracting system is upper-bounded by a
constant after exponential transients. In contrast with previous works on
incremental stochastic stability [5], we consider the case when the two tra-
jectories are subject to distinct and independent noises, as detailed in sec-
tion 2.2.1. This specificity enables our theory to have a number of new and
practically important applications. However, the fact that the noise does
not vanish as two trajectories get very close to each other will prevent us
from obtaining asymptotic almost-sure stability results (see section 2.3.2).
In section 3, we show that results on combinations of deterministic con-
tracting systems have simple analogues in the stochastic case. These combi-
nation properties allow one to build by recursion stochastically contracting
systems of arbitrary size.
Finally, as illustrations of our results, we study in section 4 several ex-
amples, including contracting observers with noisy measurements, stochas-
tic composite variables and synchronization phenomena in networks of noisy
dynamical systems.
2 Main results
2.1 Background
2.1.1 Nonlinear contraction theory
Contraction theory [24] provides a set of tools to analyze the incremental
exponential stability of nonlinear systems, and has been applied notably
to observer design [24, 25, 1, 21, 36], synchronization analysis [35, 28] and
systems neuroscience modelling [15]. Nonlinear contracting systems enjoy
desirable aggregation properties, in that contraction is preserved under many
types of system combinations given suitable simple conditions [24].
While we shall derive global properties of nonlinear systems, many of our
results can be expressed in terms of eigenvalues of symmetric matrices [19].
Given a square matrix A, the symmetric part of A is denoted by As. The
smallest and largest eigenvalues of As are denoted by λmin(A) and λmax(A).
Given these notations, the matrix A is positive definite (denoted A > 0) if
λmin(A) > 0, and it is uniformly positive definite if
∃β > 0 ∀x, t λmin(A(x, t)) ≥ β
The basic theorem of contraction analysis, derived in [24], can be stated
as follows
Theorem 1 (Contraction) Consider, in Rn, the deterministic system
ẋ = f(x, t) (2.1)
where f is a smooth nonlinear function. Denote the Jacobian matrix of f
with respect to its first variable by ∂f
. If there exists a square matrix Θ(x, t)
such that M(x, t) = Θ(x, t)TΘ(x, t) is uniformly positive definite and the
matrix
F(x, t) =
Θ(x, t) +Θ(x, t)
Θ−1(x, t)
is uniformly negative definite, then all system trajectories converge exponen-
tially to a single trajectory, with convergence rate | sup
x,t λmax(F)| = λ > 0.
The system is said to be contracting, F is called its generalized Jacobian,
M(x, t) its contraction metric and λ its contraction rate.
2.1.2 Standard stochastic stability
In this section, we present very informally the basic ideas of standard stochas-
tic stability (for a rigourous treatment, the reader is referred to e.g. [22]).
This will set the context to understand the forthcoming difficulties and dif-
ferences associated with incremental stochastic stability.
For simplicity, we consider the special case of global exponential stability.
Let x(t) be a Markov stochastic process and assume that there exists a non-
negative function V (V (x) may represent e.g. the squared distance of x
from the origin) such that
∀x ∈ Rn ÃV (x) ≤ −λV (x) (2.2)
where λ is a positive real number and Ã is the infinitesimal operator of the
process x(t). The operator Ã is the stochastic analogue of the deterministic
differentiation operator. In the case that x(t) is an Itô process, Ã corre-
sponds to the widely-used [27, 34, 8] differential generator L (for a proof of
this fact, see [22], p. 15 or [3], p. 42).
For x0 ∈ Rn, let Ex0(·) = E(·|x(0) = x0). Then by Dynkin’s formula
([22], p. 10), one has
∀t ≥ 0 Ex0V (x(t))− V (x0) = Ex0
ÃV (x(s))ds
≤ −λEx0
V (x(s))ds = −λ
Ex0V (x(s))ds
Applying the Gronwall’s lemma to the deterministic real-valued function
t → Ex0V (x(t)) yields
∀t ≥ 0 Ex0V (x(t)) ≤ V (x0)e−λt
If we assume furthermore that Ex0V (x(t)) < ∞ for all t, then the above
implies that V (x(t)) is a supermartingale (see lemma 3 in the Appendix for
details), which yields, by the supermartingale inequality
T≤t<∞
V (x(t)) ≥ A
≤ Ex0V (x(T ))
≤ V (x0)e
(2.3)
Thus, one obtains an almost-sure stability result, in the sense that
∀A > 0 lim
T≤t<∞
V (x(t)) ≥ A
= 0 (2.4)
2.2 The stochastic contraction theorem
2.2.1 Settings
Consider a noisy system described by an Itô stochastic differential equation
da = f(a, t)dt+ σ(a, t)dW d (2.5)
where f is a Rn × R+ → Rn function, σ is a Rn × R+ → Rnd matrix-valued
function and W d is a standard d-dimensional Wiener process.
To ensure existence and uniqueness of solutions to equation (2.5), we
assume, here and in the remainder of the paper, the following standard
conditions on f and σ
Lipschitz condition: There exists a constant K1 > 0 such that
∀t ≥ 0, a,b ∈ Rn ‖f(a, t)− f(b, t)|+ ‖σ(a, t) − σ(b, t)‖ ≤ K1‖a− b‖
Restriction on growth: There exists a constant K2 > 0
∀t ≥ 0, a ∈ Rn ‖f(a, t)‖2 + ‖σ(a, t)‖2 ≤ K2(1 + ‖a‖2)
Under these conditions, one can show ([3], p. 105) that equation (2.5)
has on [0,∞[ a unique Rn-valued solution a(t), continuous with probability
one, and satisfying the initial condition a(0) = a0, with a0 ∈ Rn.
In order to investigate the incremental stability properties of system (2.5),
consider now two system trajectories a(t) and b(t). Our goal will consist of
studying the trajectories a(t) and b(t) with respect to each other. For this,
we consider the augmented system x(t) = (a(t),b(t))T , which follows the
equation
f(a, t)
f(b, t)
σ(a, t) 0
0 σ(b, t)
dW d1
dW d2
f(x, t)dt+
σ(x, t)dW 2d (2.6)
Important remark As stated in the introduction, the systems a and
b are driven by distinct and independent Wiener processes W d1 and W
This makes our approach considerably different from [5], where the authors
studied two trajectories driven by the same Wiener process.
Our approach enables us to study the stability of the system with respect
to variations in initial conditions and to random perturbations: indeed, two
trajectories of any real-life system are typically affected by distinct “real-
izations” of the noise. In addition, it leads very naturally to nice results on
the comparison of noisy and noise-free trajectories (cf. section 2.4), which
are particularly useful in applications (cf. section 4).
However, because of the very fact that the two trajectories are driven
by distinct Wiener processes, we cannot expect the influence of the noise
to vanish when the two trajectories get very close to each other. This con-
strasts with [5], and more generally, with the standard stochastic stability
case, where the noise vanishes near the origin (cf. section 2.1.2). The con-
sequences of this will be discussed in detail in section 2.3.2.
2.2.2 The basic stochastic contraction theorem
We introduce two hypotheses
(H1) f(a, t) is contracting in the identity metric, with contraction rate λ,
(i.e. ∀a, t λmax
≤ −λ)
(H2) tr
σ(a, t)Tσ(a, t)
is uniformly upper-bounded by a constant C (i.e.
∀a, t tr
σ(a, t)Tσ(a, t)
In other words, (H1) says that the noise-free system is contracting, while
(H2) says that the variance of the noise is upper-bounded by a constant.
Definition 1 A system that verifies (H1) and (H2) is said to be stochas-
tically contracting in the identity metric, with rate λ and bound C.
Consider now the Lyapunov-like function V (x) = ‖a−b‖2 = (a−b)T (a−
b). Using (H1) and (H2), we derive below an inequality on ÃV (x), similar
to equation (2.2) in section 2.1.2.
Lemma 1 Under (H1) and (H2), one has the inequality
ÃV (x) ≤ −2λV (x) + 2C (2.7)
Proof Since x(t) is an Itô process, Ã is given by the differential operator
L of the process [22, 3]. Thus, by the Itô formula
ÃV (x) = L V (x) =
∂V (x)
f(x, t) +
σ(x, t)T
∂2V (x)
σ(x, t)
1≤i≤2n
f(x, t)i +
1≤i,j,k≤2n
σ(x, t)ij
∂xi∂xk
σ(x, t)kj
1≤i≤n
f(a, t)i +
1≤i≤n
f(b, t)i
1≤i,j,k≤n
σ(a, t)ij
∂ai∂ak
σ(a, t)kj
1≤i,j,k≤n
σ(b, t)ij
∂bi∂bk
σ(b, t)kj
= 2(a− b)T (f(a, t) − f(b, t))
+tr(σ(a, t)Tσ(a, t)) + tr(σ(b, t)Tσ(b, t))
Fix t ≥ 0 and, as in [10], consider the real-valued function
r(µ) = (a− b)T (f(µa+ (1− µ)b, t)− f(b, t))
Since f is C1, r is C1 over [0, 1]. By the mean value theorem, there exists
µ0 ∈]0, 1[ such that
r′(µ0) = r(1)− r(0) = (a− b)T (f(a)− f(b))
On the other hand, one obtains by differentiating r
r′(µ0) = (a− b)T
(µ0a+ (1− µ0)b, t)
(a− b)
Thus, one has
(a− b)T (f(a)− f(b)) = (a− b)T
(µ0a+ (1− µ0)b, t)
(a− b)
≤ −λ(a− b)T (a− b) = −2λV (x) (2.8)
where the inequality is obtained by using (H1).
Finally,
ÃV (x) = 2(a− b)T (f(a)− f(b)) + tr(σ(a, t)Tσ(a, t)) + tr(σ(b, t)Tσ(b, t))
≤ −2λV (x) + 2C
where the inequality is obtained by using (H2). �
We are now in a position to prove our main theorem on stochastic incre-
mental stability.
Theorem 2 (Stochastic contraction) Assume that system (2.5) verifies
(H1) and (H2). Let a(t) and b(t) be two trajectories whose initial condi-
tions are given by a probability distribution p(x(0)) = p(a(0),b(0)). Then
∀t ≥ 0 E
‖a(t)− b(t)‖2
+ e−2λt
‖a0 − b0‖2 −
dp(a0,b0)
(2.9)
where [·]+ = max(0, ·). This implies in particular
∀t ≥ 0 E
‖a(t)− b(t)‖2
‖a(0) − b(0)‖2
e−2λt (2.10)
Proof Let x0 = (a0,b0) ∈ R2n. By Dynkin’s formula ([22], p. 10)
Ex0V (x(t)) − V (x0) = Ex0
ÃV (x(s))ds
Thus one has ∀u, t 0 ≤ u ≤ t < ∞
Ex0V (x(t)) − Ex0V (x(u)) = Ex0
ÃV (x(s))ds
≤ Ex0
(−2λV (x(s)) + 2C)ds (2.11)
(−2λEx0V (x(s)) + 2C)ds
where inequality (2.11) is obtained by using lemma 1.
Denote by g(t) the deterministic quantity Ex0V (x(t)). Clearly, g(t) is a
continuous function of t since x(t) is a continuous process. The function g
then satisfies the conditions of the Gronwall-type lemma 4 in the Appendix,
and as a consequence
∀t ≥ 0 Ex0V (x(t)) ≤
V (x0)−
e−2λT
Integrating the above inequality with respect to x0 yields the desired
result (2.9). Next, inequality (2.10) follows from (2.9) by remarking that
‖a0 − b0‖2 −
dp(a0,b0) ≤
‖a0 − b0‖2dp(a0,b0)
‖a(0) − b(0)‖2
(2.12)
Remark Let ǫ > 0 and Tǫ =
E(‖a0−b0‖2)
. Then inequal-
ity (2.10) and Jensen’s inequality [30] imply
∀t ≥ Tǫ E(‖a(t)− b(t)‖) ≤
C/λ+ ǫ (2.13)
Since ‖a(t)−b(t)‖ is non-negative, (2.13) together with Markov inequal-
ity [11] allow one to obtain the following probabilistic bound on the distance
between a(t) and b(t)
∀A > 0 ∀t ≥ Tǫ P (‖a(t)− b(t)‖ ≥ A) ≤
C/λ+ ǫ
Note however that this bound is much weaker than the asymptotic
almost-sure bound (2.4).
2.2.3 Generalization to time-varying metrics
Theorem 2 can be vastly generalized by considering general time-dependent
metrics (the case of state-dependent metrics is not considered in this article
and will be the subject of a future work). Specifically, let us replace (H1)
and (H2) by the following hypotheses
(H1’) There exists a uniformly positive definite metric M(t) = Θ(t)TΘ(t),
with the lower-bound β > 0 (i.e. ∀x, t xTM(t)x ≥ β‖x‖2) and f(a, t)
is contracting in that metric, with contraction rate λ, i.e.
Θ(t) +Θ(t)
Θ−1(t)
≤ −λ uniformly
or equivalently
M(t) +
M(t) ≤ −2λM(t) uniformly
(H2’) tr
σ(a, t)TM(t)σ(a, t)
is uniformly upper-bounded by a constant C
Definition 2 A system that verifies (H1’) and (H2’) is said to be stochas-
tically contracting in the metric M(t), with rate λ and bound C.
Consider now the generalized Lyapunov-like function V1(x, t) = (a −
b)TM(t)(a − b). Lemma 1 can then be generalized as follows.
Lemma 2 Under (H1’) and (H2’), one has the inequality
ÃV1(x, t) ≤ −2λV1(x, t) + 2C (2.14)
Proof Let us compute first ÃV1
ÃV1(x, t) =
f(x, t) +
σ(x, t)T
σ(x, t)
= (a− b)T
(a− b) + 2(a− b)TM(t)(f(a, t) − f(b, t))
+tr(σ(a, t)TM(t)σ(a, t)) + tr(σ(b, t)TM(t)σ(b, t))
Fix t > 0 and consider the real-valued function
r(µ) = (a− b)TM(t)(f(µa+ (1− µ)b, t)− f(b, t))
Since f is C1, r is C1 over [0, 1]. By the mean value theorem, there exists
µ0 ∈]0, 1[ such that
r′(µ0) = r(1)− r(0) = (a− b)TM(t)(f(a) − f(b))
On the other hand, one obtains by differentiating r
r′(µ0) = (a− b)TM(t)
(µ0a+ (1− µ0)b, t)
(a− b)
Thus, letting c = µ0a+ (1− µ0)b, one has
(a− b)T
(a− b) + 2(a− b)TM(t)(f(a) − f(b))
= (a− b)T
(a− b) + 2(a− b)TM(t)
(c, t)
(a− b)
= (a− b)T
M(t) +M(t)
(c, t)
(c, t)
(a− b)
≤ −2λ(a− b)TM(t)(a− b) = −2λV1(x) (2.15)
where the inequality is obtained by using (H1’).
Finally, combining equation (2.15) with (H2’) allows to obtain the de-
sired result. �
We can now state the generalized stochastic contraction theorem
Theorem 3 (Generalized stochastic contraction) Assume that system
(2.5) verifies (H1’) and (H2’). Let a(t) and b(t) be two trajectories whose
initial conditions are given by a probability distribution p(x(0)) = p(a(0),b(0)).
∀t ≥ 0 E
(a(t)− b(t))TM(t)(a(t) − b(t))
+ e−2λt
(a0 − b0)TM(0)(a0 − b0)−
dp(a0,b0) (2.16)
In particular,
∀t ≥ 0 E
‖a(t)− b(t)‖2
‖a(0) − b(0)‖2
e−2λt
(2.17)
Proof Following the same reasoning as in the proof of theorem 2, one
obtains
∀t ≥ 0 Ex0V1(x(t)) ≤
V1(x0)−
e−2λt
which leads to (2.16) by integrating with respect to (a0,b0). Next, observing
‖a(t)− b(t)‖2 ≤ 1
(a(t)− b(t))TM(t)(a(t) − b(t)) = 1
EV1(x(t))
and using the same bounding as in (2.12) lead to (2.17). �
2.3 Strength of the stochastic contraction theorem
2.3.1 “Optimality” of the mean square bound
Consider the following linear dynamical system, known as the Ornstein-
Uhlenbeck (colored noise) process
da = −λadt+ σdW (2.18)
Clearly, the noise-free system is contracting with rate λ and the trace of
the noise matrix is upper-bounded by σ2. Let a(t) and b(t) be two system
trajectories starting respectively at a0 and b0 (deterministic initial condi-
tions). Then by theorem 2, we have
∀t ≥ 0 E
(a(t)− b(t))2
(a0 − b0)2 −
e−2λt (2.19)
Let us verify this result by solving directly equation (2.18). The solution
of equation (2.18) is ([3], p. 134)
a(t) = a0e
−λt + σ
eλ(s−t)dW (s) (2.20)
Next, let us compute the mean square distance between the two trajec-
tories a(t) and b(t)
E((a(t)− b(t))2) = (a0 − b0)2e−2λt +
((∫ t
eλ(s−t)dW1(s)
((∫ t
eλ(u−t)dW2(u)
= (a0 − b0)2e−2λt +
(1− e−2λt)
(a0 − b0)2 −
e−2λt
The last inequality is in fact an equality when (a0 − b0)2 ≥ σ
. Thus,
this calculation shows that the upper-bound (2.19) given by theorem 2 is
optimal, in the sense that it can be attained.
2.3.2 No asymptotic almost-sure stability
From the explicit form (2.20) of the solutions, one can deduce that the
distributions of a(t) and b(t) converge to the normal distribution N
([3], p. 135). Since a(t) and b(t) are independent, the distribution of the
difference a(t)−b(t) will then converge to N
. This observation shows
that, contrary to the case of standard stochastic stability (cf. section 2.1.2),
one cannot – in general – obtain asymptotic almost-sure incremental stability
results (which would imply that the distribution of the difference converges
instead to the constant 0).
Compare indeed equations (2.2) (the condition for standard stability, sec-
tion 2.1.2) and (2.7) (the condition for incremental stability, section 2.2.2).
The difference lies in the term 2C, which stems from the fact that the influ-
ence of the noise does not vanish when two trajectories get very close to each
other (cf. section 2.2.1). The presence of this extra term prevents ÃV (x(t))
from being always non-positive, and as a result, it prevents V (x(t)) from be-
ing always “non-increasing”. As a consequence, V (x(t)) is not – in general –
a supermartingale, and one cannot then use the supermartingale inequality
to obtain asymptotic almost-sure bounds, as in equation (2.3).
Remark If one is interested in finite time bounds then the supermartin-
gale inequality is still applicable, see ([22], p. 86) for details.
2.4 Noisy and noise-free trajectories
Consider the following augmented system
f(a, t)
f(b, t)
0 σ(b, t)
f(x, t)dt +
σ(x, t)dW2d
(2.21)
This equation is the same as equation (2.6) except that the a-system is
not perturbed by noise. Thus V (x) = ‖a − b‖2 will represent the distance
between a noise-free trajectory and a noisy one. All the calculations will be
the same as in the previous development, with C being replaced by C/2.
One can then derive the following corollary
Corollary 1 Assume that system (2.5) verifies (H1’) and (H2’). Let a(t)
be a noise-free trajectory starting at a0 and b(t) a noisy trajectory whose
initial condition is given by a probability distribution p(b(0)). Then
∀t ≥ 0 E
‖a(t)− b(t)‖2
‖a0 − b(0)‖2
e−2λt
(2.22)
Remarks
• One can note here that the derivation of corollary 1 is only permitted
by our initial choice of considering distinct driving Wiener process for
the a- and b-systems (cf. section 2.2.1).
• Corollary 1 provides a robustness result for contracting systems, in the
sense that any contracting system is automatically protected against
noise, as quantified by (2.22). This robustness could be related to the
exponential nature of contraction stability.
3 Combinations of contracting stochastic systems
Stochastic contraction inherits naturally from deterministic contraction [24]
its convenient combination properties. Because contraction is a state-space
concept, such properties can be expressed in more general forms than input-
output analogues such as passivity-based combinations [29]. The following
combination properties allow one to build by recursion stochastically con-
tracting systems of arbitrary size.
Parallel combination Consider two stochastic systems of the same di-
mension {
dx1 = f1(x1, t)dt+ σ1(x1, t)dW1
dx2 = f2(x2, t)dt+ σ2(x2, t)dW2
Assume that both systems are stochastically contracting in the same
constant metric M, with rates λ1 and λ2 and with bounds C1 and C2.
Consider a uniformly positive bounded superposition
α1(t)x1 + α2(t)x2
where ∀t ≥ 0, li ≤ αi(t) ≤ mi for some li,mi > 0, i = 1, 2.
Clearly, this superposition is stochastically contracting in the metric M,
with rate l1λ1 + l2λ2 and bound m1C1 +m2C2.
Negative feedback combination In this and the following paragraphs,
we describe combinations properties for contracting systems in constant met-
rics M. The case of time-varying metrics can be easily adapted from this
development but is skipped here for the sake of clarity.
Consider two coupled stochastic systems
dx1 = f1(x1,x2, t)dt+ σ1(x1, t)dW1
dx2 = f2(x1,x2, t)dt+ σ2(x2, t)dW2
Assume that system i (i = 1, 2) is stochastically contracting with respect to
Mi = Θ
i Θi, with rate λi and bound Ci.
Assume furthermore that the two systems are connected by negative
feedback [33]. More precisely, the Jacobian matrices of the couplings are of
the form Θ1J12Θ
2 = −kΘ2JT21Θ
1 , with k a positive constant. Hence,
the Jacobian matrix of the augmented system is given by
J1 −kΘ−11 Θ2JT21Θ
J21 J2
Consider a coordinate transform Θ =
associated with
the metric M = ΘTΘ > 0. After some calculations, one has
ΘJΘ−1
Θ1J1Θ
Θ2J2Θ
≤ max(−λ1,−λ2)I uniformly (3.1)
The augmented system is thus stochastically contracting in the metric
M, with rate min(λ1, λ2) and bound C1 + kC2.
Hierarchical combination We first recall a standard result in matrix
analysis [19]. Let A be a symmetric matrix in the formA =
A21 A2
Assume that A1 and A2 are definite positive. Then A is definite positive if
sing2(A21) < λmin(A1)λmin(A2) where sing(A21) denotes the largest singu-
lar value of A21. In this case, the smallest eigenvalue of A satisfies
λmin(A) ≥
λmin(A1) + λmin(A2)
λmin(A1)− λmin(A2)
+ sing2(A21)
Consider now the same set-up as in the previous paragraph, except that
the connection is now hierarchical and upper-bounded. More precisely, the
Jacobians of the couplings verify J12 = 0 and sing
2(Θ2J21Θ
1 ) ≤ K. The
Jacobian matrix of the augmented system is then given by J =
J21 J2
Consider a coordinate transform Θǫ =
0 ǫΘ2
associated with the
metric Mǫ = Θ
ǫ Θǫ > 0. After some calculations, one has
ΘJΘ−1
Θ1J1Θ
ǫ(Θ2J21Θ
ǫΘ2J21Θ
Θ2J2Θ
Set now ǫ =
2λ1λ2
. The augmented system is then stochastically con-
tracting in the metric Mǫ, with rate
(λ1 + λ2 −
λ21 + λ
2)) and bound
2C2λ1λ2
Small gains In this paragraph, we require no specific assumption on the
form of the couplings. Consider the coordinate transformΘ =
associated with the metric Mk = Θ
Θk > 0. Aftersome calculations, one
Θ1J1Θ
Θ2J2Θ
where Bk =
kΘ2J21Θ
Θ1J12Θ
Following the matrix analysis result stated at the beginning of the pre-
vious paragraph, if infk>0 sing
2(Bk) < λ1λ2 then the augmented system is
stochastically contracting in the metric Mk, with bound C1 + kC2 and rate
λ verifying
λ ≥ λ1 + λ2
λ1 − λ2
+ inf
sing2(Bk) (3.2)
4 Some examples
4.1 Effect of measurement noise on contracting observers
Consider a nonlinear dynamical system
ẋ = f(x, t) (4.1)
If a measurement y = y(x) is available, then it may be possible to choose
an output injection matrix K(t) such that the dynamics
˙̂x = f(x̂, t) +K(t)(ŷ − y) (4.2)
is contracting, with ŷ = y(x̂). Since the actual state x is a particular
solution of (4.2), any solution x̂ of (4.2) will then converge towards x expo-
nentially.
Assume now that the measurements are corrupted by additive “white
noise”. In the case of linear measurement, the measurement equation be-
comes y = H(t)x+Σ(t)ξ(t) where ξ(t) is a multidimensional “white noise”
and Σ(t) is the matrix of measurement noise intensities.
The observer equation is now given by the following Itô stochastic dif-
ferential equation (using the formal rule dW = ξdt)
dx̂ = (f(x̂, t) +K(t)(H(t)x −H(t)x̂))dt+K(t)Σ(t)dW (4.3)
Next, remark that the solution x of system (4.1) is a also a solution of
the noise-free version of system (4.3). By corollary 1, one then has, for any
solution x̂ of system (4.3)
∀t ≥ 0 E
‖x̂(t)− x(t)‖2
+ ‖x̂0 − x0‖2e−2λt (4.4)
where
λ = inf
∣∣∣∣λmax
∂f(x, t)
−K(t)H(t)
)∣∣∣∣
C = sup
Σ(t)TK(t)TK(t)Σ(t)
Remark The choice of the injection gain K(t) is governed by a trade-
off between convergence speed (λ) and noise sensitivity (C/λ) as quantified
by (4.4). More generally, the explicit computation of the bound on the
expected quadratic estimation error given by (4.4) may open the possibility
of measurement selection in a way similar to the linear case. If several
possible measurements or sets of measurements can be performed, one may
try at each instant (or at each step, in a discrete version) to select the
most relevant, i.e., the measurement or set of measurements which will best
contribute to improving the state estimate. Similarly to the Kalman filters
used in [9] for linear systems, this can be achieved by computing, along
with the state estimate itself, the corresponding bounds on the expected
quadratic estimation error, and then selecting accordingly the measurement
which will minimize it.
4.2 Estimation of velocity using composite variables
In this section, we present a very simple example that hopefully suggests the
many possibilities that could stem from the combination of our stochastic
stability analysis with the composite variables framework [31].
Let x be the position of a mobile subject to a sinusoidal forcing
ẍ = −U1ω2 sin(ωt) + 2U2
where U1 and ω are known parameters. We would like to compute good
approximations of the mobile’s velocity v and acceleration a using only mea-
surements of x and without using any filter. For this, construct the following
observer
−αv 1
−αa 0
(αa − α2v)x
−αaαvx− U1ω3 cos(ωt)
3 cos(ωt)
(4.5)
and introduce the composite variables v̂ = v + αvx and â = a + αax. By
construction, these variables follow the equation
v̂ − v
−U1ω3 cos(ωt)
(4.6)
and therefore, a particular solution of (v̂, â) is clearly (v, a). Choose now
αa = α
v = α
2 and let Mα =
α2 −α/2
−α/2 1
. One can then show that
system (4.6) is contracting with rate λα = α/2 in the metric Mα. Thus, by
the basic contraction theorem [24], (v̂, â) converges exponentially to (v, a)
with rate λα in the metric Mα. Also note that the β-bound corresponding
to the metric Mα is given by βα =
1+α2−
α4−α2+1
Next, assume that the measurements of x are corrupted by additive
“white noise”, so that xmeasured = x+ σξ. Equation (4.5) then becomes an
Itô stochastic differential equation
3 cos(ωt)
By definition of B, the variance of the noise in the metric Mα is upper-
bounded by α
. Thus, using again corollary 1, one obtains (see Figure 1
for a numerical simulation)
∀t ≥ 0 E
‖v̂(t)− v(t)‖2 + ‖â(t)− a(t)‖2
‖v̂0 − v0‖2 + ‖â0 − a0‖2
4.3 Stochastic synchronization
Consider a network of n dynamical elements coupled through diffusive con-
nections
dxi =
f(xi, t) +
j 6=i
Kij(xj − xi)
 dt+ σi(xi, t)dW di i = 1, . . . , n
(4.7)
x, t) =
f(x1, t)
f(xn, t)
 , ⌢σ(⌢x, t) =
σ1(x1, t) 0 0
. . . 0
0 0 σn(xn, t)
The global state
x then follows the equation
x, t)− L⌢x
x, t)dW nd (4.8)
 0  2  4  6  8  10
 0  2  4  6  8  10
Figure 1: Estimation of the velocity of a mobile using noisy measurements
of its position. The simulation was performed using the Euler-Maruyama
algorithm [18] with the following parameters: U1 = 10, U2 = 2, ω = 3,
σ = 10 and α = 1. Left plot: simulation for one trial. The plot shows
the measured position (red), the actual velocity (blue) and the estimate of
the velocity using the measured position (green). Right plot: the average
over 1000 trials of the squared error ‖v̂ − v‖2 + ‖â − a‖2 (green) and the
asymptotic bound
= 200
given by our approach (red).
In the sequel, we follow the reasoning of [28], which starts by defining an
appropriate orthonormal matrix V describing the synchronization subspace
(V represents the state projection on the subspace M⊥, orthogonal to the
synchronization subspace M = {(x1, . . . ,xn)T : x1 = . . . = xn}, see [28] for
details). Denote by
y the state of the projected system,
y = V
x. Since the
mapping is linear, Itô differentiation rule simply yields
y = Vd
x, t) −VL⌢x
x, t)dW nd
y, t) −VLVT⌢y
y, t)dW nd (4.9)
Assume now that ∂f
is uniformly upper-bounded. Then for strong
enough coupling strength, A = V ∂f
VT − VLVT will be uniformly neg-
ative definite. Let λ = |λmax(A)| > 0. System (4.9) then verifies condition
(H1) with rate λ. Assume furthermore that each noise intensity σi is upper-
bounded by a constant Ci (i.e. supx,t tr(σi(x, t)
Tσi(x, t)) ≤ Ci). Condition
(H2) will then be satisfied with the bound C =
i Ci.
Next, consider a noise-free trajectory
yu(t) of system (4.9). By theo-
rem 3 of [28], we know that
yu(t) converges exponentially to zero. Thus,
by corollary 1, one can conclude that, after exponential transients of rate λ,
E(‖⌢y(t)‖2) ≤ C
On the other hand, one can show that
‖⌢y(t)‖2 = 1
‖xi − xj‖2
Thus, after exponential transients of rate λ, we have
‖xi − xj‖2 ≤
Remarks
• The above development is fully compatible with the concurrent syn-
chronization framework [28]. It can also be easily generalized to the
case of time-varying metrics by combining theorem 3 of this paper and
corollary 1 of [28].
• The synchronization of Itô dynamical systems has been investigated
in [4]. However, the systems considered by the authors of that article
were dissipative. Here, we make a less restrictive assumption, namely,
we only require ∂f
to be uniformly upper-bounded. This enables us
to study the synchronization of a broader class of dynamical systems,
which can include nonlinear oscillators or even chaotic systems.
Example As illustration of the above development, we provide here a de-
tailed analysis for the synchronization of noisy FitzHugh-Nagumo oscillators
(see [35] for the references). The dynamics of two diffusively-coupled noisy
FitzHugh-Nagumo oscillators can be described by
dvi = (c(vi + wi − 13v
i + Ii) + k(v0 − vi))dt+ σdWi
dwi = −1c (vi − a+ bwi)dt
i = 1, 2
Let x = (v1, w1, v2, w2)
T and V = 1√
1 0 −1 0
0 1 0 −1
. The Jaco-
bian matrix of the projected noise-free system is then given by
c− c(v
− k c
−1/c −b/c
Thus, if the coupling strength verifies k > c then the projected system
will be stochastically contracting in the diagonal metric M = diag(1, c) with
rate min(k − c, b/c) and bound σ2. Hence, the average absolute difference
between the two “membrane potentials” |v1 − v2| will be upper-bounded by
min(1, c)min(k − c, b/c) (see Figure 2 for a numerical simulation).
Acknowledgments We are grateful to Dr S. Darses, Prof D. Bennequin
and Prof M. Yor for stimulating discussions, and to the Associate Editor
and the reviewers for their helpful comments.
 0  2  4  6  8  10
Figure 2: Synchronization of two noisy FitzHugh-Nagumo oscillators. The
simulation was performed using the Euler-Maruyama algorithm [18] with
the following parameters: a = 0.3, b = 0.2, c = 30, k = 40 and σ = 1. The
plot shows the “membrane potentials” of the two oscillators.
A Appendix
A.1 Proof of the supermartingale property
Lemma 3 Consider a Markov stochastic process x(t) and a non-negative
function V such that ∀t ≥ 0 EV (x(t)) < ∞ and
∀x ∈ Rn ÃV (x) ≤ −λV (x) (A.1)
where λ is a non-negative real number and Ã is the infinitesimal operator
of the process x(t). Then V (x(t)) is a supermartingale with respect to the
canonical filtration Ft = {x(s), s ≤ t}.
We need to show that for all s ≥ t, one has E(V (x(s))|Ft) ≤ V (x(t)).
Since x(t) is a Markov process, it suffices to show that
∀x0 ∈ Rn E(V (x(t))|x(0) = x0) ≤ V (x0)
By Dynkin’s formula, one has for all x0 ∈ Rn
Ex0V (x(t)) = V (x0) + Ex0
ÃV (x(s))ds
≤ V (x0)− λEx0
V (x(s))ds ≤ V (x0)
where Ex0(·) = E(·|x(0) = x0).
A.2 A variation of Gronwall’s lemma
Lemma 4 Let g : [0,∞[→ R be a continuous function, C a real number
and λ a strictly positive real number. Assume that
∀u, t 0 ≤ u ≤ t g(t)− g(u) ≤
−λg(s) + Cds (A.2)
∀t ≥ 0 g(t) ≤ C
g(0) − C
e−λt (A.3)
where [·]+ = max(0, ·).
Proof Case 1 : C = 0, g(0) > 0.
Define h(t) by
∀t ≥ 0 h(t) = g(0)e−λt
Remark that h is positive with h(0) = g(0), and satisfies (A.2) where the
inequality has been replaced by an equality
∀u, t 0 ≤ u ≤ t h(t)− h(u) = −
λh(s)ds
Consider now the set S = {t ≥ 0 | g(t) > h(t)}. If S = ∅ then the
lemma holds true. Assume by contradiction that S 6= ∅. In this case, let
m = inf S < ∞. By continuity of g and h and by the fact that g(0) = h(0),
one has g(m) = h(m) and there exists ǫ > 0 such that
∀t ∈]m,m+ ǫ[ g(t) > h(t) (A.4)
Consider now φ(t) = g(m)− λ
g(s)ds. Equation (A.2) implies that
∀t ≥ m g(t) ≤ φ(t)
In order to compare φ(t) and h(t) for t ∈]m,m+ ǫ[, let us differentiate the
ratio φ(t)/h(t).
φ′h− h′φ
−λgh+ λhφ
λh(φ − g)
Thus φ(t)/h(t) is increasing for t ∈]m,m + ǫ[. Since φ(m)/h(m) = 1, one
can conclude that
∀t ∈]m,m+ ǫ[ φ(t) ≥ h(t)
which implies, by definition of φ and h, that
∀t ∈]m,m+ ǫ[
g(s)ds ≤
h(s)ds (A.5)
Choose now t0 such that m < t0 < m+ ǫ, then one has by (A.4)
g(s)ds >
h(s)ds
which clearly contradicts (A.5).
Case 2 : C = 0, g(0) ≤ 0
Consider the set S = {t ≥ 0 | g(t) > 0}. If S = ∅ then the lemma holds
true. Assume by contradiction that S 6= ∅. In this case, let m = inf S < ∞.
By continuity of g and by the fact that g(0) ≤ 0, one has g(m) = 0 and
there exists ǫ such that
∀t ∈]m,m+ ǫ[ g(t) > 0 (A.6)
Let t0 ∈]m,m+ ǫ[. Equation (A.2) implies that
g(t0) ≤ −λ
g(s)ds
which clearly contradicts (A.6).
Case 3 : C 6= 0
Define ĝ = g − C/λ. One has
∀u, t 0 ≤ u ≤ t ĝ(t)−ĝ(u) = g(t)−g(u) ≤
−λg(s)+Cds = −
λĝ(s)ds
Thus ĝ satisfies the conditions of Case 1 or Case 2, and as a consequence
∀t ≥ 0 ĝ(t) ≤ [ĝ(0)]+e−λt
The conclusion of the lemma follows by replacing ĝ by g−C/λ in the above
equation. �
References
[1] N. Aghannan and P. Rouchon. An intrinsic observer for a class of la-
grangian systems. IEEE Transactions on Automatic Control, 48, 2003.
[2] D. Angeli. A lyapunov approach to incremental stability properties.
IEEE Transactions on Automatic Control, 47:410–422, 2002.
[3] L. Arnold. Stochastic Differential Equations : Theory and Applications.
Wiley, 1974.
[4] T. Caraballo and P. Kloeden. The persistence of synchronization under
environmental noise. Proceedings of the Royal Society A, 461:2257–
2267, 2005.
[5] T. Caraballo, P. Kloeden, and B. Schmalfuss. Exponentially stable
stationary solutions for stochastic evolution equations and their pertu-
bation. Applied Mathematics and Optimization, 50:183–207, 2004.
[6] L. d’Alto and M. Corless. Incremental quadratic stability. In Proceed-
ings of the IEEE Conference on Decision and Control, 2005.
[7] B. Demidovich. Dissipativity of a nonlinear system of differential equa-
tions. Ser. Mat. Mekh., 1961.
[8] H. Deng, M. Krstic, and R. Williams. Stabilization of stochastic nonlin-
ear systems driven by noise of unknown covariance. IEEE Transactions
on Automatic Control, 46, 2001.
[9] E. Dickmanns. Dynamic Vision for Intelligent Vehicles. Course Notes,
MIT EECS dept, 1998.
[10] K. El Rifai and J.-J. Slotine. Contraction and incremental stability.
Technical report, MIT NSL Report, 2006.
[11] W. Feller. An Introduction to Probability Theory and Its Applications.
Wiley, 1968.
[12] P. Florchinger. Lyapunov-like techniques for stochastic stability. SIAM
Journal of Control and Optimization, 33:1151–1169, 1995.
[13] P. Florchinger. Feedback stabilization of affine in the control stochastic
differential systems by the control lyapunov function method. SIAM
Journal of Control and Optimization, 35, 1997.
[14] V. Fromion. Some results on the behavior of lipschitz continuous sys-
tems. In Proceedings of the European Control Conference, 1997.
[15] B. Girard, N. Tabareau, Q.-C. Pham, A. Berthoz, and J.-J. Slotine.
Where neuroscience and dynamic system theory meet autonomous
robotics: a contracting basal ganglia model for action selection. Neural
Networks, 2008.
[16] P. Hartmann. Ordinary differential equations. Wiley, 1964.
[17] R. Has’minskii. Stochastic Stability of Differential Equations. Sijthoff
and Nordhoff, Rockville, 1980.
[18] D. Higham. An algorithmic introduction to numerical simulation of
stochastic differential equations. SIAM Review, 43:525–546, 2001.
[19] R. Horn and C. Johnson. Matrix Analysis. Cambridge University Press,
1985.
[20] J. Jouffroy. Some ancestors of contraction analysis. In Proceedings of
the IEEE Conference on Decision and Control, 2005.
[21] J. Jouffroy and T. Fossen. On the combination of nonlinear contracting
observers and uges controllers for output feedback. In Proceedings of
the IEEE Conference on Decision and Control, 2004.
[22] H. Kushner. Stochastic Stability and Control. Academic Press, 1967.
[23] D. Lewis. Metric properties of differential equations. American Journal
of Mathematics, 71:294–312, 1949.
[24] W. Lohmiller and J.-J. Slotine. On contraction analysis for nonlinear
systems. Automatica, 34:671–682, 1998.
[25] W. Lohmiller and J.-J. Slotine. Nonlinear process control using con-
traction theory. A.I.Ch.E. Journal, 2000.
[26] W. Lohmiller and J.-J. Slotine. Contraction analysis of nonlinear dis-
tributed systems. International Journal of Control, 78, 2005.
[27] X. Mao. Stability of Stochastic Differential Equations with Respect to
Semimartingales. Longman, White Plains, NY, 1991.
[28] Q.-C. Pham and J.-J. Slotine. Stable concurrent synchronization in
dynamic system networks. Neural Netw, 20(1):62–77, Jan. 2007.
[29] V. Popov. Hyperstability of Control Systems. Springer-Verlag, 1973.
[30] W. Rudin. Real and complex analysis. McGraw-Hill, 1987.
[31] J.-J. Slotine and W. Li. Applied Nonlinear Control. Prentice-Hall, 1991.
[32] E. Sontag and Y. Wang. Output-to-state stability and detectability of
nonlinear systems. Systems and Control Letters, 29:279–290, 1997.
[33] N. Tabareau and J.-J. Slotine. Notes on contraction theory. Technical
report, MIT NSL Report, 2005.
[34] J. Tsinias. The concept of “exponential iss” for stochastic systems and
applications to feedback stabilization. Systems and Control Letters,
36:221–229, 1999.
[35] W. Wang and J.-J. E. Slotine. On partial contraction analysis for cou-
pled nonlinear oscillators. Biol Cybern, 92(1):38–53, Jan. 2005.
[36] Y. Zhao and J.-J. Slotine. Discrete nonlinear observers for inertial
navigation. Systems and Control Letters, 54, 2005.
	Introduction
	Main results
	Background
	Nonlinear contraction theory
	Standard stochastic stability
	The stochastic contraction theorem
	Settings
	The basic stochastic contraction theorem
	Generalization to time-varying metrics
	Strength of the stochastic contraction theorem
	``Optimality'' of the mean square bound
	No asymptotic almost-sure stability
	Noisy and noise-free trajectories
	Combinations of contracting stochastic systems
	Some examples
	Effect of measurement noise on contracting observers
	Estimation of velocity using composite variables
	Stochastic synchronization
	Appendix
	Proof of the supermartingale property
	A variation of Gronwall's lemma
