SARS传播模型

成员：王策 刘丽英 汪先鹏

摘要：本文讨论了SARS疫情的传播规律和对经济方面的影响。首先本文对题中的早期模型进行了评价，认为其最大的优点是：能较好的描述疫情早期的发展情况，并在理论上可大致预测出疫情的爆发点以便卫生部门及早控制；最大的不足是：原模型在求解过程中参数Ｋ经过多次手工调整，而且L取为一个定值，此做法主观性太强，缺乏普适性。针对上述模型的不足，本文在原模型基础上进行了改进，得到了一个更为优化的Logistic SARS模型，并且它与实际情况符合得较好。同时，我们也通过Matlab语言对北京、山西等的计算值和实际数据进行了拟合，进而验证了这个模型的可靠性。当然，要建立一个最优模型还需要考虑更多因素，在考虑了传播途径及易感染人群等因素后，也可以建立一个最优的SEIRQ模型。但这样考虑就需要大量的数据采集整理工作，但在实际中不易实现。对卫生部所采取部分措施的评析中，我们引入了小世界网络模型，对政府措施做出定量的评论，并用图形直观的表示出来。我们还分析了Logistic SARS模型特点，并对其改进与应用做出了展望。

本文还通过对北京市1997～2002年各月接待的海外旅游人数的分析并建立了时间序列模型，“预测”出2003年疫情期间本应接待的人数，对比实际接待人数，计算出在非典期间少接待的游客人数约为115万人。最后我们给当地报刊写了一篇短文，说明了建立传染病数学模型的重要性。

关键词：Logistic SARS模型，SEIRQ模型，小世界网络模型，自回归。

１

一、问题重述

SARS（Severe Acute Respiratory Syndrome，严重急性呼吸道综合症, 俗称：非典型肺炎）是21世纪第一个在世界范围内传播的传染病。SARS的爆发和蔓延给我国的经济发展和人民生活带来了很大影响，我们从中得到了许多重要的经验和教训，认识到定量地研究传染病的传播规律、为预测和控制传染病蔓延创造条件的重要性。请你们对SARS 的传播建立数学模型，具体要求如下： （1）对附件1所提供的一个早期的模型，评价其合理性和实用性。 （2）建立你们自己的模型，说明为什么优于附件1中的模型；特别要说明怎样才能建立一个真正能够预测以及能为预防和控制提供可靠、足够的信息的模型，这样做的困难在哪里？对于卫生部门所采取的措施做出评论，如：提前或延后5天采取严格的隔离措施，对疫情传播所造成的影响做出估计。附件2提供的数据供参考。

（3）收集SARS对经济某个方面影响的数据，建立相应的数学模型并进行预测。附件3提供的数据供参考。

（4）给当地报刊写一篇通俗短文，说明建立传染病数学模型的重要性。

三、合理假设及说

1、假设所研究的人口为理想下的人群，对该病普遍易感染，每个发病病人单位时间内传染的易感染者人数与未被感染的人数是成正比，隔离或预防意识增强可在一定程度上影响病人单位时间内传染易感染者人数的比例。 2、不考虑温度、气压等自然因素对SARS发病的影响。 3、假设预测地区足够大，患病人数足够多。

4、在整个过程中不考虑由人口流动因素、人口自然出生率与死亡率。

２

二、问题分析

（一）通过对早期模型和实际情况的分析，我们认为影响SARS传播因素众多，大致可分为时域因素和地域因素。

(1)时域因素

a．媒体宣传：初期疫情较轻，媒体宣传强度很弱，导致民众的自我保护意识不足，容易感染；后期疫情较重，媒体宣传强度很大，民众的自我保护意识大大加强。

b．政府干预：初期疫情较轻，政府介入不足，后期疫情较重，政府加强干预（如：强行隔离，公共场所消毒等行为）。

c．认识程度：当一种新的传染病出现时，初期由于人们的认识程度不足，无法采取有效的预防和治疗措施，但随着研究的深入，认识程度会越来越高。

(2)地域因素

a．经济水平和医学水平：经济水平和医学水平高的地区的疫情控制情况会明显比水平低的地方好。

b．人口密度和人口流动：人口密度和人口流动大的城市若爆发传染病，疫情程度会比人口密度和人口流动小的城市大。

c．气候：SARS适合在春秋两季传播，且各城市的气候会疫情程度。

（二）我们认为在SARS疫情期间考察的人群大致可分为三类：健康人群，感染人群，治愈人群。而感染人群又可分为非传染源和传染源两类。 （三）一个较好的传染病传播模型因该具有如下功能： a．能较好的描述疫情的大致走势。

b．能较精确的给出关键时间（初期爆发时刻；中期稳定时刻；高峰期；0病例增长的时刻），以便政府和卫生部门针对不同情况作出及时而正确的措施。 c．能给出描述疫情的指标，以便政府和卫生部门决定其各项工作的力度。

四、变量说明

