基于HMM的食品安全风险预测方法研究

马永军，刘宇姗，侯阳阳，刘 洋

（天津科技大学计算机科学与信息工程学院，天津 300222）

基于HMM的食品安全风险预测方法研究

马永军，刘宇姗，侯阳阳，刘 洋

（天津科技大学计算机科学与信息工程学院，天津 300222）

以乳品为例，提出了一种基于隐马尔可夫模型（HMM）的食品安全风险预测方法.依据HACCP对乳品供应链进行分析，找出供应链各个环节的关键控制点，并将其危害影响因素作为HMM的量化指标.分析时将供应链分为4个环节，每个环节最终的安全等级概率作为下一环节的初始状态概率分布，利用HMM预测供应链风险等级和计算风险值，从而为供应链风险评价提供依据.

食品安全；HMM模型；乳品；供应链；风险预测

食品安全是关系国计民生的重大问题，近年来食品安全问题的预警和治理已受到高度关注.其中，风险预测作为风险预警和治理的基础受到格外重视.风险预测是利用适当的风险预测工具和方法来确定安全风险等级和优先控制顺序的过程，从而为实施有效的风险管理措施提供决策支持.

目前常用的食品安全风险预测方法主要分为定性和定量两大类，其中定性的方法主要有灰色系统法、德尔菲法、层次分析法等.灰色系统法［1］是对不确定系统的研究方法；德尔菲法［2］能够充分发挥各位专家的不同意见，评定过程公正；层次分析法［3］是食品安全风险分析与预测中应用最为广泛的方法.但是它们的共性缺点就是在运算过程中受主观因素影响较大，使得结果的精确度不够，并且不能根据系统状态的变化实时改变风险预测结果［4］.定量的分析方法主要有贝叶斯网络法、人工神经网络和支持向量机（SVM）等.贝叶斯网络法［5］是一种基于概率推理的图形化概率网络，能在有限的不确定信息条件下进行学习和推理，但是贝叶斯网络构造繁琐，实际应用时还需反复交叉不断完善，易用性不好.人工神经网络［6］具有高度的自适应能力，但是学习速度较慢，而且算法容易陷入局部极值.SVM［7］本质上是将问题转换为凸优化问题，可以保证找到全局极值，但在将低维空间向高维空间映射时，又存在计算耗时问题，因此不适合对于大样本数据的分析处理.

食品安全分析和食品供应链密不可分，食品供应链主要由原料生产、生产加工、储运和消费4个环节组成，食品安全风险存在于食品供应链的各个环节中，而且存在食品安全风险迁移问题，即一个环节的风险会向下一个环节传递.在目前的方法中，对风险的判断一般是对4个环节的风险进行简单平均或加权平均，未充分考虑食品安全风险的形成机理，割裂了4个环节之间的风险传递关系.

HMM模型是序列数据处理和统计学习的一种重要概率模型，具有建模简单、数据计算量小、运行速度快等特点［8］，其在风险预测分析领域的应用较广泛［9］.HMM模型被用来对网络系统安全进行实时风险预测，表现出了很好的适用性和扩展性［10］.此外，HMM模型中的初始状态分布矢量参数，对应初始状态时各风险等级出现的概率，该参数可以有效表示食品供应链各环节对下一环节的影响程度.因此本文以乳品为例，利用HMM模型进行食品安全风险分析与预测.

1 HMM模型

HMM模型是在马尔可夫链的基础上发展起来的，在HMM中观察到的事件是状态的随机函数，该模型是一个双重随机过程，即一个观察状态，一个隐藏状态.

HMM的数学表达式为

其中：N为HMM中状态的个数；M为HMM中对应的观测值个数；π是初始状态的分布矢量；A为状态转移概率矩阵；B是观察向量的概率矩阵.

初始状态的分布矢量可表示为

也可以简单表示为

状态转移概率矩阵可表示为

元素aij表示从状态i转移到状态j的概率，即

观察向量的概率矩阵可表示为

元素bjk表示在j状态的情况下，观察状态k出现的概率，即

HMM需要解决3个问题：

（1）评估问题：给定五元组模型λ=（N， M，π，A， B）和观察序列O=（O1， O2，…，OT），计算这个观察序列出现的概率.

