第十七章马氏链模型 §1随机过程的概念 个随机试验的结果有多种可能性,在数学上用一个随机变量(或随机向量)来描 述。在许多情况下,人们不仅需要对随机现象进行一次观测,而且要进行多次,甚至接 连不断地观测它的变化过程。这就要研究无限多个,即一族随机变量。随机过程理论就 是硏究随机现象变化过程的概率规律性的 定义1设{1,t∈T}是一族随机变量,T是一个实数集合,若对任意实数t∈T5 是一个随机变量,则称{1I∈T}为随机过程 T称为参数集合,参数t可以看作时间。5的每一个可能取值称为随机过程的一个 状态。其全体可能取值所构成的集合称为状态空间,记作E。当参数集合T为非负整 数集时,随机过程又称随机序列。本章要介绍的马尔可夫链就是一类特殊的随机序列。 例1在一条自动生产线上检验产品质量,每次取一个,“废品”记为1,“合格品” 记为0。以5n表示第n次检验结果,则n是一个随机变量。不断检验,得到一列随机 变量51252…,记为{n,n=1,2,…}。它是一个随机序列,其状态空间E={0,1} 例2在m个商店联营出租照相机的业务中(顾客从其中一个商店租出,可以到m 个商店中的任意一个归还),规定一天为一个时间单位,“51=j”表示“第1天开始营 业时照相机在第j个商店”j=1,2,…,m。则{n,n=12,…}是一个随机序列,其状 态空间E={1,2…,m}。 例3统计某种商品在t时刻的库存量,对于不同的t,得到一族随机变量, {5,t∈[0,+∞)}是一个随机过程,状态空间E=[0,R],其中R为最大库存量 我们用一族分布函数来描述随机过程的统计规律。一般地,一个随机过程 51,t∈T},对于任意正整数n及T中任意n个元素,…,n相应的随机变量54,…,5 的联合分布函数记为 F.(x1,…,xn)=P{5,≤x1…,.≤xn} (1) 由于n及t(i=1,…,n)的任意性,(1)式给出了一族分布函数。记为 {F….(x1,…,xn),1∈T,i=1,…,n,n=12… 称它为随机过程{,t∈T}的有穷维分布函数族。它完整地描述了这一随机过程的统计 规律性。 §2马尔可夫链 2.1马尔可夫链的定义 现实世界中有很多这样的现象:某一系统在已知现在情况的条件下,系统未来时刻 的情况只与现在有关,而与过去的历史无直接关系。比如,研究一个商店的累计销售额 如果现在时刻的累计销售额已知,则未来某一时刻的累计销售额与现在时刻以前的任一 时刻累计销售额无关。上节中的几个例子也均属此类。描述这类随机现象的数学模型称 为马氏模型。 定义2设{n,n=1,2,…}是一个随机序列,状态空间E为有限或可列集,对于任 意的正整数m,n,若,j∈E(k=1,…,n-1),有
-207- 第十七章 马氏链模型 §1 随机过程的概念 一个随机试验的结果有多种可能性，在数学上用一个随机变量（或随机向量）来描 述。在许多情况下，人们不仅需要对随机现象进行一次观测，而且要进行多次，甚至接 连不断地观测它的变化过程。这就要研究无限多个，即一族随机变量。随机过程理论就 是研究随机现象变化过程的概率规律性的。 定义 1 设{ ,t T} ξ t ∈ 是一族随机变量，T 是一个实数集合，若对任意实数 T t t ∈ ,ξ 是一个随机变量，则称{ ,t T} ξ t ∈ 为随机过程。 T 称为参数集合，参数t 可以看作时间。ξ t 的每一个可能取值称为随机过程的一个 状态。其全体可能取值所构成的集合称为状态空间，记作 E 。当参数集合T 为非负整 数集时，随机过程又称随机序列。本章要介绍的马尔可夫链就是一类特殊的随机序列。 例 1 在一条自动生产线上检验产品质量，每次取一个，“废品”记为 1，“合格品” 记为 0。以ξ n 表示第n 次检验结果，则ξ n 是一个随机变量。不断检验，得到一列随机 变量ξ 1 ,ξ 2 ,L，记为{ ,n = 1,2,L} ξ n 。它是一个随机序列，其状态空间 E = {0,1}。 例 2 在m 个商店联营出租照相机的业务中（顾客从其中一个商店租出，可以到m 个商店中的任意一个归还），规定一天为一个时间单位，“ j ξ t = ”表示“第t 天开始营 业时照相机在第 j 个商店”， j = 1,2,L,m 。则{ ,n = 1,2,L} ξ n 是一个随机序列，其状 态空间 E = {1,2,L,m}。 例 3 统计某种商品在 t 时刻的库存量，对于不同的 t ，得到一族随机变量， { ,t ∈[0,+∞)} ξ t 是一个随机过程，状态空间 E = [0,R]，其中 R 为最大库存量。 我们用一族分布函数来描述随机过程的统计规律。一般地，一个随机过程 { ,t T} ξ t ∈ ，对于任意正整数n 及T 中任意n 个元素 n t , ,t 1 L 相应的随机变量 n ξ t ξ t , , 1 L 的联合分布函数记为 ( , , ) { , , } t1 t 1 n t1 1 t n F x x P x x n n L L = ξ ≤ L ξ ≤ （1） 由于n 及t (i 1, ,n) i = L 的任意性，（1）式给出了一族分布函数。记为 { ( , , ), , 1, , ; 1,2, } Ft1Ltn x1 L xn ti ∈T i = L n n = L 称它为随机过程{ ,t T} ξ t ∈ 的有穷维分布函数族。它完整地描述了这一随机过程的统计 规律性。 §2 马尔可夫链 2.1 马尔可夫链的定义 现实世界中有很多这样的现象：某一系统在已知现在情况的条件下，系统未来时刻 的情况只与现在有关，而与过去的历史无直接关系。比如，研究一个商店的累计销售额， 如果现在时刻的累计销售额已知，则未来某一时刻的累计销售额与现在时刻以前的任一 时刻累计销售额无关。上节中的几个例子也均属此类。描述这类随机现象的数学模型称 为马氏模型。 定义 2 设{ ,n = 1,2,L} ξ n 是一个随机序列，状态空间 E 为有限或可列集，对于任 意的正整数m,n ，若i, j,i ∈ E(k = 1, ,n −1) k L ，有
P{nm=儿n=1,n1=n1…,51=1}=P{5n+m=n=l} (2) 则称{n,n=1,2,…}为一个马尔可夫链(简称马氏链),(2)式称为马氏性。 事实上,可以证明若等式(2)对于m=1成立,则它对于任意的正整数m也成立。 因此,只要当m=1时(2)式成立,就可以称随机序列{n,n=1,2,…具有马氏性 即{n,n=1,2,…}是一个马尔可夫链 定义3设{n,n=12,…}是一个马氏链。如果等式(2)右边的条件概率与n无 关,即 PIMm=jIm=1)=Pu,(m) (3) 则称{5n,n=1,2…为时齐的马氏链。称P(m)为系统由状态i经过m个时间间隔(或 m步)转移到状态j的转移概率。(3)称为时齐性。它的含义是:系统由状态i到状态 j的转移概率只依赖于时间间隔的长短,与起始的时刻无关。本章介绍的马氏链假定都 是时齐的,因此省略“时齐”二字 22转移概率矩阵及柯尔莫哥洛夫定理 对于一个马尔可夫链{n,n=12…},称以m步转移概率P2(m)为元素的矩阵 P(m)=(P2(m))为马尔可夫链的m步转移矩阵。当m=1时,记P()=P称为马尔可 夫链的一步转移矩阵,或简称转移矩阵。它们具有下列三个基本性质: (i)对一切,j∈E,0≤P2(m)≤1 (i)对一切i∈E,∑P2(m)=1 (i)对一切1∈E,P=9={当=时 0,当≠j时 当实际问题可以用马尔可夫链来描述时,首先要确定它的状态空间及参数集合,然 后确定它的一步转移概率。关于这一概率的确定,可以由问题的内在规律得到,也可以 由过去经验给出,还可以根据观测数据来估计。 例4某计算机机房的一台计算机经常出故障,研究者每隔15分钟观察一次计算 机的运行状态,收集了24小时的数据(共作97次观察)。用1表示正常状态,用0表 示不正常状态,所得的数据序列如下 111001001111111001111011111100111111111000110110 111011011010111101110111101111110011011111100111 解设X(n=1…97)为第n个时段的计算机状态,可以认为它是一个时齐马氏 链,状态空间E={0,1},编写如下 Matlab程序 a1=11110010011111110011110111111001111111110001101101; a2=11110110110101111011101111011111100110111111001111 a=[a1a2] f00=length(findstr(00,a) f0l=length(findstr('01,a) f10=length(findstr(10,a) fll=length(findstr('ll,a) 或者把上述数据序列保存到纯文本文件data1.txt中,存放在Mat1ab下的work 子目录中,编写程序如下: C⊥c,c1e
-208- { | , , , } { | } 1 1 1 1 P j i i i P j i ξ n+m = ξ n = ξ n− = n− L ξ = = ξ n+m = ξ n = （2） 则称{ ,n = 1,2,L} ξ n 为一个马尔可夫链（简称马氏链），（2）式称为马氏性。 事实上，可以证明若等式（2）对于m = 1成立，则它对于任意的正整数m 也成立。 因此，只要当 m = 1时（2）式成立，就可以称随机序列{ ,n = 1,2,L} ξ n 具有马氏性， 即{ ,n = 1,2,L} ξ n 是一个马尔可夫链。 定义 3 设{ ,n = 1,2,L} ξ n 是一个马氏链。如果等式（2）右边的条件概率与n 无 关，即 P{ j | i} p (m) ξ n+m = ξ n = = ij （3） 则称{ ,n = 1,2,L} ξ n 为时齐的马氏链。称 p (m) ij 为系统由状态i 经过m 个时间间隔（或 m 步）转移到状态 j 的转移概率。（3）称为时齐性。它的含义是：系统由状态i 到状态 j 的转移概率只依赖于时间间隔的长短，与起始的时刻无关。本章介绍的马氏链假定都 是时齐的，因此省略“时齐”二字。 2.2 转移概率矩阵及柯尔莫哥洛夫定理 对于一个马尔可夫链{ ,n = 1,2,L} ξ n ，称以 m 步转移概率 p (m) ij 为元素的矩阵 P(m) ( p (m)) = ij 为马尔可夫链的m 步转移矩阵。当m = 1时，记 P(1) = P 称为马尔可 夫链的一步转移矩阵，或简称转移矩阵。它们具有下列三个基本性质： （i）对一切i, j ∈ E ，0 ≤ pij(m) ≤ 1； （ii）对一切i ∈ E ， ∑∈ = j E pij(m) 1； （iii）对一切i, j ∈ E ， ⎩ ⎨ ⎧ ≠ = = = 当 时 当 时 ， i j i j pij ij 0 1, (0) δ 。 当实际问题可以用马尔可夫链来描述时，首先要确定它的状态空间及参数集合，然 后确定它的一步转移概率。关于这一概率的确定，可以由问题的内在规律得到，也可以 由过去经验给出，还可以根据观测数据来估计。 例 4 某计算机机房的一台计算机经常出故障，研究者每隔 15 分钟观察一次计算 机的运行状态，收集了 24 小时的数据（共作 97 次观察）。用 1 表示正常状态，用 0 表 示不正常状态，所得的数据序列如下： 1110010011111110011110111111001111111110001101101 111011011010111101110111101111110011011111100111 解 设 X (n = 1,L,97) n 为第n 个时段的计算机状态，可以认为它是一个时齐马氏 链，状态空间 E = {0,1}，编写如下 Matlab 程序： a1='1110010011111110011110111111001111111110001101101'; a2='111011011010111101110111101111110011011111100111'; a=[a1 a2]; f00=length(findstr('00',a)) f01=length(findstr('01',a)) f10=length(findstr('10',a)) f11=length(findstr('11',a)) 或者把上述数据序列保存到纯文本文件data1.txt中，存放在Matlab下的work 子目录中，编写程序如下： clc,clear
format rat fid=fopen('datal. txt,'r)i while (feof(fid) a=[a fgetl(fid)] end for 1=0: 1 for j=0: 1 s=[int2str(),int2str (j)] f(i+1,j+1)=length(findstr(s, a))i end fs=sum(f)i f(i,:)=f(i,:)/fs(i) 求得96次状态转移的情况是 0→0,8次;0→1,18次 1→0,18次;1→1,52次 因此,一步转移概率可用频率近似地表示为 Poo= P(Xn=0 Xn=0 a 8+1813 o=P{Xn+1=1|Xn=0}≈ P1o=P(m=01,=l-18 73 8+ 9 18+5235 P1=PXn+=1|Xn=1≈-5226 18+5235 例5设一随机系统状态空间E={12,3,4},记录观测系统所处状态如下: 2232 14311 2143 若该系统可用马氏模型描述,估计转移概率P 解首先将不同类型的转移数nn统计出来分类记入下表 i→j转移数n 4121 行和n 2324 11 各类转移总和∑∑n等于观测数据中马氏链处于各种状态次数总和减1,而行和几是
-209- format rat fid=fopen('data1.txt','r'); a=[]; while (~feof(fid)) a=[a fgetl(fid)]; end for i=0:1 for j=0:1 s=[int2str(i),int2str(j)]; f(i+1,j+1)=length(findstr(s,a)); end end fs=sum(f'); for i=1:2 f(i,:)=f(i,:)/fs(i); end f 求得 96 次状态转移的情况是： 0 → 0，8 次； 0 →1，18 次； 1→ 0 ，18 次； 1→1，52 次， 因此，一步转移概率可用频率近似地表示为 13 4 8 18 8 { 0 | 0} 00 1 = + p = P Xn+ = Xn = ≈ 13 9 8 18 18 { 1| 0} 01 1 = + p = P Xn+ = Xn = ≈ 35 9 18 52 18 { 0 | 1} 10 1 = + p = P Xn+ = Xn = ≈ 35 26 18 52 52 { 1| 1} 11 1 = + p = P Xn+ = Xn = ≈ 例 5 设一随机系统状态空间 E = {1,2,3,4}，记录观测系统所处状态如下： 4 3 2 1 4 3 1 1 2 3 2 1 2 3 4 4 3 3 1 1 1 3 3 2 1 2 2 2 4 4 2 3 2 3 1 1 2 4 3 1 若该系统可用马氏模型描述，估计转移概率 pij 。 解 首先将不同类型的转移数nij 统计出来分类记入下表 i → j 转移数nij 1 2 3 4 行和ni 1 2 3 4 4 4 1 1 3 2 4 2 4 4 2 1 0 1 4 2 10 11 11 7 各类转移总和 ∑∑ i j nij 等于观测数据中马氏链处于各种状态次数总和减 1，而行和ni 是
系统从状态i转移到其它状态的次数,n是由状态i到状态j的转移次数,则pn的估 值P 算得 2/52/51/101/10 3/112/114/112/11 4/114/112/111/11 01/74/72/7 Matlab计算程序如下: format rat =[432143 2123443311 1332122244 for 1=1: 4 for j=l: 4 f(i,3)=length(findstr([i 3l, a)) ni=(sum(f)) for i=1: 4 p(i,:)=f(i,:)/ni(i) end 例6(带有反射壁的随机徘徊)如果在原点右边距离原点一个单位及距原点S(s>1) 个单位处各立一个弹性壁。一个质点在数轴右半部从距原点两个单位处开始随机徘徊。 每次分别以概率p(0<p<1)和q(q=1-p)向右和向左移动一个单位;若在+1处,则 以概率P反射到2,以概率q停在原处;在s处,则以概率q反射到s-1,以概率P停 在原处。设n表示徘徊n步后的质点位置。{n,n=1,2,…}是一个马尔可夫链,其状 态空间E={1,2,…,S},写出转移矩阵P。 当i=2时 0,当≠2时 q,当=1时 P={P,当=2时 0.其它 p,当=s时 Py={q,当 0,其它
-210- 系统从状态i 转移到其它状态的次数， nij 是由状态i 到状态 j 的转移次数，则 pij 的估 计值 i ij ij n n p = 。计算得 ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ = 0 1/ 7 4 / 7 2 / 7 4 /11 4 /11 2 /11 1/11 3/11 2 /11 4 /11 2 /11 2 / 5 2 / 5 1/10 1/10 Pˆ Matlab 计算程序如下： format rat clc a=[4 3 2 1 4 3 1 1 2 3 ... 2 1 2 3 4 4 3 3 1 1 ... 1 3 3 2 1 2 2 2 4 4 ... 2 3 2 3 1 1 2 4 3 1]; for i=1:4 for j=1:4 f(i,j)=length(findstr([i j],a)); end end f ni=(sum(f'))' for i=1:4 p(i,:)=f(i,:)/ni(i); end p 例 6（带有反射壁的随机徘徊）如果在原点右边距离原点一个单位及距原点 s(s > 1) 个单位处各立一个弹性壁。一个质点在数轴右半部从距原点两个单位处开始随机徘徊。 每次分别以概率 p(0 < p < 1) 和 q(q = 1− p) 向右和向左移动一个单位；若在+1 处，则 以概率 p 反射到 2，以概率q 停在原处；在 s 处，则以概率 q 反射到 s −1，以概率 p 停 在原处。设ξ n 表示徘徊 n 步后的质点位置。{ ,n = 1,2,L} ξ n 是一个马尔可夫链，其状 态空间 E = {1,2,L,s}，写出转移矩阵 P 。 解 ⎩ ⎨ ⎧ ≠ = = = 当 时 当 时 0， 2 1, 2 { } 0 i i P ξ i ⎪ ⎩ ⎪ ⎨ ⎧ = = = 其它 当 时 当 时 0, , 2 , 1 1 p j q j p j ⎪ ⎩ ⎪ ⎨ ⎧ = − = = 其它 当 时 当 时 0, , 1 , q j s p j s psj
p,当-1=1时 P={q,当-i=-1时(=23…,s-1) 0,其它 因此,P为一个s阶方阵,即 P0g 000 0 0 000 0 p 0000qp 定理1(柯尔莫哥洛夫一开普曼定理)设{n,n=1,2,…}是一个马尔可夫链,其 状态空间E={1,2,…},则对任意正整数m,n有 Pu, (n+m)=2P&(n)Pi, (m) 其中的i,∈E。 定理2设P是一个马氏链转移矩阵(P的行向量是概率向量),PO是初始分布 向量,则第n步的概率分布为 p(n)= p(o)pi 例7若顾客的购买是无记忆的,即已知现在顾客购买情况,未来顾客的购买情况 不受过去购买历史的影响,而只与现在购买情况有关。现在市场上供应A、B、C三个 不同厂家生产的50克袋状味精,用“5n=1”、“5n=2”、“5n=3”分别表示“顾客 第n次购买A、B、C厂的味精”。显然,{5n,n=12,…}是一个马氏链。若已知第 次顾客购买三个厂味精的概率依次为0.2,04,0.4。又知道一般顾客购买的倾向由表2 给出。求顾客第四次购买各家味精的概率。 表2 A B 上次 ABC 0.8 0 0.1 0.1 购买 0.5 0.3 0.2 解第一次购买的概率分布为 02040 0.80.10.1 转移矩阵P=0.50.104 0.50.30.2 则顾客第四次购买各家味精的概率为 P+=PP=[070040.3601636] 23转移概率的渐近性质一极限概率分布
-211- ⎪ ⎩ ⎪ ⎨ ⎧ − = − = − − = = 其它 当 时 当 时 0, , 1 ( 2,3, , 1) , 1 q j i i s p j i pij L 因此， P 为一个 s 阶方阵，即 ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ = q p q p q q p q p P 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 L L L L L L L L L 。 定理 1 （柯尔莫哥洛夫—开普曼定理）设{ ,n = 1,2,L} ξ n 是一个马尔可夫链，其 状态空间 E = {1,2,L}，则对任意正整数 m,n 有 ∑∈ + = k E pij(n m) pik (n) pkj(m) 其中的i, j ∈ E 。 定理 2 设 P 是一个马氏链转移矩阵（ P 的行向量是概率向量）， (0) P 是初始分布 行向量，则第n 步的概率分布为 n n P P P ( ) (0) = 。 例 7 若顾客的购买是无记忆的，即已知现在顾客购买情况，未来顾客的购买情况 不受过去购买历史的影响，而只与现在购买情况有关。现在市场上供应 A、B、C 三个 不同厂家生产的 50 克袋状味精，用“ξ n = 1”、“ξ n = 2 ”、“ξ n = 3”分别表示“顾客 第n 次购买 A、B、C 厂的味精”。显然，{ ,n = 1,2,L} ξ n 是一个马氏链。若已知第一 次顾客购买三个厂味精的概率依次为 0.2，0.4，0.4。又知道一般顾客购买的倾向由表 2 给出。求顾客第四次购买各家味精的概率。 表 2 下 次 购 买 A B C 上次 购买 A B C 0.8 0.5 0.5 0.1 0.1 0.3 0.1 0.4 0.2 解 第一次购买的概率分布为 [ ] 0.2 0.4 0.4 (1) P = 转移矩阵 ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎣ ⎡ = 0.5 0.3 0.2 0.5 0.1 0.4 0.8 0.1 0.1 P 则顾客第四次购买各家味精的概率为 [0.7004 0.136 0.1636] (4) (1) 3 P = P P = 。 2.3 转移概率的渐近性质—极限概率分布