N: 某地累计发病人数 t: 计算病例的初始时间 N0: 时的累计发病人数

r: 发病率

k: 预防效果指数

Mmax: 理论预计累计发病最多人

r2: Logistis模型的决定系数 T2: 发病高峰时间

３

五、早期模型评价

（一）早期模型重述

假定初始时刻的病例数为N0，平均每病人每天可传染K个人，平均每个病人可以直接感染他人的时间为L天。则在L天之内，累计病例数N(t)随时间t（单位天）的关系是：

N(t)?N0(1?K)t （1-1）

如果不考虑传染期的限制，病例数将按指数规律增长，考虑传染期的限制后，则采用半模拟循环计算的办法，把达到L天的病例从可以引发直接传染的基数中去掉。然后假定从开始至高峰期之间均采用同样的K值（从拟合这一阶段的数据得出），到达高峰期后，在10天的范围内逐步调整K值至较小，然后保持不变，拟合其后在控制阶段的全部数据（认为社会在短期剧烈调整后，进入对疫情控制较好的常态）。

（二）早期模型的合理性和实用性的评价 A．早期模型的优点

dN?K?(N(t)?N(t?L))在（0～L）区间内的1、模型（1-1）实际上是微分方程dt特解。其中N(t)表示t时刻的累计病例数，则(N(t)-N(t-L))表示传染源数量，为病例总数减去失去传染能力的病例数。

2、参数K和L是描述SARS传播的两个重要参数，并且具有实际的意义：L可理解为平均每个病人在被发现前后可以造成直接传染的期限，在此期限后失去传染作用，可能原因是被隔离、病愈或死去等等。K表示某种社会条件下平均每病人每天传播的人数（但并非文中所述的一个病人的感染他人的平均概率）。

3、通过对模型的分析，可得到一个预测疫情发展的参数R＝K×L，R表示平均每个病人在其传染期内感染的总人数。若R＞1，说明社会上现存传染源人数在上升，疫情将因失控而爆发；若R＜1，说明社会上现存传染源人数在下降，疫情将得到控制；若R＝1，说明社会上现存传染源人数不变，疫情将持续下去。

4、从拟合的图形来看，此模型对各城市早期的SARS疫情描述的较好，具有一定的通用性。其实际意义就是此模型可大致预测出疫情的爆发点和发展趋势。因为通过对早期数据的拟合确定参数，得出形如（1-1）的一个指数形式的模型，而指数函数的曲线初期增长较慢，后期增长急速，必可大致找到一个“转折点”，而“转折点”所对应的时间便是预测的疫情爆发点。这一数据对于卫生部门十分重要，因为控制疫情的最好时间是在疫情爆发之前。

B．早期模型的不足之处

1、首先模型并未给出拟合程度的参数，而当我们试图通过计算得到该模型的拟合程度参数时发现无法进行。原因是原模型求解过程的中间阶段参数K多次手工调整，而且模型中并未给出调整的标准和相关理论，所以我们很难重复该求解过程。由此我们得出结论：该模型的参数取值主观性太强，此用法给阅读者运用并改进模型带来了极大的困难，所以此模型的普适性较差。

2、在数据不足的情况下因无法进行手工调整，所以该模型用香港后期拟合的K值去预测北京后期疫情的发展趋势不太准确。但如问题分析所述，地域因素会造成不同地区的K值不同（如人口密度和人口流动大的城市若爆发传染病，初期

４

的K值会比人口密度和人口流动小的城市大，等等），而很难找到地域因素几乎相同的两城市。所以此作法可能导致预测结果相差较大。图1为用此方法预测的北京后期疫情情况与实际情况的对比图

3500300025002000150010005000010203040506070 图1 对该模型的评价：

从图中可以看出，预测值与真实值偏差越来越大。该模型具有较好的实际意义，能比较合理的反映非典的传播情况和发展趋势，模型中的参数设置也较科学，能较好地反映非典传播的影响因素，但模型的求解过程主观性太强，很难重现，而且参数的取值也很主观，没有取值的标准和理论。

综上所述，我们得到结论：该模型较合理，但不实用。

六、SARS传播模型的建立---Logistic SARS模型

1、模型建立及其对题中模型的优越性

题中模型只考虑了传染期限和传染率的问题，涉及的参数及考虑因素存在上述的不足。而实际情况中，的发病规律并不为我们所熟知，目前也没有治疗的有效方法，那么，以最原始的预防手段—隔离预防是最为有效的。而且，经实践证明，隔离防治是最为有效的。而且，经实践证明，隔离防治也确切在控制疾病的蔓延上起到至关重要的作用。于是，我们引入了预防效果指数k，用来反映疾病控制程度，它直接影响SARS流行趋势、发病时间、发病高峰出现时间累积发病人数。有因为SARS发病传染迅猛，为了描述这个特征，我们又引入了参数r，用来表示发病率。应用回归研究各地SARS的疫情资料，其流行趋势可用式1表示。 ?dN?rN?kN2? ?dt?N|t?t?N00?（1）

