黃春賢，周效良

(1. 閩南師范大學 數(shù)學與統(tǒng)計學院，福建 漳州 363000; 2. 嶺南師范學院 數(shù)學與統(tǒng)計學院，廣東 湛江 524048)

傳染病嚴重威脅人類的生存和發(fā)展。14世紀中葉，黑死病經(jīng)過幾次爆發(fā)，由亞洲傳到歐洲，三分之一的歐洲人因此而死亡。1720年至1722年間爆發(fā)在法國的腹股溝腺炎瘟疫令馬賽失去了一半的人口[1]。在現(xiàn)代，人們衛(wèi)生意識和社會醫(yī)療水平已今非昔比，但過度藥物治療使病毒產(chǎn)生抗藥性等問題，往往導致新型病毒的產(chǎn)生，這使得對傳染病的防治仍需嚴陣以待。因此，傳染疾病的控制仍然是現(xiàn)代社會亟待解決的問題[2-3]。

研究傳染病模型能為傳染疾病的控制提供有效的措施[4]。Kermack等[5]于1927年應用動力系統(tǒng)原理在假設人口不同類成員之間的流率相對簡單的基礎上建立經(jīng)典的SIR傳染病模型；隨后一些學者在此基礎上考慮了SIS模型[6-8]。SIR模型和SIS模型均是在極端情況下建立的，SIR模型僅考慮患病個體經(jīng)過治療必定康復成為具有免疫能力的個體[9]；SIS模型只考慮患病個體完全沒有獲取免疫能力。Gomes等[3]就該問題建立了一類非理想免疫的傳染病模型，通過分析平衡點的性質，給出一些降低傳染病流行的疫苗接種措施；Gomes等[10]又考慮一類具有非理想免疫和非線性發(fā)生率的傳染病模型，得到該模型參數(shù)穿越對應臨界值時會有后向分岔、振動和Bogdanov-Takens分岔。為了更好地描述傳染病的傳播，馬知恩等[11]建立了一類SIRS模型，該模型考慮患者經(jīng)過治療，可能獲取完全的免疫能力或者沒有獲得免疫能力而再次成為易感者。在某些情況下，傳染病的發(fā)生率并不是雙線性的，因此，Lu等[12]考慮了一類具有廣義非單調和飽和發(fā)生率的SIRS模型，得到了該模型隨著參數(shù)的變化會發(fā)生鞍結分岔、Bogdanov-Takens分岔、Hopf分岔等。

White等[13]建立了一類海洛因毒品傳播模型，并給出該模型無海洛因吸食患者平衡點的局部穩(wěn)定性分析以及有吸食患者平衡點的存在性；Mulone等[14]通過應用特征方程和龐加萊-本迪克松定理對文獻[13]的模型的全局穩(wěn)定性進行了分析；海洛因吸食個體通過治療成為康復個體需要一定的時間周期，因此Liu等[15]以及Huang等[16]都考慮了分布式的時滯項，前者得到了無海洛因吸食患者平衡點的全局穩(wěn)定和有吸食患者平衡點的局部穩(wěn)定條件，后者得到有吸食患者平衡點的全局穩(wěn)定性的條件。Muroya等[17]對輕度海洛因傳播的情況提出如下帶有等級治療率以及不完全康復率SIRS模型

(1)

得到模型(1)的基本再生數(shù)為

(2)

并給出了無海洛因吸食患者平衡點的全局穩(wěn)定和有吸食患者平衡點的全局穩(wěn)定條件分別為R0<1和R0>1，式中：S、I、R分別是總人群中易感個體的比例、患病個體的比例和康復個體的比例；d是新出生的人口比例，假定出生人口比例等于死亡人口的比例；β是易感個體與患病個體接觸成為新的患病個體的概率；ε是等級治療率；δ表示康復個體因為沒有獲得抗體而轉化成易感個體的比率；γ是患病個體的成為康復個體的比率；σ代表康復個體因為不完全治療而再次成為患病個體的比率。該系統(tǒng)所有參數(shù)均是非負常數(shù)。

2016年，史學偉等[18]在模型(1)基礎上加入信息變量，得到了系統(tǒng)存在Hopf分岔；2017年，Ma等[19]在模型(1)的基礎上，進一步考慮非線性發(fā)生率以及海洛因吸食者不能通過自身意志克服毒癮再次成為易感者，對一類具有非線性發(fā)生率的海洛因模型進行了分岔分析，得到系統(tǒng)存在Bogdanov-Takens分岔。

本文研究模型(1)的分岔性質。在系統(tǒng)(1)的極限集S+I+R=1上考慮系統(tǒng)(1)，可以得到系統(tǒng)(1)的約化系統(tǒng)[20-21]如下：

