具有非線性發(fā)生率的諾如病毒傳播動(dòng)力學(xué)模型分析

秦海燕 侯強(qiáng)

摘 要：為了減少因諾如病毒感染引起的感染性腹瀉對(duì)人們身體健康造成的危害，在明確諾如病毒傳播特征的基礎(chǔ)上，研究了諾如病毒的傳播動(dòng)力學(xué)行為，考慮感染諾如病毒的潛伏者也傳染疾病的特性，建立具有非線性發(fā)生率的諾如病毒傳播動(dòng)力學(xué)模型，在計(jì)算模型的基本再生數(shù)R0的基礎(chǔ)上，利用Lyapunov函數(shù)和幾何方法證明了無(wú)病平衡點(diǎn)和地方病平衡點(diǎn)的穩(wěn)定性，并進(jìn)行了數(shù)值模擬。結(jié)果表明，當(dāng)R0≤1時(shí)，無(wú)病平衡點(diǎn)全局漸近穩(wěn)定，疾病消失;當(dāng)R0>1時(shí)，在一定條件下，地方病平衡點(diǎn)全局漸近穩(wěn)定。數(shù)值模擬驗(yàn)證了理論結(jié)果的正確性。研究結(jié)果豐富了感染性病毒傳播理論，對(duì)進(jìn)一步研究病毒的傳播機(jī)理具有借鑒意義。

關(guān)鍵詞：穩(wěn)定性理論;諾如病毒;傳染病模型;非線性發(fā)生率;基本再生數(shù)

中圖分類(lèi)號(hào)：O175.1 文獻(xiàn)標(biāo)識(shí)碼：A

doi：10.7535/hbkd.2021yx02005

Dynamic model analysis of Norovirus transmissionwith nonlinear incidence

QIN Haiyan，HOU Qiang

（School of Science，North University of China，Taiyuan，Shanxi 030051，China）

Abstract：In order to reduce the great harm of infectious diarrhoeal disease caused by Norovirus infection to human health，based on the transmission characteristics of Norovirus，the transmission dynamics behavior of Norovirus was studied.Taking into account the characteristics that the latent infected with Norovirus can also transmit the disease，a dynamic model of Norovirus transmission with nonlinear incidence was established.The basic reproduction number R0of the model was calculated and then the stability of the disease-free equilibrium point and the endemic equilibrium point were proved by using the Lyapunov function and the geometric method.The results show that when R0≤1，the disease-free equilibrium point is globally asymptotically stable and the disease disappears; when R0>1，under certain conditions，the endemic equilibrium point is global asymptotically stable.The theoretical results are verified by numerical simulation.The research results have enriched the theory of infectious virus transmission and provide a reference for the study of virus transmission mechanism.

Keywords：

stability theory;Norovirus;infectious disease model;nonlinear incidence;basic reproduction number

諾如病毒是全球急性腸胃炎的主要致病原，致人感染發(fā)病的主要表現(xiàn)為腹瀉和嘔吐，其傳播途徑包括人與人之間的傳播和食源性傳播。人傳人通過(guò)糞口途徑傳播（包括吸入糞便或嘔吐物產(chǎn)生的氣溶膠），食源性傳播通過(guò)食用被諾如病毒污染的食物和水進(jìn)行傳播。諾如病毒變異很快，對(duì)環(huán)境的抵抗力也很強(qiáng)，感染后潛伏期較短，且在潛伏期便可排出病毒，排毒時(shí)間長(zhǎng)，免疫保護(hù)時(shí)間卻較短[1-9]。自2013年以來(lái)，感染性腹瀉病暴發(fā)大多以諾如病毒感染為主，尤其是自2014年冬季起，諾如病毒的感染暴發(fā)疫情有較大幅度增加，顯著高于歷年水平[10-11]。