对式（1）求解得式（2）

1N? (1)

k1k?rt?(?)erN0r其中，Nmax?r/k为预期传染发病总人数，即理论上最多累积发病人数，

５

t?r/2k所对应的时间为发病高峰时间。

依据题中的数据并按所建模型拟合，可得表1中北京地区的参数估计值 N0 k r 地区 339 0.000063202 0.16152 北京 根据北京地区确诊病例累计，经过Matlab编程拟合，得到疫情走势及预测与真实值的比较，见图1

图1 北京地区疫情走势及预测与真实值的比较

为了检测我们的模型是否能较好的反映各地区实际情况、具有普遍意义，我们又以山西为例，对模型进行了验证。山西的疫情数据见附表2，参数估计值。（见表2）

表2 logistic回归对山西疫情的参数估计值

地区 山西 拟合（同上）结果见图2。

N0 416 k 0.000303971 r 0.13894

６

图2 山西区疫情走势及预测与真实值的比较

以上两个地区的模型拟合结果与实际非常接近。为了更进一步证明我们的模型对疫情发病情况的拟合程度，我们引入决定系数R作为检测标准。决定系数(R)=1-残差平方和。经过计算，两地区模型决定系数R均高于0.99，预测值与真实值非常接近，拟合优度检验无显著的差异，说明 回归模型较好的描述 的发病、流行情况，适合于发病拟合及流行特征研究。在模型中，r表示发病增加速率，r越大，疾病发展变化越快，反映是初级段发病人数增长速度快，高峰到来时间越早，持续时间越短，r越小，高峰到来时间越晚，持续时间越长。另外，从医学的角度来讲，对SARS采取相应的预防措施后，病例数有所下降说明预防措施与该疾病发生发展密切相关，在我们所建立的 Logistic 模型中引入了预防指数k，恰能更贴切的反映实际情况。而模型1并未对该因素给予考虑，这是它一个欠佳的方面。

1、 建立更优模型及困难所在

SARS是流行传染病，对于传统的流行数学模型通常假设：平均每个传染者在单位时间内可与βN个种群中的其他成员进行接触(其中N表示人群的总规模，β为接触率)。在对病愈后不具免疫力传染病模型的研究中，以前疾病的潜伏区都被忽视，假设易感染者一旦被感染立即变成了感染者，即为SIS模型。但是对于SARS来说，在易感染者被感染成为感染者之前，存在一段时间的潜伏期，为了掌握具有潜伏期的传染病的传播规律，建立一个真正能够预测以及能为预防和控制提供可靠、足够信息的模型，就应该进行全面考虑。SARS的传播基本可以用下图描述：

７

∏(1-p) 易感者S β ∏ 潜伏者E ∏p

γ1 k

隔离着Q Qμ γ2 μ1 感染者I

σ2 d2 σ1 d1

Rμ 死亡者D 恢复者R

其中，∏为初始人数，p为感染率，μ为治愈率，γ为隔离率，σ为隔离治愈率，d为死亡率，k为潜伏发病率。于是，可以建立如下模型——SEIRQ模型：

S(I??qE???Q?)?dS??(1-P)-?S??dtN??dE??p?S(I??qE???Q?)?(??k??)E1?dtN?dI???Ek?(?1?k??)I ?dt?dQ??1E??1I?(?1?k??)Q?dt??dR?dt??1I??2Q?R????N?N(t)?S(t)?E(t)?Q(t)?I(t)?R(t)由于这种模型是建立在对发病后调查分析的基础上，所要算的系数的精确度与调查的数据密切相关。但是对于数据的采集，通常只包括存活者，而对那些已死的病人，或对病程短，已经病愈的病例以及对轻型不典型病例或隐伏性的病例，我们通常很难调查。此外，某些病人在患病后，可能会改变原来的暴漏状况，如生活习惯的改变等。这样使病例对照研究后横断面研究所采取的病例类型，会与队列研究或实验所获得新病例不同。这就是现患病例—新病例偏奇。而SARS刚刚出现时，由于人们茫然，没有引起足够的重视，致使其在很短的时间内就扩散到全世界32个国家和地区。要统计较为准确的数据，因受人为、自然等因素的影响，其难度可想而知。并且各个地区采取的措施不同，人口流动不定，所以，其预测会受到限制。

８

3.、对卫生部所采取部分措施的评价 在SARS流行期间卫生部所采取的措施（见附录3）主要有:卫生部门控制人们之间的密切联系、控制传染期时间、引入反馈机制（如：政府强行措施）、加强疾病危险性的宣传教育、加强信息的透明等方面。为了能定量的评价这些措施的得力性，我们拟用小世界网络模型模拟卫生部针对SARS病毒的传播所采取的这些措施对疫情传播造成的影响。

