一类具有疫苗接种的双菌株流感模型的动力学分析

王晓静,梁宇,郭松柏,陈靖宜,李佳慧,郭德玉

(北京建筑大学 理学院,北京 102616)

甲型H1N1流感(简称甲流)为急性呼吸道传染病,其病原体是一种新型甲型H1N1流感病毒,人群对甲型H1N1流感病毒普遍易感,而且可以人传染人.感染甲流后,人们表现出的早期症状与普通流感相似,包括发热、咳嗽、喉痛、身体疼痛、头痛、发冷和疲劳等,有些还会出现腹泻或呕吐、肌肉痛或疲倦、眼睛发红等.2009年开始,甲型H1N1流感在全球范围内大规模流行.2010年8月,世界卫生组织宣布甲型H1N1流感大流行期已经结束[1].流感病毒很容易产生耐药毒株,耐药性的出现极大增加了流感治疗的难度,甚至使药物对病毒的治疗失效.因此接种疫苗成为控制流感的有效方法.在中国,甲流疫苗接种的成功率很高,但是任何一种疫苗都只能起到预防作用,不可能完全阻止[2-3].

QIU等[1]探讨了疫苗接种和抗病毒治疗对具有耐药性菌株流感模型的影响,提出了两种菌株共存的阈值,研究结果表明,在过度治疗的情况下,会促进耐药性发展,从而增加感染规模,因此应适当实施抗病毒治疗.WANG 等[4]建立了两个H1N1流感的动力学模型,研究了抗病毒治疗对耐药性菌株的影响,证明了无病平衡点的全局稳定性,并且推导了基本和控制再生数.ROBINSON等[5]研究了由药物治疗引起的耐药性对流感传播的影响,并且证明了当R1>1且R1>R2时药物的敏感性菌株和耐药菌株是共存的.BABA等[6]考虑了敏感性和耐药性两种菌株的流感模型,将双线性和饱和发生率分别用于刻画敏感性和耐药性菌株的传染率,证明了边界平衡点和地方病平衡点的全局稳定性.文献[7]在文献[1]的基础上加入了潜伏期感染者仓室,分别分析了敏感性和抗药性个体的再生数,得到了两种菌株共存的阈值,通过数值模拟表明应该进行适当的抗病毒治疗.PEDR等[8]建立了一个包含治疗和疫苗接种的流感模型.对模型参数值进行了不确定性分析和敏感性分析,结果表明治疗措施和流行病规模之间呈非线性关系,并且治疗率存在一个临界阈值.

基于已有的研究成果,本文在模型中考虑了疫苗接种、耐药性菌株以及潜伏期感染者可以自愈等因素对流感传播的影响,以便更准确地刻画流感传播动态.

1 模型的建立

将总人口分为易感者(S)、潜伏期感染者(E)、注射疫苗的群体(V)、敏感性菌株感染者(IS)、耐药性菌株感染者(IR)、恢复者(R).基于文献[4,9]建立如下的流感模型:

(1)

流感的传播流程图见图1,模型(1)中相关参数见表1.

2 模型的稳定性分析

令模型(1)右端的方程等于零,

(2)

控制再生数Rc是反映传染病传播能力的一种指标,表明在实施一定的防控措施下一个感染者在其感染期内产生的继发染病者的平均数量.当Rc<1时,疫情会得到有效的控制,当Rc>1时,会形成地方病.通过下一代矩阵的方法[11]可计算出控制再生数.

其生物意义如下:RS和RR分别表示一个敏感性菌株和耐药性菌株的感染者在其感染期内平均感染的人数.每个敏感性菌株的感染者以f的概率接受治疗,以c的比例转化为耐药性菌株的感染者.并且ε+με+μ+η表示易感者的比例,1f+γS+μ为敏感性菌株的感染期,1γR+μ为耐药性菌株的感染期.

定理1当Rc>1且RS>RR时,模型(1)存在唯一的正平衡点

P*=(S*,V*,E*,I*S,I*R,R*).

表1 模型(1)的参数变量描述及取值

因为模型(1)的前5个方程和最后一个方程是相互独立的,因此只需考虑如下子系统:

(3)

令系统(3)的右端都等于0,可得子系统的无病平衡点E1=(S0,V0,0,0,0).

下面对子系统(3)的无病平衡点的全局稳定性进行证明,为此首先引入如下引理.

引理1[12]如果一个模型可以表示为如下形式:

(4)

其中X∈Rm表示未感染类,Z∈Rn表示感染类,包括潜伏期感染者.E2=(X*,0)为系统(4)的无病平衡点,满足如下两个条件:

则当Rc<1时,无病平衡点E2=(X*,0)是全局渐近稳定的.

定理2当Rc<1时,系统(3)的无病平衡点E1是全局渐近稳定的.

在平衡点E1=(X*,0)处,模型简化为:

(5)

其特征方程为:

λ2+(ε+2μ+η)λ+μ(ε+μ+η)=0.

由λ1λ2>0,λ1+λ2<0,知方程有两个负的特征根,故X*是局部渐近稳定的.

