王明齋,余珊珊,芮佳,楊蒙,王瑤,王琦琦,陳田木,鄭蓉蓉?,張小芬?

結(jié)核病是由結(jié)核分枝桿菌引起的一種呼吸道疾病[1]。近年來我國結(jié)核病疫情形勢雖有所下降，但是結(jié)核病特別是耐藥結(jié)核病負擔仍很嚴重[2]，與世界衛(wèi)生組織提出的消除結(jié)核病目標仍有較大差距[3]。由于結(jié)核病在3 ～6 月存在一個明顯的發(fā)病高峰[4，5]，疾病防控人員可以在流行高峰前提前采取干預措施以控制疾病流行。因此，對各地結(jié)核病流行高峰進行早期識別和預警具有重要意義。近年來，已經(jīng)越來越多的數(shù)學模型運用到結(jié)核病預測預警研究中，常見的基于數(shù)理統(tǒng)計方法的預測模型有灰色GM(1，1)模型[6]、自回歸移動平均模型[7]、Azuma 模型[8]等。基于微分方程的模型有SEIT 模型[9]等。

Logistic微分方程模型(以下簡稱LDE模型)早期主要用于描述人口的增長模式，后被陸續(xù)應用在醫(yī)學等多個領(lǐng)域。該模型由于實操性強、操作簡便，在手足口病[10]、流感[11，12]、流行性腮腺炎[13]、感染性腹瀉疾病[14]等急性傳染病中已經(jīng)有一定的預測應用，但在慢性傳染病中的研究尚未見報道。本研究嘗試利用該模型對結(jié)核病進行模擬研究，為開展類似研究和指導防控提供參考。

原理與方法

1模型原理

模型模擬兩條曲線，一條為累計病例曲線(圖1中的虛線)，另一條為新發(fā)病例曲線(圖1中的實線)。其中，模型描述的新發(fā)病例曲線出現(xiàn)先增加，達到高峰后逐漸下降；而疾病累計發(fā)病數(shù)隨著時間推移，呈現(xiàn)出“慢—快—慢”的過程，整體曲線呈現(xiàn)出“S”型(圖1)。

圖1 Logistic微分方程模型累積和新發(fā)結(jié)核病病例數(shù)曲線(虛擬數(shù)據(jù))

LDE模型可以精確計算疫情加速時間，進而判斷流行周期的早期、中期和晚期節(jié)點，這有助于實施針對性防控措施。因此，對于疾控部門應對結(jié)核病等公共衛(wèi)生事件中具有重要參考意義。Logistic微分方程模型如下:

其中，dy/dt為t時刻傳染病累計病例y的變化速率。y為累計病例數(shù)，x為病例變化速率參數(shù)，Y為累計病例上限。方程(1)的通解為:

其中，a為微分方程求解過程產(chǎn)生的常數(shù)。將方程(2)分別求二階、三階導數(shù)，可以獲得模型的3個拐點，其橫坐標分別為:

結(jié)核病的疫情發(fā)展被t1、t2、t3這三個節(jié)點分為漸增期、快增期、緩增期3個階段。在漸增期采取相應措施防控效果最優(yōu)。而越往后的時期不僅成本高昂，而且防控難度大、效果不佳。因此，Logistic微分方程模型能夠分期的這一優(yōu)勢無疑具有重要的流行病學意義。t1為“疫情加速時間”，在該時間點結(jié)核病流行周期內(nèi)的疫情發(fā)展速度由慢變快。鑒于落實衛(wèi)生決策和實施干預措施需要時間，若放任疫情發(fā)展到“疫情加速時間”時才進行預警會導致滯后。因而，鑒于以往研究和專家意見，將“疫情加速時間”提前1-2個標準差作為“建議預警時間”，如圖2 所示。該模型詳細推導過程參照文獻[10-14]。

圖2流行時期內(nèi)結(jié)核病預警時間的確定方法示意圖(s為標準差，設(shè)為1個月)

2建模步驟

Logistic微分方程模型建模步驟包括:確定研究目的、分析流行特征、建立模型、收集數(shù)據(jù)、參數(shù)估計與初始值設(shè)置、模型擬合優(yōu)度檢驗、疫情加速和建議預警時間計算。

結(jié)果

1確定研究目的

以廈門市2005-2018年結(jié)核病月發(fā)病例數(shù)為例進行模型研究。計算廈門市2005-2018年各年間結(jié)核病疫情加速時間和建議預警時間。

