中國H5N1禽流感病毒傳播的動力學研究

高 瑞，張 娟

(1.山西大學 復雜系統研究所， 太原 030006;2.疾病防控的數學技術與大數據分析山西省重點實驗室， 太原 030006)

1997年，人感染H5N1禽流感病例首次在中國香港地區出現，之后又有多個國家相繼出現[1-2]。據世界衛生組織公布數據知，截止2017年12月，我國實驗室確診感染H5N1禽流感病毒累計病例數為53例，死亡人數為31例，死亡率高達58%[3-5]。據獸醫公報公布，截止2017年11月，H5N1高致病性禽流感病毒已蔓延至18個省份[6-7]，感染家禽宿主包括：雞、鴨、鵝、鵪鶉等，野生宿主包括：江湖海周圍的水鳥以及各種遷徙野鳥等[8-10]。由于H5N1高致病性禽流感病毒的傳播，已經有近690萬家禽被捕殺，給養殖業帶來毀滅性的打擊[11-12]。現階段暫無證據表明H5N1高致病性禽流感病毒能在人與人之間進行有效傳播，但隨著H5N1高致病性禽流感病毒的持續流行與進化，該病毒對人類仍存在嚴重威脅[13]。

目前，研究禽流感傳播與控制的動力學模型大多數以禽個體、人為研究對象，且種間只考慮了禽對人的單向傳播。而中國H5N1禽流感病毒的傳播特征為：傳播速度快，一只家禽攜帶病毒很快引起整個場群的感染；我國目前以疫點為單位上報政府并進行相應的捕殺；病毒在場群間主要通過養殖業相關人員在場群間的活動進行傳播，例如養殖業相關人員的衣著、攜帶的交通工具等，在禽的調運與銷售過程中，將會攜帶病毒，并在場群之間進行傳播與擴散。基于此，本文以場群為單位，以禽養殖與銷售場所、養殖業相關人員、普通人群為研究對象，既考慮禽養殖與銷售場所對人的感染，又考慮養殖業相關人員攜帶的病毒對禽養殖與銷售場所的感染，建立中國H5N1禽流感場群間傳播動力學模型，挖掘其內在傳播機制，評估預防控制其流行的有效措施。

1 動力學模型的建立

本文研究對象分為3類：全國禽養殖場與銷售場所、普通人群和養殖業相關人員。其中，禽養殖場與銷售場所指禽養殖場、活禽交易市場和零售市場；養殖業相關人員指養殖場與銷售場所的工作人員，即養殖、運輸以及銷售活禽的高危人群；普通人群指前往活禽交易市場或者零售市場購買禽的普通民眾。根據染病和免疫狀態，禽養殖場與銷售場所FN分為易感、染病和免疫3種狀態，分別記為FS、FI、FV，其中FN=FS+FI+FV；禽養殖場/銷售場所中一旦有一只禽染病，該場將認為處于染病狀態。根據染病狀態，普通人群HN分為易感和染病2種狀態，分別記為HS和HI，其中HN=HS+HI；根據衣著、攜帶工具、交通工具等是否攜帶病毒，將養殖業相關人員MN分為無攜帶病毒者、攜帶病毒者，分別記為MS和MI，其中MN=MS+MI。基于H5N1禽流感病毒的傳播特點，該疫病在各倉室之間的轉化見流程圖1，對應的動力學模型見式(1)。模型中涉及參數的具體意義見表1，其中所有參數都是非負的，且0<α<1， 0<γ<1。

(1)

圖1 H5N1禽流感病毒在人和禽養殖場與銷售場所之間傳播流程

表1 模型(1)中有關參數的說明

觀察模型(1)，可得

從而可以得出可行域

Γ={(HS,HI,MS,MI,FS,FI,FV)|HS,

HI,MS,MI,FS,FI,FV≥0,

2 模型的動力學分析

由于模型(1)的后5個方程與前2個方程無關，因此只分析后5個方程即可，組成系統(2)如下：

(2)

2.1 平衡點的存在性

令系統(2)方程右邊為0，可以解出模型的平衡點。

其中：

其中,

根據上述討論得出：當R0>1時，系統(2)存在正平衡點E*。

2.2 無病平衡點的局部穩定性

定理1若R0<1，那么系統(2)的無病平衡點E0是局部漸近穩定的。

證明系統(2)在無病平衡點E0處的雅克比矩陣為

可以計算出其特征多項式：

其中，

A=2δ+3d2+d+λ+φ,

(2δ+d+d2)(λ+φ+2d2)+

很明顯，-d是p(λ1)的一個根。

可以得出：

Δ1=A=2δ+3d2+d+λ+φ>0

Δ2=AB-C=(2δ+d+d2) [(2δ+d+d2)×

(2d2+λ+φ)+(δ+d)(δ+d2)×

Δ3=ABC-C2-A2D=

(2d2+λ+φ)(2δ+d+d2)×

(2δ+d+3d2+λ+φ)+(δ+d2)×

(2δ+d)2]}>0