（2）解码问题：给定五元组模型λ=（N， M，π，A， B）和观察序列O=（O1， O2，…，OT），求可能性最大的隐藏状态序列.

（3）学习问题：也叫训练问题.给定观察值序列O=（O1， O2，…，OT），据此确定HMM模型λ=（N， M，π，A， B）.即如何调整λ=（N， M，π，A， B），使得P（ O| λ）最大.

上述的3个问题分别对应3个经典的解决算法［11-12］，分别是前向-后向算法、Viterbi算法、Baum-Welch算法.

2 食品安全风险预测方法

食品供应链主要包括原料生产、生产加工、储藏运输和消费4个环节，每个环节都会受很多因素的影响.假设影响因素为

其中XN（1≤N≤4）为风险预测的观测值，XNi为第N个环节的第i个影响因素，通过一定的方法反映出风险预测值和其影响因素之间的关系，数学模型表达式为

HMM模型有2个随机过程：可见的状态序列，本文指供应链各个环节检测到的数据；系统的真实状态，本文指供应链各环节的真实风险值.供应链每个环节的HMM三元组参数模型包括状态转移概率矩阵A，观测概率矩阵B和最初状态分布π，字母表示为λ=（π， A， B）.除了第一个环节外，其他3个环节的参数模型都受上一环节的影响，本环节最终的风险概率直接影响下一环节的初始状态概率分布.例如，原料生产环节通过HMM模型计算，最终得出各个安全状态出现的概率分布为δ=｛δ1， δ2，…，δi｝（i为真实状态的个数），那么下一环节，即生产加工环节的初始状态分布概率π=δ.

设有A1—A5共5个状态，表示供应链风险等级，其中A1表示正常安全状态，供应链受威胁的概率为P=0；A2、A3、A4表示危险等级逐级加深，对应供应链受威胁的概率分别为0＜P≤0.2（低风险）、0.2＜P≤0.5（中级风险）、0.5＜P≤0.8（中高级风险）；A5表示重大危险状况，供应链受威胁的概率为0.8＜P≤1.

这里的HMM为离散HMM，在建立离散的HMM模型系统时，可以假定：系统状态包含系统的所有信息，并且在当前状态下的观察是独立的.则HMM模型的状态转移模型可以表示为图1.从一个节点移动到另一个节点，表示系统在源节点显示的状态，会转变成目的节点的状态，该模型图是一个完全连接图，表明任意安全状态都有转变为其他任意安全状态的可能.

图1 HMM状态转移模型Fig. 1 Transfer model of HMM state

据以上分析，提出图2所示的风险预测流程.

图2 风险预测流程Fig. 2 Risk assessment process

在用HMM模型进行风险预测时，首先对数据进行预处理，采用的归一化映射为

其中x， y∈Rn.归一化后的原始数据在［0，1］范围内，即yi∈［0，1］，i=1，2，…，n .然后，将归一化的数据通过DBSCAN聚类去除噪音点，得到初始数据.DBSCAN算法不需要事先知道要形成的簇类的数量，可以发现任意形状的簇类，并且能够识别出噪声点，适用于类别数未知的聚类问题［13］.把聚类之后的数据通过Baum-Welch算法训练得到对应的HMM三元组模型后，即可通过Viterbi算法对供应链进行相应的风险预测.

3 实验分析

在食品安全方面，乳品一直受高度关注.所以本文以乳品为例进行分析.乳品在从农场到餐桌的过程经历了4个环节：原奶生产、乳品加工、储藏运输和消费.HACCP［14］（hazard analysis critical control point）表示危害分析和关键控制点，它是对食品加工过程中可能发生的食源性危害进行识别和评估，进而采取一定控制的食品安全控制体系.利用HACCP确定关键控制点的判断流程如图3所示.

图3 确定关键控制点的判断树Fig. 3 Critical control point decision tree

根据HACCP体系并结合乳品供应链的4个环节，对各个环节中的生物性、物理性和化学性的潜在危害进行分析，确定了关键控制点.以第一环节为例，原奶生产环节的关键点选择见表1.

