2013年“深圳杯”数学建模A题论文

基于多目标规划的食品安全抽检模型

摘 要

本文通过对深圳的食品质量安全的影响因素以及食品抽检数据的分析，建立了基于多目标规划的食品安全抽检模型，优化了食品抽检的方法。

针对问题一，首先将所给的数据进行统计，再用Matlab软件描绘出其三年中各影响因素的变化趋势，最后用模糊集对分析方法与密切值法评价其变化趋势，可知微生物超标和重金属引起的食品不合格的比重有变小的趋势，说明微生物和重金属对食品安全的影响会越来越小。食品添加剂在食品不合格的因素中所占的有变大的趋势，说明食品添加剂会越来越多的影响食品质量。

针对问题二，仅考虑食品质量与季节之间的变化规律，建立起初等数学模型，并用Matlab软件进行数据拟合，得出食品质量与抽检期数之间的关系

S??0.0010x4?0.0154x3?0.1089x2?0.x3100?0.。 7199针对问题三，基于多目标规划建立了一个以检测误差、检测时间、检测成本

i?n,j?p?qi,j?UQi,j???minf1?i?1,j?1?n?p?i?n,j?p?为目标函数的食品抽检模型?minf2??tij，从而得出最佳解决办法，

i?1,j?1?i?n,j?n??minf3??eij?i?1,j?1??既节约成本费用又有较高的检测效果。

关键词：模糊集对分析；多目标规划 ；食品抽检

1

一、问题重述

食品质量安全状况是一个国家经济发展水平和人民生活质量的重要标志。我国改革开放三十年来，人民生活水平在不断地提高，食品安全和卫生问题越来越受到人们的关注。近年来，先后出现苏丹红、瘦肉精、三聚氰胺等事件，以及各种不利于健康的食品添加剂、强化剂问题的出现，食品安全和卫生的检测已成为全社会重点关注的问题之一。食品的质量和卫生问题涉及到原材料的使用、生产加工、运输与贮存、流通与销售等环节，在每一个环节上出现差错，都会导致食品出现安全和卫生问题，食品质量和卫生的检测工作在实际显得非常重要。但是，由于食品的种类、品牌和批次繁多，从生产加工到销售食用中间环节复杂，质检部门不可能对所有食品做到全面的质量检测，一般只能做一定的抽检。对食品进行抽检也需要一定人力、物力等成本费用，抽检的越多检测效果就越好，但需要的时间就越长，其成本费用也就越高。为此，合理抽检，既能保证较好的检测效果，又能节省时间和成本费用。

针对问题一，只考虑微生物、重金属、食品添加剂对食品质量的影响。从2010年到2012年，根据所给出的数据，在不断变化的抽检结果中找出微生物、重金属、食品添加剂的变化趋势，用模糊集对分析方法加以评价。

针对问题二，考虑到食品质量可能还与生产地点、抽检地点和季节等因素有关，我们从每年抽检的数据中总结出相应的规律，并用Matlab软件进行数据拟合。

针对问题三，食品抽检需要一定人力、物力等成本费用，抽检的越多检测效果就越好，但需要的时间就越长，其成本费用也就越高。故本文基于多目标规划建立了一个以检测误差、检测时间、检测成本为目标函数的食品抽检模型，从而得出最佳解决办法，既节约成本费用又有较高的检测效果。

二、模型假设

1、假设主要食品分为蔬菜、鱼类、鸡鸭，主要影响因素为微生物、重金属、添加剂，其他因素造成的影响忽略不计；

2、不考虑食品因保质期的问题被抽检出不合格；

2

3、考虑季节因素对食品质量的影响时，控制抽检地点等因素不变。

三、符号说明

本文中所用的符号和说明见表1 表1：符号说明

S 测点集(期数) 第k个测点 指标集 第r个指标 决策矩阵 为测点sk关于指标er的属性值 sk E er D dkr ur、vr 分别为指标er的最优值和最劣值 Ai 第i种食品 wi 第i种食品的安全指标 四、模型的建立与求解

