http://www.ma-xy.com
目录
第一章 随机微分方程 1
1.1 符号注记 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.2 随机微分方程建模 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2
1.2.1 人口增长模型 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2
1.3 随机微分方程基本理论 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.3.1 随机微分方程的一般形式 . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.3.2
随机过程简介
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
1.3.3 随机积分的初步尝试 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.3.4 随机积分 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
1.3.5 Ito 过程与 Ito 公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
1.3.6 解的存在唯一性及解析计算公式 . . . . . . . . . . . . . . . . . . . . . . . . 22
1.4 随机微分方程数值方法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
1.4.1 随机泰勒展开 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
1.4.2 Euler 方法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
1.4.3 解的收敛性与稳定性 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
1.4.4 数值方法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
1.4.5 数值实验 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
1.5 SDE 解存在唯一条件的进一步讨论 . . . . . . . . . . . . . . . . . . . . . . . . . . 38
1.5.1 非全局 Lipschitz 条件下解的存在唯一性 . . . . . . . . . . . . . . . . . . . 38
1.5.2 非全局 Lipschitz 条件下数值方法 . . . . . . . . . . . . . . . . . . . . . . . 40
1.6 SDE 结构的进一步讨论 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
1.6.1 SDE 基本理论 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
1.6.2 含 Markov 切换的 SDE . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
1.6.3 含 Poisson 跳的 SDE . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
1.6.4 随机延迟微分方程 SDDE . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
1.6.5 中立型 SDE . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
1.6.6 分段连续型 SDE . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
1.6.7 SDE 参数估计问题 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
1.6.8 随机偏微分方程 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
1
http://www.ma-xy.com
目录 目录
1.7 MATLAB 与随机微分方程 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
1.7.1 MATLAB 自带工具箱 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
1.7.2 SDEtoolbox 工具箱 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
http://www.ma-xy.com 2 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程
1.1 符号注记
R : R = R
1
，实数集
R
+
: [0, ∞]
N = N ：整数集
R
d
：d 维欧几里得空间 R
d×m
：d × m 维矩阵空间
a.s：依概率收敛
L
1
(0, T ; R
n
)：所有实值可测 F
t
适应的且
T
0
|f
t
|ds < ∞, a.s 的 n 维随机过程集合。
L
2
(0, T ; R
n
)：所有实值可测 F
t
适应的且
T
0
|f
t
|
2
ds < ∞, a.s 的 n 维随机过程集合。
C([a, b]; R
n
)：在 [a, b] 上连续的函数的集合。
C
b
F
t
([a, b]; R
n
)：有界，F
t
可测的 C([a, b]; R
n
) 的随机变量集合。
C(R
d
× R
+
; R)：R
d
× R
+
上关于 x2 次连续可微，关于 t 一次可微的所有连续非负函数
v(x, t) 的集合。
a ∨ b：max{a, b}
a ∧ b：min{a, b}
下面，给出统计中的 4 中收敛的定义：
1. 几乎处处收敛 (依概率 1 收敛)
lim
k→∞
X
k
= X a.s
记为 X
n
→ X a.s。
2. 依概率收敛。∀ε > 0，当 k → ∞ 时，有
P {|X
k
− X| > ε} → 0
称 X
k
随机收敛到 X，记为 X
n
P
−→ X。
3.L
p
上 p 阶矩收敛。
E|X
k
− X|
p
→ 0
4. 依分布收敛。
lim
k→∞
E
g
(X
k
) = E
g
(X)
1
http://www.ma-xy.com
1.2 随机微分方程建模 第一章 随机微分方程
记为 X
n
L
−→ X。
上面 4 中统计下的收敛之间的关系是：
几乎必然收敛 ⇒ 依概率收敛 ⇒ 依分布收敛
p 阶矩收敛 ⇒ 依概率收敛 ⇒ 依分布收敛
1.2 随机微分方程建模
1.2.1 人口增长模型
在前面的常微分方程章节中，我们已经讨论过基础的人口增长模型
˙x
t
=
dx
t
dt
= ax
t
(1.1)
其中：x
t
表示 t 时刻人口数量， ˙x
t
表示 t 时刻人口增长，a 为人口增长率。
上述模型假设了人口增长率 a 为常值，且 t 时刻的人口增长 ˙x
t
与人口基数 x
t
成正比，由此
建立了一阶常系数线性微分方程，其解 x
t
所表现出来的特点是“光滑”。如果假设 a 是时变的，
即 a
t
，不同时刻 t 的增长率 a
t
是不同的，进而可以将原 ODE 模型变为一阶变系数线性微分方
程。下面，我们再对 a
t
做一个详细的分析。考虑 t 时刻的人口增长率 a
t
，a
t
是由一些确定因素
和不确定因素组成的，我们假设确定因素 (确定增长率) 为 r
t
，不确定因素 (战争、干涸等) 为 ε
t
，
并假设二者之间是加法关系，而非乘法等其他关系，于是有
a
t
= r
t
+ µε
t
其中：µ 是参数，为了简单，设置 µ = 1；ε
t
是一干扰，可以是多种不确定因素的共同作用，是
一个随机变量，记为 ε
t
∼ f (θ)。
由于 ε
t
是一随机变量，导致 t 时刻的增长率 a
t
也变为随机变量。将 a
t
带入原人口增长的
常微分方程 (1.1) 当中，有
dX
t
dt
= (r
t
+ ε
t
)X
t
(1.2)
其中：为了区别变量 x
t
和随机变量 x
t
，我们用大写的 X
t
来表示随机变量。
由于 ε
t
是一随机变量，导致 t 时刻的人口增长
dx
t
dt
也变为一随机变量，从而 t + 1 时刻的
人口数量 x
t+1
亦为一个随机变量，从而 {x
t
}
T
t=1
是一个随机过程。并且注意到时间 t 的连续性，
所以 {x
t
}
T
t=1
是一个连续随机过程。注意到，我们常用的马氏链、Poisson 过程和 Brown 运动是
在离散情况下讨论的，因此，有必要将离散随机过程扩展到连续随机过程。将上式 (1.2) 的 dt 右
移，有
dX
t
= r
t
X
t
dt + ε
t
X
t
dt
其中：r
t
dt 可以视为确定性增长率的微分，不妨记为 dx
t
= r
t
dt；ε
t
dt 可以视为随机性增长率的
微分，不妨记为 dW
t
= ε
t
dt。
http://www.ma-xy.com 2 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
但是，如果 dW
t
不可微，即 dW
t
在古典微分中不可微，那么，我们就不能将其写为 dW
t
= ε
t
dt，
而是要记回 dW
t
，即
dX
t
= r
t
X
t
dt + X
t
dW
t
(1.3)
显然，对于随机过程 W
t
而言，W
t
不一定是可微的，这一点我们在后面会用 Kolmogorov 定理
来进行说明。我们用式 (1.3) 来表示随机微分方程 (SDE)，并给其一个初始条件：t
0
= 0, X
t
0
= X
0
，
形成完整的随机微分方程模型
dX
t
= r
t
X
t
dt + X
t
dW
t
t ∈ [t
0
, T ]
X
t
0
= X
0
(1.4)
下面，我们来讨论 SDE 的解 X
t
。如果一个随机过程 X
t
满足上面的 SDE，则 X
t
为 SDE
问题的解。像 ODE 那样，我们尝试对方程两边求积分，有
X
t
= X
0
+
t
0
r
s
X
s
ds +
t
0
X
s
dW
s
称上式为随机微分方程对应的随机积分方程。然而，在 W
t
不可微的情况下，积分
t
0
X
s
dW
s
又
是什么呢？
上述讨论遗留了许多问题，下面，我们将人口增长的随机微分方程模型 (1.3) 推广到一般的
形式，并在 SDE 的一般形式下对其进行讨论。下面的内容包括：1、SDE 的一般形式。我们将人
口增长模型推广到高维情况，然后正则化，标准化，并讨论 SDE 的变量域。2、介绍概率空间和
随机过程。3、Ito 积分、Ito 积分的性质、Ito 过程和 SDE 的解 X
t
。4、解的存在唯一性。5、解
的稳定性和收敛性。
1.3 随机微分方程基本理论
1.3.1 随机微分方程的一般形式
SDE 起源于马尔科夫过程、Kolmogorov 分析法和 Feller 半群方法。1902 年，Gibbs 和
Boltzmann 等人在研究统计力学的积分问题时，引入了随机微分方程问题。1923 年，Wiener 在
数学上给出了 Brown(布朗) 运动严格的理论研究成果，所以 Brown 运动也称为 Wiener 过程。由
于
Brown
运动的随机性，
Ito
在
1944
给出与
Wiener
过程相应的随机系统的
Ito
积分。
1951
年，
Ito 提出 Ito 公式，也即随机链式法则。1951 年，Ito 著作《Stochastic Dierential Equation》中
详细讨论了如下 Ito 型随机微分方程
dX
t
= b(X
t
, t)dt + σ(X
t
, t)dW
t
X
t
0
= X
0
1961 年，Ito 发表了《论随机微分方程》，至此，标志着随机微分方程基本理论的成熟。
http://www.ma-xy.com 3 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
下面，我们将前面引入的人口增长的随机微分方程 (1.4) 规范一下：将 dt 和 dW
t
前面的系
数项规范化，写为一般形式，有
d
X
t
=
f
(
X
t
, t
)
d
t
+
g
(
X
t
, t
)
d
W
t
t
∈
[
t
0
, T
]
X
t
0
= X
0
(1.5)
其中：f, g 是 [t
0
, T ] 上的 Borel 连续可测函数，E[|X
0
|
2
] < ∞，W
t
是某一随机过程。
注：Borel 连续可测函数 (可测函数)：设 (Ω, F) 为一可测空间，R 为实数域，
¯
R = R ∪{−∞, ∞}，
分别用 B(R) 和 B(
¯
R) 表示 R 和
¯
R 上的 Borelσ 代数，令 f 为 Ω 到
¯
R 中的一个映射，如果
f
−1
[B(
¯
R)] ∈ F，则称 f 为 Borelσ 代数，简称为可测函数，若进一步 f 只取实值，则称 f 为实
值可测函数。Borel 可测包含于勒贝格可测，勒贝格可测将在后面介绍。
我们称
f
为漂移函数，
g
为扩散函数，
f
(
·
)
d
t
为漂移项，
g
(
·
)
d
W
t
为扩散项。上式是一维下
的 SDE 的一般形式，实际上，模型中有两个随机过程 {X
t
} 和 {W
t
}，不考虑其他问题，定义一
个广义映射 f，有
˙
X
t
= f (X
t
, W
t
, t) t ∈ [t
0
, T ]
X
t
0
= X
0
(1.6)
式 (1.6) 定义的一般形式的 SDE 是一维情况下的，{X
t
} 是一个随机过程，而不是多个随机过程。
我们易将式 (1.6) 推广到高维—随机微分方程系统
˙
X
t
= f (X
t
, W
t
, t) t ∈ [t
0
, T ]
X
t
0
= X
0
其中：X
t
∈ R
n
是 n 维随机过程，n 为方程中变量的个数，或者说是待求随机过程 X
t
的个数；
W
t
∈ R
m
是 m 维随机过程，m 为随机过程 W
t
的个数。
下面，我们在式 (1.5) 定义的一维 SDE 下进行讨论，并且假设随机过程 W
t
是一维 Wiener
过程。
1.3.2 随机过程简介
随机过程
随机过程是概率空间 (Ω, F, P ) 上的一族随机变量 {X
t
}
T
t=t
0
，对每个时刻 t，X
t
是一个随机
变量，因此，在 [t
0
, T ] 内 {X
t
} 为一随机过程。X
t
表示随机过程 (随机系统) 在时刻 t 的状态。
当 [t
0
, T ] ∈ R × R
+
时，称 {X
t
} 为连续型随机过程；当时间 [t
0
, T ] 取值如 1, 2, . . . , T 时，
称 {X
t
} 为离散型随机过程。为了便于书写，我们一般默认了 X
t
就是一个随机过程。注意：一
般随机过程教材上，讨论的随机过程 X
t
都是离散型随机过程，而我们式 (1.5) 中待求的随机过
程 X
t
是连续型随机过程。
注：概率空间 (Ω, F, P ) 的简单注释：Ω 为给定集合，w ∈ Ω 为样本点；F 为 Ω 子集的 σ 代数，
事件 A ∈ F；P 为 F → [0, 1] 的概率测度，P (A) 为事件 A 的概率。
http://www.ma-xy.com 4 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
随机过程的基本性质
上面介绍了随机过程的基本定义，下面，我们给出随机过程的一些基本性质。如果 {X
t
} 为
一随机过程，则
µ
x
(t) = E[X
t
]
为过程 {X
t
} 的均值函数。因为我们在每一个时刻 t 上求均值 E[X
t
]，从而，均值函数是时间 t
的函数，所以写为 µ
x
(t)。
如果 ∀t ∈ [t
0
, T ]，有
E[X
2
t
] < ∞
即任意时刻 t，随机变量 X
t
平方的期望 E[X
2
t
] 存在 (有限)，则称 {X
t
} 为二阶矩过程。
如果 {X
t
} 为二阶矩过程，则称
r(t
1
, r
2
) = E[(X
t
1
− µ
x
(t
1
))(X
t
2
− µ
x
(t
2
))]
为 {X
t
} 的协方差。称
V ar[X
t
] = r(t, t)
为 {X
t
} 的方差函数。称
R
x
(s, t) = E[X
s
X
t
]
为 {X
t
} 的自相关函数。
定义 (宽平稳过程) 若 {X
t
} 为二阶矩过程，且随机过程的均值函数为常值 E[X
t
] = µ，协
方差 r(t, s) 只与时间间隔 (t − s) 有关，则 {X
t
} 为宽平稳过程
¬
。
平稳性是时间序列和随机过程的一个要点，实验中，我们只有随机过程 {X
t
} 的一次样本
{X
t
}
T
t=1
(因为时间不可逆)。我们想到，如果要求某一时刻 t 的均值 E[X
t
]，那么，能否用样本
{X
t
}
T
t=1
的均值来代替 t 时刻的均值呢？即能否用时间均值来替代时刻均值呢？我们知道，如果
{X
t
} 每一时刻的随机变量 X
t
相同，那这种替代就是可行的，这也是我们讨论 {X
t
} 的平稳性的
原因之一。
在式 (1.5) 中我们假设 W
t
是 Brown 运动，下面，我们就来介绍一下这个随机过程。
定义 (Brown 运动) 如果随机过程 {X
t
} 满足
(1)X(0) = 0；
(2){X
t
} 是平稳独立增量 (t ⩾ 0)；
(3)∀t > 0, X
t
∼ N (0, σ
2
t)。
则称随机过程 {X
t
} 为 Brown 运动。其中：N (·) 为正态分布，平稳独立增量是指：增量
X(t
i+1
) − X(t
i
) 平稳且相互独立。
¬
注：严平稳过程的有限维分布在时间上是不变的，可以想象为各时刻 t 上的分布相同。
http://www.ma-xy.com 5 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
Brown 运动的基本性质
上面给出了 Brown 运动的定义，下面简单讨论一下 Brown 运动的性质。设 {W
t
} 是一个
Brown 运动，Brown 运动具有如下一些性质：
1). 正态增量。∀s, 0 ⩽ s ⩽ t，有
W
t
− W
s
∼ N (0, t − s)
即增量 (W
t
− W
s
) 服从正态分布。
2). 独立增量。∀s, 0 ⩽ s ⩽ t，有 W
t
− W
s
独立于过程 {W
t
} 的过去状态 W
u
,0 ⩽ u ⩽ s。
3). 路径连续性。W
t
是时间 t 的连续函数。
样本路径是指随机过程 {W
t
} 具体的一次样本实现，如果我们采用数值模拟的方法，则每次
模拟的过程数据都是该过程的路径。不妨记 W (·, w) 为过程 {W
t
} 的一个样本路径，则 W (·, w)
有如下性质：
1. W (·, w) 是时间 t 的连续函数；
2. W (·, w) 在任意区间上不单调；
3. W (·, w) 在任意点处不可微；
4. W (·, w) 在任意区间上有无限变差；
5. ∀t，W (·, w) 在 [0, t] 上的二次变差为 t；
Brown 运动又称为 Wiener 过程，一般简写为 W
t
或者 B
t
。根据前面的分析，我们知道，式
(1.5) 定义的 SDE 中的 X
t
, W
t
为随机过程，而且该随机过程为连续随机过程。但由于 W
t
是
Brown 运动，所以 W
t
几乎处处不可微，这导致 dW
t
不能像一般的函数那样微分，不仅是因为
dW
t
不易处理，随机积分
t
0
g(·)dW
s
也是不易处理的。由于 W
t
不可导，我们不妨直接使用 dW
t
来标记，只要我们能解决
t
0
g(·)dW
s
，
那么问题也就差不多了。关于 W
t
的几乎处处不可微问题，我们有下面的 Kolmogonov 定理
定理 (Kolmogonov 定理) 如果 {X
t
} 几乎处处的样本路径是连续的，并且 ∃α > 0, β > 0，
有
E[|X
t
− X
s
|
β
] ⩽ c|t − s|
1+α
则对 0 < r <
α
β
，T > 0，以及几乎处处的样本 w 路径 X(·, w) 在 [0, T ] 上一致 γ −Holder 连续。
我们将 Kolmogonov 定理应用到 Brown 运动中，有
E[|W
t
− W
s
|
4
] ⩽ n(n + 1)|t − s|
2
http://www.ma-xy.com 6 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
这里取 α = 1, β = 4, c = n(n + 1)，于是，有
1. 对 γ <
1
4
，以及几乎处处的样本路径 W (·, w) 是一致 γ − Holder 连续的。
2. 对
1
4
< γ < 1，以及几乎处处的样本路径 W (·, w) 处处不可能 γ − Holder 连续，从而几
乎处处不可微。
1.3.3 随机积分的初步尝试
下面，我们来尝试着求解随机积分
t
0
g(·)dW
t
。一维随机微分方程的一般形式为
dX
t
= f (X
t
, t)dt + g(X
t
, t)dW
t
t ∈ [t
0
, T ]
X
t
0
= X
0
(1.7)
如果一个随机过程 X
t
满足上面的方程，则称 X
t
为 SDE 的解。对上式 (1.7) 两边积分，有
X
t
= X
0
+
t
0
f(X
s
, s)ds +
t
0
g(X
s
, s)dW
s
上式中的 X
t
只是前部分 X
s
(0 ⩽ s ⩽ t) 的积分，我们称上式为 (1.7) 的随机积分方程，其解也
为 X
t
。上式中的第一个积分是普通意义下的 Rieman-Stireltjes 积分并无特别之处，但其第二个
积分则是其独特之处，这是一个随机积分，是随机过程 W
t
的积分，由于 W
t
的不可微性，其积
分形式应该与普通的黎曼-斯蒂尔切斯积分不同，我们不妨用 S-T 积分对其试算一下。
T
0
W
t
dW
t
t ∈ [0, T ]
用普通积分方法求上述随机积分，有
Step1. 将 [0, T ] 离散为 n 段
0 = t
0
< ··· < t
k
< t
k+1
< ··· < t
n
= T
为了简便，我们将其简单的均分成 n 段。
Step2. 取点 τ
k
。不妨记小区间 [t
k
, t
k+1
] 上的某一时刻点为 τ
k
τ
k
= (1 − λ)t
k
+ λt
k+1
则 W (τ
k
) 为随机过程 W
t
在 τ
k
处的状态值。
Step3. 计算 dW
t
，有
dW
t
= W (t
k+1
) − W (t
k
)
则整体的积分可以表示为
T
0
W
t
dW
t
≈
n−1
k=0
W (τ
k
)[W (t
k+1
) − W (t
k
)]
http://www.ma-xy.com 7 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
Step4.W (τ
k
) 的替代。我们用
V (t) = (1 − λ)W (t
k
) + λW (t
k+1
)
来代替 W (τ
k
)。
经过上面的替换，我们有
T
0
W
t
dW
t
≈
n−1
k=0
W (τ
k
)[W (t
k+1
) − W (t
k
)]
≈
n−1
k=0
V (t)[W (t
k+1
) − W (t
k
)]
= (1 − λ)
n−1
k=0
W (t
k
)[W (t
k+1
) − W (t
k
)] + λ
n−1
k=0
W (t
k+1
)[W (t
k+1
) − W (t
k
)]
= (1 − λ)
1
2
W (t
n
)
2
−
1
2
W (t
0
)
2
−
1
2
n−1
k=0
(W (t
k+1
) − W (t
k
))
2
+ λ
1
2
W (t
n
)
2
−
1
2
W (t
0
)
2
+
1
2
n−1
k=0
(W (t
k+1
) − W (t
k
))
2
=
1
2
W (t
n
)
2
−
1
2
W (t
0
)
2
+
1
2
(2λ − 1)
n−1
k=0
(W (t
k+1
) − W (t
k
))
2
所以，当 n → ∞ 时，有
T
0
W
t
dW
t
= lim
n→∞
1
2
W (w
n
)
2
−
1
2
W (t
0
)
2
+ (λ −
1
2
)
n
k=0
(W (t
k+1
) − W (t
k
))
2
=
1
2
W (T )
2
−
1
2
W (0)
2
+ (λ −
1
2
)T (1.8)
我们对上式中划横线的部分进行证明。先引入随机过程的 (二次) 变差的概念，然后再对其进行
证明。
定义 (二次变差) 当 {t
i
}
n
i=0
为时间区间 [0, T ] 上的分割 (分割为 n 段)，且当
δ
n
= max
0⩽k⩽n−1
{t
k+1
− t
k
} → 0
时，[W, W ](T ) 依概率收敛于
[
W, W
](
T
) =
lim
δ
n
→0
n−1
k
=0
|
W
(
t
k+1
)
−
W
(
t
k
)
|
2
则称 [W, W ](T ) 为二次变差。
http://www.ma-xy.com 8 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
注：如果是 n 等分割，则 δ
n
=
T −0
n
；如果是不等分割，则将
[0
, T
]
分割为
n
段，令最大段长度
为 δ
n
。
根据上面二次变差的定义，证明上式 (1.8) 划横线的部分，即证明：
lim
n→∞
n−1
k=0
|W (t
k+1
) − W (t
k
)|
2
= T (1.9)
即 [W, W ](T ) = T 。
证明 证明 [W, W ](T ) = T ，即是说对 ∀ T ，在 [0, T ] 上，随机过程 W
t
的二次变差为 T 。为
便于书写，我们记
S
n
=
n−1
k=0
[W (t
k+1
) − W (t
k
)]
2
则其均值为
E(S
n
) =
n−1
k=0
E([W (t
k+1
) − W (t
k
)]
2
)
=
n−1
k=0
(t
k+1
− t
k
)
= t
注意到 Brown 运动正态增量性质：W (t) − W (s) ∼ N 。由标准正态分布的四阶矩公式，我们可
以计算其方差，有
V ar(S
n
) = V ar
n−1
k=0
[W (t
k+1
) − W (t
k
)]
2
=
n−1
k=0
V ar{[W (t
k+1
) − W (t
k
)]
2
}
=
n−1
k=0
3(t
k+1
− t
k
)
2
⩽ 3 max{t
k+1
− t
k
}t
= 3tδ
n
所以
∞
n=1
V ar(S
n
)∞
由单调收敛定理，得
E
∞
n=1
(S
n
− T )
2
< ∞
http://www.ma-xy.com 9 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
因此，
∞
n=1
(S
n
− T )
2
< ∞, a.s
于是，S
n
− T → 0, a.s(依概率收敛)，故 [W, W ](T ) = T 。2
上面的单调收敛定理后面介绍。由式 (1.8) 和上述证明我们发现，当离散段数 n → ∞ 时，
随机积分
T
0
W
t
dW
t
并不是固定的，它的取值与参数 λ 有关，导致这种现象的根本原因是：W
t
的方差是时变的，如果过程 (函数)w
t
不是随机的，那么它的方差为 0，这样式 (1.8) 划横线的
部分为 0，就不存在参数 λ 了。随机过程 W
t
的方差不仅是时变的，而且还随着时间增大，有
V ar(W
t
) = t，因此，我们不能对
T
0
W
t
dW
t
使用黎曼-斯蒂尔切斯积分，需要定义一个新的积分
形式，这个积分形式是定义在随机过程上的。
下面，我们先给出两类常用的随机积分：Ito 积分和 Stratonovich 积分。对于式 (1.8)
T
0
W
t
dW
t
=
1
2
W (T )
2
−
1
2
W (0)
2
+ (λ −
1
2
)T
当 λ = 0 时，上述积分为 Ito 积分，即
b
a
W
t
dW
t
=
1
2
(W (b)
2
− W (a)
2
) −
1
2
(b − a)
当 λ =
1
2
时，上述积分为 Stratonovich 积分，即
b
a
W
t
dW
t
=
1
2
(W (b)
2
− W (a)
2
)
由于有不同的随机积分定义，我们不得不对随机微分方程加以区分，将式 (1.7) 形式的方程
称为 Ito 型随机微分方程，其解 X
t
称为 Ito 型随机过程，其对应的随机积分方程
X
t
= X
0
+
t
0
f(s, X
s
)ds +
t
0
g(s, X
s
)dW
s
为 Ito 型随机积分方程。
记
dX
t
= f (t, X
t
)dt + g(t, X
t
) ◦ dW
t
X(t
0
) = X
0
为 Stratonovich 型随机微分方程，其解 X
t
为 S 型随机过程，其对应的随机积分方程
X
t
= X
0
+
t
0
f(s, X
s
)ds +
t
0
g(s, X
s
) ◦ dW
s
为 S 型随机积分。
下面，我们将导出 Ito 积分的定义、给出 Ito 积分的性质、给出 Ito 过程的定义并给出 Ito
积分的计算公式
(
就像黎曼积分有牛顿—莱布尼茨公式那样
)
。
http://www.ma-xy.com 10 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
1.3.4 随机积分
前面我们提到过，不能对
T
0
W
t
dW
t
使用黎曼-斯蒂尔切斯积分，我们需要定义一个新的积
分形式。下面，就来一步步引入随机积分—Ito 积分。
黎曼积分
黎曼积分也称 R 积分。如果你不太熟悉定积分，不妨去查阅一下《数学分析讲义 (第五
版)》上册刘玉琏等著 P369-P410，书中的定积分讲的深入浅出，是非常好的。我们的目标是定
义
T
0
X
t
dW
t
，并计算它，让我们从最基本的定积分开始，一步步扩展到随机积分。黎曼积分的
定义如下：设 x
t
是一个普通的函数，而非随机过程，后面我们慢慢将其引入到随机过程 X
t
当
中。考虑 x
t
在区间 [0, T ] 上的积分
T
0
x
t
dt
Step1. 划分 n 段。我们将区间 [0, T ] 离散化为 n 段，并且记分割方法为 A = {t
k
}
n
，有
0 = t
0
< ··· < t
k
< t
k+1
< ··· < t
n
= T
记第 k 段的长度为 ∆t
k
= t
k+1
− t
k
, k = 0, 1, . . . , n −1，当然，也可以写为 ∆t
k
= t
k
− t
k−1
, k =
1, 2, . . . , n，不过一定要注意下标别混乱了。记 n 个小段的最大长度为 δ
δ = max
0⩽k⩽n−1
{∆t
k
}
为了简便，我们将其简单的均分 n 段。
Step2. 取点 τ
k
。在第 k 段 [t
k
, t
k+1
] 中任取一点 τ
k
∈ [t
k
, t
k+1
]，则该点的函数值为 x(τ
k
)。
Step3. 计算黎曼积分
T
0
x
t
dt = lim
δ →0
n−1
k=0
x(τ
k
)∆t
k
Step4.R 可积。若在 δ → 0 时，极限
lim
δ →0
n−1
k=0
x(τ
k
)∆t
k
= I < ∞
存在，且 I 与分割方法 A 无关，与取点 {τ
k
} 点无关，则称 x
t
在区间 [0, T ] 上黎曼可积，即
∀ε > 0, ∃N, ∀A : δ < N, ∀{τ
k
}，有
n−1
k=0
x(τ
k
)∆t
k
− I
⩽ ε
则称 x
t
在 [0, T ] 上黎曼可积，I 为其黎曼积分。
对于黎曼积分，我们有如下可积准则 (黎曼可积的充分条件)：如果 x
t
在区间 [0, T ] 上连续，
则 x
t
可积；如果 x
t
在 [0, T ] 上有界且有有限个间断点，则 x
t
可积；如果 x
t
在 [0, T ] 上单调，
则 x
t
可积。
定理 (勒贝格定理 (黎曼可积的充分必要条件)) 有界函数 x
t
在 [0, T ] 上可积充分必要条件
为 x
t
在 [0, T ] 内的所有间断点是零集。
http://www.ma-xy.com 11 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
注：零集是指可被长度总和任意小的可列个开区间所覆盖的集合。1. 定义一个事物 2. 然后讨论
其充分必要条件 3. 然后再讨论事物的性质 (例如：定义解，解的充分必要条件，解的存在唯一性，
稳定性)。在程序中，1. 定义一个类 2. 然后讨论类的属性 3. 然后讨论类的方法。这些做法都是
很自然的，也是人类认识事物 (类) 的一个过程。
黎曼积分的性质：上述普通定积分具有数乘性和可加性等性质，这里不做介绍。像一般情况
那样，我们不能通过黎曼积分的定义来计算积分，而是需要根据它的充要条件和性质来计算积分。
下面，给出普通定积分的计算公式：牛顿—莱布尼茨公式。若 x
t
在 [0, T ] 连续，且 φ(t) 是 x
t
的
原函数，则
T
0
x
t
dt = φ(t)
T
0
= φ(T ) − φ(0)
勒贝格积分
法国数学家 Cauchy(1789-1857) 对连续函数定义了积分，之后德国数学家黎曼 (1826-1866)
对有界函数定义了黎曼积分，也称 R 积分。但 R 积分对被积函数有很高的要求，他要求函数 f
几乎处处连续，并且 R 积分的“积分号下取极限”运算要求函数列 f
n
一致收敛于 f。1904 年，
昂利. 勒贝格引入了勒贝格积分的概念，也称为 L 积分，解决了一些 R 积分不能积的情况。一般
书上用 f (x), g(x) 或 u(x) 表示被积函数，自变量为 x，由于我们是在时间 t 上进行研究的，所以
我们用 x
t
或者 x(t) 表示一般被积函数，用 X
t
表示随机过程 (函数)，自变量为 t。L 积分和 R
积分对比示意图如图 (1.1)
(a) L 积分示意图 (b) R 积分示意图
图 1.1: L 积分和 R 积分对比示意图
设 E = [0, T ] ⊂ R 是可测集，且 E 的长度 m(E) 有界，m(E) < ∞，x
t
是 E 上的可测函数
(或者有界函数)，记函数 x
t
在域 E 上的界 x
t
(E) 为 M, N，即 M ⩽ x
t
(E) ⩽ N，下面，我们在
域 E 上介绍函数 x
t
的勒贝格积分。
Step1. 划分 n 段。我们将界 [M, N] 离散化为 n 段，有
M = x
0
< ··· < x
k
< x
k+1
< ··· < x
n
= N
并且记分割方法为 A = {x
k
}
n
，记第 k 段的长度为 ∆x
k
= x
k+1
− x
k
, k = 0, 1, . . . , n − 1，当然，
也可以写为 ∆x
k
= x
k
− x
k−1
, k = 1, 2, . . . , n，不过一定要注意下标别混乱了。记 n 个小段的最
大长度为 δ
δ = max
0⩽k⩽n−1
{∆x
k
}
http://www.ma-xy.com 12 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
为了简便，我们将其简单的均分 n 段。
Step2. 取点 ξ
k
。在第 k 段 [x
k
, x
k
+1
] 中任取一点 ξ
k
∈ [x
k
, x
k
+1
]，则其对应的自变量集为 E
k
，
E
k
= {t|x
k
⩽ x
t
⩽ x
k+1
}，且 E
k
的长度为 m(E
k
)。
Step3. 计算勒贝格积分
T
0
x
t
dt = lim
δ →0
n−1
k=0
ξ
k
m(E
k
)
Step4.L 可积。若在 δ → 0 时，极限
lim
δ
→
0
n−1
k=0
ξ
k
m(E
k
) = L < ∞
存在，且 L 与分割方法 A 无关，与取点 {ξ
k
} 点无关，则称 x
t
在区间 E = [0, T ] 上勒贝格可积，
L 为其黎曼积分。
注：L
p
(Ω) 为 Ω 上 L 可积的函数集；C
∞
0
为 L 可积且光滑的函数集。勒贝格积分的性质如下：
1).R 可积必 L 可积，且二者积分值相等；
2). 线性性和可加性；
3). 绝对值不等式
E
fdm
⩽
E
|f|dm
积分号下取极限
在前面的随机过程中，我们提到过样本路径 X(t, w)(或者写为 X
t
(w)) 是随机过程 X
t
的一
次样本实现。现在，不妨想象我们有 n 次样本实现，于是，我们就有 n 个样本路径 X
n
(t, w)，当
然，我们可以将 n 次样本实现绘成图像。再考虑一个问题，我们“放弃”随机，不妨考虑普通函
数 x
t
，如果用 x
n
(t) 来做函数 x
t
的逼近 (这可以是函数估计，函数逼近和 ODE 等等问题)，我
们会使函数列 {x
n
(t)} 尽可能逼近目标函数 x
t
，如图 (1.2) 所示
图 1.2: 函数列逼近函数示意图
在用函数列 {x
n
(t)} 逼近 x
t
时，我们对 x
t
求积分，如果
lim
n→∞
x
n
(t) = x
t
那么，我们很希望有下面的情况 (积分与极限可以互换位置)
lim
n→∞
T
0
x
n
(t)dt =
T
0
lim
n→∞
x
n
(t)dt =
T
0
x
t
dt
http://www.ma-xy.com 13 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
然而，上面的积分与极限互换位置不是随便就能发生的。在 R 积分中，只有当函数列 {x
n
(t)}
一致收敛于 x
t
时，我们才能将积分与极限互换位置，这对函数列 {x
n
(t)} 提出了很高的要求；在
L 积分中，这种要求相对弱一些，有如下的勒贝格控制收敛定理
定理 (勒贝格控制收敛定理) E 为可测集，m(E) < ∞，{x
n
(t)} 是 E 上有界可测函数列，并
且满足：
1.{x
n
(t)} 在 E 上几乎处处收敛于 x
t
lim
n→∞
x
n
(t) = x
t
a.e于E
2. 存在 E 上 L 可积函数 g(t)，在 E 上有 |x
n
(t)| ⩽ g(t),a.e 于 E。这里 g(t) 是控制函数，可以
是非负可测函数。
则 x
t
在 E 上 L 可积，并且
E
x
t
d
t
=
E
lim
n→∞
x
n
(t)dt = lim
n→∞
E
x
n
(t)dt
用数学语言描述勒贝格控制收敛定理为：L 为可测集 E 上 L 可积函数集 (有界可测)。设
x
n
(t) ∈ L, x
n
(t)
a.e
−−→ x
t
，若 ∃g ∈ L, ∀n ⩾ 1，有 |x
n
(t) ⩽ g|, a.e，则 x
t
∈ L 且
lim
n→∞
E
x
n
(t)dt =
E
x
t
d
t
引理 (Fatou 引理) {x
n
(t)} 为可测集 E 上的非负可测函数列，且 x
n
(t)
a.e
−−→ x
t
，则有
E
x
t
dt ⩽ sup
n⩾1
E
x
n
(t)dt
引理 (Levi 引理) 设 x
1
(t) ⩽ x
2
(t) ⩽ ··· ⩽ x
n
(t) ⩽ . . . 是定义在可测集 E 上的非负可测函
数列，且 x
n
(t)
a.e
−−→ x
t
，则有
E
x
t
dt = lim
n→∞
E
x
n
(t)dt =
E
lim
n→∞
x
n
(t)dt
注：上面的积分中要求 E 有界，x
t
有界，后续可推广到 ∞。上述定理更准确的数学定义参考《测
度论讲义》严加安 P50-P54。
黎曼-斯蒂尔切斯积分
黎曼-斯蒂尔切斯积分 (Riemann-Stieltjes 积分，简称 R-S 积分)，其实在前面试解
T
0
W
t
dW
t
时，我们已经初见 R-S 的身影。R-S 积分是将 R 积分在 Stielejes 下进行推广的，将 R 积分
fdx =
f∆x 推广到
fdg(x) =
f∆g(x)，例如
T
0
x
t
dw
t
=
n−1
k=0
x(τ
k
)[w(t
k+1
) − w(t
k
)]
http://www.ma-xy.com 14 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
其中：x
t
, w
t
是普通函数，而非随机过程。
Step1. 划分 n 段。我们将区间 [0, T ] 离散化为 n 段，有
0 = t
0
< ··· < t
k
< t
k+1
< ··· < t
n
= T
并且记分割方法为 A = {t
k
}
n
，记第 k 段的长度为 ∆x
k
= x
k+1
−x
k
, k = 0, 1, . . . , n −1，记 n 个
小段的最大长度为 δ
δ = max
0⩽k⩽n−1
{∆x
k
}
Step2. 取点 τ
k
。在第 k 段 [t
k
, t
k+1
] 中任取一点 τ
k
∈ [t
k
, t
k+1
]，则该点的函数值为 x(τ
k
)。
Step3. 计算 ∆w(x
k
) 。
∆w(x
k
) = w(t
k+1
) − w(t
k
)
Step4. 计算 R-S 积分
T
0
x
t
dw
t
= lim
δ →0
n−1
k=0
x(τ
k
)∆w(x
k
)
Step5.R-S 可积。若在 δ → 0 时，极限
lim
δ →0
n−1
k=0
x(τ
k
)∆w(x
k
) = R < ∞
存在，且极限值 R 与分割方法 A 无关，与 {τ
k
} 无关，则称 x
t
关于 w
t
在区间 [0, T ] 上的 R-S
积分存在且 R 为其积分值。
R-S 积分存在的充分条件为：函数 x
t
连续，w
t
单调。上面定义的积分界 [0, T ] 有界，我们
可以将它推广到无限区间 [0, ∞] 上
∞
0
x
t
dw
t
= lim
T →∞
T
0
x
t
dw
t
此外，R-S 积分也有一些性质，这里不做讨论。
Ito 积分
前面讨论的 R 积分、L 积分和 R-S 积分中的函数 x
t
, w
t
都是普通函数，而非随机过程。下
面将函数引入到随机过程当中，仍然用大写的
X
t
, W
t
来表示随机过程。考虑如下积分形式
T
0
X
t
dW
t
¬先假设 X
t
为普通函数 x
t
，W
t
为 Brown 运动，且 x
t
不依赖与 W
t
，于是，经过划分、取
点、求积等过程，我们可以定义其积分形式为
T
0
x
t
dW
t
=
n−1
k=0
x(τ
k
)[W (t
k+1
) − W (t
k
)]
http://www.ma-xy.com 15 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
可以定义 ∆W
k
= W (t
k+1
) − W (t
k
)。由于 Brown 运动的正态增量的性质 W (t) − W (s) ∼
N(0, t − s)，所以上式定义的积分是一个服从高斯分布的随机变量，且其均值为
E
T
0
x
t
dW
t
= E
n−1
k=0
x(τ
k
)[W (t
k+1
) − W (t
k
)]
= 0
其方差为
V ar
T
0
x
t
dW
t
= V ar
n−1
k=0
x(τ
k
)[W (t
k+1
) − W (t
k
)]
=
n−1
k=0
x(τ
k
)
2
(t
k+1
− t
k
)
先将 X
t
视为简单的随机过程 (例如：X
t
在各个时间段 [t
k
, t
k+1
] 为随机变量 ξ
k
，并且各
个时刻之间的随机变量毫不相关)，经过划分、取点、求积等过程，我们可以定义其积分形式为
T
0
X
t
dW
t
=
n−1
k=0
X(τ
k
)[W (t
k+1
) − W (t
k
)]
因为 X
t
是一个简单的随机过程，则 X(τ
k
) 可以记为一个简单的随机变量 ξ
k
，ξ
k
依赖于 W
t
(t ⩽ t
k
)，
但不依赖于 W
t
(t > t
k
)，由此，其积分形式还可以写为
T
0
X
t
dW
t
=
n−1
k
=0
ξ
k
[
W
(
t
k+1
)
−
W
(
t
k
)]
上面定义的这个积分具有以下性质：
1. 线性性。设 X
t
, Y
t
为简单的随机过程，有
T
0
[αX
t
+ βY
t
]dW
t
= α
T
0
X
t
dW
t
+ β
T
0
Y
t
dW
t
2. 零均值性。如果 E[ξ
2
i
] < ∞, (k = 0, 1, . . . , n − 1)，则
E
T
0
X
t
dW
t
= 0
3. 等距性。如果 E[ξ
2
i
] < ∞, (k = 0, 1, . . . , n − 1)，则
E
T
0
X
t
dW
t
2
=
T
0
E
X
2
t
dt
®现在将 X
t
扩展到更广泛的随机过程 (注意：X
t
即为随机微分方程的解)。为了方便讨论，
我们给 X
t
规定一个随机过程类 V ，在定义随机过程类之前，我们先给出适应过程的定义：
http://www.ma-xy.com 16 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
定义 (适应过程) 设 {X
t
} 是一个随机过程，{F
t
} 是 σ 代数流。若对 ∀t，x
t
是 F
t
可测的，
则称 {X
t
} 是 {F
t
} 适应过程。
下面，我们来定义一个随机过程类 V
V =
X
t
X
t
是[0, T ]上可测适应过程，满足E
T
0
X
2
t
< ∞
进一步，定义
V
∗
=
X
t
X
t
是[0, T ]上可测适应过程，且∀T > 0, 满足
T
0
X
2
t
< ∞, a.s
注：a.s 表示几乎必然收敛 (以概率 1 收敛)，英文为 almost sure。
定义 (Ito 积分) 设 X
t
∈ V (0, T )，X
t
的积分可以定义为
T
0
X(t, w)dW
t
= lim
n→∞
T
0
ϕ
n
(t, w)dW
t
∈ L
2
(P )
其中：L
2
(P ) 表示 P 上平方 L 可积函数集；{ϕ
n
(t)} 是初等随机过程序列，且当 n → ∞ 时，有
E
T
0
[X(t, w) − ϕ
n
(t, w)]
2
dt
→ 0
上面定义的积分形式即为 Ito 积分。
关于 Ito 积分中的初等随机序列 {ϕ
n
} 的构造，可以查阅《随机微分方程及其在数理金融中
的应用》蒲兴成 P16-P17。在实际问题中，我们遇到的过程常常不在随机过程类 V 中，为此，我
们可以将 Ito 积分推广到 V
∗
函数类。
Ito 积分的性质
上面给出了 Ito 积分的定义，下面，我们来介绍 Ito 积分的一些性质：1. 可加性；2. 线性性；
3. 对任意固定的 T 0，随机过程 X
t
的 Ito 积分
T
0
X
t
dW
t
是一个随机变量 (其实，不仅仅是 Ito
积分，任何带有随机性 dW
t
的积分都是一个随机变量)；4. 零均值性。
E
T
0
X
t
dW
t
= 0
引理 (Ito 引理) ∀X
t
∈ V (0, T )，有
E
T
0
X(t, w)dW
t
2
= E
T
0
X(t, w)
2
dt
下面来看一下随机过程的逼近情况：若 X(t, w) ∈ V, X
n
(t, w) ∈ V ，且当 n → ∞ 时，有
E
T
0
(X
n
(t, w) − X(t, w))
2
dt
→ 0
http://www.ma-xy.com 17 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
则有
T
0
X
n
(t, w)dW
t
→
T
0
X(t, w)dW
t
∈ L
2
(P )
在上面的性质中，我们知道，对任意固定的 T 0，随机过程 X
t
的 Ito 积分
T
0
X
t
dW
t
是一个
随机变量，由随机变量，我们就可以再构造一个随机过程。我们定义一个新的随机过程
Y
t
=
t
0
X
s
dW
s
其中：t ∈ [0, T ]，s ∈ [0, t]。我们称随机过程 Y
t
是 Ito 积分过程，可以发现，前面的随机微分方程
X
t
= X
0
+
t
0
f(s, X
s
)ds +
t
0
g(s, X
s
)dW
s
的解 X
t
其实就类似一个 Ito 积分过程。可以证明，Ito 积分过程 Y
t
存在连续的样本路径，下面，
我们来讨论一下随机过程 Y
t
的随机性质。
定理 (鞅性质) 设 X
t
∈ V
∗
，且
∞
0
E[X
2
t
]dt < ∞，则
Y
t
=
t
0
X
s
dW
s
0 ⩽ t ⩽ T
是一个零均值的连续平方可积鞅 (过程)。
定义 (Gauss 过程) 如果 X
t
非随机，且
T
0
X
2
t
dt < ∞，则对 ∀t
Y
t
=
t
0
X
s
dW
s
是一个服从正态分布的随机变量。如果 t ∈ [0, T ]，，则随机过程 {Y
t
} 是一个 Gauss 过程，其均
值函数为常值 0，协方差函数为
Cov[Y
t
, Y
t+u
] =
t
0
X
2
s
ds, u ⩾ 0
1.3.5 Ito 过程与 Ito 公式
Ito 过程
在前面随机积分的初步尝试部分，我们讨论了 Brown 运动的二次变差 (1.9)。Brown 运动
W
t
在时间区间 [0, t] 上的二次变差为 t，即
[W, W ](t) = lim
δ
n
→0
n−1
k=0
|W (t
k+1
) − W (t
k
)|
2
= t
其中：{t
k
} 为 [0, t] 上的 n 段分割，δ
n
= max
0⩽k⩽n−1
{t
k+1
− t
k
}。上式从形式上可以表示为
t
0
[dW
s
]
2
=
t
0
ds = t
http://www.ma-xy.com 18 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
或者
[dW
t
]
2
= dt
更常见的写法是
dW
t
≈
√
dt (1.10)
更一般的，我们有下面定理
定理 设 x 是有界连续函数，{t
k
} 为 [0, t] 上的 n 段分割，则 ∀θ
k
∈ [W (t
k+1
), W (t
k
)]，依概
率收敛意义下的极限
lim
δ
n
→0
n−1
k=0
x(θ
k
)[W (t
k+1
) − W (t
k
)]
2
=
t
0
x[W
s
]ds
前面，我们简单的介绍了 Ito 积分过程，下面，我们给出 Ito 过程的定义：
定义 (Ito 过程) 设 W
t
表示一维 Brown 运动，如果一个一维的过程 X
t
可以表示成一个随
机积分
X
t
= X
0
+
t
0
µ(s)ds +
t
0
σ(s)dW
s
或者随机微分方程
dX
t
= µ(t)dt + σ(t)dW
t
X(t
0
) = X
0
并且过程 µ(t), σ(t) 满足：
1.µ(t) 是适应的，并且满足
T
0
|µ(t)|dt < ∞, a.s
2. σ(t) ∈ V
∗
；
则称 X
t
是一个 Ito 过程。
由 Ito 过程的定义可以看出，前面要求解的人口增长随机微分方程的解 X
t
是一个 Ito 过程。
值得一提的是，并不是所有的随机微分方程都是
Ito
过程，对许许多多的
SDE
，我们都无法求出
其解析解。
Ito 公式
上面的内容介绍了 Ito 积分的定义、Ito 积分的性质并介绍了 Ito 过程。但是如果我们仅根
据定义来计算 Ito 积分的话，那就太麻烦了，那么 Ito 积分有没有积分公式呢？像 R 积分中有牛
顿—莱布尼茨共识那样，有的。下面，我们就来介绍 Ito 积分的计算公式—Ito 公式 (非线性链式
法则)
http://www.ma-xy.com 19 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
定理 (Ito 公式) 设 X
t
是一个 Ito 过程
dX
t
= µ(t)dt + σ(t)dW
t
g(t, x) ∈ C
2
[R
+
× R](即 g 是 [0, +∞) × R 上的二次连续可微函数)，则
Y
t
= g(t, X
t
)
仍是一个 Ito 过程，且
dY
t
=
∂g
∂t
(t, X
t
)dt +
∂g
∂x
(t, X
t
)dX
t
+
1
2
∂
2
g
∂x
2
(t, X
t
) · d(X
t
)
2
(1.11)
其中：d(X
t
)
2
= d(X
t
) · d(X
t
) 是按照如下规则计算的
dt · dt = dt · dW
t
= dW
t
· dt = 0
dW
t
· dW
t
= dt
(1.12)
证明 利用式 (1.12) 有
g(t, X
t
) = g(0 + X
0
) +
t
0
∂g
∂s
(s, X
s
) + µ
s
∂g
∂x
(s, X
s
) +
1
2
σ
2
s
·
∂
2
g
∂x
2
(s, X
s
)
ds
+
t
0
σ
s
·
∂g
∂x
(s, X
s
)dW
s
其中：µ
s
= µ(s, w), σ
s
= σ(s, w)。我们可以知道 g(t, X
t
) 是一个随机过程。为简单，设 g,
∂g
∂t
,
∂g
∂x
,
∂
2
g
∂x
2
是有界的，能够证明，在此情况下可以在 C
2
中构建一个函数列 g
n
，使得当 n → ∞ 时，g
n
,
∂g
n
∂t
,
∂g
n
∂x
,
∂
2
g
n
∂x
2
在 [0, ∞) × R 的紧子集中分别一致趋近于 g,
∂g
∂t
,
∂g
∂x
,
∂
2
g
∂x
2
。此外，我们可以假定 µ(t, w), v(t, w) 是
初等函数。利用泰勒展开，有
g(t, X
t
) = g(0, X
0
) +
j
∆g(t
j
, X
j
)
= g(0, X
0
) +
j
∂g
∂t
∆t
j
+
j
∂g
∂x
∆X
j
+
1
2
j
∂
2
g
∂t
2
(∆t
j
)
2
+
j
∂
2
g
∂t∂x
(∆t
j
)(∆X
j
) +
1
2
∂
2
g
∂x
2
(∆X
j
)
2
+
R
j
其中：
∆t
j
= t
j+1
− t
j
, ∆X
j
= X
t
j+1
− X
t
j
,
∆g(t
j
, X
j
) = g(t
j+1
, X
t
j+1
) − g(t
j
, X
t
j
), R
j
= O(|∆t
j
|
2
+ |∆X
j
|
2
)
则当 ∆t
j
→ 0 时，有
j
∂g
∂t
∆t
j
=
j
∂g
∂t
(t
j
, X
j
)∆t
j
→
t
0
∂g
∂t
(s, X
s
)ds
j
∂g
∂x
∆X
j
=
j
∂g
∂x
(t
j
, X
j
)∆X
j
→
t
0
∂g
∂x
(s, X
s
)dX
s
http://www.ma-xy.com 20 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
此外，因为 µ, σ 是初等函数，所以有
j
∂
2
g
∂x
2
(∆X
j
)
2
=
j
∂
2
g
∂x
2
µ
2
j
(∆t
j
)
2
+ 2
j
∂
2
g
∂x
2
µ
j
σ
j
(∆t
j
)(∆W
j
)
+
j
∂
2
g
∂x
2
σ
2
j
· (∆W
j
)
2
(1.13)
其中：µ
j
= u(t
j
, w), σ
j
= σ(t
j
, w)。易证在 ∆t
j
→ 0 时，上式 (1.13) 中的前两项都趋近于 0，例
如
lim
∆t
j
→0
E
j
∂
2
g
∂x
2
µ
j
σ
j
(∆t
j
)(∆W
j
)
2
= lim
∆t
j
→0
E
j
∂
2
g
∂x
2
µ
j
σ
j
2
(∆t
j
)
3
= 0
下面证明，当 ∆t
j
→ 0 时，在 L
2
(P ) 中式 (1.13) 的最后一项趋于
t
0
∂
2
g
∂x
2
σ
2
ds。为证明这一点，记
a(t) =
∂
2
g
∂x
2
(t, X
t
)σ
2
(t, w), a
j
= a(t
j
)
因为
E
j
a
j
(∆W
j
)
2
−
j
a
j
∆t
j
2
=
i,j
E[a
i
a
j
((∆W
j
)
2
− ∆t
j
)]
且当 i < j 时，a
i
a
j
((∆W
j
)
2
− ∆t
j
) 与 (∆W
j
)
2
− ∆t
j
是独立的，因而此时这项可以忽略。类似
的，当 i > j 时，有
lim
∆t
j
→0
j
E[(a
j
)
2
((∆W
j
)
2
− ∆t
j
)
2
]
= lim
∆t
j
→0
j
E[a
2
j
] · E[(∆W
j
)
4
− 2(∆W
j
)
2
∆t
j
+ (∆t
j
)
2
]
= lim
∆t
j
→0
j
E[a
2
j
] · [3(∆t
j
)
2
− 2(∆W
j
)
2
∆t
j
+ (∆t
j
)
2
]
= lim
∆t
j
→0
j
E[a
2
j
] · (∆t
j
)
2
= 0
即在 L
2
(P ) 中，有
lim
∆t
j
→0
E
j
a
j
(∆W
j
)
2
=
t
0
a(s)ds
上式通常写为 (dW
t
)
2
= dt。类似上面的讨论，同样可以证明：lim
∆t
j
→0
R
j
= 0。综上所述，Ito
公式成立。□
http://www.ma-xy.com 21 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
上面介绍了 Ito 公式，下面我们就用 Ito 公式来求解前面建立的人口增长的随机微分方程模
型 (1.4)，由
dX
t
X
t
= rdt + dW
t
对上式积分，有
t
0
dX
s
X
s
= rt + W
t
(W
0
= 0)
为计算左边积分，对函数 g(t, x) = ln x(x > 0) 运用 Ito 公式，有
d(ln X
t
) =
1
X
t
· dX
t
+
1
2
−
1
X
2
t
(dX
t
)
2
=
dX
t
X
t
−
1
2X
t
· X
2
t
dt
=
dX
t
X
t
−
1
2
dt
所以
dX
t
X
t
= d(ln X
t
) +
1
2
dt
代入所求方程，有
ln
X
t
X
0
=
r −
1
2
t + W
t
或者
X
t
= X
0
exp
r −
1
2
t + W
t
1.3.6 解的存在唯一性及解析计算公式
解的存在唯一性
现在，回到标准随机微分方程模型当中 (1.5)，对于一个 SDE 问题，我们基本的目标就是求
解方程的解析解或者数值解，在有解之后再来了解一下解的性质。对于 SDE 的解，我们已经知
道其解是一个随机过程。在对前面 ODE 和 PDE 解的分析中，我们讨论了解的定义、解的存在
唯一性、解的解析求法、解的数值求法和解的稳定性，那么，对 SDE 的解，我们仍旧要讨论上
面的问题 (前面我们已经给出 SDE 解的定义)。首先，我们来讨论解的存在唯一性，我们并不打
算直接在 SDE 一般形式下讨论
dX
t
= f (X
t
, W
t
, t)
那样太过于泛化。我们来讨论 Ito 型随机微分方程的解的存在唯一性 (Stratonovich 型随机微分
方程等其他形式这里不做讨论)。前面已经介绍过 Ito 型随机微分方程 (1.7)，这里再重写一下
dX
t
= f (X
t
, t)dt + g(X
t
, t)dW
t
t ∈ [t
0
, T ]
X
t
0
= X
0
http://www.ma-xy.com 22 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.3 随机微分方程基本理论
Ito 型随机积分方程为
X
t
= X
0
+
t
0
f(s, X
s
)ds +
t
0
g(s, X
s
)dW
s
下面规范一下上式的 1. 量空间、2. 映射空间、3. 解域。t ∈ [0, T ]；X
0
∈ R 可以是常量，也
可以是一随机变量，但要求其平方可积 E[|X
0
|
2
] < ∞；W
t
∈ R 是一个 Brown 运动。既然有随
机，我们定义一个概率空间 (Ω, F, P )，W
t
是概率空间 (Ω, F, P ) 上的 Brown 运动，{F
t
} 是 W
t
产生的 σ 代数流；f(X
t
, t) : R × [0, T ] → R；g(X
t
, t) : R × [0, T ] → R。
上面只给出了映射 f, g 的量变换，并没有给出 f, g 的类型，我们要求 f, g 是 [0, T ] 上的
(Borel) 连续可测函数，并且要求 f, g 满足 Ito 过程的要求，即
1).f(X
t
, t) 是适应的并且
T
0
|f(X
t
, t)|dt < ∞, a.s；
2).g(X
t
, t) ∈ V
∗
。
上述 Ito 型随机微分方程可能存在两种形式的解：¬给定 f, g，解 (过程)X
t
依赖于时间 t 与
W
t
在
[0
, t
]
的值有关，显然
X
t
是强吻合上述方程的，故称
X
t
为强解；
给定
f, g
，
X
t
不依赖
与时间 t，仅与 W
t
的前值有关，这样的解 X
t
称为弱解，不妨记为
˜
X
t
，
˜
X
t
并不完全吻合随机
方程。下面，我们给出 Ito 型随机微分方程强解的存在唯一性定理：
定理 (强解的存在唯一性) 设 f (t, x), g(t, x) 是 [0, T ]×R 上的 (Borel) 可测函数且 |f(t, x)|
1
2
, |g(t, x)|
平方可积。如果 f, g 满足：
1.(Lipschitz 条件)∃L > 0，∀t ∈ [0.T ]，使得 ∀x, y ∈ R，有
|f(t, x) − f (t, y)| ⩽ L|x − y|
|g(t, x) − g(t, y)| ⩽ L|x − y|
2.(线性增长条件)∃K > 0，使得 ∀t ∈ [0, T ], x ∈ R，有
|f(t, x)| + |g(t, x)| ⩽ K(1 + |x|)
并且，初始值 X
0
独立于由 {W
s
, s ⩾ 0} 产生的 σ 代数 {F
t
}，且 E[|X
0
|
2
] < ∞，则 Ito 型随机
微分方程存在唯一解 X
t
，并且解的路径 X
t
(w) 连续，还有 X
t
适应于由 X
0
和 {W
s
, s ⩽ t} 产生
的流域 {F
t
}，并且 E
T
0
|X
t
|
2
dt < ∞。
注：上述证明可以参考《Stochastic Dierential Equations》B.Qksendal；或者《随机微分方程及
其在数理金融中的应用》蒲兴成 P37-P40。
解析计算公式
在 SDE 存在唯一解的情况下，针对某一特殊形式的 SDE，我们有没有相应的解析解计算公
式呢，就像 ODE 中针对不同类型的微分方程会有相应的计算公式那样？有的，我们介绍一类标
量线性 Ito 型 SDE 的解析计算公式
http://www.ma-xy.com 23 http://www.ma-xy.com
http://www.ma-xy.com
1.3 随机微分方程基本理论 第一章 随机微分方程
¬标量线性 Ito 型随机微分方程
dX
t
= (A
t
X
t
+ a
t
)dt +
m
i=1
[B
i
(t) + b
i
(t)]dW
i
(t)
X(t
0
) = X
0
t ∈ [t
0
, T ]
其解析解为
X
t
= Φ
t
X
0
+
T
t
0
Φ
−1
(s)
a(s) +
m
i=1
B
i
(s)b
i
(s)
ds
+
m
i=1
T
t
0
Φ
−1
(s)b
i
(s)dW
i
(s)
其中：
Φ
t
= exp
T
t
0
A(s) −
m
i=1
B
2
i
(s)/2ds +
m
i=1
T
t
0
B
i
(s)W
i
(s)
线性齐次 SDE
dX
t
= A
t
X
t
dt +
m
i=1
B
i
(t)X
t
dW
i
(t)
X(t
0
) = X
0
t ∈ [t
0
, T ]
其解析解为
X
t
= exp
T
t
0
A(s) −
m
i=1
B
2
i
(s)/2 +
m
i=1
T
t
0
B
i
(s)dW
i
(s)
®自治标量线性 SDE
dX
t
= (AX
t
+ a)dt +
m
i=1
(B
i
X
t
+ b
i
)dW
i
(t)
X(t
0
) = X
0
t ∈ [t
0
, T ]
其解析解为
X
t
= Φ
t
X
0
+
T
t
0
Φ
−1
(s)
a −
m
i=1
B
i
b
i
ds +
m
i=1
T
t
0
Φ
−1
(s)b
i
dW
i
(t)
其中：
Φ
t
= exp
A −
m
i=1
B
2
i
/2
(T − t
0
) +
m
i=1
B
i
(W
i
(T ) − W
i
(t
0
))
关于随机微分方程的解析理论，才刚刚开始，其余内容留日后补充。下面，我们将进入到
SDE 的数值解法。
http://www.ma-xy.com 24 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.4 随机微分方程数值方法
1.4 随机微分方程数值方法
1.4.1 随机泰勒展开
跋山涉水走到了这里。下面，我们来讨论随机微分方程的数值解法，对数值方法，我们的目
标仍然是给定 X
n
后如何求解 X
n+1
，与 ODE 不同的是，这里的 X
n
和 X
n+1
是随机变量。ODE
解法中的欧拉方法和 R-K 方法等仍然可以推广到 SDE 当中，有限元方法和谱方法也可以进行尝
试。为简便，我们在下面的自治型 SDE 进行讨论
dX
t
= f (X
t
)dt + g(X
t
)dW
t
X(t
0
) = X
0
∈ R
t ∈ [t
0
, T ]
(1.14)
SDE 常见的数值方法：Euler 法、Milstein 法和 R-K 方法等都是基于随机泰勒展开的，下面
我们先来看一下随机泰勒展开。将时间区间 [t
0
, T ] 划分为 N 段，为简便，采用 N 等均分，设各
段长度为 ∆t = h =
T −t0
N
，记第 n 段为 [t
n
, t
n+1
]，t
n
= t
0
+ nh，我们的问题是：给出 X(t
n
) 之
后，如何求解 X(t
n+1
)？对方程 (1.14) 由随机泰勒公式和 Ito 公式，我们有
X(t
n+1
) = X(t
n
) + I
0
f(X(t
n
)) + I
1
g(X(t
n
))
+ I
11
L
1
g(X(t
n
)) + I
00
L
0
f(X(t
n
)) + R
其中：R 为余项，算子 L
0
, L
1
分别为
L
0
= f (x)
∂
∂x
+
1
2
g
2
(x)
∂
2
∂x
2
L
1
= g(x)
∂
∂x
I
0
= h, I
1
= ∆W
n
, I
00
=
1
2
h
2
, I
11
=
1
2
[∆W
2
n
− h]。
上面的展开式也可以写为
X
(
t
n+1
) =
X
(
t
n
) +
hf
(
X
(
t
n
)) + ∆
W
n
g
(
X
(
t
n
)) +
1
2
[∆W
2
n
− h]g(X(t
n
))g
′
(X(t
n
))
+
1
2
h
2
f(X(t
n
))f
′
(X(t
n
)) +
1
2
g
2
f
′′
(X(t
n
))
+ R (1.15)
Euler 法和 Milstein 法皆是在 (1.15) 式上截断得到的。我们先给出显式 Eiler 方法，并用其导出
数值解收敛与稳定的概念，之后再介绍一些其它的数值方法，最后给出
MATLAB
仿真结果。
1.4.2 Euler 方法
方程 (1.14) 的显式 Euler 方法的迭代公式为
X
n+1
= X
n
+ f(X
n
)∆t + g(X
n
)∆W
n
http://www.ma-xy.com 25 http://www.ma-xy.com
http://www.ma-xy.com
1.4 随机微分方程数值方法 第一章 随机微分方程
其中：∆t = h 为离散步长，∆W
n
是第 n 步的 Brown 运动增量。可以发现
X
n+1
− X
n
= f (X
n
)∆t + g(X
n
)∆W
n
和
dX
t
= f (X
t
)dt + g(X
t
)dW
t
是近似相同的。
对于上面的 Euler 方法，你肯定会问 ∆W
n
是如何计算的？其实，我们应该写的详细一点
∆W
n
:= W
t
n+1
− W
t
n
由 (1.12)dW
t
· dW
t
= dt，或者 (1.10)dW
t
≈
√
dt，我们知道
∆W
n
≈
t
n+1
− t
n
=
h
n
注意：这里我们对 [t
0
, T ] 用了 N 等均分，所以对 ∀i, j ∈ N，h
i
= h
j
= h。
√
h
n
不具备随机性，并
且由于 Brown 运动的正态离差性，∆W
n
∼ N (0, ∆t
n
)，即 ∆W
n
∼
√
∆t
n
N(0, 1)，记 ξ ∼ N (0, 1)，
所以，我们用
∆W
n
=
h
n
ξ
来模拟 Brown 运动增量 ∆W
n
，则 Brown 运动为 W
t
=
t
n=1
∆W
n
。
我们对如下的 SDE 应用 Euler 方法进行仿真模拟
dX
t
= λX
t
dt + µX
t
dW
t
X(t
0
) = X
0
t ∈ [0, T ]
上述方程是经济学中的 Black-Scholes 方程，和人口增长的随机微分方程相似，是期权定价模型，
可以用来表示股票价格变化。我们可以看到其 f, g 满足全局 Lipschitz 条件和线性增长条件，故
其存在唯一解。用 Ito 公式可求该 SDE 的解析解为
X
t
= X
0
exp
λ −
1
2
µ
2
t + µW
t
令方程中的参数为：λ = 2，µ = 1，t
0
= 0，X
0
= 1，T = 1，N = 2
8
，h = ∆t =
T −t
0
N
。其仿真
程序为
1 %% 此 程 序 用 于SDE章 节 Euler 数 值 方 法 的 展 示
2 %−−−−−−−−−−−−−−−−−−−−−−−−−−要 求 解 的 方 程 为−−−−−−−−−−−−−−−−−−−−−−−%
3 %%%% dWt = lambda Xt dt + mu Xt dWt
4 %%%% X(0 ) = X0
5 %%%% 0<=t<=T
6 % 其 精 确 解 Xt = X0 exp ( ( lambda − 1/2 mu^2) t + mu Wt)
7 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
8 %% 精 确 解
http://www.ma-xy.com 26 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.4 随机微分方程数值方法
9 randn( ’ s t a t e ’ ,100)
10 lambda = 2 ;
11 mu = 1 ;
12 X0 = 1;
13 t0 = 0;
14 T = 1 ;
15 N = 2^8;
16 h = (T − t0 ) /N;
17 dt = h ;
18 %
19 dW = sq r t ( dt ) *randn (1 ,N) ;%Brown运 动 增 量dW = s q r t (h)* xi , xi 服 从N(0 , 1 )
20 Wt = cumsum(dW) ;%Brown模 拟 Wt
21 % 解 析 解
22 Xt = X0 * exp ( ( lambda − 0.5*mu^2) * ( [ dt : dt :T] ) + mu*Wt) ;
23 % 画 图
24 f i g u r e
25 p lo t ( [ 0 : dt :T] , [ X0, Xt ] , ’m−’ )
26 hold on
27 %% 数 值 解
28 % Xn+1 = Xn+lambda Xn del ta t + mu Xn del ta Wn
29 R = 4 ;
30 Dt = R*dt ;
31 L = N/R; % 设 置 数 值 离 散 化 步 长 ， 由 于 dt 太 小 ， 故 将 其 放 大R倍
32 Xt_hat = zero s (1 ,L) ; % Xt的 估 计 值
33 Xt_hat ( 1) = X0 ;
34 f o r n = 1:L
35 delta_W = sum(dW(R*(n−1)+1:R*n) ) ; % 计 算Brown运 动 增 量
36 Xt_hat(n+1) = Xt_hat(n) + lambda * Xt_hat(n) *Dt + mu *Xt_hat(n) *delta_W ;
37 end
38 % 绘 图
39 p lo t ( [ 0 : Dt :T] , Xt_hat , ’ r−−*’ )
40 hold o f f
41 xl a b el ({ ’ $t$ ’ } , ’ I n t e r p r e t e r ’ , ’ lat e x ’ )
42 yl a b el ({ ’ $x$ ’ } , ’ I n t e r p r e t e r ’ , ’ l a t e x ’ )
43 legend ( ’ 解 析 解 ’ , ’ 数 值 解 ’ , ’ lo c a t i o n ’ , ’ b e s t ’ )
44
上述程序运行结果如图 (1.3) 所示
http://www.ma-xy.com 27 http://www.ma-xy.com
http://www.ma-xy.com
1.4 随机微分方程数值方法 第一章 随机微分方程
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
t
0
1
2
3
4
5
6
x
图 1.3: SDE 的 Euler 方法模拟图
1.4.3 解的收敛性与稳定性
收敛性
下面，我们来讨论一下数值解的收敛性。对某一点 t
n
而言，如果 X
t
n
与 X
n
并非随机变量，
而是一个变量，我们自然希望 X
n
能够接近 X
t
n
，即估计量与真实值靠近 (或者说数值解和解析
解靠近)。但是在 SDE 中，估计量 X
n
和真实量 X
t
n
都是随机变量，这使得我们要重新找方法来
衡量二者之间的接近度。
假设我们可以重复计算 N 次 SDE 的解析解与数值解 (注意，这里的 N 是对随机过程模拟
N 次)，即对某一个 t
n
而言，进行 N 次采样 X
t
n
, X
n
，那么，就某时刻随机变量 X
n
和 X
t
n
的
接近程度可以用
E[|X
n
− X
t
n
|] ≈
1
N
|X
n
− X
t
n
|
来衡量 (当然，也可以采用其它方法)。如果我们要求在每个时间点 t
n
上 E[|X
n
− X
t
n
|] 都非常
小，那么就可以说数值解 X
n
是强收敛的。下面，给出 SDE 数值方法的收敛性和稳定性的概念。
定义 (强收敛) 在 [t
0
, T ] 上，令 h =
T −t
0
N
，这个 N 为分割段数，h 为步长。若 ∃c > 0，c 与
h 无关，∃δ > 0，h ∈ (0, δ)，使得 ∀t
n
= t
0
+ nh ∈ [t
0
, T ]，有
E[|X
n
− X
t
n
|] ⩽ ch
r
其中：X
t
n
为 t
n
时刻的解析解，X
n
为 t
n
时刻的数值解。则称该数值方法是 r 阶强收敛的。
定义 (局部收敛阶) 若 ∃c > 0，c 与 h 无关，令 δ
n
= X
n
− X
t
n
为局部截断残差，有
max
0⩽n⩽N
|E(δ
n
)| ⩽ ch
p
1
h → 0
以及
max
0⩽n⩽N
|E(δ
n
)
2
|
1
2
⩽ ch
p
2
h → 0
并且 p
1
⩾
1
2
, p
1
⩾ p
2
+
1
2
，则称 p
2
为该数值方法均方意义上的局部收敛阶。
http://www.ma-xy.com 28 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.4 随机微分方程数值方法
定义 (弱收敛) 若对适当的 2(p + 1) 次多项式 φ，∃c > 0, δ > 0, h ∈ (0, δ)，使得 ∀t
n
=
t0 + nh ∈ [t
0
, T ]，有
E[φ(X
n
)] − E[φ(X(t
n
))]
⩽ ch
p
则称该数值方法为 p 阶弱收敛。
稳定性
关于解 X
t
的稳定性，可以参考 ODE 中的李雅普诺夫稳定。强 (弱) 收敛只是在有限时间
[t
0
, T ] 内讨论 X
t
的收敛效果，而当 t → ∞ 时，SDE 的性质也是非常值得注意的。我们知道解
X
t
对初值是连续的，现在考虑初值微小波动对解的影响。
平凡解稳定性 称 X(t; X
0
, t
0
) ≡ 0 是方程 dX
t
= f (X
t
)dt + g(X
t
)dW
t
的平凡解。1980 年
Has’mrnshiz 给出了平凡解稳定性的定义：
¬若 ∀ε > 0, t
0
> 0，有
lim
X
0
→0
P
sup
t⩾t
0
|X(t; X
0
, t
0
| ⩾ ε
= 0
则称 X(t; t
0
, X
0
) ≡ 0 是随机稳定的。
若在随机稳定的前提下，有
lim
X
0
→0
P
lim
t→∞
|X(t; X
0
, t
0
| → 0
= 1
则称 X(t; t
0
, X
0
) ≡ 0 是随机渐进稳定的。
®若在随机渐进稳定的前提下，对 ∀X
0
∈ R，有
P
lim
t→∞
|X(t; X
0
, t
0
| → 0
= 1
则称 X(t; t
0
, X
0
) ≡ 0 是大范围随机渐进稳定的。
p 阶矩稳定性 1992 年，Kloeden 和 Dleten 给出了 p 阶矩稳定的定义。若 ∃ε > 0, ∃δ = δ(ε, t
0
) >
0，使得 ∀t ⩾ 0, |X
0
| < δ，有
E (|X(t; X
0
, t
0
)|
p
) < ε
即 ∃δ
0
= δ(t
0
) > 0，使得 ∀|X
0
| < δ，有
lim
t→∞
E (|X(t; X
0
, t
0
)|
p
) = 0
则称 X(t; t
0
, X
0
) ≡ 0 是 p 阶矩渐进稳定的。特别地，当 p = 2 时，称为均方稳定 (MS 稳定)。
http://www.ma-xy.com 29 http://www.ma-xy.com
http://www.ma-xy.com
1.4 随机微分方程数值方法 第一章 随机微分方程
MS 稳定区域 对于任意一种数值方法，若迭代方程写为
X
n+1
= ϕ(h, λ, µ, J
n
)X
n
= ϕ(h, λ, µ,
√
hJ)X
n
(1.16)
其中：
J
n
=
W
(
t
n+1
)
−
W
(
t
n
)
∼
N
(0
, h
)
J ∼ N(0, 1)
对上式 (1.16) 两边取 2 阶原点矩，有
E
|X
n+1
|
2
= E
ϕ
2
+ E
|X
n
|
2
我们称 R
1
(h, λ, µ) = E
ϕ
2
(h, λ, µ,
√
hJ)
为均方稳定函数。若 |R
1
(h, λ, µ)| < 1，则该数值方法
为均方稳定 (MS 稳定)。令 p = hλ, q =
√
h，则迭代方程 (1.16) 也可以写为
X
n+1
= R
2
(p, q)X
n
称 R
2
(p, q) 为均方稳定函数。若 |R
2
(p, q)| < 1，则该数值方法 MS 稳定，并称 S = {(p, q)||R
2
(p, q)| <
1} 为该数值方法的 MS 稳定区域。
T 稳定 1993 年，Saito 和 Mitsui 提出 T 稳定。对于一般的数值迭代公式 (1.16)
X
n+1
= ϕ(h, λ, µ, J
n
)X
n
= ϕ(h, λ, µ,
√
hJ)X
n
令
¯
X
n
=
n
n
i=1
X
i
¯
X
n+1
=
n
n
i=1
X
i+1
进一步得到一步差分方程
¯
X
n+1
= R
T
(h, λ, µ)
¯
X
n
其中：
R
T
(h, λ, µ) =
n
n
i=1
ϕ(h, λ, µ, J
i
)
称 R
T
(h, λ, µ) 为 T 稳定函数。
http://www.ma-xy.com 30 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.4 随机微分方程数值方法
1.4.4 数值方法
上面，我们利用 Euler 方法引出了数值解收敛性和稳定性的概念，下面再来介绍一些数值方
法，并给出这些数值方法收敛性和稳定性的数值验证。考虑如下 SDE 问题
dX
t
= f (X
t
)dt + g(X
t
)dW
t
1. 显示数值方法
1). 显示 Euler-Maruyama 方法
X
n+1
= X
n
+ f(X
n
)∆t + g(X
n
)∆W
n
该方法是强 0.5 阶收敛，弱 1 阶收敛的。
2). 显示 Milstein 方法
X
n+1
= X
n
+ f(X
n
)∆t + g(X
n
)∆W
n
+
1
2
g(X
n
)g
′
(X
n
)(∆W
2
− ∆t)
该方法是强 1 阶收敛的。
3).
显示
R-K
方法。
Platen
的一级显式
R-K
方法为
X
n+1
= X
n
+ f(X
n
)∆t + g(X
n
)∆W
n
+
1
2
√
σ
g(X
n
+ g(X
n
)
√
σ) − g(X
n
)
(∆W
2
− ∆t)
其中：σ 为 ∆W
n
∼ N (0, σ)，例如 σ = ∆t
n
。
二级显式 R-K 方法
X
n+1
= X
n
+
3
4
K
1
+
1
4
K
2
+
1 −
f
′
(X
n
)
4
∆t
−2
g(X
n
)∆W
n
K
1
=
1 −
f
′
(X
n
)
4
∆t
−1
f(X
n
)∆t
K
2
=
1 −
f
′
(X
n
)
4
∆t
−1
f(X
n
+ K
1
)∆t
选取适当的参数，此方法可达到强 1 阶收敛。
三级显式 R-K 方法
X
n+1
= X
n
+
1
6
(K
1
+ K
2
+ K
3
) + g(X
n
)∆W
n
K
1
= f (X
n
)∆t
K
2
= f (X
n
+ K
1
)∆t
K
3
= f
X
n
+
K
1
+ K
2
4
∆t
http://www.ma-xy.com 31 http://www.ma-xy.com
http://www.ma-xy.com
1.4 随机微分方程数值方法 第一章 随机微分方程
2. 隐式数值方法
1). 隐式 Euler-Maruyama 方法
X
n+1
= X
n
+ f(X
n+1
)∆t + g(X
n+1
)∆W
n
该方法是强 1 阶收敛的。
2). 隐式 Milstein 方法
X
n+1
= X
n
+ f(X
n+1
)∆t + g(X
n+1
)∆W
n
+ g(X
n+1
)g
′
(X
n+1
)(∆W
2
− ∆t)
3). 隐式 R-K 方法 (二级)
X
n+1
= X
n
+
3
4
K
1
+
1
4
K
2
+
1 −
f
′
(X
n+1
)
4
∆t
−2
g(X
n+1
)∆W
n
K
1
=
1 −
f
′
(X
n+1
)
4
∆t
−1
f(X
n+1
)∆t
K
2
=
1 −
f
′
(X
n+1
)
4
∆t
−1
f(X
n+1
+ K
1
)∆t
注意到，隐式数值方法中的 X
n+1
并不能直接计算，解决这个问题的一种可行的方法是：预
估-矫正方法。例如，我们先用显式 Euler-Maruyama 方法预估 X
n+1
为
˜
X
n+1
，然后再带入到隐
式数值方法中进行计算。
˜
X
n+1
= X
n
+ f(X
n
)∆t + g(X
n
)∆W
n
X
n+1
= X
n
+ f(
˜
X
n+1
)∆t + g(
˜
X
n+1
)∆W
n
+ g(
˜
X
n+1
)g
′
(
˜
X
n+1
)(∆W
2
− ∆t)
半隐式数值方法
以半隐式
EM(Euler-Maruyama)
方法为例，其余半隐式方法类推。半隐式
EM
方法为
X
n+1
= X
n
+ (1 − θ)f(X
n
)∆t + θf(X
n+1
)∆t + g(X
n
)∆W
n
其中：θ ∈ [0, 1]，半隐式方法的 f 是显式方法和隐式方法的结合，而扩散函数 g 是显式方法。特
别地，当 θ = 0 时，为显示 EM 方法，当 θ = 1 时，称为后退 (back)EM 方法，即
X
n+1
= X
n
+ f(X
n+1
)∆t + g(X
n
)∆W
n
注意：
1.EM 方法收敛于 Ito 型 SDE 而非 Stratonovich 型 SDE；
2. 前面讨论的 EM、Milstein 和 R-K 方法都是一步数值方法，X
n+1
的计算只需要 X
n
；
3. 前面讨论的数值方法都是在相应 SDE 下的。
http://www.ma-xy.com 32 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.4 随机微分方程数值方法
双隐式数值方法
以双隐式 Milstein 方法为例，其余双隐式方法类推。双隐式 Milstein 方法为
X
n+1
= X
n
+ (1 − θ)f(X
n
)∆t + θf(X
n+1
)∆t + g(X
n
)∆W
n
+ L
′
g(X
n
)∆W
2
n
−
1 − σ
2
L
′
g(X
n
)∆t −
σ
2
L
′
g(X
n+1
)∆t
其中：θ ∈ [0, 1]，σ ⩽ 1 ，L
′
= g
∂
∂x
。
双隐式 Milstein 方法是强 1 阶收敛的，其全局 MS 稳定的充要条件为
(2λ + µ
2
) + ∆tλ
2
(1 − 2θ) +
∆tµ
2
2
(2σλ + µ
2
) < 0
特别地，当 θ = σ = 1 时，其 MS 稳定的充要条件为
(2λ + µ
2
) − ∆tλ
2
+
∆tµ
2
2
(2λ + µ
2
) < 0
1.4.5 数值实验
前面给出一类 SDE 在全局 L(Lipschitz) 条件及线性增长条件下的解的存在唯一性、SDE 解
析解数值解及解的收敛性和稳定性。下面，我们对数值方法的收敛性和稳定性进行数值验证。仍
然考虑如下 Black-Scholes 方程 (在 Euler 方法中我们用之进行了数值实验)
dX
t
= λX
t
dt + µX
t
dW
t
X(t
0
) = X
0
t ∈ [t
0
, T ]
上述 SDE 的解析解为
X
t
= X
0
exp
λ −
1
2
µ
2
t + µW
t
下面，我们分别对 EM 法、Milstein 方法和 R-K 方法的收敛性和稳定性进行验证。仍采用之
前的参数设置，令方程中的参数为：λ = 2 ，µ = 1，t
0
= 0，X
0
= 1，T = 1，N = 2
8
，h = ∆t =
T −t
0
N
。
Euler 方法的数值实验
Euler 方法的强 0.5 阶收敛性验证 由强收敛性的定义我们知道，在任意时刻 t
n
上，E(|X
n
−
X
t
n
|) 都非常小，满足 E(|X
n
−X
t
n
|) ⩽ ch
1
2
。由于是任意时刻，这里不妨取末端时刻 T (T = Nh)(当
然可以取其他时刻)，令
e
h
= E(|X
n
− X
T
|)
假设我们取定某个步长 h，且进行 10000 次模拟，那么我们可以计算出该步长下的 e
h
。如果 EM
方法的强收敛阶为 0.5，那么应该有
e
h
⩽ ch
1
2
http://www.ma-xy.com 33 http://www.ma-xy.com
http://www.ma-xy.com
1.4 随机微分方程数值方法 第一章 随机微分方程
上式两边都是步长 h 的函数，两边同时取对数，有
log e
h
≈ log c +
1
2
log h
我们取不同的步长 h
1
, h
2
, h
3
, . . . 。如果 log e
h
和 log h 所绘成的直线斜率接近
1
2
，则说明 EM
方法的强收敛阶确实为 0.5。
作为实验，这里，我们取步长 h
1
= 2
−8
, h
2
= 2
−7
, h
3
= 2
−6
, h
4
= 2
−5
, h
5
= 2
−4
，对每个步
长 h
i
，在 [0, 1] 上计算解析解与数值解 10000 次，然后计算 e
h
i
= E(|X
n
− X
T
|)，最终画 log e
h
和 log h 图，如图 (1.4) 所示
图 1.4: EM 方法的强收敛阶验证
Euler 方法的稳定域的求取 在判定 Euler 方法的稳定性之前，我们先来求解 Euler 方法的稳定
域。
1). 显式 Euler 方法的稳定域
X
n+1
= X
n
+ λX
n
h + µX
n
∆W
n
= X
n
+ λh + µX
n
√
hξ
= (1 + λh + µ
√
hξ)X
n
令 p = λh, q = µ
√
h，有
X
n+1
= (1 + p + qξ)X
n
对上式两边求 2 阶原点矩，有
E(|X
x+1
|
2
) = E
(1 + p + qξ)
2
X
2
n
注意到 ξ ∼ N (0, 1)，我们得到 EM 方法的 MS 稳定函数为
R
2
= (1 + p)
2
+ q
2
当 R
2
< 1 时得到稳定域如图 (1.5) 所示
http://www.ma-xy.com 34 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.4 随机微分方程数值方法
1 h = e z p l o t ( ’ (1+p )^2+q^2−1 ’ , [ − 2 ,0 ] , [ 0 , 1 . 5 ] ) ;
2 se t (h , ’ LineStyle ’ , ’− ’ , ’ c o lo r ’ , ’ r ’ )
3 xl a b e l ( ’ p = lambda*h ’ )
4 yl a b e l ( ’ q = mu * s q r t (h) ’ )
5
图 1.5: 显式 EM 方法的稳定域示意图
2). 半隐式 EM 方法的稳定域。迭代公式为
X
n+1
=
1 + qξ
1 − p
X
n
其均方稳定函数为
R
2
=
1 + q
2
(1 − p)
2
3). 隐式 EM 方法的稳定域。其迭代公式为
X
n+1
=
1
1 − p − qξ + q
2
ξ
2
X
n
其均方稳定函数为
R
2
=
1
√
2π
∞
−∞
1
((qx −
1
2
)
2
+
3
4
− p)
2
e
−
x
2
2
dx
Euler 方法的稳定性验证 取不同步长 h
1
= 1, h
2
=
1
2
, h
3
=
1
4
，在 [0, 20] 区间上进行观察，对每
个步长 h
i
，用 Euler 方法迭代 10000 次，求取均值 E(|X
t
|
2
)，并进行绘图，如图 (Euler 方法稳
定性验证图)
http://www.ma-xy.com 35 http://www.ma-xy.com
http://www.ma-xy.com
1.4 随机微分方程数值方法 第一章 随机微分方程
图 1.6: Euler 方法稳定性验证图
观察上图中当 t → ∞ 时 lim
t→∞
E(|X
t
|
2
) 的性质，如果 lim
t→∞
E(|X
t
|
2
) = 0，说明 EM 方法 MS
稳定，并且发现 p, q 在稳定域中。
Milstein 方法的数值实验
Milstein 方法的稳定域的求取 Milstein 方法的强 1 阶收敛性的验证和稳定性检验可以参考
Euler 方法。下面，我们主要介绍 Milstein 方法的稳定域的求解。
1. 显式 Milstein 方法的迭代公式为
X
n+1
= (1 + p + qξ + q
2
(ξ
2
− 1)/2)X
n
其均方稳定函数为
R
2
(1 + p)
2
+ q
2
+
1
2
q
4
2. 半隐式 Milstein 方法的均方稳定函数为
R
2
=
2(1 + (1 − θ)p)
2
+ 2q
2
+ q
4
2(1 − θp)
3. 隐式 Milstein 方法的均方稳定函数为
R
2
=
1
√
2π
∞
−∞
4
(3 − 2p + q
2
− (qx − 1)
2
)
2
e
x
2
2
dx
R-K 方法的数值实验
上面提到的 EM 方法和 Milstein 方法是收敛于 Ito 型 SDE，而对 Stratonovich 型 SDE 而
言，R-K 方法是较高效的。
R-K 方法的稳定域的求取 R-K 方法的强 1 阶收敛性的验证和稳定性检验可以参考 Euler 方
法。下面，我们主要介绍 R-K 方法的稳定域的求解。
1.
显式
2
级
R-K
方法的迭代公式为
X
n+1
= X
n
1 +
(2λh − µ
2
h + 2µ
√
hξ)(4 + 2λh − µ
2
h + 2µ
√
hξ)
8
http://www.ma-xy.com 36 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.4 随机微分方程数值方法
令 p = λh, q =
√
hµ, ξ ∼ N (0, 1)，有 MS 稳定函数为
R
2
=
1
64
[64 + 16(2p − q
2
)(2p − q
2
+ 4) + 64q
2
+ (2p − q
2
)
2
(2p − q
2
+ 4)
2
+ 4q
2
(2p − q
2
)
2
+ 4q
2
(2p − q
2
+ 4)
2
+ 48q
4
+ 16q
2
(2p − q
2
)(2p − q
2
+ 4)]
2. 半隐式 2 级 R-K 方法的迭代公式为
X
n+1
= X
n
1 +
(2λh − µ
2
h + 2µ
√
hξ)(2 − 6λh + 3µ
2
h + 4µ
√
hξ)
4(2 − 2λh + µ
2
h)
2
其 MS 稳定函数为
R
2
= 1 +
4p − 12p
2
+ 12pq
2
+ 6q
2
− 3q
4
2(2 − 2p + q
2
)
2
+
(2p − q
2
)
2
(2 − 6p + 3q)
2
+ 16q
2
(2p − q
2
)
2
+ 192q
4
16(2 − 2p + q
2
)
4
+
32q
2
(2p − q
2
)(2 − 6p + 3q
2
) + 4q
2
(2 − 6p + 3q
2
)
2
16(2 − 2p + q
2
)
4
3. 隐式 2 级 R-K 方法的迭代公式为
X
n+1
= X
n
1 +
(2λh − µ
2
h + 2µ
√
hξ)(4 − 6 + 3µ
2
h − 6µ
√
hξ)
2(2 − 2λh + µ
2
h − 2µ
√
hξ)
2
其 MS 稳定函数为
R
2
=
1
√
2π
∞
−∞
1 +
(2p − q
2
+ 2qx)(4 − 6p + 3q
2
− 6qx)
2(2 − 2p + q
2
− 2qx)
2
e
−
x
2
2
dx
1. 对于上面介绍的数值方法的收敛性和稳定性，我们是直接给出其收敛阶和稳定性，然后再
用数值方法进行验证。那么，如何依据收敛性和稳定性的定义来证明某些数值方法的收敛性和稳
定性呢？2. 可以看到，解析解 X
t
n
与数值解 X
n
都是随机过程，所以对于解析解和数值解都会
有收敛性和稳定性。
注记
1992.Kloedn 和 Platen 在文献31中通过随机泰勒展开方法得到了一系列数值方法。1995.Mil-
stein 在文献49中研究了泰勒展开及随机 R-K 方法的弱收敛性和均方强收敛性。2004.Milstein 和
Tretyakov 在文献50中进一步研究了 R-K 方法的 MS 稳定性。2005.Rodkina 研究了半隐式变步
长 θ 方法的几乎渐近稳定性。2006.Buchwar 和 Winkler 研究了一类随机线性多步方法。同年，利
用 Lypronov 函数研究了两步 Maruyama 方法的均方渐进稳定性。2007.Kossler 基于彩色根树理
论研究了 R-K 方法的 2 阶弱收敛性。2008Kim 和 Stanescu 给出了低存储高效率的 1.5 阶 R-K
方法。2008.Sickenberger 研究了变步长多步方法的均方收敛性。2010.Debrabant 研究了 3 阶弱
收敛的 R-K 方法。
http://www.ma-xy.com 37 http://www.ma-xy.com
http://www.ma-xy.com
1.5 SDE 解存在唯一条件的进一步讨论 第一章 随机微分方程
1.5 SDE 解存在唯一条件的进一步讨论
考虑如下问题：CEV(Constarit Elastictiy of Variance Model)
dX
t
= µX
t
+ σX
r
t
dW
t
X(t
0
) = X
0
上面的 CEV 模型中，当 r = 1 时，模型是我们前面所讨论的 SDE 模型；当 0 < r < 1 时，模型
中的漂移项 σX
r
t
显然不满足全局 L 条件。我们说 SDE 在全局 L 条件和线性增长条件下存在唯
一解，那么，如果 f, g 不满足全局 L 条件或者线性增长条件，我们应该怎样解决？后面我们将尝
试着解决如下的问题：
(1) 在模型结构不变的情况下，即 SDE 的形式仍为
dX
t
= f (X
t
)dt + g(X
t
)dW
t
若 f, g 不满足全局 L 条件和线性增长条件，SDE 的解存在唯一吗？解析解和数值解如何求取？
解的收敛性和稳定性如何？
(2) 如果模型形状改变，如引入 Poisson 跳、Markov 切换、延迟项和分数阶等等，那么模型
的性质如何？
针对问题 (1)，我们给出 f, g 在其它条件下 (如局部 L 条件和 Khasminskii 条件) 解的存在
唯一性，并介绍 EM 方法、双隐 Milstein 方法、(半) 驯服方法和分裂步等方法用于数值计算。
1.5.1 非全局 Lipschitz 条件下解的存在唯一性
1975 年，Fridman 在文13中研究了局部 L 条件想和线性增长条件下的解的存在唯一性。1980
年，Khasminskii 在文29中研究了在不满足线性增长条件而是满足 Khasminskii 条件下，解的存
在性，并给出 Euler 方法的依概率收敛性。
非全局 L 条件
下面，我们来介绍一些非局部 L 条件。考虑如下 SDE 模型
dX
t
= f (X
t
) + g(X
t
)dW
t
X(t
0
) = X
0
t ∈ [t
0
, T ]
其中：X
0
∈ R，W
t
是一维 Brown 运动 (可以将 Brown 运动扩展到高维)。
1). 全局 L 条件。如果 f ∈ C
1
(R) 且 ∃L > 0，使得 ∀x, y ∈ R，有
|f(x) − f (y)| ⩽ L|x − y|
则称函数 f 满足全局 L 条件。
http://www.ma-xy.com 38 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.5 SDE 解存在唯一条件的进一步讨论
2). 单边 L 条件。如果 f ∈ C
1
，且 ∃L > 0，使得 ∀x, y ∈ R，有
x − y, f (x) − f (y)
⩽ L|x − y|
2
则称 f 满足单边 L 条件。
3). 局部 L 条件。如果 f ∈ C
1
，对每个常数 j，∃L
j
> 0，使得 ∀x, y ∈ R, |x| ∨ |y| < j，有
|f(x) − f (y)| ⩽ L
j
|x − y|
则称 f 满足局部 L 条件。其中：|x| ∨ |y| = max{x, y} 。
4). 线性增长条件。如果 f ∈ C
1
，并且 ∃L > 0，使得 ∀x ∈ R，有
|f(x)|
2
⩽ L(1 + |x|
2
)
则称 f 满足线性增长条件。
5).Khasminskii 型条件。∀(x, y) ∈ R
n
× R
n
，∃v ∈ C
2
(R
n
; R
+
), α > 0，使得
lim
|x|→∞
inf v(x) = ∞
和
Lv(x, y) ⩽ α(1 + v(x) + v(y))
其中：Lv : R
d
× R
d
→ R。若 v ∈ C
2,1
(R
d
; R)，则 Lv 定义为
Lv(x, y, t) = v
t
(x, t) + v
x
(x, t)f(x, y) +
1
2
trac
g
T
(x, y)v
xx
(x, t)g(x, y)
解的存在唯一性
上面介绍了一些非全局 L 条件，下面我们给出在这个非全局 L 条件下，SDE 解的存在唯一
性。
单边 L 条件下的存在唯一性 如果 f, g ∈ C
1
(R)，且 f 满足常数 µ 的单边 L 条件，g 满足常数
¯
C 的全局 L 条件，即 ∃µ,
¯
C > 0，使得 ∀x, y ∈ R，有
x − y, f (x) − f (y)
⩽ µ|x − y|
2
|g(x) − g(y)|
2
⩽
¯
C|x − y|
2
则 SDE 存在唯一解，并且有对 ∀p ⩾ 2，∃C = C(p, T ) > 0，使得
E
sup
0⩽t⩽T
|X
t
|
p
⩽ C(1 + E|X
0
|
p
)
注：上面仅仅是解存在唯一的充分条件，而非必要。
http://www.ma-xy.com 39 http://www.ma-xy.com
http://www.ma-xy.com
1.5 SDE 解存在唯一条件的进一步讨论 第一章 随机微分方程
局部 L 条件下的存在唯一性 如果 f, g ∈ C
1
(R)，并且 f, g 满足局部 L 条件，即对每个正整数
j，∃K
j
，使得 ∀x, y ∈ R, |x| ∨ |y| < j，有
|f(x) − f (y)| ⩽ K
j
|x − y|
则 SDE 存在唯一局部解。在此基础上，如果 ∃γ > 0，使得 ∀x, y ∈ R，有
2xf(x) + |g(x)|
2
< −γ|x|
2
则 SDE 存在唯一全局解，并且 X
t
指数均方稳定
E|X
t
|
2
⩽ C(X
0
)e
−γt
其中：C(X
0
) 是一个依赖于初始值 X
0
的常数。
值得一提的是，当 f 满足单边 L，g 满足全局 L 时，可以由均值定理得到 f, g 满足局部 L，
另外
x, f(x)
∨ |g(x)|
2
⩽ α + β|x|
2
其中：α =
1
2
|
f
(0)
|
2
∨ 2|g(0)|
2
，β = (µ +
1
2
)
∨
2
¯
C
。更多的内容可以参考
1997.Mao30。
1.5.2 非全局 Lipschitz 条件下数值方法
上面，我们给出了 SDE 在某些条件下解的存在唯一性，下面，我们将要讨论某些数值方法
的收敛性和稳定性。
Euler 方法 以 Euler 方法为例，2002 年 Higham, Mao 和 Stuor18给出在 f, g 满足局部 L 条
件，或者 f 满足单边 L 条件，g 满足全局 L 条件下，当单边 L 条件常数 c、Euler 方法中的 θ 和
步长 ∆ = h 满足 cθ∆ < 1 以及 SDE 解析解和数值解矩有界条件下，Euler 方法是强收敛的，即
对 ∀p > 2，有
E
sup
0⩽n⩽N
|X
n
|
p
< ∞
E
sup
0⩽n⩽N
|X
t
n
|
p
< ∞
则
lim
∆t→0
E
sup
0⩽n⩽N
|X
n
− X
t
n
|
2
= 0
Mao 和 Szpvuch 在《Strong Convergance rates for back ward》还证明了在单边 L 条件结合一类
多项式增长条件结合多项式 L 条件下，BackEuler 方法的强收敛阶为 0.5。
现在，我们介绍具体的 Euler 方法。考虑如下两类 Euler 方法
1).θ-Euler 方法
X
n+1
= X
n
+ θf(X
n+1
)∆t + (1 − θ)f(X
n
)∆t + g(X
n
)∆W
n
http://www.ma-xy.com 40 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.5 SDE 解存在唯一条件的进一步讨论
2). 分裂步-Euler 方法
X
n
= y
n
+ θf(X
n
)∆t
y
n+1
= y
n
+ f(X
n
)∆t + g(X
n
)∆W
n
n = 0, 1, . . . , N − 1
其中：X
0
= X(t
0
)，t
0
= 0，θ ∈ [0, 1]。
对于上述 Euler 方法，我们关心两个问题：1. 收敛性 2. 稳定性。对于上述两种 Euler 方法
的收敛性，我们有如下结论：
定理 (Euler 方法的收敛性) 在 SDE 存在唯一解的情况下，当 θ ∈ [0,
1
2
) 且 f 满足线性增
长条件时，或者当 θ ∈ (
1
2
, 1] 时，对 ∀p ⩾ 2，∀t
n
∈ [0, T ]，有
lim
∆t→0
E
sup
t
n
∈[0,T ]
|X
n
− X(t
n
)|
p
= 0
lim
∆t→0
E
sup
t
n
∈[0,T ]
|X
n
− y(t
n
)|
2
= 0
但上述结论并没有给出 Euler 方法的收敛阶数，对于两类 Euler 方法的强 0.5 阶收敛性，我
们有如下结论：
定理 (强 0.5 阶收敛性) 在 SDE 存在唯一解的情况下，若 ∃D, q > 0，使得 ∀x ∈ R，有
|f
′
(x)| ⩽ D(1 + |x|
q
)
且对 ∆t < ∆
，
θ
∈
[0
,
1]
，当
θ
∈
[0
,
1
2
)
时，
f
满足线性增长条件，则
∀
p
⩾
2
，有
E
sup
t
n
∈[0,T ]
|X
n
− X(t
n
)|
p
⩽ C(∆t)
p
2
E
sup
t
n
∈[0,T ]
|X
n
− y(t
n
)|
p
⩽ C(∆t)
p
2
其中：∆t = h，∀t
n
∈ [0, T ]
∆ =
∞, θ = 0
1
2θβ
, θ > 0
上述结论证明可参考75宗小峰. 华中科技大学博士论文《随机微分方程的数值分析及随机稳定
化》。
上面给出了 Euler 方法的收敛性，下面讨论 Euler 方法的稳定性。在 SDE 精确解存在唯
一且指数均方稳定条件下，我们来研究 Euler 数值方法的稳定性，并且为了保证 Euler 方法在
θ ∈ (0, 1] 时有意义，我们假设 f 满足单边 L 条件。
定理 (Euler 稳定性) 如果 f, g 满足局部 L 条件，且 ∃γ > 0，使得 ∀x, y ∈ R，有
2xf(x) + |g(x)|
2
⩽ −γ|x|
2
http://www.ma-xy.com 41 http://www.ma-xy.com
http://www.ma-xy.com
1.5 SDE 解存在唯一条件的进一步讨论 第一章 随机微分方程
则 SDE 存在唯一全局解，并且满足
E|X
t
| ⩽ C(X
0
)e
−γt
在此情况下，如果 θ ∈ [0, 1] 且在 θ ∈ (0, 1] 时，f 满足单边 L 条件，
¬在 θ ∈ [0,
1
2
] 时，f 满足常数 K 的线性增长条件，则 ∀∆t < ∆
∗
:=
γ
(1−2θ)K
∧ ∆，有
E|y
t
|
2
⩽ C(X
0
)e
−γ
∆
K∆t
E|X
t
|
2
⩽ C(X
0
)e
−γ
∆
K∆t
其中：
γ
∆
= −
1
∆t
log
1 −
γ − (1 − 2θ)∆tK
(θ∆t
√
K + 1)
2
∆t
在 θ ∈ [
1
2
, 1] 时，对 ∀∆t < ∆，有
E|y
t
|
2
⩽ C(X
0
)e
−γ
∆
K∆t
E|X
t
|
2
⩽ C(X
0
)e
−γ
∆
K∆t
其中：
γ
∆
= −
1
∆t
log
1 −
γ(2θ − 1)
2θ − 1 + γ∆tθ
2
∆t
上述定理中的 ∆ 为
∆ =
1
µθ
, θ ∈ (0, 1]且µ > 0
∞, θ = 0或µ < 0
在
θ
∈
[0
,
1
2
]
时，经常需要
f
满足单边
L
条件的同时，还满足
µ <
0
，而在
θ
∈
[
1
2
,
1]
时，
Euler
方法是无条件保持稳定的；在 θ ∈ [0,
1
2
] 时，在 f 满足线性增长条件下，两类 Euler 方法稳定，
而在 θ ∈ [
1
2
, 1] 时，则不需要 f 的线性增长条件。
Milstein 方法 考虑如下两种 Milstein 数值方法
1.θ-Milstein 方法
X
n+1
= X
n
+ f(X
n
)∆t + g(X
n
)∆W
t
+
1
2
L
′
g(X
n
)(|∆W
n
|
2
− ∆t)
其中：∆t = h 为步长，L
′
= g(x)
∂
∂x
。
2. 分裂步 Milstein 方法
X
n
= y
n
+ θf(X
n
)∆t
y
n+1
= y
n
+ f(X
n
)∆t + g(X
n
)∆W
n
+
1
2
L
′
g(X
n
)(|∆W
n
|
2
− ∆t)
k = 0, 1, . . . , N − 1
y
0
= X(0)
http://www.ma-xy.com 42 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.6 SDE 结构的进一步讨论
2013 年，Higham 等在文献17中提出双隐 Milstein 方法
X
n+1
= X
n
+ [(1 − θ)f(X
n
) + θf(X
n+1
)]∆t + g(X
n
)∆W
n
+
1
2
L
′
g(X
n
)(∆W
n
)
2
−
1 −
¯
θ
2
L
′
g(X
n
)∆t −
¯
θ
2
L
′
g(X
n+1
)∆t
并证明在单边 L 条件和单调条件下，该数值方法是强收敛的。其中：θ ⩾ 0,
¯
θ ⩽ 1, L
′
= g(x)
∂
∂x
。
驯服方法 2012 年，Hutzentler 在文38中提出驯服 (Tamed)-Euler 方法
X
n+1
= X
n
+
f
(
X
n
)∆
t
1 + |f(X
n
)|∆t
+ g(X
n
)∆W
n
并证明 taned-Euler 方法是标准 0.5 阶收敛的。2013 年，Wang 和 Gan 在文58中提出 taned-
Milstein 方法
X
n+1
= X
n
+
f(X
n
)∆t
1 + |f(X
n
)|∆t
+ g(X
n
)∆W
n
+
1
2
L
′
g(X
n
)(|∆W
n
|
2
− ∆t)
并证明该方法是强 1 阶收敛的。
1.6 SDE 结构的进一步讨论
前面讨论的 SDE 模型的结构都是
dX
t
= f (X
t
)dt + g(X
t
)dW
t
无论 f, g 是什么类型的函数，也不论 f, g 满足什么样的条件，SDE 的结构都是上面的结构。下
面，对 SDE 的结构进行如下的扩展：
1. SDE 基本理论
2. 含 Markov 切换的 SDE
3. 含 Poisson 跳的 SDE
4. 随机延迟微分方程 SDDE
5. 中立型 SDE
6. 分段连续 SDE
7. SDE 参数估计问题
8. 随机偏微分方程 SPDE
9. * 随机控制
10. * 分数阶 SDE
http://www.ma-xy.com 43 http://www.ma-xy.com
http://www.ma-xy.com
1.6 SDE 结构的进一步讨论 第一章 随机微分方程
1.6.1 SDE 基本理论
一般形式的 SDE
dX
t
= b(t, X
t
)dt + σ(t, X
t
)dW
t
其中：X ∈ R
n
，b : [0, T ] × R
n
→ R
n
，σ : [0, T ] × R
n
→ R
n×d
，W
t
是 d 为 Brown 运动。
上述 SDE 可以用于描述液体中由分子无规则碰撞的小颗粒的运动，因此，亦称 SDE 的解
X
t
为 Ito 分布。下面，我们来看一下时齐 Ito 分布的 Markov 性，时齐 Ito 分布是如下 SDE 的
唯一解
dX
t
= b(X
t
)dt + σ(X
t
)dW
t
X
s
= x
t ⩾ s
由于上述 SDE 有 X
s
= x 的初始条件，所以将其解记为 X
t
= X
s,x
t
(t ⩾ s)，若初始时刻 s = 0，
则记为 X
t
= X
0,x
t
。
注：在后面的 BSDE 部分，我们用 X
t,x
t
(t ⩽ s ⩽ T ) 表示初始时刻为 t，初始条件为 X
t
= x 的
随机过程。
引进 {X
t
}
t⩾0
的概率分布 Q
x
(x ∈ R
n
, X
0
= x)，且令 F
t
是 Brown 运动 W
t
产生的 σ 代数，
M
t
是 X
t
产生的 σ 代数。由于 X
t
关于 F
t
可测，所以 M
t
⊆ F
t
。下面，给出 X
t
的 Markov
性。
定理 (Ito 分布的 Markov 性) 设 f 为 R
n
→ R 的 Borel 函数，则对 ∀t，对 h ⩾ 0，有
E
x
[f(X
t+h
)|F
t
] = E
X
t
(w)
[f(X
h
)]
即 X
t
在 t 时刻后的行为与 X
t
在 t 时刻前的行为一样。
注：
X
t
的
Markov
性质是在
X
t
时齐性基础上得到的。
X
t
时齐性：
{
X
s,x
s+h
}
h
⩾
0
和 {X
0,x
h
}
h⩾0
具
有相同分布。
上面给出了 X
t
的 Markov 性，下面我们来介绍 X
t
的强 Markov 性。首先，给出停时 τ 的
定义
定义 (停时) 以 {N
t
} 表示一个单增的 σ 代数族。如果 ∀t ⩾ 0，有
{w|τ(w) ⩽ t} ∈ N
t
其中：τ : Ω → [0, ∞) 是一个随机函数。则称函数 τ 为关于 {N
t
} 的停时。
停时的示例 (首达时)：{X
n
} 是一个随机变量序列，A 是一个事件集，令
T (A) = inf{n|X
n
∈ A}
并且约定 T (1) = inf{n|X
n
∈ Φ} = ∞。可见 T (A) 是 {X
n
} 首次进入 A(即发生 A 中所包含事
件) 的时刻，称 T (A) 是 {X
n
} 到集合 A 的首达时。可证明 T (A) 是关于 {X
n
} 的停时。
http://www.ma-xy.com 44 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.6 SDE 结构的进一步讨论
定理 (Ito 分布的强 Markov 性) 已知 f 是 R
n
上的 Borel 函数，τ 是关于 F
t
的停时，
τ < ∞, a.s，则
E
x
[f(X
τ+h
)|F
t
] = E
X
τ
[f(X
h
)] h ⩾ 0
粗略来讲，强 Markov 性就是将 Markov 性中的 t 变为停时 τ。
定理 (Dynkin 公式) 已知 f ∈ C
2
0
(R
n
)。设 τ 是一个停时，E
x
(τ) < ∞，则
E
x
[f(X
τ
)] = f (x) + E
x
τ
0
Af(x
s
)ds
todo: 未完待续。
1.6.2 含 Markov 切换的 SDE
1961 年，Karsovskii 和 Lidslcii 第一次在线性模型中引入切换，形成线性切换模型。他们在
模型中使用连续马尔科夫链来描述这种不同模型结构之间的切换。我们将这种 Markov 切换的
思想引入到 SDE 模型当中，带 Markov 切换的 SDE 的状态空间分为离散状态和连续状态，用
Markov 链来模拟离散状态，用 SDE 来模拟连续状态。并且由于不同过程之间相互干扰，这使得
随机系统变得更加难于处理。大量的实验表明：Brown 运动、Markov 切换和 Poisson 过程等随
机因素的引入可以使一个不平稳的系统稳定化 (当然，反之亦可)，这个观点被称为随机稳定化。
下面，我们将简单介绍一下带 Markov 切换的 SDE，以及其随机稳定化。
令 r (t) ≡ r
t
为取之于有限空间 S = {1, 2, . . . , N} 的右连续 Markov 链，并且有生成元 (状
态转移矩阵)Γ = (r
ij
)
N×N
，其转移概率为
p(r(t + δ) = j|r(t) = i) =
r
ij
δ + o(δ), i = j
1 + r
ij
δ + o(δ), i = j
其中：δ > 0，r
ij
是从状态 i 转移到状态 j 的转移率。同样，我们可以给出连续 Markov 链的
C-K 方程和 Kolmogorov 微分方程等定义。
1999 年，Mao 在文62中给出了带 Markov 切换的随机微分方程
dX
t
= f (t, X
t
, r
t
)dt + g(t, X
t
, r
t
)dW
t
X(0) = x
0
0 ⩽ t ⩽ T
并给出在局部 L 条件和线性增长条件下精确解的 p 阶指数稳定和 T 稳定。2004 年，Yuan 和
Mao 在文71中讨论了带 Markov 切换的 SDE 的 EM(Euler-Maruyama) 方法，并在全局 L 条件
下，给出了 EM 方法的强收敛性。关于带 Markov 切换的随机微分方程更详细的内容可以参考专
著：2006 年 Mao,Yuan 的书籍47。
http://www.ma-xy.com 45 http://www.ma-xy.com
http://www.ma-xy.com
1.6 SDE 结构的进一步讨论 第一章 随机微分方程
1.6.3 含 Poisson 跳的 SDE
模型结构
向 SDE 中加入 Poisson 过程以形成带跳 SDE 模型。关于带跳 SDE 问题，可以参考王小捷
2012 博士论文57、李洋 2012 博士论文35和于辉 2013 博士论文70。
引例：考虑如下期权定价 Black-Scholes 模型
dX
t
= µX
t
dt + σX
t
dW
t
其中：X
t
是股票价格，µ 为收益率，σ 为波动率 (即风险)。由于金融市场的某些风险，股票价格
X
t
可能会在时刻 τ
1
, τ
2
, . . . , τ
i
, . . . 发生比例为 ξ
1
, ξ
2
, . . . , ξ
i
, . . . 的跳跃。在区间 [τ
j
, τ
j+1
)，有
dX
t
= µX
t
dt + σX
t
dW
t
设在跳跃时刻 τ
j
，股票 X
t
的跳跃量为
∆X(τ
j
) = X(τ
j
) − X(τ
−
j
) = X(τ
−
j
)ξ
j
从而，当 t ∈ [0, τ
1
) 时，有
X
t
= X
0
e
(µ−
σ
2
2
)t+σW
t
因此
X(τ
−
1
) = X
0
e
(µ−
σ
2
2
)τ
1
+σW
τ
1
X(τ
1
) = X
0
(1 + ξ
1
)e
(µ−
σ
2
2
)τ
1
+σW
τ
1
以此类推，在时间段 [τ
1
, τ
2
), [τ
2
, τ
3
), . . . ，有
X
t
= X
0
p
ϕ
(t)
i=1
ξ
i
e
(
µ
−
σ
2
2
)τ
i
+σW
τ
i
d 维带跳 SDE 的一般形式为
dX
t
= a(t, X
−
t
))dt + b(t, X
−
t
)dW
t
+
E
c(t, X
−
t
, v)p
ϕ
(dv × dt)
X(x
−
) = X(0) = x
0
∈ L
2
F
0
(Ω : R
d
)
或者写为
dX
t
= f (t, X
−
t
)dt + g(t, X
−
t
)dW
t
+ σ(t, X
−
t
)dN
t
X(0
−
) = X
0
= x
0
如果简单考虑，可以将
N
t
视为一标量
Poisson
过程，强度为
λ >
0
。带跳
SDE
的积分形式
为
X
t
= X
0
+
t
0
a(s, X
s
)ds +
t
0
b(s, X
s
)dW
s
+
t
0
E
c(s, X
−
s
, v)N (dv, ds)
http://www.ma-xy.com 46 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.6 SDE 结构的进一步讨论
其中：X
−
t
= lim
s→t
−
X
s
，模型中的过程 (函数)a, b 为 F
t
适应过程，c 为可料过程。
a : [0, ∞) × R
d
× Ω → R
d
，b : [0, ∞) × R
d
× Ω → R
d×d
，c : [0, ∞) × R
d
× Ω → R
d
。
设
(Ω
,
F
, P
)
为带流
{F
t
}
的完备概率空间，之前的
σ
代数定义为
F
t
=
{
N
∩
σ
{
W
s
|
0
⩽
s
⩽
t
}}
，
这里的 F
t
定义为
σ
(0,s]
N(dv, ds)
s ⩽ t, A ∈ B(E)
∩ σ{W
s
|s ⩽ t} ∩ N
定义 (测度论下的 Poisson 测度) 设 E ∈ R
l
/{0} 并且有 Borel 域 B(E)。在 R
+
× E 上的
一个泊松测度 N ，其强度测度为
ˆ
N(dt, dv) = dtλ(dv)。对所有满足 λ(A) < ∞ 的 A ∈ B(E) ，N
的补偿测度
˜
N((0, t] × A) = (N −
˜
N)((0, t] × A)
t⩾0
是一个鞅 (过程)，λ 是在 [E, B(E)) 上的
任意 σ 有限测度：
|v|
λ(dv) < ∞
Poisson 测度 p
ϕ
(dv × dt) 定义在 E × [0, ∞)，其中：E ⊂ R
l
/{0}，l ∈ N。补偿 Poisson 测度为
ϕ(dv)dt = λf (v)dvdt，这里要求强度 λ = ϕ(E) < ∞。
Poisson 测度是定义在概率空间 (Ω
J
, F
J
, {F
t
}
J
, P
J
) 上的；Brown 运动是定义在概率空间
(Ω
W
, F
W
, {F
t
}
W
, P
W
) 上的，因此，X
t
定义在 (Ω, F, {F
t
}, P ) 上，其中：Ω = Ω
J
× Ω
W
，F =
F
J
× F
W
，{F
t
} = {F
t
}
J
× {F
t
}
W
，P = P
J
× P
W
，且 F 包含所有令测集 N，Brown 运动
和 Poisson 测度相互独立。p
ϕ
= {p
ϕ
(t) = p
ϕ
(E × [0, t))} 表示计算到某给定时间之前跳跃次数的
随机过程。当 T 是给定的有限常数，且 λ < ∞ 时，p
ϕ
(dv × dt) 产生序列对 {(τ
i
, ξ
i
)}
p
ϕ
(T )
i=1
，其
中：{τ
i
|Ω → R
1
, i = 1, . . . , p
ϕ
(T )} 是非负随机变量序列，表示强度为 λ 的泊松过程的跳跃时刻。
{ξ
i
|Ω → E, i = 1, . . . , p
ϕ
(T )} 是独立随机变量序列，服从分布函数 ϕ(dv)/ϕ(E)。
解存在唯一性
上面介绍了带 Poisson 跳的 SDE 模型的一般形式，下面，给出其解的存在唯一性即数值解
方法。
定理 (解存在唯一性) 设 a, b, c 满足全局 L 条件和线性增长条件，即 ∃K
1
> 0, K
2
> 0，使
得 ∀x, y ∈ R
d
，有
|a(t, x) − a(t, y)|
2
∨ |b(t, x) − b(t, y)|
2
∨
E
|c(t, x, v) − c(t, y, v)|
2
ϕdv ⩽ K
1
|x − y|
2
|a(t, x)|
2
∨ |b(t, x)|
2
∨
E
|c(t, x, v)|
2
ϕ(dv) ⩽ K
2
(1 + |x|
2
)
则解 X
t
存在唯一，且
sup
t
E|X
t
|
2
< ∞
E
sup |X
s
|
2
⩽ C(1 + E|X
0
|
2
)
http://www.ma-xy.com 47 http://www.ma-xy.com
http://www.ma-xy.com
1.6 SDE 结构的进一步讨论 第一章 随机微分方程
上述定理的证明可以参考 1972.Gikhman 文献15。
上述解存在唯一性是在全局 L 条件和线性增长条件下进行讨论的，下面，我们给出非全局 L
条件下的解存在唯一性。
定理 (非全局 L 条件下解存在唯一性) 设 a, b, c 满足非全局 L 条件：∃c
k
> 0, C > 0, L > 0，
使得 ∀x, y ∈ R
n
,k ⩾ 1, |x| ∨ |y| ⩽ k，有
|a(t, x) − a(t, y)|
2
∨ |b(t, x) − b(t, y)|
2
⩽ C
k
|x − y|
2
E
|c(t, x, v) − c(t, y, v)|
2
ϕ(dv) ⩽ C|x − y|
2
2
x, a(t, x)
+ |b(x)|
2
+
E
|c(t, x, v)|
2
ϕ(dv) ⩽ L(1 + |x|
2
)
上述定理的证明可以参考 2005.Bruti-Liberati3或者 2007.Bruti-Liberati4。2002 年，Li 利用负的
Lyapunov 指数判定了带 Poisson 测度的 SDE 的稳定性，并推导出方程几乎必然渐进稳定的充
要条件。
注：关于 Markov 切换的带跳 SDE 系统的随机稳定性问题可以参考：2001.Swishchuk 的文章27和
2014. 宗小峰. 华中科技大博士论文75第七章。
数值方法
上面介绍了带跳 SDE 的结构和解的存在唯一性，下面，我们来介绍一些数值方法。考虑如
下带跳 SDE 问题
dX
t
= f (X
−
t
)dt + g(X
−
t
)dW
t
+ σ(X
−
t
)dN
t
X(0)
−
= X
0
(1.17)
上述模型的 θ-Euler 方法为
X
t
n+1
= X
t
n
+ (1 − θ)∆t
n
f(X
t
n
) + θ∆t
n
f(X
t
n+1
)
+ g(X
t
n
)∆W
t
n
+ σ(X
t
n
)∆N
t
n
定义补偿 Poisson 过程为
˜
N
t
= N
t
− λt
补偿泊松过程是一个鞅，且满足
E(
˜
N
t+s
−
˜
N
t
) = 0
E|
˜
N
t+s
−
˜
N
t
|
2
= λs, t, s ⩾ 0
现在，我们将 (1.17) 的带跳 SDE 改为如下新的 SDE
dX
t
= f
λ
(X
−
t
)dt + g(X
−
t
)dW
t
+ σ(X
−
t
)d
˜
N
t
(1.18)
http://www.ma-xy.com 48 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.6 SDE 结构的进一步讨论
其中：f
λ
(x) = f (x) + λσ(x)。注意到，新的 f
λ
(x) 仍然满足全局 L 条件和线性增长条件，只是
相应的 L 常数变为
Lλ = (λ + 1)
2
, Kλ = (λ + 1)
2
k
新的补偿带跳 SDE(1.18) 方程的 θ-Euler 方法为
X
t
n+1
= X
t
n
+ (1 − θ)∆t
n
f
λ
(X
t
n
) + θ∆t
n
f
λ
(X
t
n+1
)
+ g(X
t
n
)∆W
t
n
+ σ(X
t
n
)∆
˜
N
t
n
注记
针对满足全局 L 条件和线性增长条件的带 PoissonJump 的 SDE 模型，1982.Platen8叙述了
关于任意强收敛阶
r
∈ {
0
.
5
,
1
,
1
.
5
, . . .
}
的定理，并构造了
Jump-Adapted
数值方法。
1987.Magh-
soodi41讨论了 Euler 方法的依概率收敛性。1995.Li56研究了数值方法的几乎必然收敛性。1996.Magh-
soodi66分析了基于 Ito-Taylor 展开的 1.5 阶数值方法。1998.Maghsoodi67分析了基于 Ito-Taylor
展开的 2 阶数值方法。2000.Liu36研究了 Euler 方法的几乎必然稳定性。2004.Gardon1给出了
Ito-type 型数值方法。2005.Higham19给出了 Euler-Maruyama 方法。2006.Higham20研究了 0.5
阶隐式方法的收敛性和稳定性。2006.Gardon2给出了 1.5 阶数值方法。2007.Higham21给出了向
后 Euler 方法的收敛性。2011.Buckwar5在单边 L 条件下，研究了两类 Runge-Kutta 方法的均方
收敛性。2011.Mishura51在非全局 L 条件下，研究了 Euler 方法的收敛性。2012.Huang22在全局
L 条件和线性增长条件下，研究了二级一阶 R-K 方法，并利用树图理论推导出了该数值方法。
1.6.4 随机延迟微分方程 SDDE
像常延迟微分方程那样，许多时候，当前时刻的增量 dX
t
不仅与当前时刻的状态量 X
t
有
关，而且还和之前时刻的状态量 X
t−τ
有关。将之前的状态量引入到 SDE 当中，从而形成随机
延迟微分方程 (思考：能否将未来时刻的状态量引入到 SDE)。给出 SDDE 的一个示例：
dX
t
= (−3X
t
− X
3t
+ X
t
sin(X
t−1
))dt + X
t
sin
3
(X
t−1
)dW
t
在上面的示例中，我们称 τ = (t −1) 为延迟数。当然，更一般的延迟数 τ 可以是时间 t 的函
数，可以是固定的，如每个 X
t
都与 X
t−τ
有关；τ 也可以是时变的，即每个时刻的 X
t
与 X
t−τ(t)
有关。下面，我们给出定延迟 τ 型 SDDE 和变延迟 τ (t) 型 SDDE 的一般形式。
1. 定延迟型 SDDE
dX
t
= f (t, X
t
, X
t−τ
)dt + g(t, X
t
, X
t−τ
)dW
t
X
t
= ξ
t ∈ [−τ, 0]
(1.19)
http://www.ma-xy.com 49 http://www.ma-xy.com
http://www.ma-xy.com
1.6 SDE 结构的进一步讨论 第一章 随机微分方程
2. 变延迟型 SDDE
dX
t
= f (X
t
, X
τ(t)
)dt + g(X
t
, X
τ(t)
)dW
t
X
t
= ξ
t ∈ [−τ, 0]
(1.20)
式 (1.19) 和式 (1.20) 型式的 SDDE 解的存在唯一性要求 f, g 满足全局 L 条件和线性增长条件。
更多内容可以参考：2013. 张玲的博士论文72、2013. 于战华的博士论文70和 2013. 陈琳的博士论
文6。
注记
2003.Mao42研究了式 (1.20) 型 SDDE 在局部 L 条件和 p 阶矩条件下 EM 方法的强收敛性。
2006.Mao47研究了式 (1.19) 型 SDDE 在局部 L 条件和 Khasminskii 型条件下 EM 方法的依概
率收敛性。2011.Mao44yanjiule 式 (1.19) 型 SDDE 在局部 L 条件和一般性 Khasminskii 型条件
下 EM 方法的依概率收敛性。2004.Liu74给出了线性 (1.19) 型 SDDE 的半隐式 Euler 方法及其
收敛性和均方渐进稳定性。2007.Mao43给出了 (1.19) 型 SDDE 在数值解收敛条件下，精确解的
均方指数稳定和 EM 方法的均方指数稳定是等价的。2010.Wu28给出 (1.19) 型 SDDE 在局部 L
条件下，EM 方法几乎处处指数稳定。
1.6.5 中立型 SDE
中立型 SDE 的一般形式
我们先来看一个中立型 SDE 的示例，然后再讨论其一般形式。示例：
d[X
t
−
1
2
sin(X
t−2
)] = (−8X
t
+ sin(X
t−2
))dt + X
t−2
dW
t
X
t
= t + 1
t ∈ [−2, 0]
由上式写出中立型 SDE 的一般形式
d(X
t
− D(X
t
)) = f (t, X
t
, X
τ
)dt + g(t, X
t
, X
τ
)dW
t
(1.21)
由 τ 是否固定，可将 (1.21) 分为下面两类
d(X
t
− D(X
t−τ
)) = f (t, X
t
, X
t−τ
)dt + g(t, X
t
, X
t−τ
)dW
t
(1.22)
d
X
t
− D
X
t−τ(t)
= f
t, X
t
, X
t−τ(t)
dt + g
t, X
t
, X
t−τ(t)
dW
t
(1.23)
如果在 (1.21) 中再考虑 Markov 切换，则有
d(X
t
− D(X
t
), r
t
) = f (r
t
, X
t
, X
τ
)dt + g(r
t
, X
t
, X
τ
)dW
t
(1.24)
http://www.ma-xy.com 50 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.6 SDE 结构的进一步讨论
注记
1、方程 (1.21)
在 1995 年，Mao61研究了一类中立型随机泛函微分方程的零解的均方指数稳定性，并且到
了充分性条件。但文61中的结论对中立型随机延迟微分方程 (1.22) 不适用。为了克服文61的结论
的局限性，Mao30利用 Razumikhin 技巧研究了中立型随机泛函微分方程 (1.21) 零解的几乎必然
和均方指数稳定性，并且文30中的结论对于中立型随机延迟微分方程 (1.22) 也适用。
Mao 还在文献45利用 Razumikhin 技巧研究了中立型随机泛函微分方程 (1.21) 零解的稳定
性，并且在文献46中总结了中立型随机泛函微分方程 (1.21) 的理论研究成果，给出了解存在唯
一和零解指数稳定的充分性条件。2007 年，Randjelovic 和 Jankovic 在文54给出了中立型随机
泛函微分方程 (1.21) 的零解 p 阶矩指数稳定的评判标准。2008 年，Wu9研究了中立型随机泛
函微分方程 (1.21) 只满足局部 Lipschitz 条件下解的存在唯一性的条件，建立了方程 (1.21) 的
Khasminskii 型定理。2009 年，Jankovic25等推广了文30的结论，利用 Razumikhin 技巧给出了
中立型随机泛函微分方程 (1.21) 的零解 p 阶矩指数稳定和几乎必然指数稳定的充分性条件。
2008 年，Huang 和 Deng23建立了中立型随机泛函微分方程 (1.21) 零解 p 阶矩渐近稳定的
Razumikhin 型定理。2010 年，Xue 等在文65研究了方程 (1.21) 全局解的存在唯一性和稳定性。
2012 年，Feng 和 Shen11给出了中立型随机泛函微分方程 (1.21) 在线性增长条件和局部 Lipschitz
条件下零解 p 阶矩渐近稳定的新的判据和标准。
2008 年，Wu 和 Mao9研究了中立型随机泛函微分方程 (1.21) 的 Euler-Maruyama 数值解，
在局部 Lipschitz 条件、线性增长条件和压缩条件等条件下，证明了 Euler-Maruyama 数值解是均
方收敛的。2012 年，Zhou 和 Fang73研究了中立型随机泛函微分方程 (1.21) 的 Euler-Maruyama
方法的收敛性。在 f, g 满足局部 Lipschitz 条件、多项式增长条件和压缩条件等条件下他们证明
了方程 (1.21) 的 Euler-Maruyama 方法数值解是依概率收敛于精确解的。
2、方程 (1.22)
在 1981 年，Kolmanovskii 和 Nosov33引入和研究了具有如下形式的中立型随机延迟微分
方程 (1.22). 方程 (1.22) 是中立型随机泛函微分方程 (1.21) 的一类特殊类型。Kolmanovskii 和
Nosov 不仅研究了方程 (1.21) 解的存在性和唯一性而且研究了它的稳定性和渐近稳定性。2000
年，Mao63利用连续的半鞅收敛定理研究了中立型随机延迟微分方程 (1.22) 的稳定性。文63中不
仅给出了方程 (1.22) 的零解 p 阶矩指数稳定和几乎必然指数稳定的一些新的充条件和标准，而
且还研究了其它类型的渐近稳定性，如多项式渐近稳定等。需要强调的是，文63中推论 5.4 改进
了文30 中推论 6.1 的结论，给出了方程 (1.22) 的零解均方指数稳定和几乎必然指数稳定的新的判
据。2010 年，Shaikkhet34建立了关于方程 (1.22) 的 Lyapunov 型稳定性定理。2011 年，Gan14等
研究了中立型随机延迟微分方程 (1.22) 的随机 θ 方法的收敛性。在 f, g 满足局部 Lipschitz 和线
性增长等条件下，他们证明了方程 (1.22) 的随机 θ 方法数值解是均方收敛于精确解的，并且均方
收敛阶为 1/2。2011 年，Wang 和 Chen59研究了中立型随机延迟微分方程 (1.22) 的半隐式 Euler
方法的均方渐近稳定性。在确保精确解均方渐近稳定的充分性条件下，证明了半隐式 Euler 方法
的数值解能够保持精确解的均方渐近稳定性。
3、方程 (1.23)
http://www.ma-xy.com 51 http://www.ma-xy.com
http://www.ma-xy.com
1.6 SDE 结构的进一步讨论 第一章 随机微分方程
2006 年，Luo37等研究了如下的中立型随机变延迟微分方程 (1.23) . 在方程 (1.23) 只满
足局部 Lipschitz 条件 (不满足线性增长条件) 下，Luo 等给出了方程 (1.23) 全局解存在唯一
和 p 阶矩指数稳定的新的判据。2011 年，Milosevic39研究了中立型随机变延迟微分方程 (1.23)
的 Euler-Maruyama 方法的收敛性。在 Khasminskii 型条件下，Milosevic 证明了方程 (1.23) 的
Euler-Maruyama 方法是依概率收敛的。2012 年，Milosevic40应用离散半鞅收敛定理研究了方程
(1.23) 的 Euler-Maruyama 方法的几乎必然指数稳定性。
4、方程 (1.24)
2003 年，Kolmanovskii32等研究了如下的具有马尔科夫转换的中立型随机延迟微分方程
(1.24)。讨论了方程 (1.24) 解的存在唯一性、p 阶矩渐近有界性和指数稳定性。2008 年，Mao48等给
出了方程 (1.24) 全局解存在唯一的充分性条件，并且利用连续的半鞅收敛定理研究了方程 (1.24)
零解的几乎必然渐近稳定性，建立了 LaSalle 型定理。文48中还特别给出了方程 (1.24) 是线性方
程时其零解几乎必然指数稳定的判据。2012 年，Li 和 Mao64改进了文48中方程 (1.24) 的 LaSalle
型定理的条件，建立了新的 LaSalle 型稳定性定理。2011 年，Yin 和 Ma69研究了具有半马尔科夫
转换的中立型随机延迟微分方程 (1.24) 的半隐式 Euler 方法的收敛性，证明了在局部 Lipschitz
等条件下方程 (1.24) 半隐式 Euler 方法数值解是均方收敛于精确解的。
数值方法
以下面的中立型 SDE 为例，展示其数值解法
d(X
t
− D(X
t−τ
)) = f (t, X
t
, X
t−τ
)dt + g(t, X
t
, X
t−τ
)dW
t
X(θ) = ξ(θ)
θ ∈ [τ, 0]
其中：τ > 0 为延迟时间，ξ ∈ C
b
F
0
([−τ, 0]; R
n
)，D : R
n
→ R
n
，f : R
+
× R
n
× R
n
→ R
n
，
g : R
+
× R
n
× R
n
→ R
n
。W
t
是定义在 (Ω, F, {F
t
}, P ) 上的 F
t
可测的一维 Brown 运动。
如果 D, f, g 满足局部 L 条件，则上述中立型 SDE 存在唯一全局解 X
t
(解析解即为零解)。
下面给出零解 (解析解) 的几乎必然稳定和指数均方稳定性。
定理 (解析解的几乎必然稳定) 零解 X
t
的几乎必然稳定描述为：∀ξ ∈ C
b
F
0
([−τ, 0]; R
n
)，有
lim
t→∞
sup
1
6
log |X
t
| < 0, a.s
定理 (解析解的指数均方稳定) 零解 X
t
的指数均方稳定描述为 ∀ξ ∈ C
b
F
0
([−τ, 0]; R
n
)，∃η >
0, C > 0，有
E|X
t
|
2
⩽ CE|ξ|
2
e
−ηt
t ⩾ 0
同样可以定义数值解的稳定性概念。下面，给出零解 X
t
几乎必然稳定和指数均方稳定的充
分条件，具体内容可以参考 2000.Mao63。
http://www.ma-xy.com 52 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.6 SDE 结构的进一步讨论
定理 (几乎必然稳定和指数均方稳定的充分条件) 若 ∃λ
1
, λ
2
, λ
3
λ
4
∈ R
+
，对 ∀t ⩾ 0, x, y ∈
R
n
，有
2(x − D(y))
T
f(t, x, y) ⩽ −λ
1
|x|
2
+ λ
2
|y|
2
|g(t, x, y)|
2
⩽ λ
3
|x|
2
+ λ
4
|y|
2
且 λ
1
− λ
3
> λ
2
+ λ
4
，则上述中立型 SDE 的零解 X
t
是几乎必然稳定和指数均方稳定的。
上面简单介绍了中立型 SDE 的一般形式，解析解的稳定性，下面，我们给出中立型 SDE 的
Euler 数值方法。
X
t
n+1
− D(X
t
n+1
) = X
t
n
− D(X
t
n
−τ
) + f(t
n
, X
t
n
, X
t
n
−τ
)∆t
n
+g(t
n
, X
t
n
, X
t
n
−τ
)∆W
t
n
X
t
n
= ξ(∆t
n
)
1.6.6 分段连续型 SDE
插播：分段连续性微分方程
近来，分段连续型微分方程 (EPCAs) 受到了很多专家学者的关注。分段连续性微分方程的
基本形式如下
x
′
(t) = f (t, x(t), x(h(t))) t ⩾ 0
x(0) = x
0
(1.25)
其中：h(t) 是常数区间，例如 h(t) = [t], h(t) = [t−n], h(t) = t−n[t], 其中 n 是正整数，[t] 是不超
过 t 最大的取整函数。分段连续型微分方程具有本身特点，该方程是一种特殊的延迟微分方程，
又有着与微分方程不同的性质，在某些区间上自变量是常数，而且该方程的解是一个连续的、局
部光滑的函数，初值是由有限集来确定，而不是一个初始函数。因此，对该方程的研究是很有必
要的。Cooke 等人7给出了该方程在生物模型和控制系统中的应用。Wiener24给出了关于该方程
的一般理论与基础的结果，如该方程解的存在唯一性、稳定性、振动性和周期性。自从 2003 年,
刘明珠教授带领的科研团队开始对分段连续型微分方程进行了数值解的研究。下面介绍该方程的
一些研究结果。
1997 年，Wang 等人60给出了分段连续型微分方程
x
′
(t) + px(t − l) + qx([t − k]) = 0 t ⩾ 0
x(0) = x
0
平衡解的全局吸引的充要条件为该方程的特征方程在 |λ| ⩾ 1 上无根。
2004 年，刘明珠等人74给出了分段连续型微分方程
x
′
(t) = ax(tl) + a
0
x([t]) t ⩾ 0
x(0) = x
0
(1.26)
的 Runge-Kutta 方法的数值解的渐近稳定性。
http://www.ma-xy.com 53 http://www.ma-xy.com
http://www.ma-xy.com
1.6 SDE 结构的进一步讨论 第一章 随机微分方程
分段连续随机微分方程的数值方法
我们知道 ODE 在随机分析上的推广形成 SDE，那么，将上面介绍的分段连续型 ODE 在随
机分析上推广就可以形成如下分段连续型 SDE
dX
t
= f (X
t
, X
[t]
)dt + g(X
t
, X
[t]
)dW
t
X(0) = x
0
其中：f : R
n
× R
n
→ R
n
，g : R
n
× R
n
→ R
n×d
，x
0
∈ R
n
，[t] 表示不超过时间 t 的最大整数。
上述模型的积分形式为
X
t
= X
0
+
t
0
f
X
s
, X
[t]
ds +
t
0
g
X
s
, X
[t]
dW
s
下面，我们来介绍上述分段连续 SDE 的解存在唯一性及其数值方法。为了讨论方便，我们
假设函数 f, g 是充分光滑的。
定理 (解存在唯一性) 若 1、X
t
是 [0, ∞] 上连续且 F
t
适应的；2、f
X
s
, X
[t]
∈ L
2
([0, ∞); R
n
)，
并且 g
X
s
, X
[t]
∈ L
2
([0, ∞); R
n×d
)；3、X
t
在 [n, n + 1) ∈ [0, ∞] 上的每个整数区间上几乎处处
满足方程。则上述分段连续 SDE 存在唯一解 X
t
。
上述分段连续 SDE 的 Euler 方法为
X
t
n+1
= X
t
n
+ f(X
t
n
, X
[t
n
]
)∆t
n
+ g(X
t
n
, X
[t
n
]
)∆W
t
n
1.6.7 SDE 参数估计问题
上面介绍的各种 SDE 模型模拟都是在模型中参数给定情况下进行的，现在，我们考虑当模
型中的参数未知时，如何根据观测到的数据来估计模型中的参数。
示例 1：考虑如下时齐股票波动模型
dX
t
= µ(X
t
)dt + σ(X
t
)dW
t
X
0
= x
当我们给出具体的 µ, σ 时，就可以对模型进行模拟 (比如用 Euler 方法)。现在，考虑 µ, σ 是具
有未知参数 α, β 的
dX
t
= µ(X
t
; β)dt + σ(X
t
|α)dW
t
X
0
= x
一般而言，求模型中参数的常用的方法有：参数估计法、非参数估计法和半参数估计法。上
述模型中参数 α, β 的最常用估计方法是极大似然估计，可以证明，对某些遍历扩散过程而言，参
数的极大似然估计是相合且渐进有效的。
1993.Florens12使用了非参数方法对扩散过程的扩散项进行了估计。1997.Stanton55使用非
参数核回归方法估计了利率期限结构回报率和波动率。1998.Fan10将局部线性回归模型应用到
SDE 参数估计当中。2003.Fan 利用局部线性回归模型估计了利率期限结构的回报率和波动率。
http://www.ma-xy.com 54 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.6 SDE 结构的进一步讨论
示例 2：考虑如下时变股票波动模型
dX
t
= µ(t, X
t
; β)dt + σ(t, X
t
|α)dW
t
X
0
= x
上述模型中的收益率函数 µ 和波动率函数 σ 是时变的。2002.Gobet16和 2003.Fan68讨论了上面
的时变股票波动模型。
1.6.8 随机偏微分方程
SPDE 的一般形式及解的存在唯一性
前面的所有问题/所有模型都是在随机 (常) 微分方程上进行的。既然 ODE 可以推广到高维
PDE，那么 SDE(SODE) 也会有相应的 SPDE(随机偏微分方程)，下面，我们给出 SPDE 的一
般形式、解的存在唯一性，然后再一维抛物型 SPDE 上讨论数值计算方法：Euler 谱方法、指数
Milstein 谱方法和指数 Runge-Kutta 谱方法。考虑如下一般形式的 SPDE
dX
t
= (AX
t
+ F (X
t
))dt + G(X
t
)dW
t
X
0
= ξ ∈ H
t ∈ [0, T ]
其中：H 是 Hilbert 空间；A 为线性算子：D(A) → H，是解析半群 S(t) = e
At
的无穷小生成元
(t ⩾ 0)；ξ 是一个 H
−
取值与 W
t
独立的随机变量；F : H → H；G : H → L(V, H) 的满足一定正
则性的映射；W
t
: [0, T ] ×Ω → V 是一个标准的关于 {F
t
} 的 Q-Wiener 过程；Q ∈ L(V ) 是一个
对称的，非负的迹类算子：Tr(Q) < ∞；(Ω, F, P ) 是带流 {F
t
} 的完备概率空间；(H, ⟨·, ·⟩
H
, ||·||
H
)
和 (V,
·, ·
V
, ||·||
V
) 是两个可分的 Hilbert 空间；L(V, H) 是从 V 到 H 的线性有界算子空间，将
L(H, H) 记为 L(H)；L
(2)
(V, H) 是从 V × V 到 H 的双线性有界算子空间。
下面给出 Q − W iener 过程的表示定理：
定理 (表示定理) 设 J 是一个有限或者可数集。(η
j
)
j∈J
∈ V 是 Q : V → V 的特征函数，使
得 Q
ij
= µ
j
η
j
(j ∈ J)，且 (η
j
)
j∈J
构成空间 V 的一组正交基底。我们可以将 W
t
表示为
W
t
(w) =
j∈J
µ
j
=0
√
µ
j
β
j
t
(w)η
j
其中：(β
j
)
j∈J,µ
j
=0
是 (Ω, F, {F
t
}, P ) 上定义的独立的实值 Brown 运动。
上面给出了 SPDE 的一般形式。和 PDE 一样，SPDE 的解亦有强解弱解适应解之分，下面，
给出 SPDE 适应解的定义以及在什么条件下 SPDE 存在唯一适应解。强解和弱解的定义可以参
考 1992.Prato 和 Zabozyk52 和 2007.Prewt53
定义 (适应解) 一个 H
−
取值的可料随机过程 X
t
，∀t ∈ [0, T ]，有
X
t
= S(t)ξ +
t
0
S(t − s)F (X
s
)ds +
t
0
S(t − s)G(X
s
)dW
s
http://www.ma-xy.com 55 http://www.ma-xy.com
http://www.ma-xy.com
1.6 SDE 结构的进一步讨论 第一章 随机微分方程
这里，上式两个积分必须是适定的。
定理 (适应解的存在唯一性) 对 F, G，若 F, G 满足
1. 全局 L 条件：∀x, y ∈ H，∃L > 0，有
||F (x) − F (y)||
2
H
∨ ||G(x) − G(y)||
2
V,H
⩽ L||x − y||
2
H
2. 线性增长条件：∀x ∈ H，∃K > 0，有
||F (x)||
2
H
∨ ||G(x)||
2
V,H
⩽ K(1 + ||x||
2
H
)
则 SPDE 存在唯一适应解 X
t
。
上述存在唯一性证明可以参考 1992.Prato52。
SPDE 的数值方法
上面介绍了 SPDE 的一般形式和解的存在唯一性，下面，我们在抛物型 SPDE 的特例上介
绍 SPDE 的数值方法。
抛物型 SPDE 令 d = {1, 2, 3}，H = v = L
2
((0, 1)
d
, R) 为 Hilbert 空间，在 H 上定义如下抛
物型 SPDE
dX
t
(x) =
k
d
j=1
∂
2
∂x
2
j
X
t
(x) + f(x, X
t
(x))
dt + g(x, X
t
(x))dW
t
X
t
∂(0,1)
d
= 0
X
0
(x) = x
0
(x) x ∈ (0, 1)
d
t ∈ [0, T ]
特别地，当 d = 1 时，有
dX
t
(x) =
k
∂
2
∂x
2
j
X
t
(x) + f(x, X
t
(x))
dt + g(x, X
t
(x))dW
t
X
t
(0) = x
0
(x) = 0
X
0
(x) = ξ(x) x ∈ (0, 1)
t ∈ [0, T ]
(1.27)
下面，我们给出 SPDE(1.27) 的数值方法：线性隐式 Euler 谱方法、指数 Milstein 谱方法和
指数 Runge-Kutta 谱方法。
线性隐式 Euler 方法 设 ( I
N
)
N∈N
, (J
K
)
K∈N
分别是 I 和 J 的有限子集序列。对 N ∈ N，我们
定义线性投影算子 P
N
: H → H 为
P
N
(v) =
i∈I
N
e
i
, v
H
e
i
v ∈ H
http://www.ma-xy.com 56 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.6 SDE 结构的进一步讨论
其中：
e
i
, v
H
=
(0,1)
e
i
(x) × v(x)dx
此外，对 ∀K ∈ N，定义截断 Q-Wiener 过程 W
K
: [0, T ] × Ω → V
0
为
W
K
t
(w) =
j∈J
K
µ
j
=0
√
µ
j
β
j
t
(w)η
j
t ∈ [0, T ], w ∈ Ω
其中：(β
j
)
j∈J,µ
j
=0
是 (Ω, F, {F
t
}, P ) 上相互独立的实取值 Brown 运动。对 w ∈ Ω, m = 0, 1, . . . , M −
1，记
∆β
j
m
(w) = β
j
t
m+1
(w) − β
j
t
m
(w)
和
∆W
M,K
m
(w) = W
K
(m+1)T /M
(w) − W
K
mT /M
(w)
=
j∈J
K
µ
j
=0
√
µ
j
∆β
j
m
(w)η
j
上式是在 [0, T ] 区间上 M 等分。线性隐式 Euler 谱方法为：对
¯
Y
N,M,K
0
= P
N
(ξ)，当 m =
0, 1, . . . , M − 1 时
¯
Y
N,M,K
m+1
= P
N
(I − hA)
−1
¯
Y
N,M,K
m
+ hf(·, Y
N,M,K
m
)
+ g(·,
¯
Y
N,M,K
m
) × ∆W
M,K
m
这里的 N, M, K 皆为划分方式。其中：∀x ∈ (0, 1) 和函数 v, w : (0, 1) → R，φ : (0, 1) × R → R。
定义 φ(·, v) : (0, 1) → R
(φ(·, v))(x) = φ(x, v(x))
定义 v × w : (0, 1) → R 为
(v × w)(x) = v(x) × w(x)
指数 Milstein 谱方法 2015.Jentzen,Rockner 在文26中构造了 SPDE 的 Milstein 方法。对
¯
Y
N,M,K
0
= P
N
(ξ)，当 m = 0, 1, . . . , M − 1 时
¯
Y
N,M,K
m+1
= P
N
S(h)
¯
Y
N,M,K
m
+ hf(·,
¯
Y
N,M,K
m
) + g(·,
¯
Y
N,M,K
m
) × ∆W
M,K
m
+
1
2
∂
∂y
g
(·,
¯
Y
N,M,K
m
) × g(·,
¯
Y
N,M,K
m
) ×
(∆W
M,K
m
)
2
− h
K
j=1
µ
j
(η
j
)
2
http://www.ma-xy.com 57 http://www.ma-xy.com
http://www.ma-xy.com
1.6 SDE 结构的进一步讨论 第一章 随机微分方程
其中：∀x ∈ (0, 1)，对函数 v, w : (0, 1) → R，定义 v
2
: (0, 1) → R 为
(v
2
)(x) = (v(x))
2
同时，定义 (
∂
∂y
g)(·, v) : (0, 1) → R 为
∂
∂y
g
(·, v)(x) =
∂
∂y
g(x, v(x))
对于连续可微函数 g : (0, 1) × R → R，
∂
∂y
g : (0, 1) × R → R 是 g 关于第二个变量的偏导数。
指数 Runge-Kutta 谱方法 将 Milstein 方法中的偏导数替换为差分，得到一个无偏导数的
R-K 方法：对
¯
Y
N,M,K
0
= P
N
(ξ)，当 m = 0, 1, . . . , M − 1 时
¯
Y
N,M,K
m+1
= P
N
S(h)
¯
Y
N,M,K
m
+ hf(·, Y
N,M,K
m
) + g(·,
¯
Y
N,M,K
m
) × ∆W
M,K
m
+
1
2
√
h
g
·,
¯
Y
N,M,K
m
+
√
hg(·,
¯
Y
N,M,K
m
)
− g(·,
¯
Y
N,M,K
m
)
×
(∆W
M,K
m
)
2
− h
K
j=1
µ
j
(η
j
)
2
其中：对所有函数 v : (0, 1) → R，g : (0, 1) ×R → R，定义 g(·, v +
√
hg(·, v)) −g(·, v) : (0, 1) → R
为
g
·, v +
√
hg(·, v)
− g(·, v)
(x)
= g
x, v(x) +
√
hg(x, v(x))
− g(x, v(x))
更详细的内容可以参考 2012. 王小捷57第六章。
数值实验 考虑一维抛物型 SPDE(1.27)，令 ξ(x) = 0，x ∈ (0, 1)，T = 1，k =
1
100
，且 f (x, y) = 1−y，
g(x, y) =
1 − y
1 + y
2
在令
µ
j
=
1
j
2
η
j
(x) = e
j
(x) =
√
2 sin(jπx) x ∈ (0, 1), y ∈ R, j ∈ N
具体的 SPDE 为
dX
t
(x) =
1
100
∂
2
∂x
2
j
X
t
(x) + 1 − X
t
(x)
dt +
1 − X
t
(x)
1 + X
t
(x)
2
dW
t
X
t
(0) = x
0
(x) = 0
X
0
(x) = ξ(x) = 0 x ∈ (0, 1)
t ∈ [0, 1]
¬Euler 方法的 Matlab 程序：
http://www.ma-xy.com 58 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.6 SDE 结构的进一步讨论
1 %% 此 程 序 是SPDE模 型 的 Euler 方 法
2 N = 128;
3 M = Nˆ3;
4 A = −pi ˆ2 * ( 1 :N) . ˆ2 /100;
5 Y = z er os (1 ,N) ;
6 mu = (1 :N) . ˆ −2; %\mu_j( j \ i n \mathbb{N})
7 f = @(x) 1−x ;
8 g = @(x) (1−x ) ./(1+ x . ˆ2 ) ;
9 f o r m = 1 :M
10 y = dst (Y) * s q r t ( 2 ) ;
11 dW = dst ( randn (1 ,N) . * s qr t (mu*2/M) ) ;
12 y = y + f ( y) /M + g(y ) .*dW;
13 Y = id s t ( y ) / s q r t (2) . / (1 − A/M ) ;
14 end
15 plo t ( ( 0 :N+1)/(N+1) , [ 0 , dst (Y) * sq rt ( 2 ) , 0 ] ) ;
16
Milstein 方法的 Matlab 程序：
1 %% 此 程 序 是SPDE的 M il s t e i n 方 法
2 N = 128 ;
3 M = Nˆ2 ;
4 A = −pi ˆ2 * ( 1 :N) . ˆ2/100 ;
5 Y = z er os (1 ,N) ;
6 mu = ( 1 :N) . ˆ −2 ;
7 f = @(x) 1−x ;
8 g = @(x) (1−x ) ./(1+ x . ˆ 2 ) ;
9 bb = @(x) (1−x) . * ( x . ˆ2−2*x−1)/2./(1+x . ˆ2 ) . ˆ3 ; %如 公 式 所 述
10 et a = ze ro s (1 ,N) ;
11 f o r n=1:N
12 eta = eta+2*s in (n * ( 1 :N) /(N+1)* pi ) . ˆ2*mu(n) /M ;
13 end
14 f o r m = 1 :M
15 y = dst (Y) * s qr t ( 2 ) ;
16 dW = dst ( randn (1 ,N) .* sq r t (mu*2/M) ) ;
17 y = y + f ( y ) /M + g(y ) . *dW + bb( y) . * (dW. ˆ2 − eta ) ;
18 Y = exp (A/M ) . * i d s t ( y) / s qr t ( 2 ) ;
19 end
20 plo t ( (0 :N+1)/(N+1) , [ 0 , dst (Y) * s q r t ( 2 ) , 0 ] ) ;
21
®R-K 方法的 Matlab 程序：
1 %% 此 程 序 是SPDE的R−K方 法
2 N = 128 ;
3 M = Nˆ2 ;
4 A = −pi ˆ2 * ( 1 :N) . ˆ2/100 ;
5 Y = z er os (1 ,N) ;
6 mu = (1 :N) . ˆ −2 ;
7 f = @(x) 1−x ;
8 g = @(x) (1−x ) ./(1+ x . ˆ 2 ) ;
9 et a = ze ro s (1 ,N) ;
10 SqrM = sq rt (M) ;%h
11 f o r n=1:N
http://www.ma-xy.com 59 http://www.ma-xy.com
http://www.ma-xy.com
1.6 SDE 结构的进一步讨论 第一章 随机微分方程
12 eta = eta+2*s in (n * ( 1 :N) /(N+1)* pi ) . ˆ2*mu(n) /M ;
13 end
14 f o r m = 1 :M
15 y = dst (Y) * s qr t ( 2 ) ;
16 dW = dst ( randn (1 ,N) .* sq r t (mu*2/M) ) ;
17 g_eva = g ( y) ;
18 y = y + f ( y ) /M + g(y ) . *dW + 0.5*SqrM*( g ( y+g ( y ) /SqrM) − g ( y ) ) . * (dW. ˆ2 − eta ) ;
19 Y = exp (A/M ) . * i d s t ( y) / s qr t ( 2 ) ;
20 end
21 plo t ( (0 :N+1)/(N+1) , [ 0 , dst (Y) * s q r t ( 2 ) , 0 ] ) ;
22
上面给出了 Euler、Milstein 和 R-K 方法的 matlab 仿真程序，下面，我们对上述算法进行
说明。以 R-K 为示例：
1. 首先，对 m = 0, 1, . . . , M − 1，x ∈ (0, 1)，定义 ξ
m
: (0, 1) → R 为
ξ
m
(x) =
¯
Y
N,M,K
m
(x) + hf(x, Y
N,M,K
m
(x)) + g(x,
¯
Y
N,M,K
m
(x)) × ∆W
M,K
m
(x)
+
1
2
√
h
g
x,
¯
Y
N,M,K
m
(x) +
√
hg(x,
¯
Y
N,M,K
m
(x))
− g(x,
¯
Y
N,M,K
m
(x))
×
(∆W
M,K
m
(x))
2
− h
K
j=1
µ
j
(η
j
(x))
2
2. 然后，对 m = 0, 1, . . . , M − 1，R-K 算法可写为
¯
Y
N,M,K
m+1
= P
N
(S(h), ξ
m
)
=
N
j=1
e
Ah
ξ
m
, e
j
H
e
j
=
N
j=1
e
−λ
j
h
ξ
m
, e
j
H
e
j
其初始值为
¯
Y
N,M,K
0
= P
N
(ξ) =
N
j=1
ξ, e
j
H
e
j
其中：H = L
2
((0, 1); R)，A 是乘以常数 K > 0 的零边值条件拉普拉斯算子。所以 A 的特征函
数和特征值分别为
e
j
(x) =
√
2 sin(jπx)
λ
j
= kπ
2
j
2
, x ∈ (0, 1), j ∈ N
对于 Fourier 分量 j = 1 , 2, . . . , N ，有
¯
Y
N,M,K
m+1
, e
j
H
= e
−λ
j
h
ξ
m
, e
j
H
=
√
2e
−λ
j
h
1
0
ξ
m
(x) sin(jπx)dx m = 0, 1, . . . , M − 1
http://www.ma-xy.com 60 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.7 MATLAB 与随机微分方程
和
¯
Y
N,M,K
0
, e
j
H
=
ξ, e
j
H
=
√
2
1
0
ξ(x) sin(jπx)dx j = 1, 2, . . . , N
因此，R-K 的实现步骤为：
Step1. 给定
¯
Y
N,M,K
m
，计算 ξ(m)。
Step2.
使用某种数值积分方法
(
复合梯形公式
)
计算内积
¯
Y
N,M,K
m+1
, e
j
H
, j
= 1
,
2
, . . . , N
。
Step3. 根据
ξ, e
j
H
, j = 1, 2, . . . , N ，求取
¯
Y
N,M,K
m+1
。
因为特征函数 e
j
(x) =
√
2 sin(jπx) 是正弦函数，可以用 MATLAB 命令 dst 和 idst 简化上
述过程。
1.dst 是离散正弦变换。将 N 个实数 z(k) 变为 N 个上述 y(j)(k, j = 1, . . . , N)
y(j) =
N
k=1
z(k) sin
jπ
k
N + 1
j = 1 , . . . , N
2.idst 是反离散正弦变换。将 N 个实数 z(k) 变为 N 个上述 y(j)
y(j) =
2
N + 1
N
k=1
z(k) sin
jπ
k
N + 1
j = 1 , . . . , N
在上式中令 z(k) = ξ
m
k
N+1
, k = 1, . . . , N ，则 y(j) 其实是一个数值求解
√
2⟨ξ
m
, e
j
⟩ 的复合梯
形公式，所以
¯
Y
N,M,K
m+1
, e
j
H
≈ e
−λ
j
h
× y(j)/
√
2 j = 1, . . . , N
在算法中：
Step1. 计算网格点
k
N+1
(k = 1, . . . , N) 上
¯
Y
N,M,K
m
(x)，的 N 个函数值，和 ξ
m
(x) 的函数值。
Step2. 然后，用 idest 近似计算
¯
Y
N,M,K
m+1
, e
j
H
, (j = 1, , . . . , N )。
Step3. 最后，调用 dst 来计算
¯
Y
N,M,K
m+1
的 N 个函数值。
Step4. 重复 Step1 到 Step3，直到得到
¯
Y
N,M,K
m
。
1.7 MATLAB 与随机微分方程
1.7.1 MATLAB 自带工具箱
MATLAB 中已经自带了 SDE 工具箱，这个工具箱只能模拟部分 (普通)SDE 模型，并且不
能进行参数估计。下面，给出其应用实例。
1 % 此 文 件 用 于 matlab 自 带 的SDE工 具 箱 的 用 法
2 % https :// cn . mathworks . com/help/ f in an ce / examples/ pri ci n g−american−basket−options−by−monte
−carl o−simula t i on . html
3 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
http://www.ma-xy.com 61 http://www.ma-xy.com
http://www.ma-xy.com
1.7 MATLAB 与随机微分方程 第一章 随机微分方程
4 % sde 的 一 般 形 式 ， 连 续 ： dXt=F( t , Xt) dt+G( t , Xt)dWt
5 % sde 的 一 般 形 式 ， 离 散 ： Xt+Δt=Xt+F( t , Xt)Δt+G( t , Xt)GΔtZ( t , Xt)
6 % 一 些 sde 的 特 殊 形 式 ：BM、GBM、CEV、CIR、HWV和 Heston 模 型 。
7 %
8 % 下 面 介 绍BM、GBM、CEV、CIR、HWV和 Heston 模 型 的 一 般 形 式
9 % Brownian Motion (BM) ：dXt=A( t ) dt+V( t )dWt
10 %
11 % Geometric Brownian Motion (GBM) ： dXt=B( t ) Xtdt+V( t )XtdWt
12 %
13 % Constant E l a s t i c i t y of Variance (CEV) ： dXt=B( t )Xtdt+V( t )X� ( t )tdWt
14 %
15 % Cox−I n g e r s o l l −Ross (CIR) ： dXt=S( t ) (L( t ) ?Xt) dt+V( t )X12tdWt
16 %
17 % Hull−White/ Vasicek (HWV) ：dXt=S( t ) (L( t ) ?Xt) dt+V( t )dWt
18 %
19 % Heston ：
20 % dX1t=B( t )X1tdt+GX2tX1tdW1t
21 % dX2t=S( t ) [ L( t ) ?X2t ] dt+V( t )GX2tdW2t
22 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
23 % in t e r p o l a t e%随 机 插 值 方 法 ( 即 Brownian 桥 )
24 %% 示 例 1 ： 模 拟 股 票 价 格
25 %1 、 加 载 matlab 的 股 票 数 据
26 load Data_GlobalIdx2
27 p r i c e s = [ Dataset .TSX Dataset .CAC Dataset .DAX . . .
28 Dataset .NIK Dataset .FTSE Dataset .SP ] ;
29 %2 、 转 化 为 日 收 益 率
30 r etu rns = ti c k 2 r e t ( p r i c e s ) ;%另 注 ： p r i c e 2 r e t ( data ) 表 示 将 日 指 数 转 化 为 日 收 益 率%r e t 2 p r i c e (
data ) 表 示 将 日 收 益 率 转 化 为 日 指 数 : )
31 % 3 、 计 算 收 益 率 信 息
32 nVariables = si z e ( returns , 2) ;
33 expReturn = mean( re t ur n s ) ;
34 sigma = std ( r etu r ns ) ;
35 c o r r e l a t i o n = co r r c o e f ( retu r ns ) ;
36 t = 0 ;
37 X = 100;
38 X = X( ones ( nVariables , 1 ) ) ;
39 % 4 、 创 建 sde 模 型
40 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
41 %高 版 本 的 matlab 将 许 多 模 型 都 模 模 式 化 了 ， 要 先 声 明 ( 定 义 ) 一 个 模 型 类 ， 然 后 解 模 型 。
42 %下 面 的 sde 、 sdedo 、 sdeld 、 cev 、gbm等 都 可 以 用 来 定 义 我 们 要 求 解 的 模 型 ：
43 %−−−−−−−−−−−−−−−−−−−SDE模 型 ： dXt=�Xtdt+D(Xt)�dWt−−−−−−−−−−−−−−−−%
44 %sdeld 、 cev 、gbm都 是 用 于 解 决 某 一 类 型 的SDE， 例 如 ： sde l d 用 来 求 解 线 性 漂 移 模 型 (SDE)
45 % 4 .1 使 用 sde 定 义 一 个 一 般 的 ( 广 泛 形 式 的 )SDE模 型
46 F = @( t ,X) diag ( expReturn ) * X;
47 G = @( t ,X) diag (X) * diag ( sigma ) ;
48 SDE = sde (F, G, ’ Corr ela tio n ’ , c o r r e l a t io n , ’ Start Sta t e ’ , X) ;
49 % 4 .2 使 用 sdedo 定 义 一 个 一 般 的 ( 广 泛 形 式 的 )SDE模 型
50 F = d r i f t ( zer o s ( nVariables , 1 ) , diag (expReturn ) ) ;
51 G = d i f f u s i o n ( ones ( nVariables , 1 ) , diag ( sigma ) ) ;
52 SDEDDO = sdeddo (F, G, ’ Cor r el a ti o n ’ , c o r r e l a t i o n , . . .
53 ’ Star t Sta te ’ , 100) ;
54 % 4 .3 使 用 sdeld 定 义 一 个 线 性 漂 移 模 型
55 SDELD = sdeld ( z er os ( nVariables , 1 ) , diag ( expReturn) , . . .
http://www.ma-xy.com 62 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.7 MATLAB 与随机微分方程
56 ones ( nVariables , 1 ) , diag ( sigma ) , ’ C o rr e la t io n ’ , . . .
57 c o r r e l a t i o n , ’ S t ar t Sta t e ’ , X) ;
58 % 4 .4 使 用 cev 定 义 一 个CEV模 型
59 CEV = cev ( diag ( expReturn ) , ones ( nVariables , 1 ) , . . .
60 diag ( sigma ) , ’ Corr e la t io n ’ , co r r e l a t i o n , . . .
61 ’ Star t Sta te ’ , X) ;
62 % 4 .5 使 用gbm定 义 一 个gbm模 型
63 GBM = gbm( diag (expReturn ) , diag ( sigma ) , ’ Co r re l at i on ’ , . . .
64 c o r r e l a t i o n , ’ S t ar t Sta t e ’ , X) ;
65 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
66 % 5 、 模 拟 股 票 指 数 价 格 的 单 一 路 径
67 nPeriods = 249; %一 个 日 历 年 ( 定 义 为 约 250 个 交 易 日 )
68 dt = 1 ; % 时 间 增 量 = 1 day
69 rng (142857 , ’ t w i s t e r ’ )
70 [ S ,T] = simulate (SDE, nPeriods , ’ DeltaTime ’ , dt ) ;%用 默 认 的 模 拟 方 法 对SDE进 行 模 拟
71 % [ S ,T] = simByEuler (SDE, nPeriods , ’ DeltaTime ’ , dt );%用 Euler 方 法 进 行 模 拟
72 f i g u r e
73 p lo t (T, S) , x l a be l ( ’ 交 易 日 ’ ) , y l a b e l ( ’ 价 格 ’ )
74 t i t l e ( ’ 多 维 市 场 模 型 的 单 路 径 ’ )
75 legend ({ ’Canada ’ ’ France ’ ’Germany ’ ’ Japan ’ ’UK’ ’US ’ } , ’ Location ’ , ’ Best ’ )
76 [X,T] = simulate (GBM, nPeriods , ’ DeltaTime ’ , dt , . . .
77 ’ nTri a ls ’ , 10000) ;%用 默 认 的 模 拟 方 法 对GBM进 行 模 , 为 了 说 明 重 载 欧 拉 逼 近 法 的 性 能 优 势 ， 增 加
试 验 次 数 为 10000
78 f i g u r e
79 histogram ( squeeze (X( end , 1 , : ) ) , 30) , x l a b e l ( ’ P r i ce ’ ) , y la b e l ( ’ Frequency ’ )
80 t i t l e ( ’ Histogram of P ri ce s a f t e r One Year : Canada (TSX Composite ) ’ )
81 % 对SDE和GBM各 模 拟 10次 路 径 ， 绘 制 第 一 次 的 路 径 ， 并 绘 制 二 者 的 第 一 个 股 票 价 格 估 计 差
82 rng (142857 , ’ t w i s t e r ’ )
83 [ S ,T] = simulate (SDE, nPeriods , ’ DeltaTime ’ , dt , ’ n T ria ls ’ , 10) ;
84 rng (142857 , ’ t w i s t e r ’ )
85 [X,T] = simBySolution (GBM, nPeriods , . . .
86 ’ DeltaTime ’ , dt , ’ nTr i al s ’ , 10) ;%解 析 解 模 拟
87 %
88 subplot (2 , 1 ,1 )
89 p lo t (T, S ( : , : , 1 ) ) , x l a be l ( ’ 交 易 日 ’ ) , y l a b e l ( ’ P r i ce ’ )
90 t i t l e ( ’ 多 个 股 票 的 第 一 次 的 模 拟 路 径−−欧 拉 逼 近 ’ )
91 subplot (2 , 1 ,2 )
92 p lo t (T, X( : , : , 1 ) ) , x l a b e l ( ’ 交 易 日 ’ ) , y la b e l ( ’ P r ice ’ )
93 t i t l e ( ’ 多 个 股 票 的 第 一 次 的 模 拟 路 径−−解 析 解 ’ )
94 %
95 f i g u r e
96 p lo t (T, S ( : , 1 , 1 ) − X( : , 1 , 1 ) , ’ blue ’ ) , g r id ( ’ on ’ )
97 xl a b el ( ’ 交 易 日 ’ ) , y l ab e l ( ’ 模 拟 误 差 ’ )
98 t i t l e ( ’ 第 一 支 股 票 第 一 次 模 拟 的 Euler 方 法 模 拟 误 差 ’ )
99 %% 示 例 2 ： 模 拟 股 票 价 格−基 于 诱 导 路 径 和 相 关 性
100 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
101 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−dXt=�Xtdt+�XtdWt−−−−−−−−−−−−−−−−−−%
102 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−
%
103 %1 、 加 载 matlab 的 股 票 数 据
104 load Data_GlobalIdx2
105 p r i c e s = [ Dataset .TSX Dataset .CAC Dataset .DAX . . .
106 Dataset .NIK Dataset .FTSE Dataset .SP ] ;
107 %2 、 转 化 为 日 收 益 率
http://www.ma-xy.com 63 http://www.ma-xy.com
http://www.ma-xy.com
1.7 MATLAB 与随机微分方程 第一章 随机微分方程
108 r etu rns = ti c k 2 r e t ( p r i c e s ) ;%另 注 ： p r i c e 2 r e t ( data ) 表 示 将 日 指 数 转 化 为 日 收 益 率%r e t 2 p r i c e (
data ) 表 示 将 日 收 益 率 转 化 为 日 指 数 : )
109 %3 、 使 用 第 一 种 方 法 指 定 西 格 玛 和 相 关 性
110 expReturn = diag (mean( retu r ns ) ) ;
111 sigma = diag ( std ( ret urn s ) ) ;%设 置 西 格 玛 作 为 一 个 对 角 矩 阵 的 资 产 回 报 标 准 偏 差
112 c o r r e l a t i o n = co r r c o e f ( retu r ns ) ;%设 置 相 关 作 为 这 些 回 报 的 样 本 相 关 矩 阵 。 在 这 种 情 况 下 ， 布 朗 运
动 的 组 件 是 依 赖
113 GBM1 = gbm( expReturn , sigma , ’ Cor r el a ti o n ’ , c o r r e l a t i o n ) ;
114 %4 、 使 用 第 二 种 技 术 指 定 西 格 玛 和 相 关 性 ：
115 cov a r i a nce = cov ( retu r ns ) ;
116 sigma = cholcov ( covari a n c e ) ’ ;%指 定 西 格 玛 的 资 产 收 益 率 的 协 方 差 矩 阵 的 Cholesky 因 子
117 GBM2 = gbm( expReturn , sigma ) ;
118 % 注 ： 在 这 里 ， 西 格 玛 捕 获 的 资 产 回 报 的 不 确 定 性 的 相 关 性 和 幅 度 。 与 第 一 种 技 术 相 比， 布 朗 运 动 的
组 件 是 独 立 的 。 此 外 ， 这 种 技 术 接 受 的 身 份 矩 阵 的 默 认 分 配 的 相 关 性 ， 并 更 简 单 。
119 %5 、 用 欧 拉 逼 近 法 对GBM1和GBM2进 行 路 径 模 拟
120 rng (22814 , ’ t w i s t e r ’ )
121 [ X1,T] = simByEuler (GBM1, 1 00 0 ) ; % co r r e l a t e d Brownian motion
122 rng (22814 , ’ t w i s t e r ’ )
123 [ X2,T] = simByEuler (GBM2, 1 00 0 ) ; % standard Brownian motion
124 subplot (2 , 1 ,1 )
125 p lo t (T, X1)
126 t i t l e ( ’ Sample Paths from Correlated Brownian Motion ’ )
127 yl a b el ( ’ Asset Pri c e ’ )
128 subplot (2 , 1 ,2 )
129 p lo t (T, X2)
130 t i t l e ( ’ Sample Paths from Standard Brownian Motion ’ )
131 xl a b el ( ’ Trading Day ’ )
132 yl a b el ( ’ Asset Pri c e ’ )
133 %% 示 例 3 ： 模 拟 股 票 价 格−基 于 动 态 模 拟 ( 系 数 是 时 变 的 )
134 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
135 %−−−−−−−−−−−−−−−−−−−−−−−−−dXt=� ( t )Xtdt+� ( t )XtdWt−−−−−−−−−−−−−−−−%
136 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
137 load Data_GlobalIdx2
138 dt = 1/250;
139 r etu rns = ti c k 2 r e t ( Dataset .CAC) ;
140 sigma = std ( r etu r ns )* sq r t (250) ;
141 y i e l d s = Dataset .EB3M;
142 y i e l d s = 360* lo g (1 + y i e l d s ) ;
143 nPeriods = length ( y i e l d s ) ;
144 rng (5713 , ’ t w i s t e r ’ )
145 obj = gbm(mean( y i e l d s ) , diag ( sigma ) , ’ Star t Sta te ’ ,100) ;%第 一 种 方 法 指 定 风 险 中 性 的 收 益 为
Euribor 收 益 样 本 平 均 值 ， 因 此 假 定 一 个 恒 定 的 ( 非 动 态 ) 无 风 险 收 益 率 ：
146 [ X1,T] = simulate ( obj , nPeriods , ’ DeltaTime ’ , dt ) ;
147 r = ts2func ( yie l d s , ’ Times ’ , ( 0 : nPeriods − 1) ’) ;
148 rng (5713 , ’ t w i s t e r ’ )
149 obj = gbm( r , diag ( sigma ) , ’ Start Sta t e ’ ,100) ;
150 X2 = simulate ( obj , nPeriods , ’ DeltaTime ’ , dt ) ;%模 拟 时 使 用 相 同 的 初 始 随 机 数
151 subplot (2 , 1 ,1 )
152 p lo t ( dates ,100* y i e l d s )
153 d a t e t i c k ( ’ x ’ )
154 xl a b el ( ’ Date ’ )
155 yl a b el ( ’ Annualized Yield (%) ’ )
156 t i t l e ( ’ Risk Free Rate(3−Mo Euribor Continuously−Compounded) ’ )
http://www.ma-xy.com 64 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.7 MATLAB 与随机微分方程
157 subplot (2 , 1 ,2 )
158 p lo t (T, X1, ’ red ’ ,T, X2, ’ blue ’ )
159 xl a b el ( ’Time ( Years ) ’ )
160 yl a b el ( ’ Index Level ’ )
161 t i t l e ( ’ Constant vs . Dynamic Rate o f Return : CAC 40 ’ )
162 legend ({ ’ Constant I n t e r e s t Rates ’ ’ Dynamic I n t e r e s t Rates ’ } , . . .
163 ’ Location ’ , ’ Best ’ )
164 f = Example_BlackScholes ( nPeriods , n T ria ls )
165 rng (88161 , ’ t w i s t e r ’ )
166 X = simBySolution ( obj , nPeriods , ’ DeltaTime ’ , dt , . . .
167 ’ nTri a ls ’ , nTrials , ’ Proces s es ’ , f . BlackScholes ) ;
168 c a l l = mean( exp(−r a t e *T) *max( squeeze (X( end , : , : ) ) − s t r i k e , 0) )
169 put = mean(exp(−rate *T) *max( s t r i k e − squeeze (X( end , : , : ) ) , 0) )
170 f . C al l Pr ic e ( s t r i k e , ra te )
171 f . PutPrice ( s t r i k e , r at e )
172
1.7.2 SDEtoolbox 工具箱
MATLAB 自带的 SDE 工具箱仅能进行 SDE 的模拟，不能进行 SDE 的参数估计。SDEtool-
box 工具箱可以进行 SDE 的模拟也可以进行参数估计，其应用实例如下
1 %% 文 件 展 示 了SDE Toolbox( by Umberto P ic ch in i ) 的 使 用
2 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
3 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−说 明 ( 中 文 )−−−−−−−−−−−−−−−−−−−−−−−−−%
4 % 0 、SDE Toolbox 作 者 ：
5 % umberto . picchini@biomatematica . i t
6 % http : //www. biomatematica . i t /pages/ p i c c h i n i . html
7 % 1 、SDE工 具 箱 可 以 进 行 I to 型 和 Stratonovich 型 随 机 微 分 方 程 的 解 的 样 本 路 径 的 模 拟 ，
8 % 以 及 随 机 微 分 方 程 参 数 估 计
9 % 2 、 可 在 http : // sdetoolbox . s o u r c e f o r g e . net 下 载
10 % 3 、 matlab 版 本 要 求 ：
11 % MATLAB 6 x 以 上 。 操 作 系 统 ： 任 何 Matlab 支 持 的 平 台 。
12 % 这 个 版 本 的SDE工 具 箱 已 经 与 Matlab 6.5 和 7 . 3 测 试 。
13 % 4 、 可 以 在 下 面 的 网 址 获 取 工 具 箱 的 使 用 指 导
14 % http : // sdetoolb ox . s ou rc ef o r g e . net /manual . pdf
15 % −−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
16 % −−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−文 件 说 明−−−−−−−−−−−−−−−−−−−%
17 % −−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
18 % CHANGES SDE工 具 箱 的 改 变 ( 历 史 信 息 )
19 % demo_sdefile SDE模 型 的 定 义 f o r sde_demo .m
20 % fminsearchbnd .m Nelder−Mead优 化 方 法 ( by John D’ Errico , woodchips@rochester .
rr .com)
21 % LICENSE 许 可 证 文 件
22 % manual . pdf 手 册 使 用 指 南
23 % README 自 述 文 件
24 %
25 % m1_sdefile .m 模 型M1的 定 义 . SDE Toolbox 共 定 义 了 10类SDE问 题 ( 当 然 你 可 以 定 义
自 己 的 问 题 )
26 % m2_sdefile .m
27 % m3_sdefile .m
28 % m4_sdefile .m
http://www.ma-xy.com 65 http://www.ma-xy.com
http://www.ma-xy.com
1.7 MATLAB 与随机微分方程 第一章 随机微分方程
29 % m5_sdefile .m
30 % m6_sdefile .m
31 % m7_sdefile .m
32 % m8_sdefile .m
33 % m9_sdefile .m
34 % m10_sdefile .m
35 %
36 % p e r c t i l e .m 百 分 位 数 计 算 器 (by Peter J . Acklam , pjacklam@online . no )
37 % sde_compute_hessian .m 计 算 参 数 Hessian 的 中 心 近 似 估 计 与 参 数 边 界 检 查 ( o r i g i n a l l y
writte n by Andrea De Gaetano , www. biomatematica . i t )
38 % sde_demo .m 演 示 文 件
39 % sde_demo_interactive_settings .m 要 求 用 户 提 供 交 互 sde_demo仿 真 ( 模 拟 ) 设 置
40 %
41 % sde_euler .m SDE数 值 积 分 ( Euler Maruyama方 法 )
42 % sde_euler_demo .m SDE数 值 积 分 ( 欧 拉 方 法 ) ， 创 建 临 时 的sde_demo .m
43 % sde_getdata .m 获 取 数 据 ( ASCII )
44 % sde_graph .m 绘 图 ： 数 据 曲 线 ( 如 果 有 ) 、 经 验 平 均 、 置 信 区 间 、 分 位 数 和 直 方 图
45 % s d e _ i ntegrator .m SDE数 值 积 分 。
46 % sde_kurtosis .m 计 算 样 本 峰 度
47 % sde_library_optsetup .m 优 化 程 序 的 设 置
48 % sde_library_run .m 模 型 l i b r a r y 的 主 程 序 。
49 % sde_library_setup .m l i b r a r y 的 设 置 ： 定 义 状 态 变 量 的 数 目 和 模 型 参 数
50 % sde_makelabel .m 查 看 到 描 述 部 分 为 sde_makelabel .m
51 % sde_milstein .m SDE数 值 积 分 ( 米 尔 斯 坦 方 法 )
52 % sde_milstein_demo .m SDE数 值 积 分 ( 米 尔 斯 坦 方 法 ) ， 创 建 临 时 的sde_demo .m和sde_demo2 .
m
53 % sde_model_description .m SDE模 型 描 述
54 % sde_moment .m 计 算 采 样 时 刻
55 % sde_npsml .m 非 参 数 参 数 估 计 方 法 ( Ito 和 Stra t onovich SDEs)
56 % sde_npsml_euler .m Euler−Maruyama积 分 可 用 于 sde_npsml .m
57 % sde_npsml_milstein .m 米 尔 斯 坦 积 分 可 用 于 sde_npsml .m
58 % sde_param_mask.m 返 回 自 由 参 数 向 量 ( o r i g i n a l l y w r i tten by Andrea De Gaetano ,
www. biomatematica . i t )
59 % sde_param_unmask .m 返 回 完 整 的 ( 免 费+常 数 ) 的 参 数 集 ( 最 初 是 由 Andrea De Gaetano 、
www. biomatematica 。 它 )
60 % sde_parcheck .m 检 查 自 由 参 数 是 否 可 以 接 受 ( 最 初 是 由 Andrea De Gaetano 、www.
biomatematica 。 它 )
61 % sde_parconfint .m 计 算 渐 近 95%的 置 信 区 间 的 近 似 参 数 的 极 大 似 然 估 计
62 % sde_predict .m 模 拟 轨 迹 的 线 性 插 值 用 于 生 成 数 据：
63 % sde_psml .m 参 数 的 最 大 似 然 值 I t o SDEs
64 % sde_psml_euler .m 过 渡 密 度 近 似 使 用 欧 拉 积 分 ( 可 用 于sde_psml 。m)
65 % sampledata1 . dat ASCII 制 表 符 分 隔 的 文 件 包 含 的 数 据 从 一 维 sde 生 成 ( 阅 读 手 册 的 内 容
)
66 % sampledata2 . dat ASCII 制 表 符 分 隔 的 文 件 包 含 的 数 据 从 二 维 sde 生 成 ( 阅 读 手 册 的 内 容
)
67 % sampledata3 . dat ASCII 制 表 符 分 隔 的 文 件 包 含 的 数 据 从 一 维 sde 生 成 ( 阅 读 手 册 的 内 容
)
68 % sde_skewness .m 计 算 样 本 偏 度
69 % sde_split_sdeinput .m 把 s d e f i l e 输 入 状 态 变 量 分 离
70 % sde_stats .m 在 过 程 结 束 时 , 计 算 描 述 性 统 计
71 % sde_xobsmatrix .m 查 看 文 件 sde_xobsmatrix .m的 描 述 部 分
72 %−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−%
73 c lc , c l e a r
http://www.ma-xy.com 66 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.7 MATLAB 与随机微分方程
74 %% 示 例 1 ： i t o 型 随 机 微 分 方 程 模 拟 ( 非 估 计 )
75 % dXt = aXtdt + �dWt , X0 = x0
76 % 其 中 ： t � [ 0 , 1 0 ] , 使 用 Euler−Maruyama方 法 ， 步 长h = 0.01 ， x0 = 1 ，
77 % ( a , � ) = ( 0 . 5 , 0 . 2 ) , 使 用 0 作 为 种 子 的 伪 随 机 维 纳 增 量 的 产 生 。
78 x0 = 1 ; %SDE初 始 条 件
79 a = 0 . 5 ; %参 数a
80 sigma = 0 . 2 ; % 参 数 � = ’ sigma ’
81 problem = ’M1 ’ ; %该 工 具 箱 一 共 有 10 个 可 解 决 的 模 型 ， 需 要 对 照 manual . pdf 进 行 查 看
82 t0 = 0; %初 始 时 间
83 T = 10; %最 终 时 间
84 h = 0. 0 1 ; %用 于 数 值 积 分 的 步 长
85 numsim = 500; %模 拟 次 数
86 sdetype = ’ Ito ’ ; % 在 使 用EM方 法 时 ， sde 的 类 型 type 必 须 是 ’ Ito ’
87 randseed = 0; % 这 是 一 个 产 生 伪 随 机 维 纳 增 量 的 种 子
88 in t e g r a t o r = ’EM’ ; % must be ’EM’ ( Euler−Maruyama, only fo r I t o SDEs) or ’ Mil ’
89 model = ’M1a ’ ;%该 工 具 箱 一 共 有 10 个 可 解 决 的 模 型 ， 需 要 对 照 manual . pdf 进 行 查 看
90 numdepvars = 1 ; %SDE变 量 数
91 yesdata = 0 ; % can be 0 (no raw data a v a i l a b l e ) or 1 ( data a v a i l a b l e )
92 xhat = SDE_euler ( [ x0 , a , sigma ] , problem , [ t0 : h :T] , numdepvars , numsim, sdetype , . . .
93 randseed ) ;% Euler 方 法 模 拟 轨 迹 ‘ xhat ’
94 SDE_graph ( [ x0 , a , sigma ] , xhat , yesdata , problem , sdetype , i n t eg ra to r , numdepvars , . . .
95 [ t0 : h :T] , model , numsim , [ ] , [ ] , randseed )% 绘 制 路 径
96 % 下 面 用 蒙 特 卡 罗 模 拟 计 算 过 程 的 近 似 均 值 、 方 差 、 偏 度 、 峰 度 ，X在T = 10 的 基 础 上 的 模 拟 轨 迹
97 % 例 如 ： 我 们 得 到E(XT) 近 似 为 6.654*10−3, Var(XT ) 近 似 为 0 . 04 3 , Skewness (XT ) 近 似 为 0 and
98 % Kurtosis (XT ) 近 似 为 2 . 4 3 3 .
99 SDE_stats ( [ x0 , a , sigma ] , xhat , problem , [ t0 : h :T] , numdepvars , numsim, sdetype , . . .
100 i n t eg ra to r , randseed )
101 %Ito 型 随 机 微 分 方 程 ，
102 in t e g r a t o r = ’ Mil ’ ;
103 xhat = SDE_milstein ( [ x0 , a , sigma ] , problem , [ t0 : h :T] , numdepvars , numsim, sdetype , . . .
104 randseed ) ;%mi l s t e i n 法 模 拟 轨 迹 ( 路 径 ) 数 据
105 SDE_graph ( [ x0 , a , sigma ] , xhat , yesdata , problem , sdetype , i n t eg ra to r , numdepvars , . . .
106 [ t0 : h :T] , model , numsim , [ ] , [ ] , randseed )
107 %S t ratonovich 型 随 机 微 分 方 程
108 sdetype = ’ S t rat ’ ;
109 xhat = SDE_milstein ( [ x0 , a , sigma ] , problem , [ t0 : h :T] , numdepvars , numsim, sdetype , . . .
110 randseed ) ;
111 SDE_graph ( [ x0 , a , sigma ] , xhat , yesdata , problem , sdetype , i n t eg ra to r , numdepvars , . . .
112 [ t0 : h :T] , model , numsim , [ ] , [ ] , randseed )
113 %% 示 例 2 ： 随 机 微 分 方 程 组 模 拟 ( 非 估 计 )
114 % 我 们 对 下 面 的 2 变 量 It ? SDE随 机 微 分 方 程 组 模 拟 1000次 ,
115 % dX1t = [ � 1 1 � 1 + � 1 2 � 2 −�11Xt 1 −�12Xt 2 ] dt + �1dWt1
116 % dX2t = [ � 2 1 � 1 + � 2 2 � 2 −�21Xt 1 −�22Xt 2 ] dt + �2dWt2
117 % 使 用 Euler−Maruyama method， 并 进 行 如 下 设 置 :
118 % 初 始 条 件 (X01, X02) = (1 , 5) ， 参 数 ( �1 , �2 , �11 , �12 , �21 , � 2 2 ) = ( . 2 , . 5 , . 1 , . 2 , . 1 ,
. 4 , . 3 , . 2 ) , [ t0 , T] = [ 1 , 5 ] , h = 0. 00 5 .
119 c lc , c l e a r
120 parameters = [ 1 , 5 , . 2 , . 5 , . 1 , . 2 , . 1 , . 4 , . 3 , . 2 ] ;
121 problem = ’M9 ’ ;
122 t0 = 1;
123 T = 5 ;
124 h = 0. 0 0 5 ;
125 numsim = 1000 ;
http://www.ma-xy.com 67 http://www.ma-xy.com
http://www.ma-xy.com
1.7 MATLAB 与随机微分方程 第一章 随机微分方程
126 sdetype = ’ Ito ’ ;
127 randseed = 0;
128 in t e g r a t o r = ’EM’ ;
129 model = ’M9a ’ ;%该 工 具 箱 一 共 有 10 个 可 解 决 的 模 型 ， 需 要 对 照 manual . pdf 进 行 查 看
130 numdepvars = 2 ;
131 yesdata = 0 ;
132 xhat = SDE_euler( parameters , problem , [ t0 : h :T] , numdepvars , numsim, sdetype , . . .
133 randseed ) ;
134 SDE_graph( parameters , xhat , yesdata , problem , sdetype , int eg ra t o r , numdepvars , . . .
135 [ t0 : h :T] , model , numsim , [ ] , [ ] , randseed )
136 %% 示 例 3 ： 数 值 模 拟
137 x0 = 5 ;
138 a = 2 ;
139 sigma = 0 . 5 ;
140 problem = ’M1 ’ ;
141 t0 = 0;
142 T = 5 ;
143 h = 0. 0 1 ;
144 numsim = 1000 ;
145 sdetype = ’ Ito ’ ;
146 randseed = 0;
147 in t e g r a t o r = ’EM’ ;
148 model = ’M1a ’ ;
149 numdepvars = 1 ;
150 yesdata = 1 ;
151 xhat = SDE_euler ( [ x0 , a , sigma ] , problem , [ t0 : h :T] , numdepvars , numsim, sdetype , . . .
152 randseed ) ;
153 [ xobs , time ] = SDE_getdata ( ’ sampledata1 ’ ) ;%加 载 数 据 ’ sampledata1 . dat ’
154 SDE_graph ( [ x0 , a , sigma ] , xhat , yesdata , problem , sdetype , i n t eg ra to r , numdepvars , . . .
155 [ t0 : h :T] , model , numsim , time , xobs , randseed )
156 %% 示 例 4 ： 参 数 估 计
157 data = load ( ’ sampledata3 . dat ’ ) ; % load the data
158 time = data ( : , 1 ) ; %时 间 t
159 xobs = data ( : , 2 ) ; %xt
160 v r bl = data ( : , 3 ) ; %标 签 变 量
161 h = 0. 0 0 1 ;
162 owntime = [ time (1 ) : h : time ( end ) ] ; %模 拟 的 时 间 跨 度
163 x0 = xobs ( 1) ; %初 始 值
164 a = 3 ; %
165 b = 0 . 5 ; %
166 sigma = 0 . 1 ; %
167 f r e e p a r s t a r t = [ a , b , sigma ] ; % 自 由 参 数 初 始 估 计 值
168 freeparmin = [1 e −6,1e−6,1e −6]; %要 估 计 的 参 数 的 上 限
169 freeparmax = [ 3 , 1 , 1 ] ; %要 估 计 的 参 数 的 上 限
170 totparmin = [ x0 , freeparmin ] ; %
171 totparmax = [ x0 , freeparmax ] ; % the upper bounds f o r the f i x e d ( x0 ) and the f r e e to vary
parameters
172 parmask = [ 0 , 1 , 1 , 1 ] ; % 用0 标 记 常 量 参 数 x0 ， 用 1 标 记 其 他 情 况 ( 待 求 )
173 parbase = [ x0 , f r e e p a r s t a r t ] ; %常 数 和 自 由 参 数 的 起 始 值 的 数 组
174 problem = ’M2 ’ ;
175 numsim = 2000 ;%模 拟 次 数
176 sdetype = ’ Ito ’ ;
177 in t e g r a t o r = ’EM’ ;
http://www.ma-xy.com 68 http://www.ma-xy.com
http://www.ma-xy.com
第一章 随机微分方程 1.7 MATLAB 与随机微分方程
178 numdepvars = 1 ;
179 randseed = 0;
180 myopt = optimset ( ’ fminsearch ’ ) ; % o p t imizat i o n s e t t in g s , type ’ help optimset ’ f o r
d e t a i l s
181 myopt = optimset (myopt , ’ MaxFunEvals ’ ,20000 , ’ MaxIter ’ ,5000 , ’ TolFun ’ , 1 . e−4, ’ TolX ’ , 1 . e −4 , . . .
182 ’ Display ’ , ’ i t e r ’ ) ;
183 %
184 f r e e p a r e s t = fminsearchbnd ( ’SDE_NPSML’ , fr e e p a r s t a r t , freeparmin , freeparmax , myopt , owntime
, . . .
185 time , vrbl , xobs , problem , numsim, sdetype , parbase , totparmin , totparmax , parmask , . . .
186 i n t eg ra to r , numdepvars , randseed ) ;% ’ f r e e p a r e s t ’ 包 含 近 似 参 数 的 最 大 似 然 估 计
187 %
188 totparam = SDE_param_unmask( fr e ep a r e s t , parmask , parbase ) ;% ’ totparam ’ 包 含 固 定 的 数 组 和 参
数
189 % 计 算 自 由 参 数 的 95%置 信 区 间
190 SDE_ParConfInt( ’SDE_NPSML’ , f re e p a r e s t , owntime , time , vrbl , xobs , problem , numsim, sdetype , . . .
191 parbase , totparmin , totparmax , parmask , in te g r a t o r , numdepvars , randseed ) ;
192 yesdata = 1 ;
193 SDE_graph(totparam , [ ] , yesdata , problem , sdetype , i nt eg ra to r , numdepvars , owntime , [ ] , numsim , . . .
194 time , xobs , 0 ) ;
195 SDE_stats ( totparam , [ ] , problem , owntime , numdepvars , numsim, sdetype , i nt eg ra to r , 0 ) ;
196
http://www.ma-xy.com 69 http://www.ma-xy.com
http://www.ma-xy.com
参考文献
[1] Gardon´ A. The order of approximation for solutions of itoˆ-type stochastic dif- ferential
equations with jumps. 2004.
[2] Gardon´ A. The order 1.5 approximation for solutions of jump-diusion equa- tions. 2006.
[3] Platen E Bruti-Liberati N. On the strong approximation of jump-diusion processes. 2005.
[4] Platen E Bruti-Liberati N. Approximation of jump diusions in nance and economics.
2007.
[5] Riedler M G Buckwar E. Runge-kutta methods for jump-diusion dierential equations.
2011.
[6] 2013. 陈琳. 随机延迟微分方程的数值解稳定性. 华中科技大博文.
[7] Wiener J Cooke K L. Retarded dierential equations with piecewise cinstant delays. 1984.
[8] Platen E. An approximation method for a class of itoˆ processes with jump com- ponent.
1982.
[9] Wu F. Khasminskii-type theorems for neutral stochastic functional dierential equations.
2008.
[10] Gijbels Fan. local polynomial modeling and its application. 1998.
[11] Shen Y Feng L. A new criteria on pth moment asymptotic stability on a class of neutral
stochastic functional dierential equations. 2012.
[12] Florens-Emirou. On estimating the diusion coecient from discrete observations. 1993.
[13] Fridman. Stochastic dierential equations and applications. 1975.
[14] Zhang H Gan S, Schurz H. Mean square convergence of stochastic �-methods for nonlinear
neutral stochastic dierential delay equations. 2011.
[15] Skorokhd Gikhaman. Stochastic dierential equation. 1972.
71
http://www.ma-xy.com
参考文献 参考文献
[16] Reib Gobet, Hoann. Nonparametric estimation of scalor diusion based on low frequency
data is ill-posed. 2003.
[17] MaoXuerong Szpruch Higham, J Desmond. Convergence non-negativity and stability of a
new milstein scheme with applications to nance. 2013.
[18] Stuor Higham, Mao. Strong convergence of euler type methods for nonlinear stochastic
dierential equations. 2002.
[19] Kloeden P E Higham D J. Numerical methods for nonlinear stochastic dierential equations
with jumps. 2005.
[20] Kloeden P E Higham D J. Convergence and stability of implicit methods for jump- diusion
systems. 2006.
[21] Kloeden P E Higham D J. Strong convergence rates for backward euler on a class of nonlinear
jump-diusion problems. 2006.
[22] Zhang Z M Huang C, Shi J. Two-stage stochastic runge-kutta methods for stochastic dif-
ferential equations with jump-diusion. 2012.
[23] Deng F Huang L. Razumikhin-type theorems on stability of neutral stochastic functional
dierential equations. 2008.
[24] Wiener J. Generalized solutions of functional dierential equations. 1993.
[25] Jovanovic M Jankovic S, Randjelovic J. Razumikhin-type exponential stabili- ty criteria of
neutral stochastic functional dierential equations. 2009.
[26] Rockner Jentzen. A milstein scheme for stochastic partical dierential equations. 2015.
[27] AV Swishchuk，YI Kazmerchu. Stability of stochastic dierential delay ito’s equations with
poisson jumps and with markovian switchings. 2001.
[28] Wu Xiu ke. Almost sure exoinential stability of numerical solutions for stochastic delay
dierential equations. 2010.
[29] Khasminskii. Stochastic stability of dierential equations. 1980.
[30] Khasminskii. Stochastic dierential equations and their applications. 1997.
[31] Platen Kloedn. Numerical solution of stochastic dierential equation. 1992.
[32] Maizenberg T Mao X Matasov A Kolmanovskii V, Koroleva N. Neutral stochastic dierential
delay equation with markovian switching. 2003.
http://www.ma-xy.com 72 http://www.ma-xy.com
http://www.ma-xy.com
参考文献 参考文献
[33] Nosov V Kolmanovskii V. Stability and periodic modes of control systems with aftereect.
1981.
[34] Shaikkhet L. Some new aspects of lyapunov-type theorems for stochastic dif- ferential
equations of neutral type. 2010.
[35] 2012. 李洋. 正倒向随机微分方程的高精度数值方法及误差估计. 山大博文.
[36] Li C W Liu X Q. Almost sure convergence of the numerical discretisation of stochastic jump
diusions. 2000.
[37] Shen Y Luo Q, Mao X. New criteria on exponential stability of neutral stochastic dierential
delay equations. 2006.
[38] Hutzenthaler M. Strong convergence of an explicit numerical method for sdes with nonglob-
ally lipschitz continuous coecients. 2012.
[39] Milosevic M. Highly nonlinear neutral stochastic dierential equations with time-dependent
delay and the euler-maruyama methods. 2011.
[40] Milosevic M. Almost sure exponential stability of solutions to highly nonlin- ear neutral
stochastic dierential equations with time-dependent delay and the euler-maruyama approx-
imation. 2013.
[41] Harris C J Maghsoodi Y. In probability approximation and simulation of non- linear jump-
diusion stochastic dierential equations. 1987.
[42] Mao. Numerical solutions of stochastic dierential delay equations under local lipschitz
condition. 2003.
[43] Mao. Exponential stability of equdistent euler-maruyama approximations of sddes. 2007.
[44] Mao. Numerical solutions of stochastic dierential delay equations under the generalized
khasmihskii-type condition. 2011.
[45] Rodkina A Mao X. Razumikhin-type theorems for neutral stochastic functional dierential
equations. 1997.
[46] Rodkina A Mao X. Stochastic dierential equations and their applications. 1997.
[47] Yuan C Mao X. Stochastic dierential equations with markovian switching. 2006.
[48] Yuan C Mao X, Shen Y. Almost surely asymptotic stability of neutral stochastic dierential
delay equations with markovian switching. 2008.
[49] Milstein. Numerical integration of stochastic dierential equation. 1995.
http://www.ma-xy.com 73 http://www.ma-xy.com
http://www.ma-xy.com
参考文献 参考文献
[50] Tretyakov Milstein. Stochastic numerics for mathematical physics. 2004.
[51] Zubchenko V P Mishura Y S. Rate of convergence in the euler scheme for s- tochastic
dierential equations with non-lipschitz diusion and poisson mea- sure. 2011.
[52] Zabozyk Prato. Stochastic equations in innite dimersions. 1992.
[53] Rockner Prewt. A concise covrse on stochastic partical dierential equations. 2007.
[54] Jankovic S Randjelovic J. On the pth moment exponential stability criteria of neutral
stochastic functional dierential equations. 2007.
[55] Stanton. A nonparametric model of term structure dynamics and the market price of interest
rate risk. 1997.
[56] Li C W. Almost sure convergence of stochastic dierential equations of jump- diusion type.
1995.
[57] 2012. 王小捷. 随机微分方程数值算法研究. 中山大学博士论文，第四章.
[58] Gan Wang. The tamed milstein method for commutative stochastic dierential equations
with non-globally lipschitz continuous coecients. 2013.
[59] Chen Y Wang W. Mean-square stability of semi-implicit euler method for nonlinear neutral
stochastic dierential delay equations. 2011.
[60] Van J R Wang Y B. Necessary and sucient condition for the global attractivity of the
trivial solution of a delay equation with continuous and piecewise constant arguments. 1997.
[61] Mao X. Exponential stability in mean square of neutral stochastic dierential functional
equations. 1995.
[62] Mao X. Stability of stochastic dierential equations with markovian switching. 1999.
[63] Mao X. Asymptotic properties of neutral stochastic dierential delay equation- s. 2000.
[64] Mao X. A note on almost sure asymptotic stability of neutral stochas- tic dierential delay
equations with markovian switching. 2012.
[65] Hu S Xue M, Zhou S. Stability of nonlinear neutral stochastic functional dif- ferential
equations. 2010.
[66] Maghsoodi Y. Mean square ecient numerical solution of jump-diusion s- tochastic dier-
ential equations. 1996.
[67] Maghsoodi Y. Exact solutions and doubly ecient approximations and simu- lation of
jump-diusion itoˆ equations. 1998.
http://www.ma-xy.com 74 http://www.ma-xy.com
http://www.ma-xy.com
参考文献 参考文献
[68] Fan Yao. Nonlinear time series:nonparametric and parametric methods. 2003.
[69] Ma Z Yin B. Convergence of the semi-implicit euler method for neutral stochas- tic delay
dierential equations with phase semi-markovian switching. 2011.
[70] 2013. 于辉. 带泊松测度随机微分方程数值解的收敛性和稳定性. 哈工大博士论文.
[71] Mao X Yuan C. onvergence of the euler–maruyama method for stochastic dierential
equations with markovian switching. 2006.
[72] 2013. 张玲. 几类随机延迟微分方程的数值分析. 哈工大博文.
[73] Fang Z Zhou S. Numerical approximation of nonlinear neutral stochastic func- tional dif-
ferential equationss. 2012.
[74] Liu Ming Zhu. Convergence and stability of the semi-implicit euler method for a linear
sddes. 2004.
[75] 2014. 宗小峰. 随机微分方程的数值分析及随机稳定化. 华中科技大学博士论文.
http://www.ma-xy.com 75 http://www.ma-xy.com