近年來(lái)，許多學(xué)者對(duì)諾如病毒的傳播動(dòng)力學(xué)進(jìn)行了研究[12-14]。史方等[15]根據(jù)諾如病毒的流行病學(xué)特點(diǎn)建立模型，利用實(shí)測(cè)數(shù)據(jù)估計(jì)了模型中的參數(shù)，分析了潛伏期內(nèi)的新發(fā)感染者、新發(fā)病人和新發(fā)無(wú)癥狀感染者之間的關(guān)系，研究了不同傳染源在疾病傳播研究過(guò)程中的貢獻(xiàn)率，以及外部干預(yù)措施對(duì)諾如病毒傳播的影響;GAYTHORPE等[16]根據(jù)年齡分層建立時(shí)間序列模型，研究發(fā)現(xiàn)，接種疫苗是預(yù)防疾病的有效策略;黃璜等[17]通過(guò)比較各模型評(píng)價(jià)隔離措施對(duì)諾如病毒在醫(yī)院內(nèi)傳播的影響，為有效防控諾如病毒在醫(yī)院內(nèi)的感染暴發(fā)提供理論依據(jù);LANE等[18]研究了諾如病毒通過(guò)人與人之間的傳播和食源性傳播導(dǎo)致人產(chǎn)生感染性腸道疾病，并進(jìn)一步了解食源性如何傳播諾如病毒。

目前的研究中，較少人考慮諾如病毒傳播與穩(wěn)定性相關(guān)的動(dòng)力學(xué)性態(tài)，因此，筆者根據(jù)諾如病毒的流行病學(xué)特性，提出了發(fā)病率形式為β1fES+β2gIS的SEIRS非線性傳播動(dòng)力學(xué)模型，并運(yùn)用微分動(dòng)力系統(tǒng)理論分析模型的穩(wěn)定性，最后通過(guò)數(shù)值模擬驗(yàn)證結(jié)果的正確性。

1 模型的建立

將人群分為4類(lèi)，即易感者S（t）、潛伏者E（t）、染病者I（t）和恢復(fù)者R（t）。用N表示總?cè)丝冢琋=S+E+I+R，易感者雖然也可通過(guò)攝取被諾如病毒污染的食物而被感染，但本文重點(diǎn)研究潛伏者、染病者對(duì)諾如病毒傳播的影響，所以忽略食源性感染。建立如下模型：

dSdt=A-β1fES-β2gIS-μS+δR，dEdt=β1fES+β2gIS-μE-εE，dIdt=εE-μI-γI，dRdt=γI-μR-δR。（1）

式中：A表示人口常數(shù)輸入率;β1表示潛伏者的傳染率;β2表示染病者的傳染率;ε表示潛伏者變?yōu)槿静≌叩谋壤?γ表示人的恢復(fù)率;δ表示恢復(fù)者失去免疫力的比率;μ表示人的自然死亡率。采用一般非線性發(fā)生率β1fES+β2gIS表示諾如病毒的傳播發(fā)生率。為了滿足模型（1）的流行病學(xué)意義，對(duì)于發(fā)生率函數(shù)，假定f0=g0=0，且f′E>0，g′I>0，f″E≤0，g″I≤0，在文獻(xiàn)[17]中，這里的發(fā)生率為雙線性發(fā)生率β1SE+β2SI。

2 平衡點(diǎn)的存在性和基本再生數(shù)

對(duì)于模型（1），其可行域?yàn)棣?S，E，I，RS≥0，E≥0，I≥0，R≥0，S+E+I+R≤Aμ。模型（1）有唯一的無(wú)病平衡點(diǎn)E0=（S0，E0，I0，R0）=Aμ，0，0，0，根據(jù)下一代矩陣法[20]，F(xiàn)=β1f′0S0β2g′0S000，V=μ+ε0-εμ+γ，則模型的基本再生數(shù)為R0=ρ（FV-1）=β1f′0S0μ+ε+β2g′0S0εμ+εμ+γ=R01+R02。

為求模型（1）地方病平衡點(diǎn)，令模型（1）的第3個(gè)和第4個(gè)等式右端為0，則有：

E=μ+γεI;R=γμ+δI。

令模型（1）的前2個(gè)等式右端為0，則有：A-μ+εE-μS+δR=0。可得：

S=Aμ-μ+εμ+γμε-δγμμ+δI。

令第1個(gè)式子右邊為0，可得：

GI=A-β1fμ+γεI+β2gI+μAμ-k1I+k2I，

其中：k1=μ+εμ+γμε-δγμμ+δ;k2=δγμ+δ。

對(duì)G（I）求導(dǎo)可得：

G′I=k2+k1β1fμ+γεI+β2gI+μ-β1μ+γεf′μ+γεI+β2g′IAμ-k1I，

則有：

G′0=k2+k1μ-β1μ+γεf′0+β2g′0Aμ=μ+εμ+γε（1-R0），

又

G″I=2k1β1μ+γεf′μ+γεI+β2g′I-β1μ+γε2f″μ+γεI+β2g″IAμ-k1I，

