灰色马尔可夫(精选8篇)
灰色马尔可夫 篇1
摘要:针对Mean Shift算法中存在的不足,提出了结合均值偏移和灰色马尔可夫预测模型的目标跟踪算法。该方法利用灰色马尔可夫模型预测目标在当前时刻的中心位置,以此点作为均值偏移算法进行目标搜索的起始位置。同时,提取目标的几何特征,根据目标的面积来改善跟踪窗口的大小,利用“目标模型”和“候选模型”之间目标特征的变化产生模型更新策略。经实验得,该方法能实时稳健地进行跟踪。
关键词:Mean Shift,灰色马尔可夫预测模型,几何特征,目标跟踪
Comaniciu等[1]提出的基于Mean Shift的跟踪算法取得了很大成功,并且吸引了越来越多学者的研究兴趣。Mean Shift算法是一种非参数密度估算法,它具有良好的实时性,且对形变、目标遮掩的稳健性良好,易于与其他算法集成,广泛地应用于目标跟踪中。但该算法也存在一定的缺陷,为了解决传统均值偏移算法不能自适应改变核函数带宽的缺陷,Collins等将Lindeberg尺度的空间理论和Mean Shift算法相结合[2],并在某种程度上改善了此缺陷;Li Jinping等提出了采用Level Set描述目标轮廓方法与Mean Shift算法组合进行跟踪,使得在光线变化、目标颜色发生改变时取得好的效果[3],改善了相似颜色干扰问题;为了解决遮挡情况下均值偏移算法的跟踪缺陷,许多学者采用了目标预测结合Mean Shift的跟踪方法,Maggio E等提出了采用粒子滤波结合Mean Shift算法的跟踪[4],取得了好的跟踪效果,但计算量大。
本文在灰色GM(1,1)模型的基础上引入马尔可夫链预测理论,建立运动目标的灰色马尔可夫GM(1,1)预测模型,利用少量的数据来预测目标的运动轨迹,并以当前时刻的目标预测位置作为Mean Shift算法进行迭代搜索的起始位置;利用提取的几何特征表示相似度函数来自适应性地更新搜索窗口,以此减少迭代次数,最后根据模型更新策略来更新目标模拟。
1 Mean Shift算法
Mean Shift算法是一个自适应地寻找概率密度局部最大值的迭代方法。在目标跟踪中,采用归一化的加权颜色直方图来描述目标模型。假定跟踪目标的中心位于x0,xi是d维Euclidean空间Rd中的一组点,用向量表示xi(i=1,2,…,n),使用带宽为h的核函数K(x)作为多变量核密度估计,则目标模型可以表示为
式中:u=1,2,…,m,为目标区域的所有特征值;Ch为归一化常数;
令候选模型的中心坐标为y,则可以描述为
pu与qu(y)的相似性用Bhattacharyya系数
从式(3)可以看出,两个模型越相似,则
式中:
2 灰色马尔可夫预测模型
灰色系统预测的概念是由邓聚龙教授首先在国内提出,之后灰色系统理论的研究得到了迅速发展。灰色预测方法的优点在于对缺少基础资料的预测能够得到较好的预测效果,与别的预测方法相比较而言它用到的样本数据较小,而预测的精度相对较高。马尔可夫概率矩阵是对随机过程每个时刻状态的描述,它是根据状态之间的转移概率来预测系统的发展,把在不相同状态范围的内在波动规律展现出来,使随机作用造成的波动得以修正。灰色马尔可夫预测模型就是将两者的优势相结合进行预测。
2.1 建立GM(1,1)模型
假设有原始序列为
x(0)={x(0)(1),x(0)(2),…,x(0)(n)} (5)
将该序列x(0)进行一次累加生成处理,得到x(1)(记作1-AGO)
x(1)={x(1)(1),x(1)(2),…,x(1)(n)} (6)
式中:
建立GM(1,1)灰微分方程
式中:a为发展系数;b为灰作用量,是待估参数。
把式(7)转化为矩阵方程
z(1)(k)=[x(1)(k)+x(1)(k-1)]/2,k=2,3,…,n (11)
由式(8)可得:
GM(1,1)模型的离散响应方程为
将
以上就是建立运动目标轨迹的灰色GM(1,1)预测模型的基础[5]。
2.2 建立转移概率矩阵
假设将目标的轨迹状态分为m个,即状态空间I={1,2,…,m},fij作为i到j状态经一步转移的数据样本数,与
由状态i经k步转到状态j的频率也可以用式(14)来表示,将状态转移概率依次排序,可得到矩阵
式中:
从以上分析可知,灰色GM(1,1)模型用来预测目标轨迹变化的总趋势,对其轨迹状态的修正则是通过马尔可夫转移概率矩阵来完成的,最终精确地预测出目标的位置,为均值偏移算法提供搜索初始位置,避免了由于遮挡造成的目标丢失。
3 Mean Shift算法和灰色马尔可夫模型的结合
跟踪的过程中,会产生一系列的时间序列{yi},i=1,2,…,表示迭代的次数,这种时间序列有平稳收敛(呈递增和递减两种形式分布)和波动收敛(呈交替状分布)两种形式。针对这两种形式,可以采用灰色马尔可夫模型进行目标中心位置的预测。在Mean Shift搜索之前先进行目标在下一帧中位置的初步预测,并用预测中心位置作为目标搜索的起始位置,之后再进行迭代搜索,不仅可以提高搜索准确性,还能够减少迭代搜索次数。灰色马尔可夫模型采用“新陈代谢”更新策略,保证用于预测的数据是最新的,一般用当前帧的前4帧数据作为基本数据,既保证了预测的准确性,也减少计算量。
如图1所示,提取的颜色特征用于传统的Mean Shift算法进行迭代计算,找出最匹配的目标位置,具体过程见本文第2节;利用文献[6]中的迭代投影算法来近视计算目标的面积,通过前后两帧面积的差值来判断是否更新跟踪窗口的大小;再设置一个阀值,根据相邻两帧的目标面积的比值是否小于阀值来判断是否更新目标模型,避免了目标模型过更新。
4 实验结果与分析
本文实验使用的硬件平台为2.53 GHz的CPU、2 Gbyte内存,在Windows7操作系统下的计算机,软件平台为MATLAB 2010a编程环境。测试视频选取一段橄榄球的运动视频序列,视频图像大小为352×288,共88帧。根据反复的实验,选取模板更新阈值为0.2。在整个视频序列中,橄榄球的运动都是随机波动的,前42帧被跟踪的橄榄球从大到小,再从小到大,其大小一直在不断变化;后40多帧发生遮挡情况。
图2为传统的均值偏移算法跟踪的仿真图,当球的大小不断变化时,跟踪窗口不能自适应地调整大小,当球快速运动时,目标窗口的中心位置严重偏移,在42帧发生遮挡时,跟踪失败,跟踪窗口停留在第42帧时的位置。本文算法的仿真图如图3所示,根据灰色马尔可夫模型和均值偏移算法相结合,在球快速运动时能准确地确定目标的中心位置,再根据提取目标的几何特征的变化来改善跟踪窗口的大小,能很好地进行实施跟踪,效果也比传统的均值偏移算法要好;在发生遮挡时,以灰色马尔可夫GM(1,1)的预测值为中心来开窗,避免丢失跟踪目标。从两个仿真图来看,本文提出的算法要优于传统的均值偏移算法,可以取得良好的实时性和稳健性。
5 结论
当在运动目标跟踪的过程中,发生遮挡、多种运动方式存在的情况下,本文提出通过灰色GM(1,1)模型对目标的运动轨迹进行总趋势预测,再用马尔可夫概率转移矩阵来对轨迹进行修正,这样能准确预测出目标在每一帧中的中心位置,并比一般的预测方法计算量小、精度高;再结合Mean Shift算法进行匹配跟踪,最终确定目标。其中通过提取的几何特征的变化来更新跟踪窗口和目标模型,这样可以适应目标的尺度变化和避免了目标模型的过更新。实验表明,该方法能够对目标进行实时有效的跟踪。
参考文献
[1]COMANICIU D,RAMESH V,MEER P.Real-time tracking of non-rigidobjects using mean Shift[C]//Proc.IEEE Conference on Computer Vi-sion and Pattern Recognition 2000.[S.l.]:IEEE Press,2000:142-149.
[2]COLLINS R T.Mean Shift blob tracking through scale space[C]//Proc.IEEE Conference on Computer Vision and Pattern Recognition 2003.[S.l.]:IEEE Press,2003:234-240.
[3]LI J P,LI Q J.Real-time tracking by combining level set and Mean Shift[J].Journal of Information&Computational Science,2008,5(2):829-836.
[4]MAGGIO E,CAVALLARO A.Hybrid particle filter and mean shifttracker with adaptive transition model[C]//Proc.ICASSP 2005.[S.l.]:IEEE Press,2005:221-224.
[5]邓聚龙.灰色系统理论教程[M].武汉:华中理工大学出版社,1992.
[6]叶有时,唐林波,赵保军,等.基于SOPC的深空目标实时跟踪系统[J].系统工程与电子技术,2009,31(12):3002-3006.
灰色马尔可夫 篇2
将Mel频率倒谱作为特征参数,实现了基于连续型隐马尔可夫模型的舰船辐射噪声目标识别.对五艘不同类别舰船的辐射噪声建立模型,并进行了识别实验.模型的状态数和高斯混合数分别确定为1和3,五艘舰船的.平均识别率达到91.6%.利用两艘不同类型舰船的噪声分析了舰船的工况对识别结果的影响.
作 者:戴卫国 张宝华 程玉胜 DAI Wei-guo ZHANG Bao-hua CHENG Yu-sheng 作者单位:戴卫国,程玉胜,DAI Wei-guo,CHENG Yu-sheng(海军潜艇学院,青岛,266071)
张宝华,ZHANG Bao-hua(解放军某部,北京,100841)
灰色马尔可夫 篇3
2009年4月以来,上海房地产市场形成供销两旺的强劲势头,商品房的开发和销售业绩屡创新高,上海房市迎来了“小阳春”。相比之下,市区居民的可支配收入和购买力水平并未随之“水涨船高”。因此,房价就成为供需双方共同关注的焦点。一方面:房地产开发的高投入、高风险,迫使开发商必须及时、准确把握市场信息、掌握房价动向,才能做出科学的决策,趋利避险;另一方面:广大购房者也只有了解未来房价的变化趋势,才会做出理性选择,避免盲目抢购和开发商恶意炒作造成无谓的经济损失。同时,稳定房地产市场更是政府义不容辞的责任。可见,进行科学的房价预测对供需双方及政府责任都具有重要理论和现实意义。
1 研究内容与方法
1.1 主要研究内容
房地产价格构成复杂,形成过程中影响因素众多。具有较大的随机波动性[1]。因此,本文以上海房地产市场2008年至2009年间二手房均价为原数据,探讨了利用灰色-马尔可夫预测模型[2]对商品房价格进行预测的可能,进而提出相关政策建议。
1.2 主要研究方法——灰色马尔可夫模型
灰色理论基础[3],其特点是“少数据建模”,灰色系统理论提供了在贫信息情况下解决系统问题的途径。灰色系统的基本模型GM(1,1)就是通过时序数据的累加,从而滤去原始序列中可能混入的随机量,从上下波动的时间序列中寻找某种隐含的规律性,得到随机性弱化而规律性强化了的新数列,然后挖掘出原始序列的内在特征。
马尔可夫链是时间离散、状态离散的随机过程。它的特点是无后效性,即在时刻t的状态只与前一时刻的状态有关,而与前面其它各时间状态无关。一个n阶马尔可夫链由n个状态集合和一组转移概率所确定。该过程的任一时刻只能处于一个状态。如果在时刻t,过程处于状态Si,在t+1时刻它将以概率P处于状态Sj。其预测是根据状态之间的转移概率来推测系统未来的发展变化[4]。
灰色预测和马尔可夫预测都可用于时间序列预测问题,但灰色预测几何曲线呈单调递增或递减趋势,能反映总体趋势走向。而马尔可夫预测则是根据状态转移概率大小来推测系统未来发展方向,转移概率反映了各种随机因素的影响,因而适合于随机波动较大的序列预测问题。因此,这两种方法具有很大的互补性,若将它们结合使用,则有较高的预测准确率。
具体模型构建方法:首先根据前n个月的数据建立GM(1,1)预测模型,作第一步预测,然后利用GM(1,1)预测值与实际能源消耗值的差值建立马尔可夫多步转移概率矩阵模型,利用前n个月的初始状态信息及相应的多步转移概率和,推测下一月所处的状态,最后根据马尔可夫模型的预测状态,对GM(1,1)模型预测值进行修正。
2 灰色-马尔可夫模型房价预测研究
已知上海市2008年1月至2009年4月的每月二手房均价,试用灰色-马尔可夫模型预测2009年5月份二手房均价。将预测值与实际房价比较,验证灰色-马尔可夫模型预测房价的可行性与准确性。
2.1 二手房价格灰色GM(1,1)预测
根据灰色GM(1,1)模型的原理,将以上二手数据代入,可以得到上海房市二手房价格变化趋势的灰色预测方程如式(1)。
式(1)中,
根据预测方程(1),计算出2008年2月至2009年4月的预测值,拟合结果验证如表1所示。
根据灰色GM(1,1)模型进行预测,得到2009年5月的二手商品住宅售价的预测值为
其实际价格为14 600 元/平米,两者偏差为预测值偏低14.626 5,相对误差0.100 2%。可见,灰色GM(1,1)模型的拟合表现已非常的好。
2.2 一次GM(1,1)预测结果的马尔可夫评价
根据马尔可夫链分析方法的应用经验和实际情况,按照房价变化趋势的增幅与灰色预测结论的比较,可以划分为5种状态。见表2。
注:误差比例为预测值的残差(预测值-实际值)占实际价格的比例。
从以上分类中可以获得2008-01~2009-04区间内的价格状态转移情况,状态1、状态5未出现,剔除,得到表3。
进而得到一步转移概率矩阵
2009年4月的灰色预测结果的相对误差为0.901 6%,处于状态3,预测准确。根据以上转移概率矩阵,2009年5月的灰色预测结果很有可能(66.7%的概率)也处于状态3,即预测准确。
实际情况,获得2009年5月价格数据后发现预测结果的相对误差为0.100 2%,确实属于状态3。可见,对灰色GM(1,1)模型预测结果的马尔可夫评价是可信的。
2.3 多次GM(1,1)预测结果的有效性评价
对GM(1,1)模型预测房价的有效性进一步进行分析。
根据2008年1月至2009年4月份二手房价准确数据,利用GM(1,1)模型进行多次预测,得到以下结果,如表4。
通过比较,我们可以发现,自第17期,即2009年5月,开始预测以来,GM(1,1) 模型模拟拟合值的准确度正急剧下降。预测一个季度后,即2009年8月,拟合值与实际值的相对误差已达到10%,拟合值已不具备参考价值。
因此,GM(1,1)模型的准确性是有严格的限定条件的:基于准确的原始数据下,其预测时效有限,仅能准确预测未来的一至两期。
2.4多次GM(1,1)预测结果的马尔可夫状态转移
根据以上分析,我们需要对GM(1,1)模型预测结果进行马尔可夫链改进。
根据马尔可夫链预测原理,可以得到原始数据之后一个季度内(2009-05~2009-07)的预测状态向量如表5所示。
由表5发现,有效预测的三个月内,预测值处于状态1(高估)和状态3(准确)的概率正在下降,而状态2(低估)的概率正在急剧上升,上升了近30个百分点。因此,有效预测的三个月内,预测值低估市场的风险越来越大。这与表3(2009-05~2009-07间数值)的事实表现是相符的。
事实上,回观2009年的上海房地产市场,自4月开始,上海房市开始进入爆发期。4月份,成交量为255.78万平方米,创下自2005年来的月度新高。之后数月,全市平均房价突破万元大关。这些表现是在政策鼓励前提下,需求和投资的全面释放导致的。
4月份是房市由复苏向爆发转变的转折点,因此,基于2009年4月之前的原始数据预测的5月及以后月份房价的预测出现严重低估的现象也就可以理解了。
3 结论
经过以上分析,可得以下小结:
(1) 将灰色-马尔可夫模型应用于商品住宅销售价格预测是切实可行和行之有效的。灰色-马尔可夫模型既考虑了从时间序列中挖掘数据的演变规律,又通过状态转移概率矩阵的变换提取数据的随机响应。因而,具有严密的科学性和适用性。
(2) 灰色-马尔可夫模型是基于对历史数据的统计分析之上的,因此历史数据越详尽,预测精度越高,预测结果越可靠。
(3) 灰色-马尔可夫模型的预测亦具有一定的局限性。该模型只是基于历史数据预测未来,未考虑未来政府政策变化、国际局势变化等市场突发现象的影响。因此,该模型只能预测近期内房价的走势。
(4) 房地产市场受各方面因素的影响很大,很容易出现突发性波动。为了提高灰色-马尔可夫模型预测的准确性,需要观察员及时更新原始数据,并关注马尔可夫状态转移概率的变化趋势。
摘要:针对房地产营销过程中,价格变化随机性的特点,将灰色GM(1,1)预测与马尔可夫预测有机结合,构成灰色—马尔可夫预测模型,对未来某一时段的房价的灰色预测结果做马尔可夫评价。结果表明,该模型是切实有效的。因此,此模型可用于指导房地产政策的制定。
关键词:灰色-马尔可夫模型,房价,预测
参考文献
[1]武永祥,王学涵.房地产开发.北京:中国建筑工业出版社,2003
[2]刘志彬,施斌.灰色马尔可夫链在深基坑沉降预测中的应用.煤田地质与勘探,2002;(6):35—37
[3]邓聚龙.灰色理论基础.武汉:华中科技大学出版社,2002
灰色马尔可夫 篇4
在股票市场中,股票价格是一个基本特征量,但是它总受政治、经济等各方面的影响,具体的影响因素的程度和信息是不完全的,所以我们可以把股市当成一个灰色系统来处理。灰色GM(1,1)预测模型是灰色系统理论的重要组成部分,主要适用于时间短、数据资料少、波动性不大的预测问题,且只需很少的几个数据即可建立模型进行预测,可以很好地解决由于数据少而导致的精确度不高的问题,但由于灰色GM(1,1)预测模型的预测曲线是一条较平滑的单调曲线,对波动性较大的股票市场中的数据列拟合较差,预测度较低。同时马尔可夫链比较适合随机波动性较大的预测问题,但是马尔可夫链要求状态无后效性,且要具有平稳过程等特点。如果灰色GM(1,1)模型对数据进行拟合,找出其变化趋势,则可以弥补马尔可夫预测的局限性,而在灰色预测基础上进行马尔可夫预测,又可弥补灰色预测对随机波动性较大的数据序列准确度低的不足,所以可将二者结合起来。形成灰色马尔可夫预测模型,对随机波动性较大的数据列进行预测,大大提高随机波动性较大数据列的预测精度,为随机波动性较大的对象的预测提供一种新的方法。
1 灰色-马尔可夫链预测模型
1.1 灰色模型GM(1,1)
1.1.1 给定原属数据列,记作x(0)=[x(0)(k)|k=1,2,…,n]
1.1.2 对原始数据进行一次累加,生成新的数据序列记x(1)=[x(1)(k)|k=1,2,…,n]=[x(1)(1),x(1)(2),…,x(1)(n)]
其中x(1)与x(0)之间满足下述关系
1.1.3 构造矩阵B和常数矩阵Yn。
1.1.4 用最小二乘法求解a,u。
1.1.5 该模型的时间响应方程为:
1.2 马尔可夫链预测模型
1.2.1 状态划分。
(1)根据GM(1,1)模型求出原始数据序列的拟合值x(t);
(2)求出残差;
(3)残差的相对值为;
(4)为了使每一状态的数据相差不多,将ε(t)的值从小到大排列,根据用户的需要和数据的多少,将状态分为自己想要的数。
1.2.2 马尔可夫链。
(1)马尔可夫链的概念。
定义1设随机过程{X(t),t∈T}的状态空间S为R中的可列集。如果对T中任意n个参数t1<t2<…<tn,以及使P{X(tn)=in|X(tn-1)=in-1,
则称{X(t),t∈T}为马尔可夫链。
定义2 若马尔可夫过程的状态空间S为R中的可列集,时间参数集T为可列离散集,则称P={Xm+1=j|Xm=i,Xm-1=im-1,…,X0=i0}=P
定义3若P{Xm+1=j|Xm=i}≡P{X1=j|X0=i}≡P(m,1)≡Pij即从i状态转移到j状态的概率与m无关,称这类马尔可夫链为齐次马尔可夫链,Pij表示由状态i经过一步转移到状态j的概率,称为一步概率转移,以Pij为元素的矩阵P=(Pij)称为状态一步转移概率矩阵,其形式为:
其中Pij叟0,(i,j∈s);∑j=sPij=1,(i∈s)。
(2)转移概率的计算。
a.前面的状态划分,将残差的相对值分为若干状态,记为E1,E2,E3,E4,残差相对值序列由状态Ei经k步转移到状态Ej的概率成为n步转移概率,记为,式中mij(k)为状态Ei经k步转移到状态Ej的次数,Mi为状态Ei出现的次数,由于数据序列的最后状态的转向不明确,故计算Mi为时要去掉数据序列中最末的k个数据;
b.当k=1时,即为一步转移概率Pij,其矩阵形式可记为:
c.P(1)=PP(0),考察P(1)中的n个值,若max Pkj=Pk1,则可以认为下一时刻系统最有可能由状态Ek转向状态E1,即下一时刻最有可能处于状态E1;
d.P(n)=PP(n-1),同理上一步可知n时刻后系统最有可能所处状态。
e.计算预测区间及预测值。确定预测对象未来的转移状态转移以后,即确定了残差预测值的变动区间Ej,我们用区间的中位数作为残差的预测值,即残差的预测值=(预测区间最大值+预测区间最小值)/2。则最终预测对象的预测值=(1+残差预测值)灰色预测值。
2 实例分析
由于收盘价是影响股票价格的最主要的因素,本文在上交所交易指数数据中选取从2008年10月到2010年5月20个月末收盘指数,其中时间序列的单位以月计。
步骤一:根据表1中已有的实际数据,用Matlab求出模型GM(1,1)中的a,u。a=0.0215,u=2160.4,则
步骤二:由残差相对值公式,求得残差相对值序列ε(t)的范围为(-25.317%,23.1901%);
步骤三:根据实际情况将残差相对值ε(t)平均分为4个状态a,b,c,d;其中a=[13%,26%);b=[0,13%);c=[-13%,0);d=[-26%,-13%)。
步骤四:残差相对值状态的一步转移矩阵为:
用2010年5月上证指数残差相对值作为初始状态,认为初始分布I(0)=(0 0 0 1)。
则下个月的绝对分布为:
所以下个月预测对象最有可能处于状态d,3,即ε(t)∈d=(13%,26%],又由,所以灰色马尔可夫预测区间x(21)∈(2691.19,3000.796),灰色马尔可夫预测值x(21)=2845.993。
步骤五:同理再进行预测可以依照上面的步骤得到以后各个月份的预测区间及预测值。预测区间及预测值详见表2。
由以上数据经灰色马尔可夫预测模型计算得,见表2。
以上计算结果表明的灰色马尔可夫预测模型较单一的灰色预测模型和马尔可夫预测模型来说精度和准确性都有了明显的提高。
3 结论
灰色马尔可夫链预测模型是根据灰色模型和马尔可夫链模型思想建立的,同时建立在历史数据的统计分析基础上,因此历史数据越准确,精度越高。该模型用灰色GM(1,1)模型预测曲线来反映上证指数的发展趋势,以此预测曲线为基准,再运用马尔可夫模型来寻找上证指数的转移状态,不仅考虑了数据序列中的演变规律,而且通过状态转移概率矩阵的变换提取数据中的随机响应,因而它将灰色模型和马尔可夫链模型的优点结合起来,克服各自的缺点提高预测精度,另外该模型的原理浅显易懂、计算过程不复杂,适用性比较强。该模型存在许多值得探讨的问题,它的预测精度与状态的划分有很大的关系,目前状态的划分没有统一的标准,所以还需要进一步的研究,另外该模型对波动性较大且有一定上升趋势的数据来说,可以取得比较好的预测效果,但并不是对任何数据都能适用。总之,灰色马尔可夫链预测模型有其应用价值,为投资者投资决策提供一定的理论依据。
摘要:用GM(1,1)预测具有良好的精确性和规律性,但对于随机波动性较大的股市行业,它的预测精度比较低,而马尔可夫模型可以克服波动性较大的局限性,弥补灰色模型的不足,因此将两者结合起来对股市进行预测将能提高预测的精度。本文依据上交所20个月末收盘指数预测后四个月的月末收盘指数范围,实证分析表明灰色马尔可夫链模型在股市预测中应用的可行性
关键词:灰色预测模型,马尔可夫模型,月末上证收盘指数,预测
参考文献
[1]牛东晓.电力负荷预测技术及其应用[M].第一版.北京:中国电力出版社,1998.
[2]刘克.实用马尔可夫决策过程[M].北京:清华大学出版社,2004.
[3]唐娜,桂预风,李宝.灰色马尔可夫模型应用于股指分析[C]//第五届中国不确定系统年会论文集,2007:195-198.
灰色马尔可夫 篇5
1 灰色模型的建立
1.1 原始数据处理
道路交通事故的原始数据具有很大的随机性, 为了找出其中的内部规律, 弱化原始数据的随机性, 增强其规律性, 把杂乱无章的统计数据列整理成有序的数列, 灰色理论中通常把原始统计数列进行累加或累减生成处理。
记X (0) 为历年来道路交通事故原始统计序列, X (1) 为X (0) 累加生成后的生成序列。则有X (0) = (x (0) (1) , x (0) (2) , x (0) (3) , …, x (0) (n) ) ,
X (1) = (x (1) (1) , x (1) (2) , x (1) (3) , …, x (1) (n) ) .其中:
1.2 GM (1, 1) 模型的建立[4,5]
作为道路交通事故的预测模型, 应该是随时间变化的单个变量的模型。因此, 可以用一阶的单个变量的微分方程来作为预测模型, 即GM (1, 1) 。对于历年交通事故统计次数的变量X (0) , 其相应GM (1, 1) 模型的微分方程为
令为待估参数向量, 即, 则在最小二乘准则下, 有
其中:
, 其中:
于是得到方程 (1) 的解为
. 从传统的GM (1, 1) 建模过程, 可以看出序列X (1) 在k点的灰导数是用差分来处理的, 即
如前所述, 灰色GM (1, 1) 模型比较适合成指数规律变化的系统预测, 灰色预测是以GM (1, 1) 模型为基础进行的预测, 几何图形的变化为指数型, 是一条较为平滑的曲线。因此, 对于时间短、数据资料少、波动性不大的预测问题, 只需要较少的几个数据就可建立模型进行预测, 而且这种短期预测精度较高。但是由于指数型变化是单调的, 所以作长期预测其预测值就会偏高或偏低。因此, 该模型不适合长期的、随机性和波动性较大的数据序列预测。
马尔可夫过程是研究事务状态及其转移的理论。它是通过对不同状态的初始概率以及状态间转移概率的研究, 来确定状态的变化趋势, 从而达到预测未来的目的。马尔可夫过程的特点是:当过程在时刻t0所处的状态为已知的情况下, 过程在时刻t (t>t0) 所处的状态与过程在时刻t0之前的状态无关, 这种特性称为无后效性[5]。因而, 将两模型结合使用, 建立灰色马尔可夫预测模型, 一方面利用GM (1, 1) 模型使数据序列满足马尔可夫的前提条件, 即无后效性和平稳过程等均值特点;另一方面, 马尔可夫模型又解决了对随机波动性较大序列的预测问题。本文尝试将两模型结合为改进灰色马尔可夫预测模型, 它的优点主要有:它既继承GM (1, 1) 模型的优点, 又对其存在的不足之处提出了改进, 并且融入了马尔可夫模型的预测方法, 提高了对随机性强的系统的预测精度。
2 灰色马尔可夫预测模型的建立[6,7,8]
2.1 状态划分
状态划分与预测正确与否有很大关系。根据数据序列不同特点采用不同方法进行划分, 主要有相对值法、平行曲线法、残差的标准化离差法及相对误差法。对于道路交通事故的预测, 采用平行曲线法比较合适。平行曲线法就是以
其中, Ai、Bi为根据待预测对象数据定的常数。由于
2.2 状态转移概率计算
由状态
2.3 预测值的计算
确定了系统的未来转移状态
3 实例分析
本实例采用1998~2007年共10年的全国道路交通事故数据 (见表1) , 来测试改进后的模型的效果, 使用改进的灰色马尔可夫模型进行预测。
3.1 建立GM (1, 1) 模型
表1为1998~2007年全国道路交通事故次数统计, 用传统GM (1, 1) 模型计算可得到 a=0.017 313, b=582 313.309 466.则
3.2状态划分
1998年~2007年的全国道路交通事故的平均值为529 552, 根据各年事故实际情况, 由于数据较为接近, 利用平行曲线法可将样本序列划分为2个状态:
⨂2:
.由此可以预测2008年的交通事故生成量最有可能处于状态⨂2, 最有可能的预测值是
. 根据公安部统计, 2008年全国共发生交通事故次数为265 204, 预测的误差为3.1%。说明使用改进的模型预测精度还是非常高的。
4 结束语
道路交通事故预测是进行公路交通网规划, 制定交通安全对策, 预防交通事故不可缺少的重要部分。传统的预测模型需要考虑的因素较多, 过程显得极为复杂, 而且预测的精度往往不尽如人意。
灰色马尔可夫预测模型综合了灰色GM (1, 1) 模型和马尔可夫模型的优点, 充分利用历史数据给予的信息。用灰色GM (1, 1) 模型预测曲线来反映系统宏观发展规律, 并以此预测线为基准, 再用马尔可夫概率矩阵来寻求系统的微观波动规律, 因而对随机波动较大的数据序列的预测具有较高的精度。本文采用灰色马尔可夫预测模型, 使用1998~2007年全国道路交通事故数据, 来预测2008年事故次数, 预测精度达到96.9%。
参考文献
[1]敖谷昌, 杨利.道路交通事故中机动车驾驶员人为因素与事故危害性关联分析[J].重庆交通大学学报, 2010 (1) :12-14.
[2]秦利燕, 邵春福, 赵亮.道路交通事故宏观预测模型[J].武汉理工大学学报:交通科学与工程版, 2010, 34 (1) :154-157.
[3]舒怀珠.基于灰色理论的道路交通事故预测模型综述[J].道路交通与安全, 2009 (6) :25-27.
[4]邓聚龙.灰色预测与决策[M].武汉:华中理工大学出版社, 1986.
[5]张大海, 江世芳, 史开泉.灰色预测公式的理论缺陷及改进[J].系统工程理论与实践, 2002 (8) :140-142.
[6]刘思峰, 党国耀, 方志耕, 等.灰色系统理论及其应用[M].北京:科学出版社, 2004.
[7]梅振宇, 王炜, 李铭, 等.高速公路交通生成预测的改进灰色马尔可夫模型[J].公路交通科技, 2004, 21 (12) :89-92.
灰色马尔可夫 篇6
关键词:灰色系统模型,马尔可夫链,预测
马尔可夫过程是一种无后效性的随机过程,马尔可夫链预测法是应用概率论中马尔可夫链的理论与方法,研究分析某些随机变化的动态系统的发展变化过程,并预测其发展变化趋势的一种预测方法,它是现代预测方法中的一种,具有较高的科学性,准确性和适应性,在现代预测方法中占有重要的地位.灰色预测就是通过原始数据的处理和灰色模型的建立,发现、掌握系统发展规律,对系统的未来状态做出科学的定量预测.本文把灰色预测和马尔可夫链预测两种方法结合起来,用GM(1,1)模型来揭示经济现象长期发展变化的某种总趋势,用马尔可夫模型来确定现象状态之间的转移,建立灰色马尔可夫预测模型,对社会物流总额进行预测研究.同时给出线性回归模型对物流总额的预测结果作为比较。
一、灰色系统GM(1,1)模型
把社会物流总额看作一个灰色系统,每年物流总额看作灰色量.选取从1995年到2004年共10年的数据作为灰色序列,即原始数据序列X(0),具体数据如下(由于数据较大,所以将数据做一处理:将小数点向后移动4位,即变为原数据的万分之一):X(0)(k)=(10.2230,11.0614,12.4138,12.9388,13.9717,17.1427,19.5422,23.3597,29.6595,38.3829).对X(0)作1-AGO得:X(1)(k)=(10.2230,21.2844,33.6982,46.6370,60.6087,77.7514,97.2936,120.6533,150.3128,188.6957).Z (1)为X (1)紧邻均值生成序列:Z(1)=(15.7537,27.4913,40.1676,53.6228,69.1800,87.5225,108.9734,135.4830,169.5042).则由
及GM(1,1)模型x(0)(k)+a z(1)(k)=b的最小二乘估计参数列满足,得到.确定模型为:,时间响应式为:.求出X(1)的模拟值并由还原求出X(0)的模拟值.记此次的拟合值为一次拟合值,它与原始数据间有残差和相对误差,具体数据见表2。
二、灰色马尔可夫链预测
马氏链预测是根据每年物流总额状态变化来确定状态转移概率矩阵从而确定预测值的方法。所以首先把10年的物流总额按相对误差划分状态。本文选择划分为3个状态,各状态划分的相对误差取值见表1。
由表1和表2计算状态转移概率矩阵P如下:
当预测对象的状态为i,若pi1,pi2,pi3中最大的为pij,j=1,2,3.则序列的下一个时期状态将处于j,但若第i行中有多个概率相近时,需要考虑P(2),P(3),….直至该行中只有一个概率最大,这样就确定了预测值变动的区间[j1,j2],得二次拟合值:
其中是原始数据的均值.1995年处于状态2,p2j中最大的是p22,所以1996年将处于状态2,该年值为:,其中为一次拟合值.依此计算下去可得二次拟合值,见表2。
三、结论
根据上述模型进行预测,结果见表3。得到2005年物流总额为47.4250,相对误差为1.60%,能够达到预测精度。从上述模型的预测结果来看灰色马尔可夫组合预测模型的收敛性比其他两种单一预测模型收敛性更好,预测精度相对较高。
参考文献
[1]刘思峰等.灰色系统理论及其应用(第三版)[M].科学出版社.2004年11月.
[2]唐启义,冯明光著.实用统计分析及其DPS数据处理系统[M].科学出版社.2002年5月.
[3]王艳玲.基于灰色马尔可夫链的农民平均收入实证研究[J].河南师范大学学报.2008年1月.
[4]中国物流年鉴2007年统计数据.
灰色马尔可夫 篇7
消费者物价指数, 俗称CPI, 是反映与居民生活有关的商品及劳务价格统计出来的物价变动指标, 常作为观察通货膨胀水平的重要指标。消费者物价指数升幅过大, 表明通胀成为经济不稳定因素, 因此过高的CPI往往不被欢迎。
我国物价自2009年7月后一直不断攀升, 足以见得此次物价上涨来势猛烈。推动本轮物价上涨主要原因主要可分为:国内因素——为应对金融危机我国实施的宽松的货币政策, 导致货币供应大于货币需求, 国内流动性充裕;国外因素——为应对经济衰退, 世界各国家为经济注入大量的流动性导致国外流动性泛滥和国际物价攀升;成本推动——原材料价格上涨与劳动力价格上涨导致成本上升, 推动了物价的上涨;供求紧张——极端的气候使各国粮食减产, 在供求关系的作用下, 各种商品价格被显著提高。
在这些因素的影响下, CPI指数不断上涨, 严重影响居民生活水平, 高涨的物价已经开始渗透到经济生活的各个角落, 老百姓可以明显感觉到物价上涨带来的负面影响。过高物价所带来的不良影响将不利于我国经济的健康与稳定, 因此对物价指数的研究一直受到各方关注。中央经济工作会议要求, “把稳定物价总水平放在更加突出的位置”逐见此物价问题的重要性。考虑到CPI具有明显的动态特征和不确定性, 符合灰色系统的特点, 本文将引入灰色—马尔可夫模型对消费物价指数进行研究预测。
灰色系统理论提供了在贫信息情况下解决系统问题的途径, 其基本模型GM (1, 1) 的基本思想是:通过时序数据累加来滤去原始数据中可能混入的随机量, 从波动的时间序列中寻找规律, 挖掘原始序列的内在特征。马尔可夫模型的最重要特点是无后效性, 其预测是根据状态之间的转移概率来推测系统未来的发展变化。此两种预测方法虽然都可用于时间序列的预测问题, 但是由于马尔可夫模型克服了物价指数随机波动性较大的局限, GM (1, 1) 预测又使预测具有精确性和规律性, 两者之间具有很大的互补性。恰恰由于以上分析得知CPI的影响因素众多, 且影响机理比较复杂, 因此将此两种预测模型合二为一, 用GM (1, 1) 模型来揭示经济现象长期发展变化的某种总趋势, 用马尔可夫模型来确定现象状态之间的转移, 建立灰色马尔可夫预测模型, 可以对消费者物价指数进行精准预测。
二、理论基础
(一) 灰色模型介绍
数列预测是对某一指标的发展变化情况所作的预测, 其预测的结果是该指标在未来各个时刻的具体数值。其基础是累加生成数列的GM (1, 1) 模型。
设x (0) (1) , x (0) (2) , …, x (0) (M) 是所要预测的某项指标的原始数据, 一般而言, 它是一个不平稳的随机数列, 非平稳数据趋势无规律可循, 也无法用回归预测法对其进行预测。
如果对原始数列作依次累加生成处理, 得到一个新的数列{x (1) (t) }tM=1。这个数列与原始数列相比较, 随机性程度大大弱化, 平稳程度大大增加。对于这样的数列, 其变化趋势可以近似用如下微分方程描述:
在 (1) 式中, a和u的值可以通过如下最小二乘法拟合得到:
在 (2) 式中, YM为列向量YM=[x (0) (2) , x (0) (3) , …, x (0) (M) ]T;
微分方程 (1) 式所对应的时间响应函数为:
由 (3) 式可以得到对一次累加生成数列的预测值x (1) (t) 可以求得原始数的还原值:
(4) 式中, t=1, 2, …, M, 并规定x (1) (0) =0。
原始数据的还原值与其观测值之间的残差值ε (0) (t) 和相对误差值q (t) 如下:
(一) 马尔可夫模型的基本概念介绍
1. 状态转移概率
事件的发展变化过程中, 从某一种状态出发, 下一时刻转移到其他状态的可能性, 称为状态转移概率。由状态Ei转为状态Ej的状态转移概率P (Ei→Ej) 就是条件概率P (Ej/Ei)
2. 状态转移概率矩阵
假定某种被预测的事件有E1, E2, …, En, 共n个可能的状态。记Pij为从状态Ei转为状态Ej的状态转移概率, 作矩阵P= (Pij) 为状态转移概率矩阵。
3. 状态转移概率矩阵的计算
计算状态转移概率矩阵P, 就是要求每个状态转移到其他任何一个状态的转移概率Pij (i, j=1, 2, …, n) 。为了求出每一个Pij, 我们采用频率近似概率的思想来加以计算。
4. 状态概率πj (k)
πj (k) 表示事件在初始 (k=0) 时状态为已知的条件下, 经过k次状态转移后, 第k个时刻 (时期) 处于状态Ej的概率。从初始状态开始, 经过k次状态转移后到达状态Ej这一状态转移过程, 可以看作是首先经过 (k-1) 次状态转移后到达状态Ei (i=1, 2, …, n) , 然后再由Ei经过一次状态转移到达状态Ej。根据马尔可夫过程的无后效性及Bayes条件概率公式, 有
若记行向量π (k) =[π1 (k) , π2 (k) , …, πn (k) ], 可得逐次计算状态概率的递推公式:
三、实证分析
(一) 模型建立过程
给出中国2010.1-2011.8共二十个月的消费者物价指数, 用灰色马尔可夫模型对2011.9-11月份的消费者物价指数进行预测, 将预测值与实际值相比较, 验证灰色马尔可夫模型预测消费者物价指数的可行性和准确性。
步骤一:从国家统计局的统计年鉴中得出数据。
步骤二:对原始数据进行一次累加, 生产新的数据列。
步骤三:构造矩阵B和常数矩阵Y。
步骤四:用最小二乘法求解, 得到a=-0.05168, b=2.5167
步骤五:将a, b代入, 得到该模型的时间相应模型然后通过累减还原, 最终得到x (1) (ˆt+) 1 (28) .2545e0.05168t。
步骤六:根据方程求出原始数据的拟合值, 拟合值如表1所示。并根据公式求出残差的相对值, 求得残差相对值序列为ε (t) ∈ (-70%, 17%) , 由于第一个月的CPI偏离常态将影响预测结果, 使其变得不准确, 在此将其舍去, 所以残差的相对值范围为ε (t) ∈ (-18%, 17%) 。
步骤七:根据实际情况将残差相对值序列平均分为五个状态划分的状态以及每个残差相对值所处的状态如表1最后一列所示。
步骤八:计算一步转移概率, 预测对象2011年4月29号属于状态1, 那么下个月的绝对分布为, 根据处在每个状态的概率大小可以得到下个预测对象最有可能处于状态2, 即ε (t) ∈E2 (-1%, -4%) , 由此可得:x (0) (21) ∈ (6.445, 6.8792) 。同理再进行10月份和11月份的预测可以依照上面的步骤依次进行, 得到的各个月份的预测区间及预测值详见表2, 其中预测值为预测区间的中点。
(二) 结果探讨
CPI数据具有滞后性, 选择合理的指标, 运用合适的方法预测CPI, 对政策制定者把握宏观经济走向, 相机抉择财政政策和货币政策起至关重要的作用, 从表2中可以看出, 灰色—马尔可夫模型的预测值与单纯的灰色模型的预测值相比要更接近实际数据, 一方面说明基于灰色—马尔可夫模型对CPI进行预测是行之有效的预测方式。另一方面说明将灰色模型与马尔可夫预测法相结合提高了数据的利用率, 因此提高了预测水平。
从CPI近期的变化趋势判断, CPI一直有呈明显增长的态势。考虑到目前影响CPI增长的因素无法立刻消除。因此, 国家此时应当制定切实可行的宏观调控政策, 应对此次物价上涨风波。采用切实可行的财政政策与稳健的货币政策, 严格控制货币供应量。加强外汇管理, 在食品流通方面增开绿色通道, 减小流通阻力带来的成本增加。同时应该建立健全的储备机制, 打击投机倒把行为。
(三) 模型的改进方向
模型虽然具有相对精准的优点, 但是可以很明显地看出其有一定的人为性, 比如在划分残差范围时没有一定的原则, 不同的划分方式将使预测结果产生一定的偏差, 因此主观性很强, 不够严谨。另外灰色—马尔可夫模型的预测值就简单设定为其对应预测区间的中点值稍微有不合理之处。因此模型还需在严谨性与客观性方面加以改进。
(四) 模型的拓展
既然前文提到灰色—马尔可夫模型适合运用于波动无规律, 影响因素较为复杂的时序数据的预测上, 而且事实也确实证明其预测比较精准, 因此除了CPI以外, 我们还可以将模型运用于其他方面, 如股价的预测、GDP的预测等, 深入扩展还可以运用到各种领域上, 可以解决现实生活中的种种问题, 方便和改善生活, 具有良好的现实意义。
参考文献
[1]简艳群.灰色-马尔可夫链模型在股市预测中的应用[J].价值工程, 2010 (24) .
[2]何俊, 王传丽.基于改进G M (1, 1) 的C P I灰色预测模型[J].商业研究.2009 (03) .
[3]严丽玉, 罗明安.厦门市GDP二次参数拟合灰色马尔科夫链预测模型[J].厦门理工学院学报, 2010 (03) .
灰色马尔可夫 篇8
高斯混合模型 (GMM) [2]是神经学中广泛使用的一种统计学模型, GMM不仅与脑部MR图像的分段常数性质一致, 而且具有较低的计算复杂度。GMM的模型参数可以利用期望最大化 (EM) 算法根据最大似然 (ML) 准则来估计, 然而, 基于EM的ML估计具有过度拟合和容易限于局部最优解的缺点。为了克服这些缺点, 使用几种全局优化技术来替代EM算法, 例如, 文献[3]在似然估计中结合了遗传算法, 提出了GA-EM算法。此外, 当先验知识可用时, 最大后验概率 (MAP) 估计是ML估计的一种常见替代方法。文献[4]提出一种MAP-MRF框架来求解图像分割问题, 通过将体素类标签建模为马尔可夫随机场 (MRF) 来表示体素空间依赖性的先验。文献[5]提出一种基于全局随机搜索的推理方法, 即马尔可夫链蒙特卡尔 (Markov chain monte carl, MCMC) 推理, 用来替代确定性程序。文献[6]受免疫机制启发, 提出一种克隆选择算法 (clonal selection algorithm, CSA) , 基于克隆选择理论, 选择能够识别抗原的抗体来进行繁殖, 繁殖的细胞会通过一个亲和力成熟过程来改进它们对抗原的亲和力, CSA模仿对抗原刺激免疫应答机制来实现全局最优。
将CSA和MCMC技术融合到隐马尔可夫随机场 (hidden Markov random field, HMRF) 模型估计中, 提出一种用于脑部MR图像分割的HMRF-CSA算法。首先, 通过MCMC方法近似最优标签配置, 然后, 由CSA算法估计HMRF模型参数。用全局随机优化技术替代确定性搜索程序, 以此提高分割算法的鲁棒性, 同时, 将MR图像建模为分段常数图像的乘法分量, 根据MCMC推断方法获得的中间分割结果来评估图像的不均匀性。通过仿真实验, 将本文HMRF-CSA算法与现有的GA-EM方法、可变形共同分割 (D-C) 算法、SPM软件包中的统一分割算法和FMRIB软件库 (FSL) 上的HMRF-EM分割算法进行比较, 结果表明该算法具有更好的分割精度。
1 相关技术
1.1 图像不均匀模式
由图像采集的不完善所导致的图像不均匀性, 或称偏场或强度非均匀性 (INU) [7]是MR图像分析的难点之一。设定y={yi;i=1, 2, …, N}表示一副脑部MR图像, 其中yi表示在体素i处的强度, N表示体素的数目, 未知偏场B={bi;i=1, 2, …, N}通常建模为y的乘法分量, 如下式所示
式 (1) 中, 是理想图像, 是附加的高斯白噪声。由于图像中偏场B变化缓慢, 所以可将它定义为在整个图像域上的一个平滑函数。采用正交多项式{Wj:j=1, 2, …, NOP}作为偏置函数来近似偏场[8]
式 (2) 中, φ={φj:j=1, 2, …, NOP}表示组合系数, NOP= (D+1) (D+2) /2是多项式的数目, D是多项式的度。
1.2 统计学模型
假设体素强度y={yj;j=1, 2, …, N}符合GMM;从有先验概率πk的高斯分布N (μkΣk) 中独立采样每个强度值yj, 观察图像的似然, 计算式如下:
通过最大化上述似然函数来估计最优GMM参数, 确定参数后, 利用贝叶斯分类器对每个体素进行分类, 以此求解脑部图像分割问题[7]。
为了将空间约束融入到这个模型中, 本文应用MRF到模型类标签x={xj;j=1, 2, …, N}中, 根据Hammersley-Clifford理论, 类标签p (x) 的先验联合分布符合Gibbs分布。在MAP-MRF框架下, 图像分割等价于通过最大化其后验概率寻找最优配置x*
式 (4) 中, Θ={μk, Σk;k=1, 2, …, K}表示模型参数, p (y|x;Θ) 是图像似然, p (x) 是空间先验。
本文将图像强度y当作在相同图像点阵中另一个随机场建模的模型, 然后将代表潜类标签的MRFx变成HMRF。在这个模型中, 将图像分割问题制定为配置x和参数Θ的最大联合概率
式 (4) 中后验与式 (5) 中联合概率之间的差是惩罚项p (Θ|y, x) , 用于检查模型参数是否与配置x给出的观察值一致。
2 HMRF-CSA算法提出
HMRF模型主要用来估计式 (5) 中的最优类标签和模型参数 (x, Θ) , 估计过程可划分成两个相互依赖的优化步骤:搜索最优配置x*和学习最匹配模型参数Θ*。使用下式三个步骤的迭代程序来实现HMRF模型估计
式 (6) 中, f (., .) 是基于观察y和分割结果x纠正偏场的函数, t∈{1, 2, …, Tmax}表示当前迭代数目。在每次迭代中, 首先采用MCMC方法实现MRF-MAP估计, 在近似的分割结果下估计偏场, 然后利用CSA学习HMRF模型参数, 当达到最大迭代数目或分割结果变成稳态时迭代停止。
2.1 MCMC体素分类
上述迭代步骤中第一步是通过MRF-MAP近似寻找最优配置x*, 使用MCMC方法求解这种优化问题, 根据式 (4) , 对于给定任意特定配置x, 假设yj相互独立且符合基于参数Θk={μk, Σk}的多元高斯分布, 则似然为
MRFx的联合分布可表示为Gibbs函数[9],
式 (8) 中, Z是规范化常量, Vc (x) 表示派系c的潜力, C是根据邻域系统确定的所有派系的集合, T是温度参数。本文使用Potts模型表示派系潜力, 运用式 (7) 和式 (8) 到式 (4) , 并对其进行负对数变换, 得到
根据模拟退火MCMC方法, 为温度参数T定义一个冷却进度表,
式 (10) 中, i=1, 2, …, I表示MCMC算法迭代的数目, C是冷却因子, 本文设置T (0) =4, C=0.97。若给定一幅脑部MR图像y和标签x (0) 的初始配置, 则可计算出模型参数。定义用来表示从x (i) 随机移动的跳跃密度Q (.|x (i) ) 符合高斯分布, 每次迭代中, 从建议密度Q (x* (i+1) |x (i) ) 提取一个候选x* (i+1) , 从均匀分布u (0, 1) 提取一个随机序列, 计算每个体素j的接受率
如果uj<αj, 接受模拟xj (i+1) =xj* (i+1) , 否则拒绝它并保持类标签与上一次迭代xj (i+1) =xj (i) 相同, 当达到最大迭代次数时停止, 如算法1所示。
算法I:体素分类的MCMC采样
2.2 偏场校正
在利用MCMC体素分类之后, 可以获得分割结果x*和最小能量E={Ekj, xj=k, j∈S}。归一化后验概率n={nkj;j=1, 2, …, N, k=1, 2, …, K}作为软分割结果
根据理想MR图像的分段常数性质, 定义软分割与对应平均μk的积作为存储的图像
本文使用奇异值分解 (SVD) 求解下列最小二乘拟合问题, 以此估计偏场。
式 (14) 中, ./表示点对点划分, 根据估计的组合系数, 获得偏场
偏场损坏的图像可恢复如下:
2.3 CSA进行参数估计
第三步是通过最大化后验概率p (Θ|y (t) , x (t) ) 来学习当前图像强度y (t) 和配置x (t) 给出的最优参数
式 (17) 中, p (Θ) 是参数的先验概率, 这个先验指的是基于马尔可夫性质信息p (Θkj) =p (xj=k|kj) 的体素, 可以通过MRF能量计算得到, 定义它作为这些项的混合来平衡参数的收敛和多样性, 对于每个体素j∈S
式 (18) 中, v是平衡常量, Πkj=πk表示分体素全球先验。给定任意具体参数集Θ, 即可计算式 (18) 中所示优化问题的目标函数。为了实现全局最优, 采用CSA[10]求解该问题, 以群体方式模拟所有可能参数。CSA是一种进化优化算法, 通过迭代生成一群编码抗体来寻找全局最优解。本文中, 抗体群np设为100, 定义每个抗体为一个候选参数集Θ, 将有特定抗原的抗体Θk的亲和力定义为后验似然p (Θ|y, x) , 迭代优化过程由下列六个主要步骤组成:
第一步:评估每个抗体的亲和力, 根据其亲和力按降序排列所有抗体;
第二步:从当前群体中选择Ns个抗体, 克隆它们形成克隆群。对于亲和力排序为j的抗体, 定义其克隆的数目正比于其亲和力排序, 如式 (19) 所示。
式 (19) 中, β是常量, 用来控制克隆率, round (.) 用来将实数变换到与其最接近的整数;
第三步:分别对概率为phm和pre的克隆群运用超突变和受体编辑操作。超突变是在动态范围±10%内随机改变抗体的值, 目的是局部搜索最优解。受体编辑是在动态范围±100%内随机改变抗体, 实现全局搜索;
第四步:评估克隆群中抗体的亲和力, 根据其亲和力按降序排列;
第五步:选择克隆群中排名靠前的抗体代替记忆细胞集中较低亲和力的40%抗体, 保证记忆细胞集保存迄今为止获得的最优解, 以便最高亲和力的抗体按代递增;
第六步:用随机生成的抗体代替剩余集中具有最低亲和力的10%抗体, 对新群体引入多样性。
重复迭代这个过程直到达到最大迭代数目, 如图1所示。
2.4 总结
给定K-平均算法产生的初始分割结果, HMRF-CSA算法迭代执行基于MCMC的体素分类、偏场校正和基于CSA的模型参数估计, 直到算法收敛。一旦达到收敛, 则获得最终分割结果、偏场和模型参数。HMRF-CSA算法的主要步骤见算法II。
算法II:HMRF-CSA脑部图像分割算法
3 实验结果
本文从Brain Web数据集[11]获取仿真TI加权脑部MR图像, 比较提出的HMRF-CSA算法与现有的e HMRF算法、GAMIXTURE包中GA-EM算法、D-C算法、FSL包中的经典HMRF-EM算法和SPM包中的统一分割程序。Brain Web数据集提供的一组仿真脑部图像, 这些图像具有各种INU和噪声级别的解剖模型仿真, 每个仿真研究的维度为181×217×181, 体素大小为1 mm×1 mm×1 mm。
图2分别显示了仿真图像中具备40%INU和7%噪声的第88个横切片, 偏长矫正图像, 估计的偏场, 使用六种算法获得的分割结果和地面实况组织图。其中图2 (a) 表示仿真图像的第88个横切片 (7%噪声和40%INU) ; (b) 表示INU校正图像; (c) 估计的INU; (d) 表示HMRF-EM算法的结果; (e) 表示D-C算法的结果; (f) 表示SPM算法的结果; (g) 表示GA-EM算法的结果; (h) 表示e HMRF算法的结果; (i) 表示HMRF-CSA算法的结果; (j) 表示地面实况。可以看出, 本文算法产生的分割结果比其他算法产生的结果更接近地面实况。
接下来, 在两组仿真MR图像上对这些算法进行进一步比较。第一组MR图像包含有20%INU和噪声级别范围从1%到7%的四个图像, 使用骰子相似度系数 (DSC) [12]定量评估每个脑部组织类型分类的性能。
式 (20) 中, Vs (k) 是分割结果中脑部组织类k的体, Vg (k) 是在地面实况上对应的体, |V|代表体V中体素的数目。通过正确分类的脑部体素百分比来计算分割精度, 并评估整体精度。图3表示六种算法获得的分割精度。
从图3可以看出, 在大部分仿真图像中, 本文算法在划分每个脑组织和分类整个脑部体方面都具有较高的精度。而且, 随着噪声和INU级别的增加, 提出算法的精度下降幅度比其他算法低, 这表明本文算法具有较强的抵制噪声和INU影响的能力。
第二个测试组包含有40%INU和噪声级别范围从1%到7%的四个图像, 六种算法获得的分割精度如图4所示。
从图4可以看出, 本文算法能在高噪声和INU级别下保持良好分割性能。
4 讨论
4.1 参数设置
本文提出的HMRF-CSA算法中, 有三组需要近似的参数, 包括MCMC推断、INU估计和基于CSA的参数近似。在INU近似中, 以正交多项式的阶来权衡考虑近似精度和计算复杂度, 由于INU变化非常慢, 对于INU近似, 10个三阶多项式已经足够。式 (19) 中权重参数v决定MRF先验, 较大的v能使MRF的作用更大, 另一方面, 小v则更支持GMM先验。CSA本身需要很多参数, 文献[13]对此进行了详细讨论。本文使用CSA程序的经验参数设置:群体大小Np=100、记忆集大小Nm=0.3Np、选定抗体的数目Ns=0.5Np、克隆率常量β=0.5, 超突变概率phm=0.8、受体编辑概率pre=0.1和最大的代Nt=20。
4.2 计算复杂度
计算机程序的性能与许多因素有关, 包括计算机处理能力、数据表示、编程语言和编码实现等[14]。本文评估了HMRF-CSA算法的计算复杂度, 本文算法在每次迭代中顺序执行MCMC推断、偏场估计和基于CSA的参数估计。设定对于有N个体素的一副图像, MCMC推断的计算复杂度为O (N) ;偏场估计仅进行一些矩阵计算;基于CSA的参数估计的复杂度为O (Np+NcK) , 其中Np表示群体大小, Nc表示总克隆数目。提出的迭代分割算法的迭代次数达到wmax后停止, 其线性整体计算复杂度O (N+Np+NcK) 。需要注意的是, MCMC方法的主要缺点是需要大量仿真图[15], 然而, 由于CSA为MCMC方法配置了一个良好的开始状态, 所以本文算法不需要许多仿真图。同时, MCMC方法的输出使CSA在有限代数之后成熟, 因此, 尽管本文HMRF-CSA算法涉及耗时的MCMC和CSA程序, 然而, 其计算复杂度只稍微高于传统分割方法。
5 结论
提出了HMRF-CSA脑部MR图像分割算法, 在基于HMRF模型分割中结合了CSA和MCMC, 本文算法能够有效的用于基于HMRF模型估计的图像分割问题。在仿真脑部MR图像上进行实验, 将本文算法与GA-EM算法、D-C算法、SPM和FSL软件包算法进行比较, 实验表明该算法获得了更好的分割精度。