（1）控制人口接触流动及隔离时间先后对SARS传播的影响

为了说明这两点，我们引入两个可调参数，在现实生情况他们分别对应W（表示人们之间联系的密切程度）和T（表示发现并隔离病原的速度）。可以预料W越大，T越延后，病毒就越容易传播：W越小，T越提前，病毒越难传播。用小世界模型拟合结果（如图1、2所示）也证明了这一点。其中，Ni为当天仍犯病人数，Nt为总患病人数。在图1中T=2，上图W=10，病毒就越难传播。用小世界模型模拟结果（如图1、如图2）也证明了这一点。其中，N 为当天仍患病人数。在图1中T=10，病毒传播自动衰减；下图W=20，病毒迅速传播。图2反映了发现并隔离病原的速度相差1天的发病变化趋势，如果拖后5天，其严重情况可想而知：若提前5天，则可使病情得到有效控制。

图1.T不变时，W对病毒传播的影响

９

图2 W不变时，T对病毒传播的影响

由此可见，不能及时发现控制病原和人们之间接触太多会非常有利于病毒的传播。初期出现病毒的爆发正式由于这两个原因，要控制病毒的蔓延应该从这两个方面入手。所以，卫生部所采取的限制人口流动，以及早发现、早诊断。早报告。早隔离、早治疗等有效措施都有效的控制W和T，从而使SARS的传播蔓延得到了有效控制。

(2)引入反馈机制后对SARS传播的影响

如果不引入其它机制，那么病毒的传播就只有两种结果，要么自动衰减，要么迅速蔓延，而实际情况中并非如此。在SARS传播过程中，人们的自觉性是一个渐变的过程，会随着疫情的变化而变化，是一个反馈过程，引入这个反馈过程同样可以减小T值和W值，从而达到抑制SARS的传播和扩散的效果（如图3所示）

图3 引入反馈机制后随时间的变化曲线

由上图我们可以看到，引入了反馈机制是人们自我隔离后，病毒的传播得到有效

１０

的控制。因此，卫生部采取的加大宣传提醒大家自觉地进行隔离、减少与别人的接触以及改变不良的卫生习惯等措施，都有效的防治了SARS传播蔓延。但是，实际中当人们发现当前患病人数N <100病且持续减小时,往往会放松警惕.而此时,SARS发病人数又会有所回升，这种趋势正如图4所示。

图4 人们会放松警惕时随时间的变化曲线

因此即使在疫情已经减轻的形势下，我们也一定不能麻痹大意，要贯彻好隔离制度，提高警惕性和自觉性，这样才能根本地战胜SARS。从这个意义讲，在实现了兵力零增长后，卫生部仍坚持通报病例的统计工作以及坚决不放松警惕的措施，都是十分正确的。

（3）信息透明度对病毒的影响 实际情况中，不是所有人都能及时获得疫情信息从而开始自我隔离的，例如在北京，知道4月20日公布了准确的患病人数才开始大规模采取自我隔离，因此这里就有一个信息透明度的问题。所以，我们也引入了一个叫信息透明度Ti=来表征这种情况，Ti的意义是知道疫情的情况从而进行自我隔离的人占总数人数的比

例。从图5我们可以看到透明度对病毒的传播也有重要的影响。

１１

图5 在有反馈的情况下不同的信息透明度对病毒传播的影响 很明显，当透明度比较高时，疫情消失需要的时间比较短，高峰期患病人数较少。因此，卫生部每天通报疫情，让人们尽早地了解，从而做好预防措施，也是控制疫情的有效方法。

七、模型的特点 1、模型的优点

首先，此模型引入了参数k（预防指数）并运用Logistic阻滞增长模型来拟合，比较符合实际情况以及相关医学解释，也即用医学与数学相结合的知识阐明了SARS发生发展规律，拟合结果与实际流行趋势贴近。

其次，用t 和N 能够预测SARS高峰期的到来时间，可能累计最大发病数，这样，人们就可以以此为参考，认为第改变一些参数或控制一些相关因素，从而达到预防疾病的传播与蔓延的效果，在现实中有实用性。如，加强消毒，控制人口流动均可以增大K，从而可以使实际的高峰期累计病例数减低。

第三，结果用图形拟合说明、验证，简介直观。 2、模型的不足

尽管我们的预测结果已相当好，但仍有一些不足。当医疗条件变化时，治疗可能成为最有效的手段，而且医源性传染的可能性也会大大降低，那么，k的意义就成为最有效的手段，而且，根据现有数据所拟合的模型都具有滞后性，原始数据的精准度也会影响到模型的效果，而下一代的SARS的模型或许应该更加复杂，涉及到更详细的空间和随机过程、作为病毒源头的动物、季节因素和多种传染模式，那么，Logistic SARS模型的应用就受到了限制。因此，它应该随条件而逐步改进。

3、模型的改进与推广