當(dāng)I∈0，Ak1μ時(shí)，G″I>0，所以G′I單增，又當(dāng)R0>1時(shí)，G′0<0，可得GI先減后增，且G0=0，GAk1μ>0，所以一定存在一個(gè)I，使得G（I）=0;當(dāng)I∈Ak1μ，Aμ時(shí)，GAk1μ>0，G′I>0，此時(shí)不存在一個(gè)I，使得G（I）=0。

因此有結(jié)論：當(dāng)R0>1時(shí)，模型（1）存在唯一正平衡點(diǎn)E*（S*，E*，I*，R*），即有：

定理1 對(duì)于模型（1）有以下結(jié)論：當(dāng)R0≤1時(shí)，模型（1）有一個(gè)無(wú)病平衡點(diǎn)E0;當(dāng)R0>1時(shí)，模型（1）存在唯一正平衡點(diǎn)E*。

3 平衡點(diǎn)的穩(wěn)定性

定理2 如果R0≤1，模型（1）的無(wú)病平衡點(diǎn)E0在Ω內(nèi)是全局漸近穩(wěn)定的;如果R0>1，無(wú)病平衡點(diǎn)是不穩(wěn)定的。

證明 矩陣在E0處的特征方程為（λ+μ）λ+μ+δ（λ2+a1λ+a2）=0。其中：a1=μ+γ+μ+ε1-R01，a2=μ+εμ+γ1-R0。

當(dāng)R0≤1時(shí)，a1>0，a2>0，又λ1=-μ<0，λ2=-μ-δ<0，所以模型（1）的無(wú)病平衡點(diǎn)是局部漸近穩(wěn)定的。當(dāng)R0>1時(shí)，有a2<0，這里存在正根，無(wú)病平衡點(diǎn)不穩(wěn)定，構(gòu)造Lyapunov函數(shù)證明無(wú)病平衡點(diǎn)的全局漸近穩(wěn)定性：因?yàn)閒x在區(qū)間0，Aμ上恒有f″x≤0，則對(duì)任意的x，y∈0，Aμ滿足x

同理可得：gII′=g′II-gII2=g′I-gIII≤0。所以，fEE和gII為減函數(shù)。

因此：β1f（E）S（μ+ε）E≤limE→0β1f′（E）S0（μ+ε）E=β1f′E（0）S0μ+εb1，β2g（I）S（μ+γ）I≤limI→0β2g（I）S0（μ+γ）I=β2g′I（0）S0μ+γb2。

定義 J=10-εμ+ε1，c1，c2=b1，b2J-1。因此，c1=R0≤1，根據(jù)（b1，b2）的正性可得（c1，c2）的正性，構(gòu)造Lyapunov函數(shù)L=R0E+c2I。

對(duì)L關(guān)于模型（1）求導(dǎo)：

dLdt=R0β1fESμ+εE，β2gISμ+γIμ+εE，μ+γIT-R0，c2Jμ+εE，μ+γIT≤

R0β1f′E0S0μ+ε，β2g′I0S0μ+γμ+εE，μ+γIT-R0，c2Jμ+εE，μ+γIT=

R0-1b1，b2μ+εE，μ+γIT。

因此，當(dāng)R0≤1時(shí)，L′≤0，而L′=0的最大不變集包含唯一的點(diǎn)E0=Aμ，0，0，0。由Lasalle不變集原理[19]可知，當(dāng)R0≤1時(shí)，無(wú)病平衡點(diǎn)在Ω內(nèi)是全局漸近穩(wěn)定的。

由定理2知，如果R0>1，無(wú)病平衡點(diǎn)是不穩(wěn)定的。根據(jù)文獻(xiàn)[20]，易證模型（1）的一致持續(xù)等價(jià)于無(wú)病平衡點(diǎn)的不穩(wěn)定性，故當(dāng)R0>1時(shí)，模型（1）是一致持續(xù)的，因此存在一個(gè)緊吸引子集K。

因此當(dāng)t>t*時(shí)，可得：

1t∫t0g1ds≤1t∫0t*g1ds-1tlnItI（t*）-b—t-t*t，

1t∫t0g2ds≤1t∫0t*g2ds-1tlnItI（t*）-b—t-t*t，

1t∫t0g3ds≤1t∫0t*g3ds-1tlnItI（t*）-b—t-t*t，