4.1 问题一： 4.1.1 问题分析：

这是个评价总结问题，要找出食品领域各因素的变化趋势，需要解决两个问题。一般来说食品安全的影响因素众多，而且随时间的波动较大，而且与抽检的地点和种类有较大关系，表现得杂乱无章。所以首先要对数据进行整理、分析，

3

用合理的方法归纳出各因素的变化趋势，然后我们需要就所得出的结论进行验证，以确定其科学性，所以需要运用模糊集对分析法建立目标函数，确定约束条件来评价该变化趋势是否可信。对于本问题而言，目标函数是变化趋势的同一隶属度、相对隶属度和相对贴近度，用模糊集对法和密切值法排序，来验证总结出的变化趋势的准确性及科学性。

4.1.2 微生物、重金属、添加剂含量的变化趋势

根据深圳市市场监督管理局所提供的数据[1]，计算出了2010年到2012年各个月份的抽检中分别由微生物、重金属和食品添加剂影响食品不合格所占的比例，所列表格如下（考虑到论文篇幅的限制，只列出微生物超标占不合格样品的比重）：

表2：微生物超标占不合格样品的比重 期数 比重（%） 期数 比重（%） 期数 比重（%） 期数 比重（%） 2 23 10 9 21 31 28 87.5 4 57 11 82 22 100 29 33.3 5 0 12 73 23 58 30 66.7 6 11 14 75 24 27 32 27.2 7 2 15 67 25 29 34 0 8 48 16 27 26 4.2 35 100 9 5 18 29 27 83.3 36 25.7 注：从2010年到2012年记为1期-36期，部分期数数据缺省。 根据各因素所占不合格样品的比例，用Matlab软件画出如下曲线：

4

微生物超标占不合格样品比重1重金属超标占不合格样品比重10.50.500102030400010203040添加剂超标占不合格样品比重10.50010203040

图1：微生物、重金属、添加剂超标分别占不合格样品比重

从以上图像可以看出，抽检结果很不稳定，随时间波动比较剧烈，但是呈现一定的趋势。微生物超标和重金属引起的食品不合格的比重有变小的趋势，说明微生物和重金属对食品安全的影响会与来越小；食品添加剂在食品不合格的因素中所占的有变大的趋势，说明食品添加剂会越来越多的影响食品质量，需要我们加大监督力量，严格控制食品添加剂的使用。

通过查阅一些资料，初步预测这三种影响食品质量的因素的变化趋势，得出了一些结论：由于食品加工和保存的时候，消毒杀菌的工作做得比较好，使微生物对食品的危害有所减小，呈现比较乐观的态势。近年来一些排放不达标的企业被勒令改进，使得重金属对食品的危害也有所减小。但是人们也越来越注重生产效率，越来越多的使用食品添加剂，对食品安全有较大的威胁，仍需我们努力改善。

4.1.3 用模糊集对的方法评价其变化趋势[2]：

设多属性决策问题Q??S,E,D?, 其中:

S??sk??k?1,2?,m?为测点集(期数) , sk为第k个测点;

5

E??er? ?r?1,2,?,n?为指标集,E有不同的类型指标, er为第r个指标; D??dkr?为决策矩阵, dkr为测点sk关于指标er的属性值。

记最优测点集 U??u1,u2,?,ur?, 最劣测点集 ??v1,v2,?,vr?, 其中: ur、vr 分别为指标er的最优值和最劣值。对于er?E1，比较区间为?ur,vr?，定义集对

?dkr,ur?同一度akr和对立隶属度ckr如下：

akr?dkr/?ur?vr? (1)

akr?urvr/??ur?vr?dkr? (2)

计算同一隶属度、对立隶属度：

ak??akr (3)

ck??ckr (4)

评价结果排序：由于ak、ck是相对确定的，分别表示对sk接近最优测点U的肯定和否定程度，那么在相对确定条件下可定义sk和U的相对贴近度为：

akr?ak/?ak?ck? (5)

表3：2010~2012三年各因素统计结果