由文献[13]中的引理5.1可知,存在一个序列{tn}使得S(tn)→S∞,S′(tn)→0,tn→∞(n→∞),再由系统(5)的第一个方程,有

另外,

3 敏感性分析

敏感性分析是通过敏感性指数来确定影响疫情传播的重要因素.敏感性指数是指测量一个状态变量在参数变化时的相对变化.变量m关于n的敏感性指数定义为:

将根据上面的公式依次求出Rc关于βS,βR,p,ε,c,η的敏感性指数.

由计算可得:

并且

同理

所以Rc随着参数βS,βR,p,ε,c的增大而增大,且对潜伏期人群转换为感染者的比例p的变化最敏感.

计算Rc对η的偏导数,有

因此,Rc随疫苗接种率η的增大而减小.

由于无法通过计算确定Rc与参数m,f之间的变化关系,因此,通过第4部分的数值模拟进行判断.

4 数值模拟

在本小节,首先对各个参数进行全局灵敏度分析并给出了参数敏感性菌株的比例m与Rc之间的关系;其次,给出了随着疫苗接种η和治疗率f的变化,累计感染人数以及耐药性菌株的累计病例数的变化规律.

灵敏度分析[14-15]是研究系统或模型的状态对系统参数变化的敏感程度,分为局部灵敏度分析和全局灵敏度分析.全局灵敏度分析可以分析多个输入参数对输出参数的影响程度,并且涵盖了参数的概率分布对结果输出的影响.将通过Sobol法进行全局灵敏度分析,Sobol法是一种基于方差的蒙特卡罗法,主要思想是方差分解,将模型的总方差分解为各个参数的方差以及各个参数之间相互作用的方差之和.

图2分别表示RS,RR,Rc关于各参数的灵敏度分析.其中ST为总灵敏度,其意义为各参数对输出变量(再生数)的影响.柱状图上的黑线是在95%的置信度下画出的误差线,刻画潜在的误差或不确定程度.

图2(a)表明敏感性菌株的比例m对敏感性菌株的控制再生数RS的影响最大,并且参数βS,p,μ,f,η的灵敏度值较大.由图2(b)可知,耐药性菌株的控制再生数RR对潜伏期人群转换为感染者的比例p最敏感,参数βR,m,μ,η的变化对RR的影响较大,由敏感性分析以及数值模拟的结果可知,随着潜伏期人群转换为感染者的比例p增大,控制再生数会增大;随着疫苗接种率η的增大,控制再生数会减小,因此为了控制疾病的暴发,需要采取干预手段以减小潜伏期人群向感染者转移的比例,并且提高疫苗接种率,形成群体免疫,减少流感感染的规模.

由图2(c)可以观察到,敏感性菌株的比例m对控制再生数Rc的影响最大,参数βR,p,μ的灵敏度值较大,并且根据图3可知Rc随着敏感性菌株的比例m的增大而减小.因此,可以加强追踪力度,合理进行药物治疗,以减少流感暴发期内的感染规模.

由图4可以发现,随着治疗率f的增加,耐药性菌株的感染人数IR逐渐增大,这表明随着抗流感病毒药物治疗使用的增多,药物治疗的副作用会越来越明显,不利于耐药性菌株传播的控制,造成这种现象的主要原因是滥用抗生素.

图5的结果表明随着疫苗接种率η的增大,累计感染人数IS+IR逐渐减少,并且当疫苗接种率达到0.3时,相比于未接种疫苗的情况,累计感染人数降低约0.47倍.这表明有效接种疫苗可以减少流感的感染人数,提高易感者的抵抗力,为控制甲型H1N1的传播起到积极作用.

图6表明随着潜伏期感染者转换为双菌株感染者比例p的增大,累计感染人数IS+IR也逐渐增加.当p=1时,累计感染人数达到最大.

5 结 论

本文在已有的研究成果基础上考虑了疫苗接种、耐药性菌株以及潜伏期感染者可以自愈等因素对甲型H1N1流感传播的影响,构建了一类SVEISIRR传染病模型,计算了控制再生数Rc,给出了在一定条件下地方病平衡点的存在性,且得到了当控制再生数Rc<1时,无病平衡点是全局渐近稳定的充分条件.对主要参数进行了敏感性分析,结果表明控制再生数Rc随参数βS,βR,p,ε,c的增大而增大,随η的增大而减小.通过Sobol法灵敏度分析证实了敏感性菌株的比例m对控制再生数Rc的影响最大,其次为耐药性菌株的传播率βR和潜伏期人群转换为感染者的比例p.

数值模拟表明,感染规模随着潜伏期感染者自愈率1-p的增大而减少,在无自愈的情况下,感染规模将达到最大.耐药性菌株的感染人数会随着治疗率f的提高而增大,因此易感人群应做好防护并积极接种疫苗,感染人群要遵循医嘱合理用药,规范用药剂量,医院也应该建立标准的药敏试验方法,通过病原性细菌耐药检查,帮助医生精准选取药物,将会有效控制流感的传播.