一种新的高阶多元马尔可夫链模型及其应用

阙素琴，张晓薇，滕忠铭

(福建农林大学 计算机与信息学院，福建 福州 350002)

关键字：马尔可夫链；多元马尔可夫链；高阶多元马尔可夫链；参数估计

马尔可夫链是人力资源[1]、金融[2]、互联网应用[3]、音乐[4]、软件测试[5]、土地覆盖变化[6]、能源消耗[7]、语音识别[8]、微生物基因[9]、DNA序列[10]、信用风险[11]等许多研究领域的重要工具。探索不同分类数据序列之间的关系，建立更准确的预测模型是一个有意义的研究课题。

分类数据序列是由相互关联的相同或相似的资源产生的，针对分类数据序列的预测，已经提出了很多不同的模型，如Ching等[12-13]于2002年提出了一阶多元马尔可夫链模型，并于2009年改进了该模型，即在正态模型的后面添加负关联部分，其特点是加快计算平稳或稳态解的收敛速度。Wang等[14]提出了一种改进的简约多元马尔可夫链模型，并发现了具有变异性的新收敛条件，可用于提高预测精度和收敛条件的最小化尺度。此外，Ching等[15]于2008年还提出了高阶多元马尔可夫链模型，并以DNA序列、销售需求为例说明了所提出模型的预测能力。

随着马尔可夫链模型的发展及其应用，序列的数目可能会越来越大，序列的增大意味着有更多的参数需要被估计，从而不可避免地导致较高的计算成本。因此，减少模型需要被估计的参数的个数对数值计算具有重要的帮助。因此，Wang等[16]于2014年提出了一种增量式的多元马尔可夫链模型，此模型探究了原始分类序列中增加一个新的序列后与原序列之间的关系，但此模型并没有考虑高阶多元马尔可夫链的情形。为此，本文提出一种新的高阶多元马尔可夫链模型，其本质是Wang的增量式模型在高阶多元马尔可夫链模型上的推广，并通过数值实例说明本文所提出模型的有效性。

1 高阶多元马尔可夫链模型

定义 1随机过程{St,t∈N}，N={0，1，2,…}，状态空间 M={0，1，2,…，m}，对于任意的 t∈N，i，j，j0，…，jt-1∈M。若满足

P{St+1=i∣S0=j0，S1=j1，…，St-1=jt-1，St=j}=P{St+1=i∣St=j}

则称此随机过程为马尔可夫链。

定义2任意i，j∈M的，条件概率Pij=P{St+1=i∣St=j}称为转移概率。

令m×n阶矩阵P=[Pij]是由Pij为元素构成的一步转移概率矩阵，其中。那么一元的马尔可夫模型可表示为xt+1=Pxt，t∈N，这里xt表示t时刻状态概率分布。当问题涉及到多个分类数据序列时，模型可以自然推广到多元马尔可夫模型。当有s个分类数据序列，多元马尔可夫链模型可表示如下形式

这里λjk是模型中需要被估计的参数。

当多元马尔可夫链中t+1时刻概率分布不仅仅依赖于时刻t的所有序列（包括它自己）的状态概率分布时，Ching等[15]提出了一种如下形式的高阶多元马尔可夫链模型

在模型（2）中，如果令

那么（2）可用矩阵的形式表示如下

注意到在模型（2）中一共有ns2个需要被估计。令X是Xt的平稳分布，即X满足表达式X=QX。通过计算每个序列中各个状态出现的比例，得到X的近似，那么参数可通过求解如下优化问题得到

这里‖·‖∞表示向量的无穷大范数。具体的估计过程可详见参考文献[15]。特别地，如果问题中增加一个分类序列，则模型（2）中已估计的ns2个需要被重新计算，从而估计参数一共为 n（s+1）2个。这意味着原先已估计的个的信息，在增加一个序列后被完全抛弃了，这给数值计算增加了较大复杂度。为此，提出了一种增量式的高阶多元马尔可夫链模型，使得已估计的信息可以得到充分利用，从而减少模型中需要估计参数的数目。

2 新高阶多元马尔可夫链模型

假设原s个分类数据序列的高阶多元马尔可夫链模型如（2）所示，当增加一个新的分类数据序列时，本文提出一种新的高阶多元马尔可夫链模型，表示如下

其中 B*（ij）的定义与式子（3）和（4）类似，只需将式子（3）和（4）中的用代替。另外，注意到当i≠j≠s+1 时，有 B*（ij）=ljB（ij）。

3 参数估计

4 数值实验

4.1 简单例子

考虑下面的 3 个分类数据序列 Y（1）={2,1,3,3,4,3,2,1,3，3,2,1}，Y（2）={2,4,4,4,4,2,3,3,1，4,3,3}，Y（3）={1,2,3,4,1,4,4,3,3，1,3,1}，其中Y（3）序列是新增加的数据序列。通过计算转移频率矩阵，并经正规化处理后得到转移概率矩阵，且分类数据序列的稳态概率分布如下

