MERS传播的数学模型的建立及分析.doc

1、2015 高教社杯全国大学生数学建模竞赛 承 诺 书 我们仔细阅读了全国大学生数学建模竞赛章程和 全国大学生数学建模竞赛 参赛规则（以下简称为“竞赛章 程和参赛规则”，可从全国大学生数学建模竞赛网站下载）。 我们完全明白，在竞赛开始后参赛队员不能以任何方式（包括电话、电子邮件、网上咨询等）与队外的任何人（包括指导教师）研究、讨论与赛题有关的问题。 我们知道，抄袭别人的成果是违反竞赛章程和参赛规则的，如果引用别人的成果或其他公开的资料（包括网上查到的资料），必须按照规定的参考文献的表述方式在正文引用处和参考文献中明确列出。 我们郑重承诺，严格遵守竞赛章程和参赛规则，以保证竞赛的公正、公平性。如有

2、违反竞赛章程和参赛规则的行为，我们将受到严肃处理。 我们授权全国 大学生数学建模竞赛组委会，可将我们的论文以任何形式进行公开展示（包括进行网上公示，在书籍、期刊和其他媒体进行正式或非正式发表等）。 我们参赛选择的题号是（从 A/B/C/D 中选择一项填写）： A 我们的参赛报名号为（如果赛区设置报名号的话）： 所属学校（请填写完整的全名）： 成都工业学院 参赛队员 (打印并签名 ) ： 1. 王 2. 卢 3. 唐 指导教师 或 指导教师组负责人 (打印并签名 )： （论文纸质 版与电子版中的以上信息必须一致，只是电子版中无需签名。以上内容请仔细核对，提交后将不再允许做任何修改。如填写错误，论

3、文可能被取消评奖资格。） 日期： 2015 年 7 月 27 日 赛区评阅编号（由赛区组委会评阅前进行编号）： 2015 高教社杯全国大学生数学建模竞赛 编 号 专 用 页 赛区评阅编号（由赛区组委会评阅前进行编号）： 赛区评阅记录（可供赛区评阅时使用）： 评 阅 人 评 分 备 注 全国统一编号（由赛区组委会送交全国前编号）： 全国评阅编号（由全国组委会评阅前进行编号）： MERS 传播的数学模型的建立与分析 摘要 本文针对 MERS 的传播建立了 传统的 SIR 仓室 数学模型。 针对问题一 ，对附件 一 提供的早期模型，认为“传染概率”的说法欠妥，传染期限 L的确定缺乏医学上的支持，使模

4、型的说服力降低 。模型中借鉴广东香港的参数来预测北京的疫情走势，不失为一种方法。但在不同国家 因政策，地域的不同，病毒 的传播和控制呈现不同的特点，使不同 国家不同 城市之间的可 比性降低。而且由于 MERS 和 SARS 的治病机理不同，传染率和死亡率不同，患病人数规模不同，所以对 MERS 的传播分析，不能简单的套用附件一所用的模型。 针对问题二 ， 我们在 WHO(World Health Orgnazation)的官方网站上查找到了韩国 MERS 疫情从 2015 年 5 月 20 日至 2015 年 7月 5日的详细数据 1。 对 MERS的传播建立传统的 SIR 仓室模型，采用最小

5、二乘法拟合参数 ，利用 MATLAB 编程求解， 画出参数感染率 的参数散点图和参数移出率 的散点图 对第三个问题，本文研究对 MERS 疫情对韩国 入境旅游收入 的影响，建立了灰色预测 GM（ 1， 1） 模型。 关键字： SIR 仓室模型 常微分方程参数拟合 灰色预测 1.问题重述 MERS（ Middle East Respiratory Syndrome）病毒是一种新型的 冠状病毒 ，这种病毒已经被命名为中东呼吸综合征冠状病毒，大多数 MERS 病毒感染病例发生在沙特。 2015 年， MERS 在韩国又有新一轮的爆发和蔓延，给韩国的经济发展和人民生活带来了较大影响。对 MERS 的传

6、播建立数学模型，具体要求如下： （ 1）对附件所提供的一个 SARS 传播的早期模型，请对其评价是否适用 MERS。 （ 2）收集 MERS 的韩国疫情数据，建立 MERS 传播的数学模型，说明优于附件1中的模型的原因；特别要说明怎 样建立一个真正能够预测以及能为预防和控制提供可靠、足够的信息的模型。对于卫生部门所采取的措施做出评论，如：提前或延后 n天采取严格的隔离措施，对疫情传播所造成的影响做出估计。 （ 3）收集 MERS对韩国旅游方面影响的数据，建立相应的数学模型并进行预测。 2.模型假设 1.假设一个 MERS 康复者不会二度感染，他们已退出传染体系，因此将其归为“退出者”； 2.模