如果能够将传染源、接触率等因素化为参数考虑进去，可以进一步完善模型，使之提供潜伏期、阀值、恢复期、最多感染系数、平均易感人数等参数。并且，若将功能基因组的分析工作扩张到冠状病毒以外与SARS相关的病原体中，则对该模型的分析还可以为疫苗和新型药物的研制提供理论支持。

八、SARS对旅游业的影响

由于SARS疫情的爆发，导致期间北京市接待海外旅游人数锐减，造成了极大的经济损失。现建立模型估计此经济损失。 建模思想

本文首先通过对北京市1997～2002年各月接待海外旅游人数的分析建立起时间序列模型，“预测”出2003年3～8月本应接待的人数，然后用3～8月本应接待的人数减去实际接待人数得到“损失人数”，最后通过分析“损失人数”大致推测出旅游业恢复时间和因疫情导致的经济损失。模型即

经济损失=（本应接待的总人数-实际接待总人数）×人均旅游消费 即 ?M?(S?S)?p

??

１２

???M为经济损失，S为本应接待的总人数，S为实际接待人数，p为人均消费。

（一）对本应接待人数的分析 1、模型的建立

首先令｛X(t)，t=1,2,3?｝为北京市1997年后各月接待海外旅游人数的随机序列（例如X(1)为1997年一月接待人数的随机变量，X(13)为1998年一月接待人数的随机变量，以后类推）。

然后分析1997～2002接待人数样本值的走势，发现X(t)具有明显的周期性（周期为12个月），其数学期望也应具有周期性，所以不能认为该序列是平稳时间序列。但可作如下处理，取

X(t)=F(t)+W(t)

其中F(t)表示X(t)中随时间变化的数学期望，此问题中可用周期函数来描述。而E(W(t))?E(X(t))?E(F(t))?0且E(W(t)W(t??t))与t不相关，则W(t)可以用平稳时间序列来拟合。

而数学期望为0且具有有理谱密度的平稳时间序列必有下三种形式： 1）．自回归模型AR(p)（以下记W(t)为Wt）

Wt??1?Wt?1??2?Wt?2??????p?Wt?p?at (t?1,2,3???)

其中at是白噪声，常数p为阶数 ,常系数?i为参数 2）．滑动平均模型MA(q)

Wt?at??1at?1??2a2?????qat?q (t?1,2,3???) 常数q为阶数 ,常系数?i为参数

3）．混合模型ARMA(p,q)

Wt??1?Wt?1??2?Wt?2??????p?Wt?p?at??1at?1??2a2?????qat?q 为了确定模型的类别和阶数，首先给出自相关函数?k和偏相关函数?kk的计算公式。

?k由下式得出：

?k??k （其中?k?E(WtWt?k)） ?0?1????k?1???k1???1???????1??????????k2???2?

???????2???????????????????????1???kk???k??kk由下述方程得出

?1???1???????k?1

于是可用?k和?kk描述时间序列各模型的性质（如下表）

AR(p) MA(q) ARMA(p,q) ?k 拖尾 截尾k=p处 截尾k=q处 拖尾 拖尾 拖尾 ?kk

如果确定了F(t) 和W(t)的解析式，就可以用其来描述X(t)；

１３

2、模型的求解

现有北京市1997～2002年各月接待人数的样本值X(t),(t?1,2???72)

1）．周期函数F(t)的确定

通过对样本值的观察，我们发现样本值具有明显的周期性，周期大致为12（月），且样本值呈逐年递增的趋势，所以我们采用傅立叶级数与线性递增函数

?2n?2n?和叠加作为周期函数F(t)。型如 F(t)??(ancos(t)?bnsin(t))?c?t?d

TTn?1通过曲线拟合确定系数得

4?2?4?2?F(t)??2.3?cos(?t)?2.7?cos(?t)?2.6?sin(?t)?3.5sin(?t)?0.2?t?16.7TTTT（2-1）

拟合图形如图8（其中实线为拟合曲线，虚线为原曲线） 35?3025201510501020304050607080 图8 拟合样本值的周期函数图形 2）．平稳时间序列来拟合W(t) 平稳时间序列样本值W(t)?X(t)?F(t),(t?1,2,???72) 首先利用样本值求得样本自协方差函数

72?k??? 72进一步可得样本自相关函数

?k???W?Wjj?1???j?k,k?0,1,2???K(K?72)，

?0而样本偏相关函数?kk可利用下面的递推公式计算：

?k???k,k?0,1,2,???K(K?72)

１４

???????111???1???kk???????????????????1?????k?1,k?1?k?1?k?1?jkj???jkj??

j?1???j?1?????????????k?1,j??kj??k?1,k?1??k,k?(j?1),j?1,2,???k?计算后画出?k和?kk的曲线图（见图9，图10） 0.30.50.20.40.30.20.100.10-0.1-0.1-0.2-0.3-0.2-0.3-0.4-0.5-0.40102030405060708002468101214161820 图9 ?k的趋势图 图10 ?kk的趋势图