由于乳品一般都难以长久保持并容易受到微生物污染，供应链上任何环节的潜在危险因素都会影响食品质量安全.乳品加工环节是整个供应链最重要的部分，负责将收集到的原奶加工完成，储运和消费环节的风险因素也基本来源于供应链的前几个阶段.因此，乳品生产供应链的各个环节紧密联系，对每一个环节的检查监管都要格外重视.最后一个环节的风险值即为整条供应链的风险大小.

表1 乳品原奶生产环节关键控制点分析Tab. 1 Critical control point analysis of raw material production

通过对天津某乳品企业的调研，采集大量数据，并建立了乳品追溯系统，利用HMM模型对乳品安全风险进行预测.

乳品供应链各环节的隐藏状态有5个：A1—A5，分别代表5个安全等级.对归一化的数据通过DBSCAN聚类进行聚类处理后，共得到4类，也就是系统的4个观察状态：V1—V4；利用HMM模型依次对4个环节进行参数评估和风险概率计算，每个环节最终的风险概率直接影响下一环节的初始状态概率分布.

在训练HMM模型时，利用Baum-Welch算法，通过让概率P（O|λ）达到局部最大值得到HMM的参数模型.

定义供应链所有特征观察值序列前向变量［11］为

由式（8）可知，供应链所有特征观察值序列概率为

通过式（9）求取HMM的三元组参数模型λ=（π， A， B），再由重估公式［15］可得

输出概率P（O|λ）会随着重估次数的增加越来越大，直至参数πi、aij、bij收敛或算法达到停止条件为止.

计算得到第一环节的初始状态分布概率为π=｛0.4612，0.1132，0.2139，0.0547，0.1570｝，状态转移矩阵为

观测向量概率矩阵为

有了HMM的三元组模型，就可以对相关数据进行解码分析.在此给定一组观察序列｛O2， O1， O3， O4，O4， O3， O2， O2， O1， O1， O4， O3｝，观察序列的状态分布如图4所示.

图4 观察状态序列Fig. 4 Sequence of observed state

继而利用Viterbi算法求解观察序列的最优路径，得到供应链最有可能的真实状态.首先将πi与所有特征观察值转移概率bi相乘，得到初始化路径δt（ i）.

将初始化的路径δt（ i）与状态转移概率矩阵的元素aij相乘，取最大值与所有特征观察值转移概率bj相乘，得到当前路径最大值δt（j）.

根据式（14）、（15）结果，可求出对应的食品安全各环节真实状态的最大概率P∗和最大概率对应状态序列qt∗分别为

因此，可以评估出该环节对应的真实状态｛T1， T2， T1， T2， T4， T4， T1， T2， T4， T2， T5， T2｝.

真实状态序列的分布如图5所示.

图5 真实状态序列Fig. 5 Sequence of true state

由结果可知，系统在t=5，6，9，11时的风险级别较大，在t=1，3，7时刻，风险级别较小.通过Viterbi算法可以计算出原奶生产环节在最后时刻处在各状态的概率为δ=｛0.1672，0.5342，0.1745，0.0723，0.0518｝，则在乳品生产环节中，初始状态分布概率就为π=｛0.1672，0.5342，0.1745，0.0723，0.0518｝.

参照对第一环节的计算过程，依次利用Baum-Welch算法对各环节进行参数计算，然后利用Viterbi算法计算系统的真实状态.

最后，计算得到最终环节——消费环节各安全等级的概率为π=｛0.1547，0.6073，0.0645，0.1043，0.0692｝.

假设在A1-A5不同安全等级的情况下，系统对应的开销，即风险对系统的影响分值为C=｛0，5，10，15，20｝.参考风险评价模型的典型代表OCTAVE［10］给出一个计算风险值的公式

其中：Rt表示在t时刻系统所处的总体风险值；αt（ i）为t时刻系统处在安全状态Ai的概率；N是安全状态的数目；c（ i）是与状态Ai关联的开销.

根据式（19），可以计算出乳品的最终风险值R= 6.049，5.

因此，可以利用上述方法对乳品各环节的不同时刻进行风险预测与开销分析，有助于生产者和决策者更全面了解乳品供应链各环节的风险态势，及时查出问题发生的缘由，根据已有流程或者规定采取相应措施，并不断更新和完善风险规避的方法，最大程度减少乳品安全对群众的影响.

