消灭埃博拉病毒美赛二等奖.doc

消灭埃博拉摘要埃博拉出血热是由埃博拉病毒引起的一种急性出血性传染病。因其极高的致死率而被世界卫生组织列为对人类危害最严重的烈性传染病之一。由于埃博拉病毒在Guinea, Liberia, Sierra Leone的感染情况最严重，本文就埃博拉病毒在这三个国家为例，探讨了埃博拉病毒的传播和药物运输系统以及药物生产的需求和速度。本文引入了经典的Susceptibles Infectives Recovered（SIR）模型来描述埃博拉病毒的传播问题，预测埃博拉病毒随时间变化的传播规律；并对一系列收集数据进行了数值模拟，得到了SIR模型中日接触率和日治愈率的时间的多项式函数；进一步，根据其拟合的和，得到感染病人数的时间的函数关系式，进行了统计检验，其准确度较高。关于药品的运输系统问题，选取了America，England，Belgium，Canada四个国家为生产地，Conakry International Airport ()，Nzrkor Airport ()，Lungi International Airport ()，Kenema Airport ()，Roberts International Airport (), James Spriggs Payne Airport()六个机场为分配地，采用了度量法对六个分配地进行中心仓库选取，将药品运输分为两地足够供应和供应不足两种情况分别进行讨论:制定了抗击埃博拉病毒的药品运输系统，为药品的分发以及药品生产地的选择提供了良好的方案，并分别得到了紧急配送方案(如下两图)。 两地足够供应的交通运输网络 两地供应不足的交通运输网络并且优化运输系统，得到了最优药物运输路线(如下图),使制药产地可以更高效率运输到各个疫情高发区。优化交通运输网络根据World Health Organization (WHO)提供的药品需求量的测算公式，得出药品需求量与患病人数之间的函数关系式，从而算得药品生产速度与时间的关系式，给出最优生产速度，求得最短生产时间。进一步，根据最短生产时间与最短运输时间求得消灭病毒所需最短时间。总之，本文针对埃博拉的传播建立了SIR模型，对所述药物的所需要的数量，可能可行递送系统，递送位置，制造药物和其它关键因素都给出了解答，并提供了一些建设性的建议用于世界医学协会。关键词: SIR，数值模拟，度量法选址，网络优化一、 问题重述埃博拉出血热是由埃博拉病毒引起的一种急性出血性传染病，有着极高的致死率。目前非洲的大部分国家都受其毒害并且已经蔓延到其他国家。而且但一旦发病，发展十分迅速，重症患者发病两三天后就可能导致死亡。世界医学组织已经宣布其新药物治疗可以阻止埃博拉并治愈非晚期病人。为了消灭埃博拉病毒，或者减小当前埃博拉病毒的压力，我们需要解决以下问题：1. 建立疾病的传播模型；2. 构建一个可行的运输系统和运送地点；3. 分析求得药物数量、生产疫苗或药物的速度以及其他因素的影响；4. 为世界医学组织准备一份1-2页非技术信函，以用于他们的宣告。二、 模型假设1. 由于国家人口基数较大，假设人口总数不变，不考虑人口的出生率和自然死亡率，即为常数；2. 将埃博拉的传播途径视为与病人的有效接触，每个健康人变成病人的途径只能是有效接触，而不是突然变异；3. 不考虑气候条件对埃博拉传播的影响；4. 假设埃博拉爆发有个初始时间；5. 疫苗和药品都可用于消灭埃博拉。三、符号说明总人口数健康人口所占比例患病人口所占比例免疫人口所占比例患病人口总数传染期接触数为药品用量为缓冲百分比药物需求量药品生产速度四、问题分析由于问题要求建立疾病的传播模型，设立一个运输系统和运送的地点，以及求出药物的需求量和药物生产速度。对于疾病的传播，可以采用经典的Susceptibles Infectives Recovered（SIR）模型，建立出新增病人与时间的关系。最后可以对模型进行检验，看其是否准确。对于运输系统的建立，我们可以找到制药产到疫情高发区作为出发地和目的地，根据其距离和各个节点赋予其权重，利用度量法评分得到最优的交通网络图。对于药品需求量可以根据患病人数求得，对于制药速度，可以先求出消灭病毒所需要的时间，然后根据时间和病人之间的关系算得其药品生产速度。五、模型求解根据维基百科中的患病人数可得折线图1图1：患病人数折线图5.1疾病传播模型5.1.1模型的建立对于疾病的传播，根据维基百科中的患病人数可看出，埃博拉病毒在潜伏期传染性极低，所以我们采用经典的SIR模型初步建立微分方程人群分为易感染者，病人和病人免疫的移除者,分别三类人在总人数N的比例，健康人口所占比例，为患病人口所占比例，免疫人口所占比例。(1)图2：病毒发展流程图记初始时刻易感染人数和患病人口的比例分别为和()此时所以，均为常数，其中病人的日接触率为，日治愈率为，传染期接触数为 。建立SIR传染模型如下：(2)利用差分法处理数据以求得和，由（2）可得：(3)(4)和分别为易感染人数占总数的比例和患病人数占总数的比例；和表示易感人群数和患病人口数；和表示相邻两天和各自增加的人数。由于无法求出和的解析式，我们将模型消去并注意到的定义，可得： (5)容易求出（5）的解为,。 从而可以求得患病人数为： (6)5.1.2模型的求解对维基百科中所查得艾埃博拉病毒在非洲的传染情况进行整理，得到了由2014年3月22日到2015年1月25日患病总人数以及三个国家分别的患病总人数。表1：患病总人数统计表Sumt020,452,970292929-58320,452,910593059-60420,452,882862760-28520,452,8591031766-23920,452,8441141170-151020,452,8161301682-281120,452,805135588-111620,452,76316328102-421920,452,7511696108-122620,452,69020940129-613020,452,6772156136-133220,452,6652205143-123920,452,64223313153-234420,452,62324310162-195320,452,6192452164-46220,452,57727025181-426620,452,51930939200-587220,452,46635445208-537520,452,35943884231-1078020,452,30247436252-578820,452,16352854337-1399420,452,09159971338-7210220,451,802759160467-28910820,451,66684485518-13610920,451,60188844539-6511420,451,46196476603-14011720,451,3481,04884632-11311920,451,2751,09345660-7312220,451,1551,201108672-12012620,450,9761,323122729-17912920,450,7621,440117826-21413120,450,5381,603163887-22413420,450,3851,711108932-15313620,450,2881,77968961-9713920,450,1671,848691449,9841,9751271449,7562,1271521,145-22814620,449,5592,2401131,229-19714820,449,2052,4732331,350-35415020,448,9862,6151421,427-21915620,448,4073,0694541,552-57916120,447,4733,7076381,848-93416420,446,9384,0012942446,4854,3663652,177-45317120,445,8074,8464802,375-67817520,445,1035,3394932,586-70417820,444,5205,7624232,746-58318220,443,8656,2635012,900-65518420,443,4116,5743113443,0616,8082343,159-35018920,442,5507,1923843,286-51119220,442,0977,4923003,439-45319620,441,1308,0335413,865-96719820,440,6548,3863533,988-47620320,439,5468,9976114,485-1,10820820,438,5249,6936964,811-1,02221020,438,1879,9642714,877-33721520,436,13411,8681,9045,026-2,05322120,435,29412,6477795434,79613,0413945,191-49822620,434,80013,2682274,960423120,433,77014,0988305,160-1,03023320,433,43814,4133155,177-33223820,432,46315,1457325,420-97524020,432,21815,3512065,459-24524520,431,40415,9355845,689-81425220,429,81317,1451,2106,070-1,59125920,428,69817,9427976,388-1,11526620,427,51018,6036616,915-1,18827320,425,94319,4978947,588-1,56728020,424,91720,2067097,905-1,02628720,424,04620,7475418,235-87129420,423,30321,2965498,429-74330120,422,66321,7244288,641-64030820,422,12622,0923688,810-537其中，表1中每列分量除以总数，即可得到和以及和的离散值。于是对和的计算，就转化为利用第天的和以及和的离散值计算出当天对应的和的离散值，同理我们可求得其他三个国家的和的离散值（具体表格见附录）。表2：和的离散值Sumtt02.0000001.0000001480.143146 0.048928 31.0169490.5084751500.083748 0.029446 40.3255810.01162 1560.188661 0.040730 50.2233010.0582521610.251956 0.079849 90.1315790.0350881640.133717 0.060235 100.2153850.0923081680.103756 0.020156 110.0814810.0444441710.139909 0.040858 160.2576690.0858901750.131860 0.039521 190.0710060.0355031780.101180 0.027768 260.2918660.1004781820.104582 0.024589 300.0604650.0325581840.069060 0.021752 320.0545450.0318181860.051410 0.017039 390.0987120.0429181890.071051 0.017659 440.0781890.03703 1920.060464 0.020422 530.0163270.0081631960.120378 0.053031 620.1555560.0629631980.056761 0.014667 660.1877020.06148 2030.123152 0.055241 720.1497180.0225992080.105437 0.033633 750.2442920.05251 2100.033822 0.006624 800.12025 0.0443042150.172986 0.012555 880.26325 0.1609852210.066419 0.004823 940.1202000.0016692240.038187 0.007975 1020.3807640.169960226-0.000301 -0.017410 1080.1611370.06042 2310.073060 0.014186 1090.0731980.0236492330.023035 0.001179 1140.1452280.0663902380.064378 0.016045 1170.1078240.0276722400.015960 0.002541 1190.0667890.0256182450.051083 0.014434 1220.0999170.0099922520.092797 0.022222 1260.1352990.04308 2590.062145 0.017724 1290.1486110.0673612660.063861 0.028329 1310.13973 0.0380542730.080371 0.034518 1340.08942 0.0263002800.050777 0.015688 1360.0545250.0163012870.041982 0.015906 1390.0654760.02813 2940.034889 0.009110 1410.0926580.0283543010.029461 0.009759 1430.1071930.0357313080.024307 0.007650 1460.0879460.037500在已经得到了和关于的离散点后 ，查找三个国家的人口数量，可以得知几内亚人口总数为10057975人，利比里亚的人口总数为3955000，塞拉利昂的人口总数为6440053人。算得三个国家总人数=20453028人。根据数据我们用MATLAB拟合，得到三个国家以及它们总的的拟合图(分别为Sum, Guinea, Liberia, Sierra Leone)：图3：的拟合图求得三个国家以及它们总体的日接触率关于时间的函数为：同样的用MATLAB拟合出三个国家以及它们总体的(分别为Sum, Guinea, Liberia, Sierra Leone)的拟合图： 图4：的拟合图求得三个国家以及它们总体的日治愈率关于时间的函数为：进一步处理数据求得健康人数与总人数的比例值，将其作为拟合对象，得到关于时间的函数：可得到关于时间 的趋势图(分别为Sum, Guinea, Liberia, Sierra Leone)： 图5：健康人口趋势图我们可以把，和带入以求得病人增长率的表达式。5.1.3模型的检验为了保证模型的准确性，我们把求得的表达式用MATLAB进行进行检验,如下图所示，可以清楚看到患病人口所占的比例的实际数据与理论曲线吻合的相当不错，也就是说所求的表达式精确度很高。图6：模型检验图5.2运输系统模型5.2.1模型的建立由于病情紧急，对于药物的供给，我们采用空运的方式，根据埃博拉病毒的疫苗生产地得知，最有可能解决埃博拉病毒的八种药物cAd3-EBO Z、VSV-EBOV、Advac /MVA-BN、EBOV GP、Nasal vaccine、Vaxart tablet、Clinical trials in West Africa、U.S. national stockpile的工厂分别位于美国 英国 比利时 加拿大四个国家，因此，我们将这四个国家作为新药物的生产地进行分析。疫情高发区的三个国家共有六个机场，其中每个国家有两个机场。如下图所示：图7：机场分布图结合疫情严重的区域我们分别选取了如下表格的六个机场。分别将三个国家的六个机场其表示为，。表3：机场一览表GuineaLiberiaSierra LeoneConakry International AirportNzrkor AirportLungi International AirportKenema AirportRoberts International AirportJames Spriggs Payne Airport 把制药国家看为一个点，分别取New York，London，Brussel，Ottawa经纬度作为其地理位置求出其到各个机场的距离。根据地球上两点的距离公式可知： (7)经过计算我们得到了4个制药点分别到6个分派点的距离，如下表所示：表4：配送点到产地距离Distance /kmAmericaEnglandBelgiumCanadaSum6841.164826.344672.347091.9523431.797376.554926.564807.557608.0624718.726946.344921.484771.347200.9123840.077165.984958.694822.77412.524359.877319.535127.684995.937573.9425017.087356.235126.894997.817608.525089.43Sum43005.7929887.6429067.6744495.86十大最具潜力的埃博拉疫苗及药物十大最具潜力的埃博拉疫苗及药物十大最具潜力的埃博拉疫苗及药物十大最具潜力的埃博拉疫苗及药物由上表可看出，England和Belgium距离西非六个机场距离最近，所以选取这两个国家为最优生产地。同时可看出Guinea ,Liberia ,Sierra Leone的最优机场分别为Conakry International Airport , Lungi International Airport , Spriggs Payne Airport。我们计算距离的权重为：(8)其中为第个影响因素的初始权数；为所有影响因素初始权数求和。由于飞行距离与成本成正比，所以根据各地之间的距离计算出由生产地到机场的权重和由一个机场到另一个机场的权重分别为：表5：权重分布表The production place to the airportThe airport to the another airportWeights0.950.05根据权重公式(9)可得如下数据计算求得弧长权重分布表：表6：各分配点到产地弧长权重分布表WeightsAmericaEnglandBelgiumCanada0.047 0.033 0.032 0.048 0.050 0.034 0.033 0.052 0.047 0.034 0.033 0.049 0.049 0.034 0.033 0.051 0.050 0.035 0.034 0.052 0.050 0.035 0.034 0.052 Sum0.2930.2050.1990.304由于各国对药品需求不同，因此在急需药品治疗的情况下，可在6个配送点之间相互补给，所以我们计算出6个机场之间的相互距离，如表7所示：表7：各个机场之间的距离Distance /kmSum574.765116.304326.999482.576516.5612017.205574.765503.016272.654283.513254.0341887.982116.304503.016236.688373.433410.0451639.486326.999272.654236.688184.78205.7111226.832482.576283.513373.433184.7844.2911368.593516.561254.034410.045205.71144.2911430.642并可根据表7所示数据分别计算个机场之间弧长权重，计算结果如下表：表8：各个机场之间弧长权重分布表Weights0.0600.0120.0340.0500.0540.0600.0530.0280.0300.0270.0120.0530.0250.0390.0430.0340.0280.0250.0190.0210.0500.0300.0390.0190.0050.0540.0270.0430.0210.005对表4中数据进行分析，由于四个机场中Belgium和England离供药地点距离较近，其配送方案为最优，考虑到其配送地点的供药能力有限，不足以供给所有病人的需求，因此我们将问题分为这两个药品产地药品供应量足够和药品供应量不足两种情况。5.2.2两个药品产地的药品供应量足够。可由England和Belgium直接供药，由于Belgium到各个机场距离均最短，所以以Belgium为优先供应区，药品不足部分由England供应。由于各国对药品需求不同，因此在急需药品治疗的情况下，可在三个国家内紧急配送。根据表6所述各机场到England和Belgium的权重，计算结果如下表：表9：各机场到两最近产地的权重分布表WeightsEnglandBelgiumSum0.033 0.032 0.0650.034 0.033 0.0670.034 0.033 0.0670.034 0.033 0.0670.035 0.034 0.0690.035 0.034 0.069将6个机场视为节点，根据表5权重分布表的数据按公式(9）式中为节点所占权重，为节点到两产地距离权重之和，为节点到另外五个节点的距离权重之和。计算可得各节点权重分布，计算结果如下表10表10：各节点所占权重分布表Weights0.1690.1720.1540.1640.1700.171由于权重最小，因此可得最优节点为，所以由到其余各节点和产地距离之和最节约成本。因此在有生产地到六个机场分配药物不能满足个别地区所需药物要求时，可将视为仓库，及时补充各地区所需药品。利用度量法评分（此时已将视为仓库，所以不考虑移走其他点对的影响）表11：度量法评分表Removed pointInfluenced pointAssigned pointsIncremental/ WeightsSum000000000.0220.02200000.0020.0210.0190.0010.0150.015由表中可以看出，将各个点移走时受影响最大的点为和，因此可将和设为仓库点。由于位于与中间，而且与点最为接近。所以选择作为其仓库点，最后得到的交通运输网络如下图所示：图8：供大于求的交通运输网络。最优化消费即最小成本公式为： (10)表示第i个需求点的需求量，表示从需求点i到配送中心仓库j的单位的运输成本（其包括时间，距离，费用等，由于疫情的紧急性，我们暂时只考虑距离）表示需求点i是否由配送中心仓库j来提供服务，表示配送中心仓库j是否被选中成为设施，如果选中则为1，未被选中则为0。式（10）应满足如下约束条件:其中表示n个需求点集合；表示m个建设配送中心仓库节点候选点集合；p表示允许投建的配送中心仓库数。5.2.3两个药品产地的药品供应量不足以Belgium为优先供应区，药品不足部分由England首先供应，次之有America和Canada供应与5.2.1相似，用同样的计算方法计算，在两个生产地药品生产供应不足的情况下，建立其最优交通运输网络。计算其弧长权重为下表所示：表12：各机场到四个国家弧长的权重WeightsAmericaEnglandBelgiumCanadaSum0.0470.0330.0320.0480.160.050.0340.0330.0520.1690.0470.0340.0330.0490.1630.0490.0340.0330.0510.1670.050.0350.0340.0520.1710.050.0350.0340.0520.171由表12中生产地到机场和由机场到机场成本的权重数据可得：表13：各节点权重分布表Weights0.1630.1700.1630.1650.1700.170由于和所占权重最小，因此其距离最短，由此可知机场最优点为和。因此在有生产地到六个机场分配药物不能满足个别地区所需药物要求时，可将和视为仓库，及时补充各地区所需药品。利用度量法评分计算：表14：度量法评分表Removed pointInfluenced pointAssigned pointsIncremental/ WeightsSum0.130.1300000.0220.02200000.0020.0210.0190.0010.0150.015由表中可得将各个点移走时受影响最大的点为和。综合以上因素可得和为最优点，因此可将和设为仓库点。由于与相对距离较小，所以由对进行药物供应；又位置较居中，所以可由对,进行药物供应。同时因为两仓库点距离较近，因此和药物可互相流通。最后我们得到最优的交通运输网络如下图所示：图9：供不应求的交通运输网络5.2.3运输网络的优化根据表8中数据可得西非三国家内药物互相补给路线，得到各节点最优分配方案，其权重如表15：表15：每个节点最优分配方案表Destination2 nodes3 nodes4 nodes5 nodes6 nodesSum0.0120.0370.0560.0610.0880.2540.0270.0320.0510.0760.0880.2740.0120.0440.0650.070.0970.2880.0190.0240.0510.1040.1160.3140.0050.0260.0510.0630.1230.2680.0050.0240.0490.0610.1210.26由于各国对药品需求不同，因此在急需药品治疗的情况下，可在三个国家内互相配送，对表15数据进行分析得到最优药物紧急配送最优路线如图10图10：优化模型路线图最优运输网络如下图图11：优化运输网络因此在有生产地到六个机场分配药物不能满足个别地区所需药物要求时，可根据图9进行紧急配送。5.3药品需求量模型根据资料显示，World Health Organization (WHO)提供的方法，药品年度需求量的测算公式为：下一年度某药品的需求量=下一年度预计正在治疗的人数*预计某药物的使用人数百分比*药物用量*365天*（1+缓冲百分比）-库存剩余量我们对其公式进行改得到其需求量随时间的函数：(11)其中为药品用量， 为缓冲百分比，根据资料显示，其值一般为25%。(12)5.4药品速度求解以一天为间隔计算药物生产速度，为生产速度，为间隔时间一天， (13)为相邻两天病人增加的人数即，即(14) 即为求得的药品生产速度。六、结论通过改进经典的SIR模型，运用MATLAB软件拟合，求得其参数值，从而使模型的预测达到了相当高的精度。十分精确的刻画出患者总数随时间变化的规律。 且模型数据充足，为我们论述的观点提供强有力的证据。但模型中的值恰当拟合，是保证预测准确的关键所在。在讨论模型细节时，将的值认为是常量，从而变为只与时间有关的变量。这一假设合理的简化了模型，但同时使模型必须在人口分布相对均匀的地区才有较高的预测精度。但是对于单个国家来说还是比较有效的。对于人口分布相对不均匀的多个国家来说，可以将模型进一步优化。考虑随机因素的影响，我们可以采用随机Susceptibles Infectives Susceptibles（SIS） 模型或者随机Susceptibles Infectives Recovered（SIR） 模型，易得到更好的结果。在运输系统的模型中，我们分供应量大于和小于需求量两种情况分别建立了最优化模型，通过综合考虑弧长和节点的权重，建立了最优化运输网络。还求得了最佳应急方案。使得各大制药公司可以最低成本的运输最大量的药物，加快了药物之间的互相补给，可以以最短时间，最高效率运输到各个疫情高发区，但其计算过程略显繁琐。参考文献1Ebola virus. Scientific & Technical Information of Gansu.2014.43(8).12Ebola virus epidemic in West Africa, /wiki/Ebola_virus_epidemic_in_West_Africa3Guinea, /wiki/Guinea4Liberia, /wiki/Liberia5Sierra Leone, /wiki/Southern_Province,_Sierra_Leone6Ebola vaccine, /wiki/ Ebola_vaccine7Pherical distance between the two places on Earth.Journal of Yulin College.1996.03.61-638Mozhi Li, Xiaodong Liu. The model of pharmaceutical distribution center warehouse location, China Economist,2014(12), 277-2789WHO. Demand forecast for antiretroviral drugs in low and middle-income countries, HIV/AIDS Programme, Strengthening health services to fight HIV/AIDS.R.2007.10Yanli Zhou, Weiguo Zhang. The model of stochastic SIS infectious diseases of nonlinear incidence rates Persistence and extinction. Journal of Shandong University(Natural Science).2013,10:68-76.11Hongjing Shi，Zhisheng Duan，Guanrong ChenAn SIS model with infective medium on complex networksJPhysicalA, 2008,387: 2133-2144.12 Linlin Meng, Sanling Yuan. The Asymptotic Behavior of an Stochastic SIR Epidemic Model. Journal of Biomathematics,2013,01:47-52.消灭埃博拉埃博拉疾病已被我们认定为21世纪人类遭遇的又一场全球性流行疾病，对人们的生活和经济发展构成了严重威胁。自2014年2月在几内亚爆发以来，已经迅速蔓延至利比里亚等其他西非国家，对人们的生活和经济发展构成了严重威胁。面对突如其来的埃博拉疫情的严峻考验，冷静客观的分析就显得尤为重要。对此我们建立SusceptiblesInfectivesRecovered（SIR）模型，对近段时间各个地区的感染情况进行分析，预测了病毒随时间的传播规律,为增加治愈率提供一些有建设性的意见。同时，对疫情重灾区进行隔离，减少日接触率，从而减少埃博拉病毒的传播。进一步，我们制定了抗击埃博拉病毒的药品运输系统：考虑到供求关系，对各种情况进行分析，制定了各种情况下的运输系统和运输路线，为药品的分发和生产地的选择提供了良好的方案。由于药物的稀缺以及疾病的大肆蔓延，各地药品需求量巨大。为了应对此严峻形势，各大制药公司可根据我们提供的路线以最低的成本运输最大量的药物，不但减少了中间环节与患病者的花销，还可以加快药物之间的互相补给，以最短时间，最经济，最高效的运输效率运送到各个疫情高发区，使病人快速得到治疗。最终各国将紧密合作，应对挑战，采取切实可行的措施，努力防治这一疾病，以更快速度消灭埃博拉，确保各国的繁荣、稳定和发展。保证各国社会安定和人民生命安全，减少埃博拉的负面影响。 附录附录1：Guinea患病总人数统计表Guinea 411,199,92886060-72511,199,9051031766-23911,199,892112970-131011,199,8721221080-201111,199,864127583-81611,199,8281512495-361911,199,8161576101-122611,199,75519740122-613011,199,7422036129-133211,199,7302085136-123911,199,70722113146-234411,199,68823110155-195311,199,6842332157-46211,199,64225825174-426611,199,60728123186-357211,199,59029110193-177511,199,51534453215-758011,199,46637228236-498811,199,41239826264-549411,199,414390-8270210211,199,35841323303-5610811,199,359408-5307110911,199,3564091309-311411,199,364406-3304811711,199,3544104310-1011911,199,3454155314-912211,199,32842712319-1712611,199,27546033339-5312911,199,25647212346-1913111,199,23148513358-2513411,199,21649510363-1513611,199,2124950367-413911,199,19550611373-1714111,199,