2分析流行特征

由圖3可知，廈門市每年結(jié)核病疫情存在一個明顯的冬春季節(jié)流行高峰。

3建立模型

如上述方程(1)所示模型。

4收集數(shù)據(jù)

搜集廈門市2005-2018年月報告發(fā)病數(shù)，數(shù)據(jù)包含14年共168個月。

5參數(shù)估計與初始值設(shè)置

采用Berkeley Madonna 軟件進行模型和數(shù)據(jù)擬合，微分方程求解方法為4階龍格庫塔法，容忍度為0.001，模型收斂方法為最小均方根(least rootmean square，LRMS)。經(jīng)模型擬合，得到相應參數(shù)x、Y、a置為1月，則冬春季建議預警時間為上年度末月或下年度首月。值。如表1所示。

圖3 Logistic微分方程模型擬合效果

表1廈門市2005-2018年結(jié)核病疫情Logistic微分方程模型擬合結(jié)果

表2廈門市2005-2018年結(jié)核病疫情預警時間

根據(jù)方程(2)，求得y值。再由方程(3)和(5)求得t1和t3值。T1為疫情加速時間，T3為疫情放緩時間，以月為單位。如2005年的T1=2，T3=9表示2005年2月廈門市結(jié)核病疫情出現(xiàn)加速而9月開始放緩，如表2所示。

6模型擬合與優(yōu)度檢驗

圖3展示了實際發(fā)病率與Logistic微分方程擬合的效果。采用SPSS 13.0軟件進行決定系數(shù)計算，結(jié)果顯示擬合效果較好(R2=0.577，P＜0.05)。

7疫情加速時間和建議預警時間計算

t1為疫情加速時期，其值為1 ～2月，因而夏門市結(jié)核病疫情加速時間為第1或第2 個月。若s設(shè)

討論

本研究采用Logistic微分方程模型，結(jié)合結(jié)核病流行特征，建立各個周期內(nèi)的Logistic微分方程模型，利用實際周發(fā)病數(shù)據(jù)計算疫情發(fā)展速度參數(shù)x和流行周期結(jié)束時的累計病例數(shù)的理論值Y，并確定了結(jié)核病的“疫情加速時間”和“建議預警時間”，這為結(jié)核病控制工作提供了一定參考。評價模型適用性時，結(jié)合既往學者研究結(jié)果推薦擬合優(yōu)度檢驗，它基于Excel 和SPSS軟件即可進行，對于基層人員而言更加容易理解、易于推廣。未來在進行實際建模研究時可以加入疾病傳播動力學特征、干預的效果、氣候因素等數(shù)據(jù)，從而為結(jié)核病的防控提供更精細的證據(jù)。

Logistic微分方程模型總體思路是:首先利用結(jié)核病發(fā)病數(shù)據(jù)分析其流行周期規(guī)律、明確流行高峰，隨后采用數(shù)學模型擬合，分析每個流行周期內(nèi)的疫情加速時間以確定建議預警時間。Logistic微分模型在現(xiàn)實防控工作中具有獨特的優(yōu)越性，主要表現(xiàn)在兩個方面:一是模型使得基層在疫情評估時更加可靠。在我國基層公共衛(wèi)生單位，受專業(yè)水平的限制，靠傳統(tǒng)的專家疫情評估方式實施有一定的難度；二是該模型簡便易操作。目前研究的數(shù)學建模方法相較均復雜，模型仿真結(jié)果需要精確判斷，且未能定量給出疫情具體的分期。但該模型也有著自身的缺陷，主要體現(xiàn)在其對數(shù)據(jù)資料要求較為嚴格，如要求疾病高峰出現(xiàn)時其左右波形要對稱分布。而在現(xiàn)實中，當結(jié)核病疫情呈現(xiàn)上升趨勢，相關(guān)部門勢必會采取干預措施，從而導致疾病發(fā)病例數(shù)比自然狀態(tài)下減少，左右兩邊會不對稱，進而擬合效果受影響。這時需要嘗試引入廣義Logistic微分方程模型進行校正，從而提高模型預測精度。

盡管如此，Logistic微分方程模型仍能在一定程度上幫助基層工作者進行結(jié)核病疫情預警工作，通過明確疫情早期速度由慢變快的時間點，為科學采取防控措施提供有益參考。