4 结 语

本文利用HMM模型对乳品质量安全风险进行分析与预测.该方法根据食品安全风险的形成机理，考虑了4个环节之间的风险传递关系以及各个环节风险的动态性，结合定性与定量评价，可以更精确地描述乳品风险的大小.该风险预测模型除了可用于对乳品进行评估外，还可以扩展应用到食品安全的其他领域.今后还需通过更多实验对各环节风险传递的影响因子进行优化，改进HMM模型的参数训练方法，以提高参数模型的训练速度和评估准确度.

［1］ 黄秋婷，丁怡，邱佩丽. 基于灰色残差修正模型的食品安全监测预警分析方法研究［J］. 现代农业科技，2013（5）：304-305.

［2］ Soon J M，Davies W P，Chadd S A，et al. A Delphi-based approach to developing and validating a farm food safety risk assessment tool by experts［J］. Expert Systems with Applications，2012，39（9）：8325-8336.

［3］ 珺安. 基于层次分析法的乳品质量安全预警系统研究［D］. 哈尔滨：东北农业大学，2012.

［4］ 孔大维. 基于供应链管理的食品可追溯系统的构建研究［D］. 成都：成都理工大学，2013.

［5］ 张丽，滕飞，王鹏. 基于贝叶斯网络的食品供应链风险评价研究［J］. 食品研究与开发，2014（18）：179-182.［6］ Xiao J，Ma Z S，Zhang D J. The research on the BP neural network in the food security risk early warning applied under supply chain environment［C］// Proceedings of Second International Symposium on Knowledge Acquisition and Modeling. Piscataway：IEEE，2009：281-284.

［7］ 曹秋勤. 基于支持向量机的蔬菜质量安全预测及溯源模型的研究与应用［D］. 广州：华南理工大学，2014.

［8］ 王志堂，蔡淋波. HMM模型（HMM）及其应用［J］. 湖南科技学院学报，2009，30（4）：42-44.

［9］ 张建华. 基于深度学习的语音识别应用研究［D］. 北京：北京邮电大学，2015.

［10］ 龙百元. 基于层次分析法和HMM模型的风险预测方法研究［D］. 长沙：湖南大学，2007.

［11］ Ieee L R R F. A tutorial on hidden Markov models and selected applications in speech recognition［J］. Proceedings of the IEEE，1989，77（2）：257-286.

［12］ 谢锦辉. 隐Markov模型（HMM）及其在语音处理中的应用［M］. 武汉：华中理工大学出版社，1995.

［13］ 慕升弟. 基于协议HMM模型的网络流量识别技术研究［D］. 重庆：重庆大学，2013.

［14］ CAC. Draft Principle and application of the hazard analysis and critical control point system：Report of twenty-fifth session of the codex committee on food hygiene［R］. Washington DC：CAC，1991.

［15］ Schulz L L，Tonsor G T. Cow-calf producer perceptions regarding individual animal traceability［J］. Journal of Agricultural and Applied Economics，2010，42（4）：659-677.

责任编辑：常涛

Food Safety Risk Prediction Based on HMM

MA Yongjun，LIU Yushan，HOU Yangyang，LIU Yang
（College of Computer Science and Information Engineering，Tianjin University of Science & Technology，Tianjin 300222，China）

Taking dairy products as the research subjects，a food safety risk assessment method based on the Hidden Markov Model（HMM） was developed.The critical control point for each sector of the supply chain was found out with the method of hazard analysis critical control point.Those critical points can be used as quantitative indicators of HMM.According to the characteristics HMM model，the supply chain was divided into four sectors.The ultimate probability of each sector is the initial state probability distribution of the next part.Finally the level of risk in the supply chain and risk level can be calculated with HMM，which can provide a reference for supply chain risk assessment.

food safety；HMM；dairy products；supply chain；risk prediction

TP311

A

1672-6510（2016）05-0069-05

10.13364/j.issn.1672-6510.20150209

2015-11-16；

2016-04-18

教育部规划项目（12YJAZH091）

马永军（1970-），男，吉林长春人，教授，yjma@tust.edu.cn.