一類潛伏期傳染且具有病毒變異的傳染病模型

馬怡婷，張?zhí)?，鄧金超，劉俊利

（1.長安大學 理學院，陜西 西安 710064；2.西安工程大學 理學院，陜西 西安 710048）

0 引言

傳染病動力學模型在理解傳染病的傳播規(guī)律和預測傳染病的傳播趨勢中起著重要作用。早在1927 年，Kermack 和Mckendrick 就建立了SIR 傳染病模型［1］，這種“倉室”建模思想一直被廣泛使用且在不斷地發(fā)展［2-7］。近年來，疾病具有潛伏期且潛伏期具有傳染性的動力學模型已引起許多學者關注［8-10］，通過對潛伏期與染病期均傳染的傳染病模型的研究，說明了忽略潛伏期會影響參數(shù)估計的準確性，進而對疾病的控制產(chǎn)生很大影響。

醫(yī)療數(shù)據(jù)統(tǒng)計發(fā)現(xiàn)，一些傳染病在其傳播過程中會發(fā)生病毒變異，使其難以控制，極易導致傳染病暴發(fā)。一旦傳染病蔓延，將對人類健康與社會發(fā)展造成嚴重威脅，例如禽流感病毒（H7N9），新型冠狀病毒（COVID-19），乙肝等疾病。特別是COVID-19，疫情暴發(fā)初期，國內(nèi)累計確診人數(shù)與死亡人數(shù)都在不斷增加，有多名醫(yī)務工作者死于COVID-19，疫情給國家造成了人力、物力和財力的巨大損失。隨著疫情的不斷發(fā)展，新冠病毒開始變異，出現(xiàn)了多種具有較強傳播力與潛伏性的變異病毒。根據(jù)文獻［11］，新冠變異病毒德爾塔毒株最早于2020 年10 月在印度被發(fā)現(xiàn)，比阿爾法毒株（英國首先報告的B.1.1.7 變異新冠病毒）更易傳播，至2021 年6 月24 日已擴散到92 個國家，更易傳播的德爾塔毒株為全球抗疫帶來新挑戰(zhàn)。因此，研究病毒變異的過程有助于了解傳染病的傳播規(guī)律和更好地控制傳染病。關于病毒變異的傳染病模型已有很多研究工作：Gao 等［12］研究了一類具有病毒變異Logistic 死亡率的SEIR 傳染病模型，證明了平衡點的全局穩(wěn)定性，并分析了模型中參數(shù)對疾病傳播的影響；李冬梅等［13］建立了具有自發(fā)病毒變異的時滯SIR 傳染病模型，給出了平衡點的存在性、局部穩(wěn)定性的充分條件，并通過數(shù)值模擬的方法分析了病毒自發(fā)變異對疾病傳播的影響；杜金姬等［14］研究了一類具有病毒變異且具有Logistic 死亡率的SEIR 傳染病模型，證明了隨機系統(tǒng)全局正解的存在唯一性，并建立了傳染病滅絕的條件。

以上模型雖然考慮到了疾病在傳播過程中發(fā)生病毒變異，但均未考慮疾病潛伏期具有傳染性，且未對模型作具體的疾病應用。因此，本文將建立一類潛伏期和染病期均傳染且具有病毒變異的SEI1I2QR 模型。在第2 節(jié)，給出模型的基本再生數(shù)?0與各類平衡點的存在性，并證明各類平衡點的全局漸近穩(wěn)定性。第3 節(jié)，根據(jù)印度COVID-19 的累計病例數(shù)對所建模型進行擬合，對疫情的發(fā)展趨勢進行回溯性研究，并對閾值?0中的部分參數(shù)進行敏感性分析，根據(jù)分析結(jié)果提出合理的建議。第4 節(jié)，對全文進行總結(jié)。

1 模型建立

從易感者到感染者會有一段潛伏期，這段時期主要表示從易感者到染病者所需要的時間，用倉室E來表示。將某一區(qū)域在t時刻的總?cè)藬?shù)N(t)分為易感者S(t)、潛伏者E(t)、病毒變異前的染病者I1(t)、病毒變異后的染病者I2(t)、隔離者Q(t)和治愈者R(t)，建立了一類潛伏期和染病期均傳染且具有病毒變異的SEI1I2QR 模型。對模型作出以下假設：