经观察我们发现?k明显拖尾，?kk截尾于k=2，所以我们采用模型AR(2) 即 Wt??1?Wt?1??2?Wt?2?at （3-1） 经过计算系数 ?1??2a????????1(1??2)1????21???-0.1787 ?2??????2??1????21?21?0.08 4 4

白噪声at的方差???0?(1??1??1??2??2)?7.6714 ，均值为0；

通过公式（3-1）可递推预测未来的Wt(t?73,74,???)例如：

W73??1?W72??2?W71?at，W71,W72为已有样本值，at为正态随机数（预测时不

???????????计其影响）。

3）．模型的结果

,(t?73,74,75) （4-1） 预测公式X(t)?F(t)?W(t)?^^其中F(t)由式（2-1）算出，W(t)由公式（3-1）递推算出。 通过对已有数据的检验，得到平均预报误差小于10%，可认为预报结果较为准确。然后预测得2003年3～8月接待量（假设3月前旅游业未受影响）为27.1 ，30.2 ， 30.4 ，29.3 ，29.7 ，32.4（单位：万人）

同时，我们观察每年同一个月的情况，发现相同月份的值有一个逐渐增大的趋势。所以我们对3到8月的数据分别进行了线性拟合，得到2003年3到8月的游客接待量分别为23.6 ，30.5 ，32.2 ，30.5 ，27.3 ，33.6（单位：万人），发现两种方法预测数据相差不大，最大相差为3.5（万人），平均相差1.7（万人）

１５

（取绝对值的均值），证明预测结果很好。

（二）对实际接待人数的分析

用“预测”所得的3～8月应接待人数减去实际接待人数，就得这6各月的“损失人数”，经过分析发现，“损失人数”的走势与传播模型中病例日增量的走势基本相同，但有大约有一个月的延迟。即病例日增量在大约5月初达到峰值，而“损失人数”在6月才达到峰值。这是因为SARS对其它行业的影响必有一定的滞后性。观察3～8月“损失人数”的散点图和病例日增量图，我们认为用二次曲线拟合“损失人数”较为适合。拟合函数为

???2.9?x2?22.2?x?14.7

如图11：

3020100-10-20-30 图11 “损失人数”拟合函数图（横坐标为月份，纵坐标为“损失人数”） 分析此图我们发现，大约到9月初“损失人数”才会降到0，即表示大约9月初疫情对旅游业的影响才会消失，恢复到预计水平。

（三）结果分析

3～8月总共“损失人数”大致为114.6万人(三月前旅游业未受影响，9月后旅游业恢复)。假定人均在京消费1035美元（2001年统计数据）[5]，计算得疫情期间总共损失约1.2亿美元。

2345678910九、控制传染病，预测靠模型

传染病是由病原体传染而引起的疾病，其传播过程有两个影响因素：社会因素、自然因素；需要三个条件：传染源、传染途径、易感人群。传染病曾经给人们带来巨大的灾难，天花、霍乱等传染病曾吞噬过无数宝贵的生命。最近的一次大规模的传染病是2002年末到2003年初在中国乃至世界爆发的SARS（Severe Acute Respiratory Syndrome，严重急性呼吸道综合症，俗称非典型肺炎）。SARS的爆发和蔓延给我国的经济发展和人民生活带来了很大影响。

古往今来，无数的科学工作者都致力于研究传染病的数学模型，来预测传染病的爆发时间、高峰期、发展趋势，为预防和控制它的传播提供可靠、足够的信息，对采取的措施提出合理的建议。

现在，有人对人们关注北京的SARS问题建立模型进行分析，计算并预测出了SARS的发展趋势，说明如果加大控制力度，对非典病人严格隔离，不久后病例的日增长数将开始下降，SARS将1个多月后日增量将逐渐减小到0，否定了当时有的人认为SARS会无限制传播的错误观点。预测结果如下图所示：

１６

日病例增长随时间变化的图

上图是以3月1日为起点，该图只反映大致趋势。从图中可以看出，该模型预测出SARS的日增病例快速上升的时间为4月15日，高峰期在4月30日左右，5月25日以后累计病例基本稳定到2500人左右，增长缓慢，非典的威胁已经逐渐消除。相信这对北京政府采取的措施有很大的帮助！北京政府是在4月20日采取措施对非典进行控制的，此时似乎有点晚了！因为在相同的控制力度下，如果提前5天采取措施，累计病例将控制在2000人以下；但值得庆幸的是，如果再延后5天，累计病历将至少达到3000多人，甚至可能超过4000人。

有人还就非典对北京旅游业的影响也建立了模型，计算知，非典期间北京的游客减少了约115万人，经济损失约1.2亿美元，约为正常情况下全年收入的33%，在经济上是个巨大的损失！我们有理由相信，如果北京在SARS刚出现时就采取相应的措施，现在这样的“灾难”完全可以避免；如果传染病的数学模型能更早建立出来，相信北京政府采取的措施会更加理智。