1t∫t0g4ds≤1t∫0t*g4ds-1tlnItI（t*）-b—t-t*t，

則有：

1t∫t0η（Q）ds≤sup1t∫0t*g1ds-1tlnItI（t*）-b—t-t*t，

1t∫0t*g2ds-1tlnItI（t*）-b—t-t*t，

1t∫0t*g3ds-1tlnItI（t*）-b—t-t*t，

1t∫0t*g4ds-1tlnItI（t*）-b—t-t*t。

從而有：

q=lim sup supt→+SymboleB@x∈K1t∫t0η（Q）ds≤ -b—<0，

即當(dāng)R0>1時(shí)，系統(tǒng)（1）的地方病平衡點(diǎn)E*全局漸近穩(wěn)定。

4 數(shù)值模擬

采用數(shù)值模擬方法證明模型（1）中非線性發(fā)生率分別為雙線性發(fā)生率和飽和發(fā)生率時(shí)，其在無(wú)病平衡點(diǎn)和地方病平衡點(diǎn)處的動(dòng)力學(xué)性質(zhì)。

圖1和圖2是非線性發(fā)生率為雙線性發(fā)生率β1SE+β2SI時(shí)，取不同參數(shù)值的數(shù)值模擬情況。

參數(shù)值如下：

1） μ=0.05，β1=0.006，β2=0.008，δ=0.01，ε=0.5，γ=0.5，A=2。I取初始值I（0）=3，I（0）=5，I（0）=7，此時(shí)R0=0.531<1。從圖1中可以看出無(wú)病平衡點(diǎn)E0在Ω內(nèi)全局漸近穩(wěn)定。

2）μ=0.05，β1=0.01，β2=0.02，δ=0.01，ε=0.5，γ=0.5，A=2。I取初始值I（0）=3，I（0）=5，I（0）=7，此時(shí)R0=1.127>1，且模型的雅可比矩陣在地方病平衡點(diǎn)處的特征方程的特征根均具有負(fù)實(shí)部，因此地方病平衡點(diǎn)E*在Ω內(nèi)局部漸近穩(wěn)定，從圖2可以看出地方病平衡點(diǎn)E*在Ω內(nèi)全局漸近穩(wěn)定。

參數(shù)值如下：

1）μ=0.05，β1=0.006，β2=0.008，δ=0.01，ε=0.5，γ=0.5，A=2，α=0.2。I取初始值I（0）=3，I（0）=5，I（0）=7，此時(shí)R0=0.531<1。從圖3可以看出無(wú)病平衡點(diǎn)E0在Ω內(nèi)全局漸近穩(wěn)定。

2）μ=0.05，β1=0.01，β2=0.02，δ=0.01，ε=0.5，γ=0.5，A=2，α=0.2。I取初始值I（0）=3，I（0）=5，I（0）=7，此時(shí)R0=1.127>1，且模型的雅可比矩陣在地方病平衡點(diǎn)處的特征方程的特征根均具有負(fù)實(shí)部，因此地方病平衡點(diǎn)E*在Ω內(nèi)局部漸近穩(wěn)定，從圖4可以看出地方病平衡點(diǎn)E*在Ω內(nèi)全局漸近穩(wěn)定。

5 結(jié) 語(yǔ)

針對(duì)諾如病毒變異快、感染性強(qiáng)、對(duì)環(huán)境抵抗力強(qiáng)的特性，在一般SEIR模型基礎(chǔ)上，提出了具有非線性發(fā)生率的諾如病毒傳播動(dòng)力學(xué)模型，在一般的生物學(xué)假設(shè)下，分析了平衡點(diǎn)的存在性和穩(wěn)定性。結(jié)果表明，如果R0<1，則無(wú)病平衡點(diǎn)全局漸近穩(wěn)定，這意味著疾病滅絕;如果R0>1，則存在唯一的地方病平衡點(diǎn)且在一定條件下全局漸近穩(wěn)定。本研究尚未證明地方病平衡點(diǎn)的局部穩(wěn)定性，這是后續(xù)要做的重要工作之一。

參考文獻(xiàn)/References：

[1] LOPMAN B，VENNEMA H，KOHLI E，et al.Increase in viral gastroenteritis outbreaks in Europe and epidemic spread of new Norovirus variant[J].The Lancet，2004，9410（363）：682-688.