（1） 易感人群常數(shù)輸入，輸入率為Λ，μ是人們的自然死亡率。

（2） 易感者與潛伏者、病毒變異前染病者和病毒變異后染病者有效接觸后都需經(jīng)過一定的潛伏期，并且潛伏期相同。

（3） 病毒進入易感者體內(nèi)經(jīng)過潛伏期后，部分易感者轉(zhuǎn)移為病毒變異前染病者，轉(zhuǎn)移率為ε，此時病毒在這過程中發(fā)生變異，部分轉(zhuǎn)化為病毒變異后染病者，變異率為σ。

其他參數(shù)是典型的傳染病動力學建模參數(shù)，其含義明確。β1表示易感者與潛伏者的有效接觸率；β2表示易感者與病毒變異前的染病者的有效接觸率；β3表示易感者與病毒變異后的染病者的有效接觸率；μ1表示潛伏者的隔離率；μ2、γ1、α1分別表示病毒變異前的染病者的隔離率、治愈率、因病死亡率；μ3、γ2、α2分別表示病毒變異后的染病者的隔離率、治愈率、因病死亡率；γ3、α3分別表示隔離者的治愈率、因病死亡率。從生物數(shù)學的角度考慮，規(guī)定所有參數(shù)都非負，所以傳染病模型倉室圖如圖1 所示。

圖1 模型的倉室結(jié)構(gòu)圖Fig.1 Compartment structure diagram of model

根據(jù)假設（1）—（3）與倉室結(jié)構(gòu)圖1，相應的傳染病模型如下：

2 平衡點的存在性與穩(wěn)定性

研究模型（1）的動力學性質(zhì)首先要確定其平衡點，且平衡點的穩(wěn)定性可以為傳染病的控制提供有力的參考價值，因此本節(jié)將討論平衡點的存在性與穩(wěn)定性。

2.1 平衡點的存在性

為方便起見，記δ1=μ+ε+μ1，δ2=μ+μ2+α1+γ1+σ，δ3=μ+μ3+α2+γ2，δ4=μ+α3+γ3。模型（1）的平衡點滿足方程組

易知模型（1）總存在無病平衡點P0(S0，0，0，0，0，0)，其中。根據(jù)文獻［15］中的再生矩陣法，可定義出模型（1）的基本再生數(shù)為，直接計算，當?0＞1 時存在唯一的地方病平衡點P*(S*，E*，I1*，I2*，Q*，R*)，其中

因此，可得如下關于平衡點存在性的定理。

定理1當?0≤1 時，模型（1）僅存在無病平衡點P0(S0，0，0，0，0，0)；當?0＞1 時，模型（1）既存在無病平衡點P0，還存在唯一的地方病平衡點P*(S*，E*，I1*，I2*，Q*，R*)。

2.2 無病平衡點的穩(wěn)定性

無病平衡點對應于疾病的消亡，無病平衡點穩(wěn)定說明易感者的人數(shù)最終會趨于一個常數(shù)，變異前后的染病者人數(shù)都將趨于零。本小節(jié)將研究無病平衡點的穩(wěn)定性。

定理2當?0＜1 時，無病平衡點P0局部漸近穩(wěn)定；?0＞1 時，無病平衡點P0不穩(wěn)定。

證明模型（1）在無病平衡P0處的Jacobi 矩陣為

顯然，J(P0)有六個特征根λi(i=1，2，…，6)，其中λ1=λ2=-μ，λ3=-δ4且λ4，λ5，λ6滿足方程λ3+a1λ2+a2λ+a3=0，這里

當?0＜1 時，

依據(jù)Routh-Hurwitz 準則，J(P0)的所有特征根均有負實部，所以無病平衡點P0局部漸近穩(wěn)定；當?0＞1 時，a3＜0，則J(P0)存在正實部特征根，因此無病平衡點P0不穩(wěn)定。

定理3當?0＜1 時，無病平衡點P0全局漸近穩(wěn)定；當?0=1 時，無病平衡點P0全局吸引。

證明構(gòu)造Lyapunov 函數(shù)

V1(t)關于模型（1）對t求導有

