重庆工商大学数学建模算法讲义第17章 马氏链模型

-207-第十七章马氏链模型§1随机过程的概念一个随机试验的结果有多种可能性，在数学上用一个随机变量（或随机向量）来描述。在许多情况下，人们不仅需要对随机现象进行一次观测，而且要进行多次，甚至接连不断地观测它的变化过程。这就要研究无限多个，即一族随机变量。随机过程理论就是研究随机现象变化过程的概率规律性的。定义1设},{Ttt∈ξ是一族随机变量，T是一个实数集合，若对任意实数tTtξ,∈是一个随机变量，则称},{Ttt∈ξ为随机过程。T称为参数集合，参数t可以看作时间。tξ的每一个可能取值称为随机过程的一个状态。其全体可能取值所构成的集合称为状态空间，记作E。当参数集合T为非负整数集时，随机过程又称随机序列。本章要介绍的马尔可夫链就是一类特殊的随机序列。例1在一条自动生产线上检验产品质量，每次取一个，“废品”记为1，“合格品”记为0。以nξ表示第n次检验结果，则nξ是一个随机变量。不断检验，得到一列随机变量L,,21ξξ，记为},2,1,{L=nnξ。它是一个随机序列，其状态空间}1,0{=E。例2在m个商店联营出租照相机的业务中（顾客从其中一个商店租出，可以到m个商店中的任意一个归还），规定一天为一个时间单位，“jt=ξ”表示“第t天开始营业时照相机在第j个商店”，mj,,2,1L=。则},2,1,{L=nnξ是一个随机序列，其状态空间},,2,1{mEL=。例3统计某种商品在t时刻的库存量，对于不同的t，得到一族随机变量，)},0[,{+∞∈ttξ是一个随机过程，状态空间],0[RE=，其中R为最大库存量。我们用一族分布函数来描述随机过程的统计规律。一般地，一个随机过程},{Ttt∈ξ，对于任意正整数n及T中任意n个元素ntt,,1L相应的随机变量nttξξ,,1L的联合分布函数记为},,{),,(1111nttnttxxPxxFnn≤≤=ξξLLL（1）由于n及),,1(nitiL=的任意性，（1）式给出了一族分布函数。记为},2,1;,,1,),,,({11LLLL==∈nniTtxxFinttn称它为随机过程},{Ttt∈ξ的有穷维分布函数族。它完整地描述了这一随机过程的统计规律性。§2马尔可夫链2.1马尔可夫链的定义现实世界中有很多这样的现象：某一系统在已知现在情况的条件下，系统未来时刻的情况只与现在有关，而与过去的历史无直接关系。比如，研究一个商店的累计销售额，如果现在时刻的累计销售额已知，则未来某一时刻的累计销售额与现在时刻以前的任一时刻累计销售额无关。上节中的几个例子也均属此类。描述这类随机现象的数学模型称为马氏模型。定义2设},2,1,{L=nnξ是一个随机序列，状态空间E为有限或可列集，对于任意的正整数nm,，若)1,,1(,,−=∈nkEijikL，有-208-}|{},,,|{1111ijPiiijPnmnnnnmn=======+−−+ξξξξξξL（2）则称},2,1,{L=nnξ为一个马尔可夫链（简称马氏链），（2）式称为马氏性。事实上，可以证明若等式（2）对于1=m成立，则它对于任意的正整数m也成立。因此，只要当1=m时（2）式成立，就可以称随机序列},2,1,{L=nnξ具有马氏性，即},2,1,{L=nnξ是一个马尔可夫链。定义3设},2,1,{L=nnξ是一个马氏链。如果等式（2）右边的条件概率与n无关，即)(}|{mpijPijnmn===+ξξ（3）则称},2,1,{L=nnξ为时齐的马氏链。称)(mpij为系统由状态i经过m个时间间隔（或m步）转移到状态j的转移概率。（3）称为时齐性。它的含义是：系统由状态i到状态j的转移概率只依赖于时间间隔的长短，与起始的时刻无关。本章介绍的马氏链假定都是时齐的，因此省略“时齐”二字。2.2转移概率矩阵及柯尔莫哥洛夫定理对于一个马尔可夫链},2,1,{L=nnξ，称以m步转移概率)(mpij为元素的矩阵))(()(mpmPij=为马尔可夫链的m步转移矩阵。当1=m时，记PP=)1(称为马尔可夫链的一步转移矩阵，或简称转移矩阵。它们具有下列三个基本性质：（i）对一切Eji∈,，1)(0≤≤mpij；（ii）对一切Ei∈，∑∈=Ejijmp1)(；（iii）对一切Eji∈,，⎩⎨⎧≠===时当时当，jijipijij0,1)0(δ。当实际问题可以用马尔可夫链来描述时，首先要确定它的状态空间及参数集合，然后确定它的一步转移概率。关于这一概率的确定，可以由问题的内在规律得到，也可以由过去经验给出，还可以根据观测数据来估计。例4某计算机机房的一台计算机经常出故障，研究者每隔15分钟观察一次计算机的运行状态，收集了24小时的数据（共作97次观察）。用1表示正常状态，用0表示不正常状态，所得的数据序列如下：1110010011111110011110111111001111111110001101101111011011010111101110111101111110011011111100111解设)97,,1(L=nXn为第n个时段的计算机状态，可以认为它是一个时齐马氏链，状态空间}1,0{=E，编写如下Matlab程序：a1='1110010011111110011110111111001111111110001101101';a2='111011011010111101110111101111110011011111100111';a=[a1a2];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-209-formatratfid=fopen('data1.txt','r');a=[];while(~feof(fid))a=[afgetl(fid)];endfori=0:1forj=0:1s=[int2str(i),int2str(j)];f(i+1,j+1)=length(findstr(s,a));endendfs=sum(f');fori=1:2f(i,:)=f(i,:)/fs(i);endf求得96次状态转移的情况是：00→，8次；10→，18次；01→，18次；11→，52次，因此，一步转移概率可用频率近似地表示为1341888}0|0{100=+≈===+nnXXPp13918818}0|1{101=+≈===+nnXXPp359521818}1|0{110=+≈===+nnXXPp3526521852}1|1{111=+≈===+nnXXPp例5设一随机系统状态空间}4,3,2,1{=E，记录观测系统所处状态如下：4321431123212344331113321222442323112431若该系统可用马氏模型描述，估计转移概率ijp。解首先将不同类型的转移数ijn统计出来分类记入下表ji→转移数ijn1234行和in123444113242442101421011117各类转移总和∑∑ijijn等于观测数据中马氏链处于各种状态次数总和减1，而行和in是-210-系统从状态i转移到其它状态的次数，ijn是由状态i到状态j的转移次数，则ijp的估计值iijijnnp=。计算得⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=7/27/47/1011/111/211/411/411/211/411/211/310/110/15/25/2ˆPMatlab计算程序如下：formatratclca=[4321431123...2123443311...1332122244...2323112431];fori=1:4forj=1:4f(i,j)=length(findstr([ij],a));endendfni=(sum(f'))'fori=1:4p(i,:)=f(i,:)/ni(i);endp例6（带有反射壁的随机徘徊）如果在原点右边距离原点一个单位及距原点)1(ss个单位处各立一个弹性壁。一个质点在数轴右半部从距原点两个单位处开始随机徘徊。每次分别以概率)10(pp和)1(pqq−=向右和向左移动一个单位；若在+1处，则以概率p反射到2，以概率q停在原处；在s处，则以概率q反射到1−s，以概率p停在原处。设nξ表示徘徊n步后的质点位置。},2,1,{L=nnξ是一个马尔可夫链，其状态空间},,2,1{sEL=，写出转移矩阵P。解⎩⎨⎧≠===时当时当，202,1}{0iiiPξ⎪⎩⎪⎨⎧===其它时当时当,02,1,1jpjqpj⎪⎩⎪⎨⎧−===其它时当时当,01,,sjqsjppsj-211-⎪⎩⎪⎨⎧−=−=−=−=其它时当时当,0)1,,3,2(1,1,siijqijppijL因此，P为一个s阶方阵，即⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎣⎡=pqpqqpqpqP000000000000000000LLLLLLLLL。定理1（柯尔莫哥洛夫—开普曼定理）设},2,1,{L=nnξ是一个马尔可夫链，其状态空间},2,1{L=E，则对任意正整数nm,有∑∈=+Ekkjikijmpnpmnp)()()(其中的Eji∈,。定理2设P是一个马氏链转移矩阵（P的行向量是概率向量），)0(P是初始分布行向量，则第n步的概率分布为nnPPP)0()(=。例7若顾客的购买是无记忆的，即已知现在顾客购买情况，未来顾客的购买情况不受过去购买历史的影响，而只与现在购买情况有关。现在市场上供应CBA、、三个不同厂家生产的50克袋状味精，用“1=nξ”、“2=nξ”、“3=nξ”分别表示“顾客第n次购买CBA、、厂的味精”。显然，},2,1,{L=nnξ是一个马氏链。若已知第一次顾客购买三个厂味精的概率依次为0.2，0.4，0.4。又知道一般顾客购买的倾向由表2给出。求顾客第四次购买各家味精的概率。表2下次购买ABC上次购买ABC0.80.50.50.10.10.30.10.40.2解第一次购买的概率分布为[]4.04.02.0)1(=P转移矩阵⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=2.03.05.04.01.05.01.01.08.0P则顾客第四次购买各家味精的概率为[]1636.0136.07004.03)1()4(==PPP。2.3转移概率的渐近性质—极限概率分布-212-现在我们考虑，随n的增大，nP是否会趋于某一固定向量？先考虑一个简单例子：转移矩阵⎥⎦⎤⎢⎣⎡=3.07.05.05.0P，当+∞→n时，⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡→125127125127nP又若取⎥⎦⎤⎢⎣⎡=125127u，则uuP=，Tu为矩阵TP的对应于特征值1=λ的特征（概率）向量，u也称为P的不动点向量。哪些转移矩阵具有不动点向量？为此我们给出正则矩阵的概念。定义４一个马氏链的转移矩阵P是正则的，当且仅当存在正整数k，使kP的每一元素都是正数。定理３若P是一个马氏链的正则阵，那么：（i）P有唯一的不动点向量W，W的每个分量为正。（ii）P的n次幂nP（n为正整数）随n的增加趋于矩阵W，W的每一行向量均等于不动点向量W。例8信息的传播一条新闻在LL,,,,21naaa等人中间传播，传播的方式是1a传给2a，2a传给3a，…如此继续下去，每次传播都是由ia传给1+ia。每次传播消息的失真概率是p，10p，即ia将消息传给1+ia时，传错的概率是p，这样经过长时间传播，第n个人得知消息时，消息的真实程度如何？设整个传播过程为随机转移过程，消息经过一次传播失真的概率为p，转移矩阵⎥⎥⎦⎤⎢⎢⎣⎡=−−ppppP11真假真假P是正则矩阵。又设V是初始分布，则消息经过n次传播后，其可靠程度的概率分布为nPV⋅。一般地，设时齐马氏链的状态空间为E，如果对于所有Eji∈,，转移概率)(npij存在极限jijnnpπ=∞→)(lim，（不依赖于i）或⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⎯⎯→⎯=∞→LLLLLLLLLLLLLLL