通过对北京SARS的分析，我们可以清楚地看到传染病数学模型在对传染病传播的预防、预测和控制方面有很大的贡献；而且，建立传染病的数学模型，也为我们以后预测和处理类似的情况积累了经验。

十、参考文献

[1] 郭东升 、吴咏时，《疫情演化方程和对SARS疫情的预测》， http://www.pkunew.pku.edu.cn ，2003.9.24 [2] 杨位钦，顾岚，《时间序列分析与动态数据建摸》，北京工业出版社1986.12 [3] 《奥运旅游瞄准》，http://www.snweb.com/gb/market_daily，2003.9.24

[4] 林国基等，用小世界网络模型研究SARS病毒的传播，http://162.105.8.101/sars/wangluomoxing.htm 2003年9月

附件1：北京市疫情的数据

日 期

4月20日

已确诊病例累计 现有疑似病例 死亡累计 治愈出院累计

339

402

18

33

１７

4月21日 4月22日 4月23日 4月24日 4月25日 4月26日 4月27日 4月28日 4月29日 4月30日 5月01日 5月02日 5月03日 5月04日 5月05日 5月06日 5月07日 5月08日 5月09日 5月10日 5月11日 5月12日 5月13日 5月14日 5月15日 5月16日 5月17日 5月18日 5月19日 5月20日 5月21日 5月22日 5月23日 5月24日 5月25日 5月26日 5月27日 482 588 693 774 877 988 1114 1199 1347 1440 1553 1636 1741 1803 1897 1960 2049 2136 2177 2227 2265 2304 2347 2370 2388 2405 2420 2434 2437 2444 2444 2456 2465 2490 2499 2504 2512 610 666 782 863 954 1093 1255 1275 1358 1408 1415 1468 1493 1537 1510 1523 1514 1486 1425 1397 1411 1378 1338 1308 1317 1265 1250 1250 1249 1225 1221 1205 1179 1134 1105 1069 1005 25 28 35 39 42 48 56 59 66 75 82 91 96 100 103 107 110 112 114 116 120 129 134 139 140 141 145 147 150 154 156 158 160 163 167 168 172 43 46 55 64 73 76 78 78 83 90 100 109 115 118 121 134 141 152 168 175 186 208 244 252 257 273 307 332 349 395 447 528 582 667 704 747 828

１８

5月28日 5月29日 5月30日 5月31日 6月01日 6月02日 6月03日 6月04日 6月05日 6月06日 6月07日 6月08日 6月09日 6月10日 6月11日 6月12日 6月13日 6月14日 6月15日 6月16日 6月17日 6月18日 6月19日 6月20日 6月21日 6月22日 6月23日 2514 2517 2520 2521 2522 2522 2522 2522 2522 2522 2523 2522 2522 2522 2523 2523 2522 2522 2522 2521 2521 2521 2521 2521 2521 2521 2521 941 803 760 747 739 734 724 718 716 713 668 550 451 351 257 155 71 4 3 3 5 4 3 3 2 2 2 175 176 177 181 181 181 181 181 181 183 183 184 184 186 186 187 187 189 189 190 190 191 191 191 191 191 191 866 928 1006 1087 1124 1157 1189 1263 1321 1403 1446 1543 1653 1747 1821 1876 1944 1994 2015 2053 2120 2154 2171 2189 2231 2257 2277 附件2山西SARS疫情发展走势（每日增量）

日期 5月11日 5月12日 5月13日 5月14日 5月15日 5月16日 5月17日 5月18日 5月19日

现患人数 311 307 293 274 263 227 220 207 200 累计确诊报告病例 416 420 431 435 442 443 442 445 444 累计治愈出院人数 87 95 119 142 160 197 203 218 232 １９

5月20日 5月21日 5月22日 5月23日 5月24日 5月25日 5月26日 5月27日 5月28日 5月29日 5月30日 5月31日 6月1日 6月2日 6月3日 6月4日 6月5日 6月6日 6月7日 6月8日 6月9日 6月10日 6月11日 6月12日 6月13日 6月14日 6月15日 6月16日 6月17日 6月18日 6月19日 6月20日 6月21日 6月22日 6月23日 6月24日 6月25日 6月26日 6月27日 6月28日 6月29日 6月30日 193 183 180 167 154 144 140 126 113 104 94 91 88 85 79 73 54 51 50 46 33 28 22 21 17 14 12 12 11 4 4 4 4 4 4 1 0 0 0 0 0 0 445 445 447 447 450 449 450 450 450 450 450 450 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 448 232 242 247 259 275 284 289 303 316 325 335 338 342 339 348 354 371 374 375 379 385 391 396 402 403 407 410 412 412 412 419 419 419 419 419 422 423 423 423 423 423 423

２０