2.3 地方病平衡點的穩(wěn)定性

地方病平衡點穩(wěn)定說明疾病會持續(xù)存在，病毒變異前后的染病者最終趨于一個常數(shù)，疾病形成一種地方病。本小節(jié)將從理論上證明地方病平衡點的穩(wěn)定性。

定理4當?0＞1 時，地方病平衡點P*局部漸近穩(wěn)定。

3 模型應用

根據(jù)模型（1）的特點，本節(jié)將對印度的COVID-19 進行研究。

3.1 數(shù)據(jù)

本文搜集了印度2020 年12 月1 日至12 月30 日的累計病例數(shù)據(jù)［16］，如圖2 所示。

圖2 2020年12月1日至12月30日的印度累計病例數(shù)Fig.2 Cumulative number of cases in India from December 1,2020 to December 30,2020

3.2 模型參數(shù)估計

印度2020 年總?cè)丝?3.8 億人［17］，出生率1.873 8%［17］，因此Λ=25 858 440。人口的平均壽命為68.89 年［17］，可取μ=3.977×10-5，新冠的潛伏期一般為7 天［18］，所以取ε=0.14。假設S(0)=13.8×108，取變異前的染病者初值為11 月30 日的累計病例數(shù)，所以I1(0)=9 463 256。

根據(jù)印度12 月1 日至12 月30 日的COVID-19 累計病例數(shù)，對模型（1）的未知參數(shù)進行估計，參數(shù)的估計值見表1。模型擬合結(jié)果如圖3 所示。根據(jù)表1 中的參數(shù)，可得到變異后的染病者的人數(shù)變化趨勢，如圖4 所示。圖4 的結(jié)果表明變異后的染病者的人數(shù)隨時間的變化在不斷增長，由于變異病毒具有更強的傳染性，印度需對其重視起來，并采取相關應對措施，防止疫情的進一步擴散。

表1 模型（1）的參數(shù)值Table 1 The parameter values of model (1)

圖3 2020年12月1日至2020年12月30日累計病例數(shù)與模型（1）擬合Fig.3 The cumulative cases from December 1, 2020 to December 30, 2020 were fitted to model (1)

圖4 2020年12月1日至2020年12月30日I2(t)變化趨勢Fig.4 Change trend of I2(t) from December 1, 2020 to December 30, 2020

3.3 新冠流行趨勢

印度后期將接種新冠疫苗，疫苗的接種會影響累計病例人數(shù)，且后期新冠康復者會感染毛霉菌病，所以模型不能進行長時間段的預測，因此僅對模型（1）關于印度未來半個月的疫情情況進行預測。對印度2020 年12 月31 日至2021 年1 月14 日累計病例數(shù)的預測如圖5 所示。從觀察到的數(shù)據(jù)可以看到，印度疫情正在不斷加重，還須保持高度警惕，加強防控措施。

圖5 2020年12月1日至2021年1月14日累計病例數(shù)與模型（1）擬合的數(shù)據(jù)對比Fig.5 Cumulative cases from December 1, 2020 to January 14, 2021 were compared with model (1) fitted data

3.4 參數(shù)敏感性分析

理論分析表明，基本再生數(shù)?0是決定傳染病流行的一個重要參數(shù)。為了更好的分析各參數(shù)與?0的關系，找到控制傳染病傳播的理論依據(jù)，利用PRCCs（偏秩相關系數(shù)） 方法對模型（1）的?0進行不確定性分析與敏感性分析，以確定各參數(shù)對閾值?0的影響大小。若PRCCs 的值為正，則表示該參數(shù)的增加會導致?0的增大；若PRCCs 的值為負，則表示該參數(shù)的增加會導致?0的減小。PRCCs 的絕對值在 0～0.2 之間表明輸入?yún)?shù)與輸出變量之間不存在顯著相關性。PRCCs 的絕對值在0.2～0.4 之間表明輸入?yún)?shù)與輸出變量之間存在中度相關性。PRCCs 的絕對值在0.4～1 之間表明輸入?yún)?shù)與輸出變量之間具有較高的相關性［19］。本文選取8 個參數(shù)β1、β2、β3、μ1、μ2、μ3、γ1、γ2進行敏感性分析，分析結(jié)果如圖6 所示。從圖6 中可以看出，在這些參數(shù)中，β1、β2、β3對?0有正影響，μ1、μ2、μ3、γ1、γ2對?0有負影響，與?0有強相關性的參數(shù)為：易感者與潛伏者的有效接觸率β1、易感者與變異前的染病者的有效接觸率β2。該結(jié)果表明β1與β2對印度疫情的控制更加有效。