2010年 2011年 2012年 微生物 38.2% 21.5% 46.2% 重金属 49.2% 32.2% 29.2% 食品添加剂 39.8% 14.9% 46.1% 依据各评价指标的性质，确定U=（21.5%，29.2%，14.9%），V=（46.2%，49.2%，46.1%），应用以上各式计算食品安全各影响因素的贴近度。

6

表4：食品质量各影响因素的平均同一隶属度、对立隶属度、相对贴近度及排序 年度 2010年 2011年 2012年 ak 0.3096 0.3412 0.6416 ck 0.5979 0.4000 0.2952 rk 0.3412 0.4603 0.6849 模糊集对法排序 3 2 1 密切值法排序 3 1 2 上表的排序结果表明我们的预测与实际情况比较符合，且模糊集对评价结果与密切值法评价结果出入不大，说明模糊集对分析模型具有较强的准确性和实用性，可用于食品安全影响因素变化趋势的综合评价。 4.2 问题二： 4.2.1 问题分析：

影响食品安全问题的因素众多，除了考虑食品加工时不同的食品产地微生物、重金属残留物与添加剂的掺入，在抽样调查时还要考虑到的抽样地点即食品的销售地点和季节等诸多因素[3]。本文只考虑季节因素，季节的变化对应着温度的变化，比如夏天，应该针对消费集中、保质期短、易发霉变质、时令性强的食品实施有效监管，重点检查食用油、肉制品、糕点等夏季高风险食品。通过三年中抽检的数据用matlab软件进行初步数据拟合，找出食品质量与季节变换之间的关系，并用于指导食品抽检。 4.2.2 数据处理：

通过对2012年26期抽检数据的分析，运用Excel软件整理出各次抽检的合格率，如下表：

表5：抽检合格率与期数的统计表 期数 合格率 期数 合格率 期数 合格率 1 0.964 9 1.000 17 0.984 2 1.000 10 0.722 18 0.875 3 0.914 11 0.800 19 0.955 4 0.993 12 0.994 20 0.760 5 0.938 13 0.964 21 0.980 6 0.988 14 1.000 22 0.980 7 0.940 15 0.875 23 0.824 8 0.578 16 0.885 24 0.768 4.2.3 模型建立与求解：

7

考虑到每一年的季节的重复，本文以2012年的数据作为基准建立模型，并用另外两年加以检验。

一般情况下食品的有效性（合格的产品）S与季节的某个特征数值c (比如温度)的成正比，S?k1c (其中k1为常数且大于0)，而每个季节的不同固然对应着不同的温度，我们假设温度随着时间x成二次曲线关系变化c?k2x2?k3x?k4,则通过比例关系推导出数学模型S?k1?k2x2?k3x?k4?，又称经验模型。 要利用上述模型计算食品质量的不合格率，必须确定其中的参数k1、k2、k3、k4。一种确定参数的办法是进行抽检得到数据或查阅先关资料，另一种办法就是把上式写成S?ax2?bx?c (其中a?k1k2、b?k1k3、c?k1k4)，只需要确定a、b、c就可以进行S与x之间的计算。显然，后一种方法比较简便。利用上述表格中所统计的食品不合格率与时间的关系，运用最小二乘法进行数据拟合，求出模型的各个参数。

用Matlab软件算出a?0.0005、b??0.0130、c?0.9773， 则食品质量与季节的变化规律的关系模型为：

S?0.0005x2?0.0130x?0.9773 (6)

图形如下：

10.950.90.85合格率0.80.750.70.650.60.550510期数152025

图2：2012年食品合格率散点图与拟合图线

8

将实际期数带入模型中去，与实际数据作比较。从上面拟合的曲线可以看出每年的夏天食品合格率与曲线出入较大。原因是该模型是通过比例关系得到的，假设条件比较粗糙，考虑比较不全面，是比较粗略的模型。

由于是比较粗略的模型，计算的结果与数据出入较大。观察

S?0.0005x2?0.0130x?0.9773，该是有一个多项式主要涉及到参数和参数的系