7、型不考虑所研究这段时间内的自然出生率和死亡率， MERS 引起的死亡人数归为“退出者”； 3.假设在疾病传播期内所考察地区总人数视为常数； 4.假设每个病人单位 时间有效接触的人数为常数； 5. 假设韩国 在 MERS 疫情流行期间和结束之后， 旅游业 数据的变化只与 MERS 疫情的影响有关，不考虑其它随机因素的影响。 3.变量说明 S ： 表示易感染人群 (susceptible)占总人数的比例； I ； 表示感染人群 (infected)占总人数的比例； R ： 表示移出人群占总人数的比例； ：表示感染者对易感染者有效感染的感染率； ：表示移出率，即移出者的增加率 01xn：表示 第 n

8、 个单位时间的旅游业的收益； 11X n 表示在 n 个单位时间时的累积旅游业的收益； ： 称为系统发展灰数； ： 称为内生控制变量 4.对早期模型的评价 附件 1 的模型主要采用“数据拟合”和“借鉴参数”的方法对北京疫情走势进行预测。 在数据拟合方面，该模型中有两个疑点： 1、感染期限 L 的确定。由于被严格隔离、治愈、死亡等原因，感染者在某一时段后不再具有对易感人群的传染力，故对病毒的传染加上感染期限是合理的。但在对该参数的确定上，作者为了较好 地拟合各阶段的数据 ,通过人为调试来确定 L的取值，缺乏医学上的支持，使模型的说服力减弱，合理性和可靠性大大降低。 2、文中认为“ K代表某种环境

9、下一个人传染他人的平均概率”。但从模型的公式中可以看出，参数 K的实际意义是一个病人平均每天传染其他人的个数。两者之间有实质的区别，文中的说法显然不妥。 从预测思想来看，该模型是借鉴先发地区 广东、香港的有关参数对北京的疫情进行预测的。由于广东、香港的疫情和控制都在北京之前，已经过了高峰期，到 5月 8日为止每日新增病例已降至 10 来例，基本处于后期控制阶段。而当时北京的疫 情刚过了高峰期，正处于社会剧烈调整时期，数据较为凌乱，略有下降趋势，但不明显。可见在当时，采取这种借鉴是无奈之举。 但是由于城市之间的政策，风俗习惯等不同，城市之间的可比性不强，借鉴存在很大的局限性。如在香港，由于对传播

10、机制认识不足，中途又出现高度感染的特殊情况。另外使用借鉴法无法对首发城市进行预测。 MERSS 和 SARS 又有许多不同 , 例如， MERS 的病死率约为 40.7%，传染性没有 SARS 强，但 SARS 的病死率为 14%-15%，低于 MERS，传染性则强于 MERS。 MERS和 SARS 治病机理不同，传染率死 亡率不同，病毒潜伏时间也不同，所以不能简单的套用附件所给模型来分析和预测 MERS 的传播。 5.模型的建立与求解 5 1 问题 2 的模型建立与求解 5.1.1 问题 2 的模型思路 对于传染病感染区，由于为了避免传染病的更大的扩散，一方面政府会对人口的流动做出限制，另

11、一方面个人由于对传染病的警惕也不会进入感染区，所以，感染区人口流动很小。即可把它看作一个封闭区，则可以建立 SIR 仓室模型 2,里面的总人数不变，里面的人分类为： 易感染者 :即正常人，但可能会被感染。 感染者：已经感染这种病的人，可以传染给周围的人。 移出者 ：包括感染者中死亡的人，感染者中自愈的人和先天对这种病毒有免疫的人，他们将不在受这种传染病的影响。 5.1.2 问题 2 的模型建立 根据传染病的感染而致的各类人口比例的变化可以建立微分方程。这里用St表示易感染人群 (susceptible)的比例关于时间 t 的函数，同样用 It表示感染人群 (infected)的比例关于时间 t

12、 的函数，用 Rt表示移出人群的比例关于时间 t 的函数。 方程一：根据感染人群比例增长相等，可建立如下方程： () ( ) ( ) ( )d I t S t I t I tdt (1) 这里 表示感染者对易感染者有效感染的感染率 ， 即单位时间内单位病人传染的人数与易感者之比值 ； 表示移出率，即 单位时间内移出者占染病者的比率 ；()It 表示感染者比例关于时间 t 的函数； ()St 表示易感染人者比例关于时间 t 的函数。 方程二：根据易感染者比例减少相等，可建立如下方程： ( ) ( )dS S t I tdt (2) 方程三：根据移出者的增长相等，建立如下方程： () ()dR t

13、 Itdt (3) 这里 ()Rt 表示 移出者比例关于时间 t 的函数。 方程四：易感染者比例，感染者比例和移出者比例之和恒为 1： ( ) ( ) ( ) 1S t I t R t (4) 综上： ()( ) ( ) ( )( ) ( )()()( ) ( ) ( ) 1d I tS t I t I tdtdSS t I tdtd R tItdtS t I t R t(5) 5.1.3 问题 2 的模型求解 5.2.1 问题 2 的模型微分方程初始值确定 我们统计的数 据 1 是 从 2015 年 5月 20 日到 2015 年 7 月 7 日的 数据（见表一） ， 共计 49 天，对于移