[2] YEN C，WIKSWO M E，LOPMAN B A，et al.Impact of an emergent Norovirus variant in 2009 on Norovirus outbreak activity in the United States[J].Clinical Infectious Diseases an Official Publication of the Infectious Diseases Society of America，2011，53（6）：568-571.

[3] GTZ H，EKDAHL K，LINDBCK J，et al.Clinical spectrum and transmission characteristics of infection with Norwalk virus：Findings from a large community outbreak in Sweden[J].Clin Infect Dis，2001，33：622-628.

[4] LEE R M，LESSLER J，LEE R A，et al.Incubation periods of viral gastroenteritis：A systematic review[J].Infect Dis，2013，13：446.

[5] ATMAR R L，OPEKUN A R，GILGER M A，et al.Norwalk virus shedding after experimental human infection[J].Emerg Infect Dis，2008，14：1553-1559.

[6] KIRBY AE，SHI J，MONTES J，et al.Disease course and viral shedding in experimental norwalk virus and snow mountain virus infection[J].J Med Virol，2014，86（12）：2055-2064.

[7] TEUNIS P F，SUKHRIE F H，VENNEMA H，et al.Shedding of norovirus in symptomatic and asymptomatic infetion[J].Epidemiol Infect，2015，143（8）：1710-1717.

[8] AOKI Y，SUTO A，MIZUTA K，et al.Duration of Norovirus excretion and the longitudinal course of viral load in Norovirus-infected elderly patients[J].Hosp Infect，2010，75（1）：42-46.

[9] THOMAS A，PARRINO M D，DAVID S，et al.Clinical immunity in acute gastroenteritis caused by norwalk agent[J].The New England Journal of Medicine，1977，297：86-89.

[10]LU Jing，SUN Limei，F(xiàn)ANG Lin，et al.Gastroenteritis outbreaks caused by Norovirus GII.17，Guangdong Province，China，2014—2015[J].Emerging Infectious Disease，2015，21（7）：1240-1241.

[11]HAN J，JI L，SHEN Y et al.Emergence and predominance of norovirus GII.17 in Huzhou，China[J].Virol Journal，2015，12（1）：139.

[12]NEWMAN K L，MOE C L，KIRBY A E，et al.Norovirus in symptomatic and asymptomatic individuals：Cytokines and viral shedding[J].Clin Exp Immunol，2016，184：347-357.

[13]VANDERPAS J，LOUIS J，REYNDERS M，et al.Mathematical model for the control of nosocomial Norovirus[J].Journal of Hospital Infection，2009，71（3）：214-222.

[14]MAESHALL J，BRUGGINK L D.The dynamics of Norovirus outbreak epidemics：Recent insights[J].International Journal of Environmental Research and Public Health，2011，8（4）：1141-1149.

[15]史方，周毅，劉茂華，等.諾如病毒感染的傳播動(dòng)力學(xué)研究[J].中國(guó)國(guó)境衛(wèi)生檢疫雜志，2015，38（6）：413-418.

SHI Fang，ZHOU Yi，LIU Maohua，et al.Study on transmission dynamics of Norovirus infection[J].Chinese Journal of Frontier Health and Quarantine，2015，38（6）：413-418.

[16]GAYTHORPE K A M，TROTTER C L，CONLAN A J K.Modelling norovirus transmission and vaccination[J].Vaccine，2018，36：5565-5571.

[17]黃璜，蘭玲，單旭征，等.SEIR模型在諾如病毒醫(yī)院感染的應(yīng)用研究[J].中國(guó)醫(yī)院感染學(xué)雜志，2018，28（14）：2205-2208.

HUANG Huang，LAN Ling，SHAN Xuzheng，et al.Application of SEIR model in nosocomial infection of Norovirus[J].Chinese Journal of Nosocomial Infectiology，2018，28（14）：2205-2208.

[18]LANE D，HUSEMANN E，HOLLAND D，et al.Understanding foodborne transmission mechanisms for Norovirus：A study for the UK's Food Standards Agency[J].European Journal of Operational Research，2019，275：721-736.

[19]LASALLE J P.The stability of dynamical systems[J].SIAM Review，1976，21（3）.doi：10.1137/1021079.

[20]MICHAEL Y L，JOHN R G，WANG Liancheng，et al.Global dynamics of a SEIR model with varying total population size[J].Math Biosci，1999，160：191-213.

[21]MICHAEL Y L，MULDOWNEY J S.A geometric approach to global stability problems[J].Siam Journal on Mathematical Analysis，1996，27：1070-1083.