数。为了对该式进行改进，我们把此式改为

S?AXm?BXn?CXr? (7)

可再利用最小二乘法进行数据拟合，求出各个参数，利用MATLAB软件进行数据拟合，得到改进后的模型：

S??0.0010x4?0.0154x3?0.1089x2?0.3100x?0.7199 (8)

图形如下：

1.0510.950.90.85合格率0.80.750.70.650.60.550510期数152025

图3:2012年食品不合格率的散点图与拟合图线

该式就是改进后的模型，从数据拟合的过程中我们可以看出该模型与实际数据有较高的符合程度，但还是有几个点拟合的效果不够好。

从中我们发现的规律是，夏季前后由于气温高的原因，消费集中、保质期短、易发霉变质、时令性强的食品不合格率比较高，应重点检查食用油、肉制品、糕点等夏季高风险食品。冬季和春季可能由于温度比较低和迫于春节压力的原因，

9

对影响食品质量的因素较少，但由于该季节期间微生物活跃度高，应该重点抽查微生物发酵的和各种微生物制品。 4.3 问题三：抽检模型及可靠性分析 4.3.1 抽检模型的建立[4]:

设共有n种食品，以字母Ai代表第i种食品。并且设每种食品都有p个安全指标，且设wi为第i种食品的安全指标，同时我们设Ai种食品现在有di个品牌，并且每个品牌批次数为Ni,j，其中i代表食品的种类，而j代表第i种食品的第j个品牌。

要建立一个数学模型来对进行合理的抽检，那么在n种食品的di个品牌里都得选取一定的样本，那么设这个样本的大小为ui,j 其中i表示的是食品的种类，而j表示的是该种类u的品牌。那么我们的模型即为求u矩阵：

?u11u12?uu22U??21?MM??un1un2u1d?Lu2d?? Ou11??Lund?L设以上矩阵为样本矩阵，那么求样本矩阵即为此模型的求解。设AU为食品的各种类的各品牌的批次矩阵，那么，U必须满足的条件为：

U?AU

设第i类食品的第j品牌的真实合格率为，另外在u中进行检测会出现一个检测的合格率矩阵，设其为UQ，那么定义检测误差矩阵为：

??|UQ?Q|/Q (9)

那么上述值的平均值，即定义为本次抽检的检测误差。 设检测成本为：

?e12?eE??12?M??e12e12e12Me12e12?Le12?? OM??Le12?L 10

而检测所需时间为：

?t11?t21T???M???tn1t12t22Mtn2t1P?Lt2P?? OM??Ltnp??L通过分析，只要求出抽检的目标函数，对于本问题是一个多目标的规划的问题，其中规划的目标函数有三个，一个是检测误差设其为f1，一个是检测时间设其为f2，最后一个是检测所花费的的费用，设其为f3。因此建立的抽检方案模型：

i?n,j?p?qi,j?UQi,j???minf1?i?1,j?1?n?p?i?n,j?p? ?minf2??tij (10)

i?1,j?1?i?n,j?n??minf3??eij?i?1,j?1??4.3.2 多目标规划模型的模拟检验

对于建立的多目标规划模型，由于没有相应的数据进行分析求解，我们根据所建立的抽检模型进行计算机模拟，以检验模型的正确性。

在抽检过程中，为了简化过程，我们以每种食品少数品牌进行计算机模拟，同时，我们设定某一种食品有10个品牌，且每个品牌有100 个批次，在模拟开始之前，随机生一个实际合格率，在一次模拟过程中，实际合格率选取为： 表6：随机产生的五个品牌实际合格率 品牌 合格率

一 0.8035

二 0.9478

三 0.9012

四 0.8107

五 0.9012

这样在完成了真合格率设定以后，就对每个品牌的100 个批次进行相应的随机处理，然后利用上述抽检模型对每个批次进行抽检，这里我们选取抽检数为17。

抽检完成以后就能得出出抽检合格率，本次随机模拟的抽检合格率为：

11

表7：第17次抽检五个品牌的合格率 品牌 合格率