圖6 部分相關系數(shù)（PRCCs）的結(jié)果Fig.6 Results of Partial Rank Correlation Coefficients (PRCCs)

因為β1與β2的PRCC 值相近且討論方式相似，因此下面僅重點討論參數(shù)β1對變異前染病者I1(t)的累計病例數(shù)的影響，結(jié)果如圖7 所示。從圖7 中可以看到，當參數(shù)β1減小0.8 倍時，隨著時間的推移，變異前染病者的累計病例數(shù)將持續(xù)減少，且效果逐漸明顯。但由于印度疫情較為嚴重，疫情的控制還需很長時間。

圖7 β1變化對I1(t)數(shù)量的影響Fig.7 Effect of β1 change on the number of I1(t)

在實際的疫情防控中，將多種有效措施結(jié)合起來會產(chǎn)生更好的效果，更加利于疫情的控制。經(jīng)過分析可以得到幾個降低與控制該傳染病的有效方法：

（1） 控制易感者與潛伏者的有效接觸率β1其降低，控制易感者與變異前的染病者的有效接觸率β2使其降低，控制易感者與變異后的染病者的有效接觸率β3使其降低，控制潛伏者、變異前的染病者與變異后的染病者的隔離率μ1、μ2、μ3使其增大，即在面對傳染病時，國家應采取嚴格的防控措施，例如；“封城”“居家隔離”等，特別是針對變異前的染病人群，要求居民自覺居家隔離，不出門，不串門，使易感人群避免接觸染病者，以切斷疫情的傳播途徑。

（2） 控制變異前的染病者與變異后的染病者的治愈率γ1、γ2使其增大，即縮短染病者的治愈周期，因此在衛(wèi)生防疫及提高醫(yī)療水平方面要有所提升。根據(jù)上面的分析，參數(shù)γ1對疫情控制更加有效，因此在醫(yī)療方面，需更加注重對變異前的染病者的治療。

（3） 控制人口的常數(shù)輸入使其減小，即降低不同區(qū)域間人口的流動性。

4 結(jié)論

本文研究了一類潛伏期具有傳染性的SEI1I2QR 傳染病模型，考慮了傳染病在傳播過程中部分染病者的病毒發(fā)生變異，成為一類新的染病者。通過研究模型平衡點的存在性與穩(wěn)定性，得到如下結(jié)論：當?0≤1 時，傳染病滅絕；當?0＞1 時，傳染病蔓延。同時結(jié)合模型特點，將印度的COVID-19 作為傳染病例子進行數(shù)值模擬，結(jié)果表明2020 年12 月1 日至12 月30 日的擬合結(jié)果與實際COVID-19 的累計病例數(shù)據(jù)基本符合，并預測了2020 年12 月31 日至2021 年1 月14 日的疫情情況，結(jié)果表明在短時間內(nèi)印度疫情會越來越嚴重，且變異的染病者會越來越多。最后對模型的基本再生數(shù)?0中的部分參數(shù)進行了敏感性分析，分析結(jié)果表明采取降低易感者與染病者的有效接觸率、加強對人口的隔離、提高染病者的治愈率、減少人口流動性等措施可以降低疫情的傳播。

年齡對于傳染病的流行規(guī)律是一個很重要的因素，因為不同年齡段的人有著不同的死亡率，年齡也影響著傳染率與恢復率等疾病傳播因素，因此考慮具有年齡結(jié)構(gòu)的模型是非常有必要的。其次，時滯在傳染病的傳播過程中也扮演著一個重要的角色，可以用時滯模擬傳染病的潛伏期、患者對疾病的感染期和恢復者對疾病的免疫期，因此考慮具有時滯的模型可以減小預測與實際之間的誤差。我們將在將來的研究工作中考慮這些問題。