已知前两个分类数据序列的高阶多元马尔可夫链模型（阶数n=2）中的为已估计参数，增加一个新的序列s3后，通过求解相应的线性规划问题估计参数lj和，从而得到新的高阶多元马尔可夫链模型为

为了说明新模型的有效性，本文探究了传统高阶多元马尔可夫链模型和新高阶多元马尔可夫链模型的CPU所用时间t（单位/s）、估计参数个数np（单位/个）、以及线性规划中的目标函数值w，其中传统高阶多元马尔可夫链模型的参数是通过使用参考文献[15]中的方法计算得到。通过MATLAB软件编程，并使用其内置函数tic和toc记录CPU所用时间t，具体的数值结果详见表1，其中新模型和传统模型估计参数的个数np分别由公式n（2s+1）+s和n（s+1）2计算得到。

表1 两种高阶多元马尔可夫链模型的w、t和np比较

由表1可以看出：传统高阶多元马尔可夫链模型和新高阶多元马尔可夫链模型在线性规划中的目标函数值w一致，且新模型CPU所用时间少于传统模型的所用时间，新模型的估计参数个数也少于传统模型的估计参数个数。因此，新模型在时间消耗和参数控制上都优于传统的模型。

4.2 销售需求预测

第二个例子以销售需求序列来说明新的高阶多元马尔可夫模型的优越性。由于市场需求波动较大，生产计划和库存控制对成本有着直接影响的关系。因此，研究库存空间需求与整体增长的销售需求之间的相互作用是公司迫切需要解决的问题。我们的目标是预测市场的销售需求，使成本最小化。假设产品被分为6种可能的状态(1,2,3,4,5,6)，例如：1=无销售量，2=非常低的销售量，3=低销售量，4=标准销售量，5=高销售量，6=非常高的销售量。从营销和生产计划的角度来看，这样的标签是有意义的。

假设已给出了产品A，产品B，产品C和产品D这四类数据的序列,数据来源于参考文献[15]。通过计算每个序列中每个状态出现的比例，可以得到四个分类数据序列的初始概率分布如下

同样的，考虑二阶多元马尔可夫模型，已知前四个分类数据序列的模型，其中的为已估计参数，此时增加一个新的序列（产品E）

利用新的高阶多元马尔可夫链模型估计得到：

根据所建的新高阶多元马尔可夫链模型得出：产品A与产品B密切相关的。特别是A产品的销售需求很大程度上取决于B产品，主要原因是A产品和B产品的化学性质是一样的，为了营销只是包装不同而已。此外，产品C和产品E的关系也很相似，原因是它们有相似的产品味道。结果表明高阶马尔可夫模型对分析销售需求关系是符合产品的实际情况，这对销售需求的预测具有重要意义。

此外，本文利用传统的和新的高阶多元马尔可夫链模型得到的预测序列Y（k）在t时刻的概率分布中取概率最大的那个状态作为它的下一时刻的状态，即

同样的，在这个例子中，比较了传统的高阶多元马尔可夫链模型和新的高阶多元马尔可夫链模型在线性规划中的目标函数值（w）、CPU所用时间（t）和参数估计数目(np)的情况，具体见表2。另外，为了比较传统和新的高阶多元马尔可夫链模型预测销售需求的效果，表2列出了它们对各数据序列下一时刻所预测的状态的概率。其中，传统高阶多元马尔可夫链模型和本文所提出的新高阶多元马尔可夫链模型对产品A、产品B、产品C、产品D和产品E的预测状态都为状态6，状态6，状态6，状态2和状态2，且它们所对应的概率分别为0.4275，0.3978，0.6270，0.3581和0.3571。由表2可以得出这两个模型在线性规划中的目标函数值w以及预测效果上基本是一致的，但是新的高阶多元马尔可夫链模型CPU所用时间（0.0320 s）少于传统模型的所用时间（0.0780 s）；新模型的估计参数个数仅仅22个，远少于传统模型的估计参数个数（50个）。为了进一步说明新模型的有效性，本文还分别分析了“产品 B、C、D、E 增加 A”、“产品 A、C、D、E 增加 B”、“产品 A、B、D、E 增加 C”以及“产品A、B、C、E增加D”的模型，得到了相似的结论。

表2 两种高阶多元马尔可夫链模型的w，t，np以及预测状态的比较

5 结论

本文提出了一种新的高阶多元马尔可夫链模型，该模型通过对原s个分类数据序列的高阶多元马尔可夫链模型的分析，建立原分类数据序列与增加一个序列后的新序列集之间的关系。正如在模型（5）的描述中所分析的那样，新的高阶多元马尔可夫链模型只需要估计n（2s+1）+s个参数，这少于传统模型的估计参数个数n（s+1）2。本文同时给出了新高阶多元马尔可夫链模型的参数估计方法，该估计方法最终将参数的计算转化为一个线性规划问题。最后通过两个具体的数值例子说明了新模型在节省计算资源方面的优越性，并且其预测性能与传统模型基本保持一致。