一 0.8824

二 0.9412

三 0.7647

四 0.9412

五 0.8824

根据以上模拟结果进行数据分析可计算得到两个出检测误差，经计算，得出检测误差为 表8：检测误差表 品牌 误差

一 0.0982

二 0.0070

三 0.1515

四 0.16090

五 0.0209

通过对检测误差算平均值，得出结果为8.77%，即建立的多目标规划抽检模型该模型的平均检测误差为8.77%。

五、模型的评价及改进

5.1 模型的优点：

（1）模型中有较多新的想法和新的思维，引用了许多经典的函数，模糊集对分析方法很好的评价了各种因素的变化趋势；

（2）模型的应用综合考虑了各种因素，比较贴近实际，在对政府抽检食品有一定的指导意义；

（3）进行了比较严格的拟合，有一定的实用价值； （4）模型带有图像，一目了然，思路清晰。 5.2 模型的缺点：

（1）虽然数据比较齐全，但部分月份数据缺省，无法让拟合的数据曲线更为准确；

（2）在问题一和问题二讨论过程中只考虑部分因素而忽略影响模型其它因素，可能对结果有一定的影响； 5.3 模型的优化与推广： 5.3.1 优化：

在考虑食品安全的影响因素时，可以再增加一些因素，也可以做更多的假设，这样可以更全面的分析食品抽检种可能遇到的问题，以便为食品监管提供依据。在模型改进方面，我们可以引入抽检地点和抽检类别等因素，并对它们进行回归

12

分析，逐步筛选出对食品安全影响较严重的因素，并以此为基础进行抽检方法的优化。 5.3.2 推广：

此模型可广泛应用于实际各种“抽样问题”，从而解决成本费用与检测可靠性之问的折中。

参考文献

[1] 深圳市市场监督管理局网站www.szaic.gov.cn

[2]姜启源，谢金星，叶俊，数学模型，北京：高等教育出版社，2004.4。 [3]关于影响食品安全因素的探讨.陈锦屏，张志国.食品科学，2005年26卷第8期

[4]郑龙，黄继达，食品安全的抽检问题研究，中国矿业大学矿业工程学院，江苏徐州 221116

附录

问题一：

x=[2 4 5 6 7 8 9 10 11 12 14 15 16 17 18 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36];

y1=[0.23 0.57 0 0.11 0.2 0.48 0.5 0.9 0.82 0.73 0.75 0.67 0.27 0.29 0.31 1 0.58 0.27 0.29 0.042 0.5 0.833 0.875 0.333 0 0.667 0 0.272 0 1 0.257];

y2=[0.076 0.06 0.5 0.56 0.14 0.54 0 0 0.16 0.11 0.5 0.13 0.47 0.14 0.56 0 0.17 0.45 0.46 0.25 0 0 0.16 0.067 0 0 1 0 0 0 0.314];

y3=[0.69 0.71 0.5 0.33 0.8 0.23 0.75 0.17 0.24 0.2 0 0.59 0.031 0.86 0.13 0 0.33 0.36 0.33 0.708 0.5 0.161 0.075 0.6 1 0.333 0 0.728 1 0 0.429]; subplot(2,2,1); plot(x,y1,'b*-');

title('微生物超标占不合格样品比重'); subplot(2,2,2); plot(x,y2,'r*-');

13

title('重金属超标占不合格样品比重'); subplot(2,2,3); plot(x,y3,'g*-');

title('添加剂超标占不合格样品比重') 问题二

x=[1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 ]

y=[0.964 1.000 0.914 0.993 0.938 0.988 0.940 0.578 1.000 0.722 0.800 0.994 0.964 1.000 0.875 0.885 0.984 0.875 0.955 0.760 0.980 0.980 0.824 ]

a=polyfit(x,y,6) z=polyval(a,x);

plot(x,y,'k+',x,z,'r') %做出数据点与拟合曲线 xlabel('期数'); ylabel('合格率');

14

本文来源：https://www.bwwdw.com/article/unp6.html