7月1日 7月2日 7月3日 7月4日0 7月5日 7月6日 7月7日 7月8日 7月9日 7月10日 7月11日 7月12日 7月13日 7月14日 0 0 0 0 0 0 0 0 0 0 0 0 0 0 448 448 448 448 448 448 448 448 448 448 448 448 448 448 423 423 423 423 423 423 423 423 423 423 423 423 423 423 附件3

中国政府所采取的9项控制措施：第一，建立一个跨部门的合作和协作机制，在各政府部门都确立了这样的机制。 第二，建立全国性的SARS监测和报告体系。第三，采取全面的防治措施来管理SARS感染病人。主要的重点就是严格控制医院内的感染和社区的传播。第四，就是加强医疗护理，并且尽可能地降低死亡率。第五，制定防治方案，组织培训项目，对于医疗人员进行培训。第六，监管和提供现场的支持。第七是在社区当中开展SARS的防治和健康教育。第八，对重点问题进行研究，也就是进行科研攻关。第九，交流信息、加强合作。

附件4：问题3的程序

%程序1 用傅立叶级数进行拟合 %问题3的数据

y0=[9.4 11.3 16.8 19.8 20.3 18.8 20.9 24.9 24.7 24.3 19.4 18.6 9.6 11.7 15.8 19.9 19.5 17.8 17.8 23.3 21.4 24.5 20.1 15.9 10.1 12.9 17.7 21.0 21.0 20.4 21.9 25.8 29.3 29.8 23.6 16.5 11.4 26.0 19.6 25.9 27.6 24.3 23.0 27.8 27.3 28.5 32.8 18.5 11.5 26.4 20.4 26.1 28.9 28.0 25.2 30.8 28.7 28.1 22.2 20.7 13.7 29.7 23.1 28.9 29.0 27.4 26.0 32.2 31.4 32.6 29.2 22.9];

%进行傅立叶级数拟合 w=2*pi/12; for t=1:72

ac(t)=cos(2*w*t);ac1(t)=cos(w*t);

２１

bc(t)=sin(2*w*t);bc1(t)=sin(w*t); end

x=[ac;ac1;bc;bc1;1:72;ones(1,72)]; b=x'\\y0';

%拟合得到的表达式的系数

%b=-2.2503 -2.7344 -2.6120 -3.5255 0.1545 16.7552

%计算出拟合出的傅立叶级数的值 t=1:72;

y=b(1)*cos(2*w*t)+b(2)*cos(w*t)+b(3)*sin(2*w*t)+b(4)*sin(w*t)+b(5)*t+b(6)*ones(1,72);

%计算出拟合出的傅立叶级数的值并计算出未来8个月的值 tt=1:80;

yy1=b(1)*cos(2*w*tt)+b(2)*cos(w*tt)+b(3)*sin(2*w*tt)+b(4)*sin(w*tt)+b(5)*tt+b(6)*ones(1,80); yy=yy1(73:80);

%程序2 用时间序列拟合

%计算真实值与傅立叶级数计算出的值的差 z=y0-y;

%计算模型的参数的值 r0=sum(z.^2)/72; t1=0;

for i=1:71

t1=t1+z(i)*z(i+1); end

r1=t1/72; t2=0;

for i=1:70

t2=t2+z(i)*z(i+2); end

r2=t2/72;

p1=r1/r0;p2=r2/r0;

phi1=p1*(1-p2)/(1-p1^2); phi2=p2*(1-p1)/(1-p1^2);

%计算方差，但此处数据较少，方差值可能不准确 segma=sqrt(r0*(1-phi1*p1-phi2*p2));

%用题目中的数据检验预测的准确性（预测效果比较好）

２２

for j=1:50 for i=3:72

d(i)=phi1*z(i-1)+phi2*z(i-2); end

y1=d+y; y0-y1;

me(j)=mean(d(3:end));

mea(j)=mean(abs(d(3:end))); end

mean(me) mean(mea)

%程序3 模式识别 r0=sum(z.^2)/72; t=zeros(1,71); for j=1:71 for i=1:72-j

t(j)=t(j)+z(i)*z(i+j); end

r(j)=t(j)/72; end

p=r/r0;

phi(1,1)=p(1); for i=2:50

for j=1:i-1

pp1(i)=0;pp2(i)=0; for jj=1:i-1

pp1(i)=pp1(i)+p(i-jj)*phi(i-1,jj); pp2(i)=pp2(i)+p(jj)*phi(i-1,jj); end

phi(i,i)=(p(i)-pp1(i))/(1-pp2(i));

phi(i,j)=phi(i-1,j)-phi(i,i)*phi(i-1,i-j); end end

%程序4 预测2003年的数据 %真实值

yn=[15.4 17.1 23.5 11.6 1.78 2.61 8.8 16.2];

d(1:2)=z(end-1:end); for i=3:10

d(i)=phi1*d(i-1)+phi2*d(i-2); end

２３

dd=d(3:10);

y1=yy+dd y1-yn;

２４