14、除者比例 ()Rt ，在病情开始没有死亡人数，也没有治愈人数，所以移除者比例 ()Rt 的初始条件 (0) 0R ，对于感染者比例 ()It ，根据数据看出开始只有 1人，所以 感染者比例 ()It 的初始条件 (0) 1 /IN ，对于易感染者 St的初始条件可有 ( ) ( ) ( ) 1S t I t R t 算出，即 (0) 1 1 /SN 。 表一 对 MERS 疫情每天数据的统计 日期 当天确诊病例人数 累计确诊病例人数 当天死亡人数 累计死亡人数 20/05/2015 3 3 0 0 21/05/2015 0 3 0 0 22/05/2015 0 3 0 0 23/05/2015

15、 0 3 0 0 24/05/2015 0 3 0 0 25/05/2015 1 4 0 0 26/05/2015 1 5 0 0 27/05/2015 0 5 0 0 28/05/2015 2 7 0 0 29/05/2015 6 13 0 0 30/05/2015 3 16 0 0 31/05/2015 8 24 0 0 01/06/2015 4 28 1 1 02/06/2015 7 35 0 1 03/06/2015 5 40 1 2 04/06/2015 13 53 1 3 05/06/2015 7 60 1 4 06/06/2015 16 76 0 4 07/06/2015 16

16、92 0 4 08/06/2015 6 98 3 7 09/06/2015 16 114 1 8 10/06/2015 12 126 3 11 11/06/2015 6 132 2 13 12/06/2015 13 145 1 14 13/06/2015 3 148 3 17 14/06/2015 5 153 4 21 15/06/2015 4 157 1 22 16/06/2015 5 162 2 24 17/06/2015 2 164 5 29 18/06/2015 3 167 0 29 19/06/2015 0 167 1 30 20/06/2015 2 169 0 30 21/06/2

17、015 3 172 0 30 22/06/2015 3 175 0 30 23/06/2015 4 179 0 30 24/06/2015 0 179 1 31 25/06/2015 1 180 1 32 26/06/2015 1 181 2 34 27/06/2015 0 181 0 34 28/06/2015 0 181 0 34 29/06/2015 0 181 0 34 30/06/2015 0 181 0 34 01/07/2015 1 182 0 34 02/07/2015 2 184 1 35 03/07/2015 0 184 0 35 04/07/2015 1 185 0 35

18、 05/07/2015 0 185 0 35 06/07/2015 0 185 0 35 07/07/2015 0 185 1 36 5.2.2 问题 2 的模型参数拟合 对于模型上的参数感染率 和移出率 受很多因素影响，很显然不是一个常数。但是对于不是常数的参数很难拟合，所以，我们把所采集的数据 （见表一）分为适当个组，把每个组的参数看作常数，这样就很方便拟合出每个组的参数，最后再综合各个组拟合的参数，把参数拟合成关于自变量时间 t 的函数 （当然这里的组数分得越多则参数拟合得更加精确） 。这样就完整拟合出了每个时刻的参数值。 对于每个组的 参数拟合， 根据参数感染率 的定义（ 单位时间内单

19、位病人传染的人数与易感者之比值 ）而易感染人数远大于单位病人传染人数，所以 1 。移出率 （ 单位时间内移出者占染病者的比率）同样可以得知 1 ，所以可以定两个参数的范围都在 0到 1 之间，这里可以用尝试参数值去逼近最佳参数 ，具体步骤如下： Step1：因为模型的参数感染率 和移出率 范围都在 0 到 1 之间，这里我们 先用 0.1 的梯度把参数分成（ 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9）九个值，两个参数共组成 81 组值。 Step2：把分成的每组参数值代入模型的微分方程组，用 matlab 中的函数可算出一个数值解。 Step3：用上面求出来的函数

20、计算出的值和实际值作差计算平方和（即最小二乘法），找出 所有组参数值中最小的二乘数，即那组的感染率 和移出率 就是在梯度 0.1 下的最佳参数值。 Step4：画出已经计算出来的函数图，与实际数据的函数图对比，如果相差太远这进行下一步，反之转到 Step7。 Step5：为了更加靠近最终的参数最佳值，这里再在上次拟合出来的参数值左右一个梯度的范围再以上次梯度的 0.1 倍进行参数再分组。例如：假设上次算出来的最佳感 染率 =0.5 移出率 =0.2，这次的参数范围就是 0.4 0.6 0.1 0.3，这次就用 0.01 的梯度把参数分为 20组参数。 Step6：回到 Step2。 Step7：到了这里我们已经算出了一组数据的最佳参数值。 然后重复以上步骤，算出每组数据的最佳参数值。 Step8： 根据各组数据的最佳参数值画出参数值得散点图，拟合成适应函数。到这里，全部参数拟合结束。 算法流程图如图一所示： 图 1 算法流程图 开始 把参数分组 将每一组参数值代入微分方程组，求方程组的数值解 由每组参数算出函数的值与实际值比较，记录最小二乘数反应的两个参数 画出函数图，与实际函数图对比 再次将参数分组，代入微分方程组，求数值解 记录这组数据的额参数值 每组数是否找到了参数 相差 不大 相差很大 否 结束