Δ4=DΔ3>0

由Routh-Hurwitz判據，J|E0所有的特征根都具有負實部，因此系統(2)的無病平衡點E0是局部漸近穩定的。

3 防控措施評估

3.1 參數估計

1) 據中國統計年鑒，2004—2017年每年出生人口數的平均值為1 637萬[14]，因此A=1.637×107。

2) 從中國畜牧業年鑒可獲得2004—2017年全國蛋雞和肉雞的總養殖場數[15]。假設每年新增禽養殖場與銷售場所為2 000個，根據A2=d2FN，可推算出每年禽養殖場與銷售場所的關閉率為6%。

3) 資料顯示，H5N1禽流感滅活疫苗的免疫保護期可達半年以上[16]，本文假設為1年，則免疫失效率φ=1。

4) 接種免疫率為接種覆蓋率和疫苗有效率的乘積。資料顯示，H5N1禽流感滅活疫苗的疫苗有效率為68.4%[17]。本文假設疫苗的接種覆蓋率為80%，所以λ=0.55。

5) 據世界衛生組織公布數據知：截至2017年12月，我國實驗室確診感染H5N1禽流感病毒累計病例數為53例，死亡人數為31例，死亡率高達58%[3]，所以1-γ=58%。

6)δ代表恢復率，一年內感染者百分之百移出，要么恢復，要么死亡，因此取δ=1。

HS(0)=1.3×109×0.7,HI(0)=2

MS(0)=1.3×109×0.3,MI(0)=1

FS(0)=1.114 927 6×108×0.7

FI(0)=50

FV(0)=1.114 927 6×108×0.3

圖2 2004—2017年禽養殖場與銷售場所染病數據

3.2 敏感性分析

利用表1的數據和基本再生數的表達式，可以估計出中國H5N1高致病禽流感在場群間傳播的基本再生數R0=0.832 4。說明在目前的免疫與撲殺作用下，該疫病隨著時間推移將會得到控制。在實際中，H5N1抗原在持續的免疫壓力下不斷進化，發生基因變異，使得疫苗對其保護效果越來越差，即接種免疫率λ會降低，因此偶爾仍會有疫病暴發。現在討論疫苗免疫率對疫病流行的影響。為了評估免疫等措施對疫病流行的影響，討論R0關于各參數的敏感性。由圖3可知：當其他參數值不變時，隨著λ降低，基本再生數R0將會增大。當λ小于0.1時，R0將大于1。

接下來實施雙參數分析，來討論C2、C3、d2和λ對R0的作用。由圖4知：隨著C2和C3的減小，基本再生數R0將會降低。隨著d2和λ的增大，也會使基本再生數R0不斷地降低。由圖4(a)～(d)知：如果按照當前每年養殖場的關閉率d2=0.06，接種免疫率λ=0.55,C2需小于1 010,C3需小于40 384,才能保證R0<1。由圖4(e)知：曲線斜率比較大，說明R0對β3較敏感，即控制染病養殖業相關人員對場群的傳染力是非常重要的。通過以上分析得到：監控養殖業相關人員對場群的感染，加強養殖業相關人員及其攜帶工具的消毒，增加養殖場接種免疫率，可有效控制傳染力度。

圖3 λ對R0的作用

圖4 雙參數敏感性分析C2、C3 、d2、λ 、β2和β3對R0的作用

4 討論

基于H5N1禽流感病毒的傳播特征，本文以場群為單位，以禽養殖場與銷售場所、養殖業相關人員、普通人群為研究對象，既考慮禽養殖場與銷售場所和養殖業相關人員之間的雙向作用，建立中國H5N1禽流感場群間傳播動力學模型，挖掘其內在傳播機制，估計了中國H5N1禽流感傳播的基本再生數R0，評估預防控制其流行的有效措施。通過分析得到，在免疫與撲殺禽的作用下，中國H5N1禽流感傳播能暫時得到控制，然而，由于病毒變異重組，接種免疫率λ會降低，會導致基本再生數R0大于1，使得疫病再次暴發。進一步分析得到，監控養殖業相關人員對場群的感染，加強養殖業相關人員及其攜帶工具的消毒，增加養殖場接種免疫率，可有效控制H5N1禽流感病毒的傳播。在未來的工作中，將研究病毒變異與疫苗替換對疫病流行的影響。