(3)

約化系統(tǒng)(3)的穩(wěn)定性與原系統(tǒng)(1)的穩(wěn)定性等價。為方便表示，令θ1=d+ε+γ，θ2=d+δ+σ。

1 平衡點分析

本章研究系統(tǒng)(3)平衡點的類型。系統(tǒng)(3)僅存在2個平衡點E0(0, 0)和E1(I*,R*)，其中：

(4)

(5)

定理1系統(tǒng)(3)的無病平衡點E0(0, 0)只為下面3種情況之一：

證明系統(tǒng)(3)在平衡點E0(0, 0)處的Jacobi矩陣為

(6)

特征方程為

λ2+a11λ+a12=0，

(7)

式中a11=θ1+θ2-β，a12=-βθ2+θ1θ2-σγ。易知

(θ1+θ2)2-2β(θ1+θ2)+β2+4βθ2-4θ1θ2+4σγ=

β2-2β(θ1-θ2)+(θ1-θ2)2+4σγ=

[β-(θ1-θ2)]2+4σγ>0。

因此，特征方程(7)存在2個不同的實特征根。根據(jù)韋達定理得：

定理2系統(tǒng)(3)的地方病平衡點E1(I*,R*)只為下面3種情況之一：

證明系統(tǒng)(3)在平衡點E1(I*，R*)處的Jacobi矩陣為

(8)

特征方程為

λ2+a21λ+a22=0。

(9)

把式(4)和式(5)代入得：

2 分岔分析

本章主要研究系統(tǒng)(3)的分岔性質。定義系統(tǒng)(3)的無病平衡點與地方病平衡點的距離為

d=‖E1-E0‖2，

(10)

(11)

引理1[22]考慮帶參數(shù)φ的常微分方程組

(12)

假定0是系統(tǒng)的平衡點，即對任意φ滿足：f(0,φ)≡0。

②對應于零特征值，A有一個非負的右特征向量ω和一個左特征向量υ，記

(13)

(14)

則當a<0，b>0時，該系統(tǒng)發(fā)生前向分岔。該分岔的分岔圖如圖1所示。

圖1 a<0,b>0時系統(tǒng)的前向分岔Fig. 1 Bifurcation diagram of forward bifurcation for a<0, b>0

(15)

可解得特征值為

(16)

因此，Jacobi矩陣J(E0,β*)有一個零特征值，且其他特征值均為負實部的。

設W=(w1,w2)T是Jacobi矩陣J(E0,β*)的右特征向量，滿足J(E0,β*)W=0，即

(17)

解得

W=(θ2,γ)T。

(18)

設V=(v1,v2)T是Jacobi矩陣J(E0,β*)的左特征向量，滿足VJ(E0,β*)T=0，解得

V=(θ2,σ)Τ。

(19)

由系統(tǒng)(3)得到：

(20)

計算該系統(tǒng)在無病平衡點E0的二次偏導數(shù)：

易知：

(21)

(22)

由引理1，定理3得證。證畢。

3 數(shù)值模擬

本章根據(jù)以上理論分析，通過數(shù)值模擬，給出系統(tǒng)(1)的相圖和狀態(tài)變量隨時間變化的曲線圖，系統(tǒng)的參數(shù)選為下面2種情形：

圖2 β<β*時系統(tǒng)(1)的相圖Fig. 2 Phase diagram of system (1) for β<β*

圖3 系統(tǒng)(1)初始值為(0.45, 0.28)的曲線Fig. 3 Curves of system for initial condition (0.45, 0.28)

圖4 β>β*時系統(tǒng)(1)的相圖Fig. 4 Phase diagram of system (1) for β>β*

圖5 系統(tǒng)(1)初始值為(0.004 5, 0.002 8)的曲線Fig. 5 Curves of system for initial condition (0.004 5, 0.002 8)

圖6 系統(tǒng)(1)初始值為(0.21, 0.13)的曲線Fig. 6 Curves of system for initial condition (0.21, 0.13)

4 結語

本文通過對一類具有等級治療率及不完全康復率的SIRS模型的平衡點類型進行分析，得到β<β*時，系統(tǒng)僅存在一個穩(wěn)定的結點E0(0, 0)；當β>β*時，系統(tǒng)僅存在一個不穩(wěn)定的鞍點E0(0, 0)和穩(wěn)定的平衡點E1(I*,R*)，發(fā)現(xiàn)該系統(tǒng)存在前向分岔，并對其進行嚴格證明；最后數(shù)值模擬驗證系統(tǒng)參數(shù)β越過臨界參數(shù)值β*時系統(tǒng)相圖的拓撲結構發(fā)生改變。