关于传染病模型.doc

关于 SARS 模型的建立与相关的预测分析 本文先根据材料提供的模型与数据较为扼要地分析了附件 1 的模型的优缺点， 摘要：全面地评价了该模型的合理性与实用性。而后在对问题进行较为全面评价的基础上引入更为全面合理的假设和建立系统分析模型。运用联立微分方程组体现疫情发展过程中各类人的内在因果联系，并在此基础上运用经典的龙格库塔微分方程求解算法结合MATLAB 编程程序在附件一拟合出与实际较为符合的曲线并进行了疫情预测。同时运用双线性函数模型对卫生部的措施进行了评价并给出建议。而后运用差分方程（程序在附件二）就 SARS 对经济（主要是旅游业）的影响进行了较为准确的分析，进而通过模型算出的理论预测数值与实际数值进行对比，以数值上的显著差异直观地表现了 SARS对经济（旅游）的影响，并对接下来的几个月进行了较为合理的预测。本文的最后，通过本次建模过程中的切身体会，以一篇短文评述去说明建立如 SARS 预测模型之类的传染病预测模型的重要意义。关键词：微分方程 龙格库塔算法 SARS 双线性函数模型 差分方程 数学模型 1一 问题的重述 SARS（Severe Acute Respiratory Syndrome，严重急性呼吸道综合症 俗称：非典型 是肺炎） 21 世纪第一个在世界范围内传播的传染病。SARS 的爆发和蔓延给我国的经济发展和人民生活带来了很大影响，我们从中得到了许多重要的经验和教训，认识到定量地研究传染病的传播规律、为预测和控制传染病蔓延创造条件的重要性。请你们对 SARS的传播建立数学模型，具体要求如下： （1）对附件 1 所提供的一个早期的模型，评价其合理性和实用性。 （2）建立你们自己的模型，说明为什么优于附件 1 中的模型；特别要说明怎样才能建立一个真正能够预测以及能为预防和控制提供可靠、足够的信息的模型，这样做的困难在哪里？对于卫生部门所采取的措施做出评论，如：提前或延后 5 天采取严格的隔离措施，对疫情传播所造成的影响做出估计。附件 2 提供的数据供参考。 （3）收集 SARS 对经济某个方面影响的数据，建立相应的数学模型并进行预测。附件 3 提供的数据供参考。 （4）给当地报刊写一篇通俗短文，说明建立传染病数学模型的重要性。（二）对附件 1 所提供的模型的评价 该模型的合理性首先体现在模型假设上：“假定初始时刻的病例数为 N0，平均每病 （K ，人每天可传染 K 个人 一般为小数） 平均每个病人可以直接感染他人的时间为 L 天。” 其一，一般来说每病人每天可传染的人数与当时的健康人数有关1，但由于北京的人数基数较大，SARS 病人数相对较少并且 SARA 持续时间不是很长，所以这样假设也是可以的。 其二，每个病人可以直接感染他人的时间是有限的，该模型考虑到了这一点，也是很合理的。 该模型的合理性还在于用数理统计的方法估计相关参数。 该模型的实用性是较好地模拟与预测了北京的 SARA 数据与发展。在传染病发病初期对疫情的预测结果还是较为理想的，这主要得益于发病初期，由于病情来得突然，有关部门没有来得及采取措施加以控制，使病情得以蔓延迅速，而且发病初期在治疗方法上不是特别有效，治愈所需的时间长，所以使用 N t N 0 1 k t 作为模型进行估计以及参数的假设均较为合理，基本上是可行的。 但是到了疫情发展中后期，由于政府部门采取强硬措施加强防治工作以及人民群众的防范意识与警觉程度上的普遍提高，加之治疗措施的改进，使得每天被传染的人数下降，并且治愈的人数在不断增加，治愈时间也在不断缩短，每天的病人数应在上一天的基础上减去治愈和死亡的人数， “ “ 并且由于采取强硬措施 L”的取值会大大的减小， K”取值也会是个变量，而不是常数。大多数疑是病人往往在早期就会被隔离，所以，基本 2上很少能转化成自由非典病人而去接触并传染别人。如果此时还是选取 N t N 0 1 k t这样的单调递增函数作为预测模型，就会有较大的误差。 该模型的另一个不足是没有考虑 SARS 的潜伏期，也没有对人群进行合理的分类（如易感染人群、病人、治愈人群等等），所以必须建立更为合理的假设与模型。（三）定义与符号说明S：表示易感染人群即健康者在人群中的比例。E：处于潜伏期人群在人群众的比率。 ，这种人暂时未发病，但他们最终将发病。I：已受感染者即病人在人群中的比例。R：移出者包括“出院者”和“死亡者”在人群中的比例。M：未被隔离的带菌者。X：疑似病人。a:每个病人每天有效接触并使之感染的平均人数（常数）。b:退出率，即 SARS 患者的每日死亡率和每日治愈率之和。c：潜伏期的病人的日发病率。d: 每个未被隔离带菌者被隔离前平均每天感染有效人数。x1: 疑似中每日被排除的人数占疑似人数的比例。x2：疑似者中每日确诊的人数占疑似人数的比例。j：每个未被隔离的带菌者转化为病人的日转化率。k：被未被隔离的带菌者有效传染的人中可以控制的比例。0:SARS 潜伏期天数 。（四）模型的假设：总体假设：材料提供的疫情统计数据真实有效。基本假设：1、潜伏期一般是 2 到 11 天，我们假设为5 天。2、据医学权威表明，传染途径主要是 SARS 患者，不包括处于潜伏期的。3、SARS 患者被治愈后不具有传染性，也不会再被传染。4、北京市的总人数可视为常数，即不考虑流入人口与流出人口的影响。也不考虑这段 （R）时间的人口出生率与自然死亡率。把由 SARS 引起的死亡人数视为“移出者” 。5、与 SARS 病人接触后都会被感染。（五）模型的分析： 根据附件三（我们自己从中国网收集得到的数据，比赛题给的数据要详细）中的数据表明，由于疫情初期政府控制力度不够，大众的对 SARS 的防范意识不强，造成病情迅速蔓延。而当政府采取有力措施，人们的防患意识增强，疫情则趋于缓和，病患者人数迅速下降。所以 SARS 传播大体上可分为两个阶段：1控制前期：即认为病毒传播方式是自然传播。2控制后期：政府强力介入之后的病毒传播模型。我们以附件三提供的北京市疫情统计数据为基础建立模型。 3 （六）模型的建立： 1、控前模型的建立： 将人群分为易感染人群、已受感染者和移出者三类。时间为 3 月 1 日到 4 月 19 日。 记 S：表示易感染人群即健康者在人群中的比例。 I：已受感染者即病人在人群中的比例。 R：移出者包括“出院者”和“死亡者”在人群中的比例。 a:每个病人每天有效接触并使之感染成为潜伏病人（在天后发病）的平均人数。 b:移出率，即 SARS 患者的每日死亡率和每日治愈率之和。 0:SARS 潜伏期天数 。 则有： dS dt aI t 0 S t 0 1 dI aI t 0 S t 0 bI t 2 dt dR bI t 3 dt S t I t Rt 1 4 初始值：S0 S1 S2 S3 S4 I0 I1 I2 I3 I4R0 该微分方程组带有时滞因素，没有解析解。因此我们考虑用差分方程的办法，因为潜 ，I（1）伏期为 5 天，所以 I（0） ，I（3） ，I（2） ，I（4）分别表示疫情开始前 5 天的病人数，从第 5 天开始有： I（5）I（4）aS（0）I（0）bI（5） 即： I（5）I（4）aS（0）I（0） /1b 同理， I（6）I（5）aS（1）I（1） /1b ： I（n）In-1aSn-5In-5/ 1b 5 对 S（t）亦有 S（5）S4-aI0S0 S（6）S5-aI1S1 ： S（n）Sn-1-aIn-5Sn-5 6 对 R（t）则有 R（n）Rn-1bIn n12 7 由此可得 SARS 的控前差分方程模型，即： In In - 1 aSn - 5In - 5/ 1 b 5 Sn Sn - 1 - aIn - 5Sn - 5 6 Rn Rn - 1 bIn 7 In Sn Rn 1 8 如果初始值给定，并将参数 a，b 确定（a：由有关数据推导得出。b：由医疗水平和有 4关数据分析得出，取其平均值）就可计算出任一天的易感染人群、已受感染者和移出者的数目， 但可惜的是，由于疫情初期政府控制力度不够，没有提供给我们真实有效数据，如北京首例 SARS 病人出现在 3 月 1 日，但只有 4 月 19 日以后的数据，所以我们只能进行模型建立和分析，而不能求解模型。 这也是建立真正有效的能预测的模型的困难之一。困难之二是这个微分方程组的求解极其困难。困难之三是我们不知道政府在何时干预及力度如何。2、控后模型的建立： 将人群分为易感染人群、已受感染者、移出者、疑似病人和未被隔离的带菌者五类。设控制开始时间为 4 月 21 日。 记 S：表示易感染人群即健康者在人群中的比例。 I：已受感染者即病人在人群中的比例。 M：未被隔离的带菌者。 X：疑似病人。 R：移出者包括“出院者”和“死亡者”在人群中的比例。 a:每个病人每天有效接触并使之感染的平均人数（常数）。 b:移出率，即 SARS 患者的每日死亡率和每日治愈率之和。 d: 每个未被隔离带菌者被隔离前平均每天感染有效人数。 x1: 疑似者中每日被排除的人数占疑似人数的比例。 x2：疑似者中每日确诊的人数占疑似人数的比例。 j：未被隔离的带菌者转化为病人的日转化率。 k：被未被隔离的带菌者有效传染的人中可以控制的比例。则有： dS dt x1 X t dM t S t 8 dI jM t bI t x 2 X t 9 dt dR bI t 10 dt dX x1 X t x 2 X t dkM t S t 11 dt dM d 1 k M t S t jM t 12 dt StItRtXtMt1 13 S0I0R0X0M0为初始值参数的确定：我们以材料提供的北京市疫情统计数据来说明参数的分析方法。（见附件三）以下全部图的坐标 0 均表示 4 月 19 日。（1）x1： x1（每天新增的疑似排除人数）/ （当天疑似病人累计人数当天移出累计人数）首先我们先直观的观察一下 X1 的变化趋势。根据材料提供的数据，用 MATLAB 来出 5x1，并画图，如图 1 所示： 图1接着用曲线拟合图 1，如图 2 所示： 图2 从上图可看出，图 2 大概有两个峰值。第一个高峰可能是疑似者中非感染者较高；第二个峰值则是因大部分真正带病的疑似者已转化为确诊后， 未带菌者相对比例增大造成的。 由此 4 阶拟合得出的曲线误差很大，为此，我们去掉几个偏差太大的点后，易看出，x1 集中分布在 0 到 0.05 之间。从图中，可以发现，最集中的数据为 0.035，这样我们就以 0.035 为 x1 的估计值。（2）x2：x2（每天新增的疑似转化为确诊的人数） / （当天疑似累计人数当天累计移出者）首先观察 x2 的变化趋势，如图 3 所示： 6 图3用 5 阶曲线拟合，如图 4 所示： 图4 从图 4 可见，x2 在疫情得到重视后一直下降。由图还可以看出 x2 的值主要分布在0 .0005 和 0.015 之间。我们根据 SARS 具有潜伏期的情况，估计 x2 分为两个阶段值：0.0223 和 0.006。 从对 x1 与 x2 的数据处理来看，我们可以将控制后期的这段时期分为两个阶段：过渡期和平稳期； 这两个阶段的产生是与非典自身的特性分不开的。 由于非典具有潜伏期，所以在控制初期，由于前一段时间对非典的控制力度不够，造成较多的人处于 SARS 潜伏期，这一部分人最终将转化为 SARS 病人；且因为他们为未被隔离带菌者，在进医院前会传染较多的人；加之各项措施从颁布到实行总会有一段反应时间，所以上述原因直接导致了过渡期的形成。（3）b：b（每天新增的出院者和死亡者的人数）/（当天病人累计人数当天累计移出者人数）首先观察 b 的变化趋势。如图 5 所示： 7 图5用 3 阶曲线拟合，如图 6 所示： 图6由图还可以看出 B 的值主要分布在 0 .005 和 0.09 之间。 也可以分为三个阶段值： b 0.0085，0.036 和 0.085。（4）j：从材料提供的数据可以估计出其值在 0.12 到 0.25 之间。现在我们令 j0.23。（5）k：根据各地的人们的意识和习惯等因素反映出来，比如在控制期对人口流动的控制严格程度，减少聚会等措施。由此我们估计 k 在 0.7 到 0.9 之间。（七）模型的解法与求解： 经过尝试，我们可以知道无法从建立的微分方程直接求出 S，I，R，X，M 的函数表达式。我们尝试编写龙格库塔微分求解算法求它们的数值解。 龙格库塔算法（求解微分方程组）： K1ftZn K2ftZnh/2K1 K3ftZnh/2K2 K4ftZnhK3 Zn1Znh/2K12K22K3K4 步长 h 取 1 8 Zn 为矩阵（S，I，R，X，M） f 为 5 个方程的形式 K（K1K4）为（4X5）的矩阵 1 由 Zn 的初值可顺序算出 K（K1，K2，K3，K4）的值，进而推出 Zn1 的值。据此 （递推可得任意的 Zn 的值。 matlab 程序附于附件一中） 我们通过求出的数值解跟根据材料算出的每一天的 S、I、R、X、M 的值作对比，并画出图像。根据图像可以看出实际数值和理论计算数值存在一定的差异。因此，我们可以通过调整由估计得来的数据来使两个图像趋于一致。 通过分析估计有： x10.035，x20.02230.006，b0.00850.1913， d0.555，j0.22，k0.71。这组数据使图像与实际最为接近。图 7 是材料提供的北京市累计病人比率图： 图7 （蓝圈表示实际的数据，红线 以下是由龙格库塔方法求得的 I 的数值解画出的图：是 I 的数值解） 图8 注：0 点表示为 4 月 19 日，从此开始往后计算天数。从上图曲线的分析可得： 91病人数在 5 月 13 号4 月 19 日往后 25 天前后达到高峰期。2病人数在 5 月 18 日左右出现最大值，且 Imax0.00015，此后开始下降。3病人数在 6 月 10 号前后趋于缓解，并不断递减，在 6 月 10 号至 7 月 10 初为缓解阶 段。这是由于政府采取了大力控制措施，使大多数带菌者被隔离，从而被传染的人 数越来越少。4病人数将在 130 天后即 8 月底左右趋于 0，即疫情的最终控制期。（八）模型的优缺点一模型的优点：1首先此模型进行了较为详细的分类，使得微分方程的联立显得较为紧凑，相关性较强。因为材料提供的数据齐全，各参数的设置合理，所以参数的设计在详细的数据的支持下和实际情况较为接近，为图形的拟合和对最终控制期的预测奠定了基础。2该模型适用范围较广，只要数据足够地精确详细，则求出估计参数便可求解。3该模型虽然是控后模型，但只要人的警惕性增强，把控制期提前到疫情初期，那么只要做适当的参数修改便可以作为疫情发展的全过程的预测模型。二模型的缺点：1首先是没有考虑年龄结构层次对疫情的影响。因为根据医学研究表明，儿童与老人极易感染非典病毒，而青壮年由于精力旺盛使的活动积极者，且由于控制后期的前期人们的活动几率仍然较高，所以接触的几率较大。而在控后期的后期由于人们的活动水平降低，使得接触的几率降低，接触几率的不断降低这一点在模型中并没有得到很好的体现。2随着医疗水平的提高与人们的意识水平、政府的严厉的控制措施，退出率实际上是不断提高的，而在我们的模型中却认为这是一个常数。（九）控后模型的改进：1在模型中引入接触率与移出率应随着病人的减少而变化，可能会是趋于 p 的变化，随着时间推移有所调整，以附和预测的需要。2在模型中引入接触感染率的概念，即体现接触不一定感染，只不过是感染率较高。3对控后模型加入潜伏期对病毒的传播造成的影响。（十）对于卫生部门采取的措施的评价 对于卫生部门提前或延后 5 天采取严格的隔离措施的影响，我们可以建立下面的模型进行辅助分析估计：1模型参数定义：Stt 时刻易感人群总数It t 时刻出现的新增患者 患者从患病起经过 时间，仍为患者的概率 患者距发病 时间，具有传染性的概率 患者与易感人群接触率 10近断时间的医学研究表明，从正式发病到治愈一般需 714 天或更长时间，假定平均治愈时间为 12 天。2基本条件假设：新患者出现的数量与现有患者的数量成正比，也与现有易感者的数量成正比，即发病率是患者人数和易感者人数的双线性函数。由基本假设条件可得：S（t1）St-It1 （1） tIt1 St I t （2） 0经整理后得： tSt1St- S t I t （3） 0 tSt1St1- I t （4） 0 tSt1/St1- I t （5） 0