中国工程科学EngineeringScienceSep.2004Vol16No19
研究报告
SARS传染病数学建模及预测预防控制机理研究
刘云忠,宣慧玉,林国玺
(西安交通大学管理学院,西安 710049)
[摘要] 首先,利用疾病传播的一般规律及人口守恒统计法则建立了四类人的SARS传染病数学模型,然后
运用数学方法对四类人的SARS传染病数学模型进行分析,得出了其生理意义和预防、控制机理。其次,利用人工神经网络理论建立了SARS的预测模型,以北京市的SARS数据为例进行了预测和分析,预测结果显示该模型简单易行,预测精度高。
[关键词] SARS;稳定性;轨线;阈值定理;人工神经网络[中图分类号]O17511 [文献标识码]A [文章编号]1009-1742(2004)09-0060-06
1 引言
传染病是当今世界最严重的问题之一。2003
年春天,非典型肺炎(SARS)传染病突然侵袭了大半个中国,由于人们对此种传染病的传播机理还不太清楚,因而一度引起了人们心理上的恐慌。现在,预防和控制SARS的研究显得极其紧迫。目前,对传染病的研究有4种方法:描述性研究,分析性研究,实验性研究和理论性研究。在理论性研究中,数学模型起着极其重要的作用。它把传染病的主要特征通过假设、参数、变量和它们之间的联系清晰地揭示出来。数学模型的分析结果能提供许多强有力的理论基础和概念。用数学模型帮助发现传染病的传播机理,预测传染病的流行趋势已成为共识。
利用非线性动力学的方法建立传染病的数学模型,来研究传染病是否会蔓延持续下去以及是否终将会被消灭具有重要的现实意义。因为这有助于人们对传染病的发展趋势进行预测,为人们预防和治疗传染病提供有益的信息和有效的措施。Kermack和Mckendrick[1~3]的论文是传染病数学模型的基石,他们首先利用非线性动力学的方法建立了传染
[收稿日期] 2003-10-28;修回日期 2004-02-26
病的数学模型,即所谓的KM模型,之后,
Cooke[4],Hethcote[5]等做了大量的工作。陈兰荪、陈健[6]详细阐述了传染病的建模思想和研究方法。笔者利用非线性动力学的方法建立了非典型肺炎(SARS)传染病四类人的数学模型,并探索了预防和控制SARS的机理。然后,利用人工神经网络理论建立了SARS的预测模型,以北京市的SARS数据为例进行了预测和分析,预测结果显示,该模型简单易行,预测精度高。
2 四类人模型[7]
SARS传染病四类人数学模型如下:
把人口分为健康人S(t)、SARS病人I(t)、病愈免疫(包括死亡)的人R(t)及SARS疑似病人P(t)四类人。
疾病传播一般服从下列法则:法则1 在所考虑的时期内,人口总数保持在固定水平N;
法则2 易受传染者S(t)人数的变化率正比于传染病患者I(t)与S(t)人数的乘积;
法则3 由I(t)向R(t)转变的速率与I(t)成正比。
[作者简介] 刘云忠(1970-),男,江西新干县人,西安交通大学博士研究生
© 1994-2006 China Academic Journal Electronic Publishing House. All rights reserved. http://www.cnki.net
第9期刘云忠等:SARS传染病数学建模及预测预防控制机理研究61
由上述疾病传播法则,不难得出SARS传染病四类人的数学模型为
dS/dt=-λISdI/dt=λIS-αI
(1)
αdR/dt=IdP/dt=βI
方程组式(1)是四维的。且初始状态为
S(0)=S0>0,I(0)=I0>0,R(0)=R0>0,P(0)=P0=0。其中常数λ,α称为传染率和移除
其中,δ>0。以下对模型式(2)进行分析,在
(S,I)相平面上考察轨线。首先由dS/dt=-λIS-δS=-(λIS+δS)<0,得到的S(t)是单调减少的,所以,对于所有t>0,有S(t) t>0成立,此时I(t)单调减少。当S0>1/σ时由S(t)单调减少可得到唯一的t1,使S(t1)=1/σ,因此有0 )为最大值。1/σ时,I(1/σβ>0。率,其值均大于零。 令σ=αλ,1σ=αλ称为相对移除率。为了讨论问题的方便,假设总体N=1。 定理1 (阈值定理)设S(t),I(t),R(t),P(t)是初值问题式(1)的解。如果σS0<1,那么,当t→+∞时,I(t)单调减少趋于零;如果σS0>1,那么,当t→+∞时,I(t)先增加达到最大值1-1σ-1σln(σS0),此时S=1σ,而后单调减少趋于零。S(t)是一个单调减少函数,并 ),是方程1-S+且其极限limS(t)=S(+∞ t→+∞ 显然方程组式(2)的轨线方程为I(t)+lnI(t)I(0) =1-S+ σ 1ln S。S0 由此得I σ1=1- σ1- 1- σ - 1ln(σS0)-ln (σ) σlnS0。 )I(1/σI0 < σ 11(lnεε)内的根。当t→+∞0)σ=0在(0,1/σ 时,R(t)→0,P(t)→0(见图1)。 上式说明采取预防措施后,可以减少得病人数,并 且I(t)的最大值小于不采取预防措施时的最大值。 令 D={(S,I)0≤S≤1,0对于方程组式(2),其轨线为: -αt S=0,I=I0e, I=0,S=S0e -δt 。 O(0,0,0,0)点是式(2)在D上唯一的平衡 图1 S(t),I(t),R(t),P(t)的变化趋势 Fig11 ThechangingtrendofS(t), I(t),R(t),P(t) 点,并且点(0,0,0,0)是局部渐近稳定的,这是因为特征根-δ及-α均小于零。 对于直线S+I=1上的所有解均有 d(S+I)dt=-αI-δS<0。 所以,在D内出发的轨线不会越出区域D。令D={(S,I)0≤S≤1,0综上所述,得出如下定理: 定理2 对于初值问题式(2),区域D是平衡点O(0,0,0,0)的渐近稳定区域(见图1)。 定理3 (阈值定理)设S(t),I(t),R(t),P(t)是初值问题式(2)的解。如果σS0<1,那 由于人类对传染病的认识提高以及现代医学水平的发展,对于许多传染病可以做到提前预防,使人群对许多种传染病具有免疫能力,例如打预防针、进行免疫接种等。在不久的将来人类将找到 SARS的免疫疫苗,因而在模型式(1)的基础上考虑这些因素,经过调整,得到如下的数学模型: dS/dt=-λIS-δS dI/dt=λIS-αIdR/dt=αIdP/dt=βI (2) © 1994-2006 China Academic Journal Electronic Publishing House. All rights reserved. http://www.cnki.net 62中国工程科学第6卷 么,当t→+∞时,I(t)单调减少趋于零;如果σS0>1,那么,当t→+∞时,I(t)先增加到达最大值I(1σ),而后单调减少趋于零。S(t)同时单调减少趋于零。当t→+∞时,R(t)→0,P(t)→0 通过建立四类人的SARS传染病数学模型并对其进行分析的基础上,可得到如下结论: 1)在不考虑自然出生和死亡的前提下,SARS传染病发生时,如果易感染人的总数小于等于该病的相对移除率,SARS不可能发生流行,将很快被消灭。如果易感染人的总数大于该病的相对移除率,SARS可能发生流行,得病的人数将猛增,当 )时,得病人数I(t)达易感染人数下降到S=(1σ)ln(σ到最大值1-1σ-(1σS0),而后得病人数逐渐减少,最终趋向于零,即SARS被消灭。在整个过程中易感染人数单调减少,最终并不是所有的易感染人都会得病。因此,疾病不是因为缺少受传染者而停止传播,而是因为没有了传染者才停止传播。 2)在SARS传染病流行之前,对易感染的人群进行有效的预防可以使易感染的人数下降,从而达到防止SARS传染病流行的目的。在SARS传染病发生之后,立即对易感染的人群进行隔离,同样可以使易感染的人数下降,从而减少得病人数。此种情况下,如果在发病初期易感染的人数S0≤(1σ),那么SARS会很快被消灭。如果在发病初期易感染人数S0>(1σ),那么得病人数先增加, )后,得病人数逐渐减少而当其达到最大值I(1σ)小后SARS被消灭。此种情况下的最大值I(1σ)ln(σ于不预防时的最大值1-1σ-(1σS0)。 和结构关系模型时,要求被预测对象必须满足预测 模型的前提条件,否则预测结果不可靠;而因果关系模型实际上是因和果之间的映射关系,几乎可以表达所有非线性关系(如人工神经网络),因而比以上二类函数的适用范围广。人工神经网络(ANN,artificialneuralnetwork)理论是在综合了众多学科理论的基础上发展起来的,其发展在预测方面为各领域的应用带来了广阔前景,特别是地质、天文、机械等理工科领域。近年来,人们才逐渐将ANN技术用于生物、医学、药学、化学等学科,但主要是应用ANN进行分类或进行静态分析,例如用蛋白质一级结构预测二级结构、药物筛选、恶性肿瘤的生存分析、疾病的诊断及预后分析等[11]。笔者拟用误差反向传播(BP,back2propagation)神经网络对北京市SARS的发病率进行分析,建立SARS发病率的ANN预测模型,旨在探讨ANN预测模型在疾病发病率或死亡率预测上的应用前景。311 三层BP神经网络的结构 BP神经网络模型是人工神经网络中应用最广泛的一种。它的功能函数为S型函数,神经元连接形式为前馈神经网络,学习方式为有监督学习。如图2和图3所示,整个学习过程的具体步骤为: 3 利用人工神经网络理论建立SARS 预测模型[11] 目前,用于预测的定量化方法尽管很多,但基本上可归纳为时间关系、结构关系和因果关系模型等三类[11]。时间关系模型的特点为被预测对象是时间变量的函数,它包括两种:一是数理统计学中的时间序列分析模型,二是用时间函数(如多项式,正余弦等)表示的趋势外推模型。结构关系模型的特点是在一定时间内被预测事件与其影响因素之间保持着某种固定的结构函数关系(如回归模型、联立方程模型、动态模型等)。因果关系模型的特点是用因果关系表达被预测事件与其影响因素之间的相互作用。一般来说,应用时间关系模型 图2 BP网络模型结构Fig12 StructureofBPANNmodelStep1初始化网络及学习参数,如设置网络初 始权矩阵、参数α等; Step2提供训练模式、训练网络、直到满足学习要求; Step3前向传播过程:对给定训练模式输入, 计算网络的输出模式,并与期望模式比较,如有误差,则执行Step4,否则返回Step2; © 1994-2006 China Academic Journal Electronic Publishing House. All rights reserved. http://www.cnki.net 第9期刘云忠等:SARS传染病数学建模及预测预防控制机理研究63 误差的变化改变学习因子η,设n为学习次数,有 (1-0.05)η En+1≤En (6)n= (1+0.05)η En+1>En BP人工神经网络模型结构和建立步骤见图2和图3。312 基于BP网络的非线性时间序列预测原理已知时间序列{Xii=1,2,…,t-1},若用过去的N(N≥1)个时刻的数值预测未来M(M≥1)个时刻的数值时,可将训练数据分为K段,长度为N+M的有一定重叠的数据段,每一段的前N个数据作为网络的输入,后M个数据作为网络的输出(见表1)。表1 训练数据的分段方法Table1 Themethodofdividingsection fortrainingdata N个输入X1,X2,…,XNX2,X3,…,XN+1X3,X4,…,XN+2 M个输出 XN+1,XN+2,…,XN+MXN+2,XN+3,…,XN+M+1XN+3,XN+4,…,XN+M+2 图3 BP神经网络模型建立步骤Fig13 StepsofSettingupBPANNModelStep4向后传播过程:计算同一层单元的误差ej;修正权值和阈值,返回Step2,直至全部m k 个学习模式对训练完毕。 各个连接权的调整量是分别与各个学习模式对的误差函数Ek成比例变化的,该方法称为标准误差反向传播算法(又称平均误差),网络的全局误差为 E= k=1 … XK,XK+1,…,XN+ K-1 … XN+K,XN+ K+1, …,XN+ K+N-1 然后应用前述的BP神经网络进行训练学习,寻找一个RN到RM的函数关系。该网络的输入层有 N(N≥1)个节点,输出层有M(M≥1)个节点。 6 m Ek= k=1t=1 66 mq ( kyt -ct)/2 2 (3) BP学习算法实质上是最小均方(LMS)算法的推广,是一种非线性梯度优化算法,因此不可避免地存在局部极小值问题。学习算法的收敛速度很慢,通常需要上千次迭代或更多。为此,采用了动态学习比率BP算法。在BP算法中,连接权值的改变规则为: Δγt=α・dk(4)t t=1,2,…,q Δθj=β・ k ej 其预测值的一般形式为: Xt=f1(X)Xt-1+f2(X)Xt-2+…+ fp(X)Xt-p +εt(7) 其中p为预测网络的输入节点数。 fi(X)=fi(Xt-1,Xt-2,…,Xt-p,εt),(i=1, 2,…,p)是一个以预测网络各输入变量为自变量的 j=1,2,…,p(5) 非线性函数。预测误差为 2 σ= 2 神经网络的训练过程就是使σ达到全局最小 其中,α,β为定义的权值步长,即学习比率,以 下均用η表示。若η选择的很小,网络的学习速度不但很慢,而且很容易使网络陷入局部最小点;而η选择的太大,又容易使网络出现振荡,使误差E始终不能达到极小点。因此,对BP网络的学习算法做以下改进: 赋予η相对较大的初始值,使网络在学习初期误差较大,以较大的步长逼近极值点,同时又容易跳开局部极小点。随着网络的学习,根据学习 6 (X^t-Xt)2 (8) 的过程,也是预测模型的建立过程。从统计学观 点看,神经网络建模的最终结果是求出非线性函数 fi(X),而fi(X)是用网络中各节点的连接权阈值 表示的。而用传统的时间序列AR(p)模型进行预测,其预测值的一般形式为 Xt=C1Xt-1+C2Xt-2+…+CpXt-p+εt(9)其中{εt}为白噪声,Ci为常系数。 © 1994-2006 China Academic Journal Electronic Publishing House. All rights reserved. http://www.cnki.net 64 313 预测模型参数的确定 中国工程科学第6卷 利用神经网络来建立SARS发病率预测模型时需要确定神经网络BP训练中的模型参数。模型参数包括模型层数、隐含层节点数和功能函数。模型层数的确定尚无统一方法。通过试验对2个隐含层的ANN与1个隐含层的ANN进行比较,前者易陷入局部极小,误差在局部极小处振荡,在同样学习次数的训练中,二者所达到的误差相近,因而选用3层BP神经网络模型。节点数包括输入输出层节点数和隐含层节点数,输入节点一般为输入变量的个数,输出层节点为SARS每天发病数。隐含层节点数目的确定难度较大,若隐含层节点太少,网络可能不能训练;若网络节点数太多,网络训练可能无穷无尽,不易收敛。在选取隐含层节点数时,根据最小上界(LUB,leastupperbound)理论[12],经大量试验取隐含层节点数为6。功能函数的确定:BP神经网络的功能函数要求处处可微分且收敛。目前的研究多用Sigmoid函数f1(x)=1(1+e-x),极少使用双曲正切函数f2(x)=(ex -e-x) (ex+e-x)作为网络功能函数。 间,由于北京市政府没有采取措施,SARS在社会上传播很快,北京市采取各种隔离和防疫措施后,疾病传播速度呈减缓趋势。研究表明,防疫措施当时就会起作用,隔离措施的效果在其后的第一、第二个周期并不明显,而在第三个周期会显现出显著作用。由于5月1日至5日大多数单位放假,外出人数明显减少,人员流动性大大降低,相当于隔离措施的进一步加强,从而使该期间的可跟踪概率变大,相当大的一部分传染源能够发现,并被隔离,这就使得5月11日至15日期间的新患病人数显著下降;5月6日上班后,人员的社会流动性增加,使可跟踪概率下降,因此,可能使5月16日至20日间的新患病人数有一定程度的反弹,或者说下降的势头与5月11日至15日相比将会减弱。如果5月16日以后,新发病人数没有反弹,说明主要的感染途径已经能够掌握,并实施了较有效的措施,北京市的SARS疫情已得到有效控制。因此,5月16日至20日是北京市SARS疫情的分水 314 预测结果 预测结果见表2。从预测结果可以看出,所建 区警报要求的数字之前,需要保持防范措施力度,立的人工神经网络模型较好地预测了SARS感染人 并进一步加强公众防范意识[13]。数,精度较高。事实上,在3月1日到4月18日 表2 SARS感染人数实际值[14]和预测值及相对误差 Table2 TherealandforecastingvaluewhichpeopleinfectSARSandtherelativeerror 日期 实际 4-194-204-214-224-234-244-254-264-274-284-294-305-015-02 10596102109113117124144126851489311383 岭[13]。 北京市防治SARS的经验以及数据分析结果表明,为阻断SARS疫情在北京传播,在新发病人数(包括疑似患者中的确诊者)在达到WHO接触疫 感染人数 预测 1109398115118112120149123811538911786 相对误差/% 4176-3112-319251504142-4127-31233147-2138-41713138-413031543161 日期 实际 5-035-045-055-065-075-085-095-105-115-125-135-145-155-16 10562946389874150383943231817 感染人数 预测 11165976092913948404142241718 相对误差/% 517141843119-417631374160-4188-410051265113-21334135-51565188 日期 实际 5-175-185-195-205-215-225-235-245-255-265-275-285-295-30 151437012925958233 感染人数 预测 141537012926958233 相对误差/% -61677114000004100000000 日期 实际 5-316-016-026-036-046-056-066-076-086-096-106-116-126-13 11000000000000 感染人数 预测 11000000000000 相对误差/% 00000000000000 © 1994-2006 China Academic Journal Electronic Publishing House. All rights reserved. http://www.cnki.net 第9期刘云忠等:SARS传染病数学建模及预测预防控制机理研究 [6][7][8][9][10][11] 65 4 结论 首先,利用疾病传播的一般规律及人口守恒统计法则建立了四类人的SARS传染病数学模型,然后运用数学方法对四类人的SARS传染病数学模型进行分析,得出了其生理意义和预防、控制机理。其次,利用人工神经网络理论建立了SARS的预测模型,以北京市的SARS数据为例进行了预测和分析,预测结果显示该模型简单易行,预测精度高。 参考文献 [1] KermackWO,Mckendrick1AG.Contributionstothemathematicaltheoryofepidemics[J].ProcRoySocA,1927,115(5):700~721 [2] KermackWO,MckendrickAG.Contributionstothemathematicaltheoryofepidemics[J].ProcRoySocA,1932,138(1):55~83[3] KermackWO,MckendrickAG.Contributionstothemathematicaltheoryofepidemics[J].ProcRoySocA,1933,141(1):94~122 [4][5] CookeKL.Stablilityanalysisforavectordiseasemode[J].RockyMountJMath1979,9(1):31~42HethcoteHW.Qualititativeanalysisofcommunicable diseasemodels[J].MathBiosci,1976,28(3):335 ~356 陈兰荪,陈 健.非线性生物动力系统[M].北京:科学出版社,1983.111~162 张运权.传染病数学模型的建立及稳定性分析[J].数理医药学杂志,1999,12(2):99~100 张芷芬.微分方程定性理论[M].北京:科学出版社,1985 张锦炎.常微分方程几何理论与分支问题[M].北京:北京大学出版社,1987 王稳地.传染病数学模型的稳定性和分枝[D].西安:西安交通大学,2002 丁守銮,王洁贞,崔希友.基于双曲正切函数 HFRS发病率的BP神经网络预测模型[J].系统工 程理论与实践,2003,23(7):126~131,136 [12]HuangSC,HuangYF.Boundsonthenumberofhiddenneuronsonmultilayerperceptions[J].IEEETransNeuralNetworks,1991,2(1):47~55 [13] 龚其国,佟仁城,吕本富,等.北京市SARS的防疫措施、疫情的发展趋势分析[J].管理评论, 2003,15(4):28~30 [14]王建锋.SARS流行预测分析[J].中国工程科学, 2003,5(8):23~29 MathematicalandPredictiveModelsofSARSEpidemic DiseaseandMechanismofPreventionandControl LiuYunzhong,XuanHuiyu,LinGuoxi (ManagementSchoolofXiπanJiaotongUniversity,Xiπan 710049,China) [Abstract] First,basedonthelawofdiseasepropagationandconservationofpopulationstatistics,thispapersetsupmathematicalmodelsoffourkindsofpersonforSARS,appliesmathematicalmethodtoanalyzethemodel,anditgetsthephysiologysensesofthemathematicalmodelandmechanismofpreventionandcontrolrespectively.Second,itsetsuppredictivemodelsofSARSepidemicdiseaseusingBPartificialneuralnetwork.AccordingtothedataofSARSinBeijing,itpredictsandanalyzesthemodel.Theexperimentalresultsshowthatthemethodissimpler,andmoreeffectivewithhighpredictionprecision。[Keywords] SARS;stability;trail;thresholdtheory;artificialneuralnetwork © 1994-2006 China Academic Journal Electronic Publishing House. All rights reserved. http://www.cnki.net 因篇幅问题不能全部显示,请点此查看更多更全内容