李貽芬

傳染病模型是研究傳染病發(fā)展規(guī)律的科學(xué)，通過數(shù)學(xué)模型可以預(yù)測傳染病的未來發(fā)展.在理論研究中，數(shù)學(xué)模型起著極其重要的作用［1］；它把傳染病的主要特征通過假設(shè)、參數(shù)、變量和它們之間的聯(lián)系清晰地揭示出來，數(shù)學(xué)模型的分析結(jié)果能提供許多強有力的理論基礎(chǔ)；用數(shù)學(xué)模型幫助發(fā)現(xiàn)傳染病的傳播機理，預(yù)測傳染病的流行趨勢已成為共識［2］.

1 常用傳染病傳播模型

近年來，越來越多的研究者開始利用隨機模型來描述環(huán)境變化下的傳染病發(fā)展過程，確定性的傳染病模型一般都是倉室模型.隨著社會流動性的增強，傳染病的空間擴展呈現(xiàn)出新的模式，倉室模型很難處理，而來源于空間生態(tài)學(xué)的復(fù)合群體一叩方法在大空間尺度的傳染病傳播建模中得到了應(yīng)用，近年來隨著復(fù)雜性科學(xué)的興起，微觀建模蓬勃發(fā)展，與社會網(wǎng)絡(luò)相結(jié)合，發(fā)展出基于網(wǎng)絡(luò)的微觀個體建模方法，為認識傳染病的傳播規(guī)律提供了新的途徑.建立傳染病傳播模型主要有三種方法.

1.1 單一群體

單一群體模型，所有疾病發(fā)生、傳播發(fā)生在每個倉室之間.最流行的單一群體模型是倉室模型，所有處于相同狀態(tài)的人構(gòu)成一個倉室，隨著狀態(tài)的變化，人員在倉室之間移動.倉室模型的基本假設(shè)是每個人與相同人群均勻混合和接觸是瞬時的，接觸沒有記憶性，每個倉室的人口數(shù)量足夠大.因此在傳染病流行過程中，感染率，恢復(fù)率是常數(shù).因此，多數(shù)模型可以直接用（偏）微分方程描述.其中非常經(jīng)典的就是1927Kermack-McKendrick模型，即“SIR”和“SIS”模型［3］.

1.2 復(fù)合群體

復(fù)合群體是指一個相對獨立地理區(qū)域內(nèi)各局域種群的集合，各局域種群通過個體遷移而連為一體，在研究傳染病傳播時，復(fù)合群體模型將人群看作由定義明確的社會單元形成的空間結(jié)構(gòu)化人口，各單元之間因個人移動產(chǎn)生聯(lián)系.每個社會單元是一個子群體，所有社會單元形成復(fù)合群體.復(fù)合群體模型多種多樣，但有兩個基本特征：一是子群體之間的藕合方式，二是子群體內(nèi)部動態(tài)的表達.

1.3 微觀個體（網(wǎng)格網(wǎng)絡(luò)模型）

單一群體模型或復(fù)合群體模型的子群體中，一般假設(shè)個體同質(zhì)，人群均勻混合，但實際上個體只能與有限個體接觸，個體接觸模式差異較大，用網(wǎng)絡(luò)描述接觸模式更符合實際.比較典型的網(wǎng)絡(luò)動力學(xué)模型有：元胞自動機、人工神經(jīng)網(wǎng)絡(luò)、隨機布爾網(wǎng)絡(luò)、無標度網(wǎng)絡(luò)等模型，其中元胞自動機、人工神經(jīng)網(wǎng)絡(luò)和無標度網(wǎng)絡(luò)是目前流行病數(shù)學(xué)模型中研究較多應(yīng)用較廣的模型.

2 改進優(yōu)化的傳染病傳播新模型

常用傳染病傳播模型關(guān)注的重點：一是傳染病在子群體之間的流轉(zhuǎn)對復(fù)合群體整體動態(tài)的影響，這類模型相對簡單，容易分析，對理解傳染病在復(fù)合群體中的傳播有一定幫助，但難以對真實系統(tǒng)進行分析和預(yù)測，二是顯式表達子群體內(nèi)部動態(tài)，多采用較簡單的模型，如果這類模型能夠刻畫真實系統(tǒng)，但有大量的微分方程，求解困難，只能對簡單結(jié)構(gòu)的復(fù)合群體模型進行解析研究，對于具有高度真實性的模型，需采用仿真方法［4］.針對常用傳染病傳播模型的不足，我們的改進點如下.

（1）在一張任意給定的地圖上構(gòu)造大范圍的城鎮(zhèn)—城鎮(zhèn)、鎮(zhèn)內(nèi)倉室—倉室的感染理論模型，并進行編程模擬分析，獲得相應(yīng)數(shù)據(jù).

（2）運用涵蓋了單一群體方法和復(fù)合群體方法，并且使用數(shù)組迭代和（偏）微分方程結(jié)合的方法，克服復(fù)合群體網(wǎng)絡(luò)中模型個體過于簡單、求解微分方程困難的問題，并且可以獲得長時間的傳染數(shù)據(jù).

（3）新引入了一個動態(tài)模型，考慮受感染后城市可能采取的隔離措施，這將導(dǎo)致人口的交流頻率減少，將城市分為兩種狀態(tài)，分別是加強隔離前后的兩個狀態(tài)，這個模型滿足馬爾科夫鏈的條件，我們構(gòu)造兩個狀態(tài)間的概率矩陣和狀態(tài)轉(zhuǎn)換矩陣，每運行一步概率矩陣都乘以狀態(tài)轉(zhuǎn)換矩陣，然后利用極大似然的思想判斷之后的狀態(tài)會怎么變化，即城市是否會采取加強隔離的措施［5］.

2.1 構(gòu)建隨機地圖

生成有限網(wǎng)絡(luò)：生成的每一個點代表一個城市，每個城市之間的人員交流用外接矩陣描述.本程序的結(jié)果是輸出一個鄰接矩陣和一個描述城市之間交流關(guān)系的三維點狀圖.鄰接矩陣是由一個低階的矩陣逐步生成的，低階矩陣有三種模式可選，一種是零矩陣，表示初始都孤立，第二個表示點都相互連接，第三個是隨機矩陣，他們對應(yīng)的鄰接矩陣分別是比較稀疏、密集、隨機.

以某方式賦予每個點一個概率值.

在生成下一個點的連接關(guān)系時，對概率進行求和，和一個隨機數(shù)比較，由大小來判斷是否連接.在輸出的三維點狀圖中，有交流的城市以直線相連，城市在圖中的位置和實際的幾何關(guān)系無關(guān)，含有連線則表示兩城市間可以流通人口，無連線表示不能直接相連（可以通過其他城市間接相連）.設(shè)定城市個數(shù)=10，初始矩陣=隨機矩陣.

三維點狀圖如圖1所示.

圖1 三維點狀圖

2.2 城市和城市間內(nèi)部模型構(gòu)建

采用與復(fù)合模型一致的方法，感染人群服從SIR模型.我們把一個城市內(nèi)的人口分為患病人群x1、易感人群x2、免疫人群x3、死亡人群x4；a、b、c分別為傳染率、治愈率、死亡率，則滿足dx1/dt=ax1x2-bx1(x2+x3)-cx1，dx2/dt=-ax1x2，dx3/dt=bx1(x2+x3)，dx4/dt=cx1.

城間傳播模型建立：考慮到現(xiàn)有模型在解決復(fù)合群體問題時微分方程復(fù)雜性過大，我們采用數(shù)列的方式來規(guī)避過大運算量，即以天為基本單位，對于相鄰的數(shù)天采取數(shù)列關(guān)系時運算，可以保證一定精度的同時，大幅減少復(fù)雜度.

考慮每一天不同城市的人口交流，人口交流的情況用鄰接矩陣刻畫，存在人口交流的城市滿足：i(t+1,m)=i(t,m)-ci(t,m)+ci(t,n)，i(t+1,n)=i(t,n)-ci(t,n)+ci(t,m)，其中，i(t,m)表示m城市第t天中某種人群的人數(shù)，c是人口遷移率，我們認為病人和健康人的遷移率不同.本模型的關(guān)鍵是將這兩種傳染的方式耦合求解，得到每個城市的實時傳染情況.

用ode45函數(shù)來解SIR模型的微分方程，輸入初始的傳染情況，就能得到一天后的傳染情況，取傳染率=0.80，治愈率=0.02，死亡率=0.01.

計算每個城市實時傳染情況：假定最后一個城市先爆發(fā)疾病，除該城市外，初始每個城市易感人群為1個單位，其余3個人群數(shù)量為0，最后一個城市有0.02個單位的人患病，0.98個單位的人為易感.先計算每一天城市內(nèi)的傳染情況，再考慮城際的人口交流，總時間為240天，最后會得到四種人群在每個城市的實時人數(shù)的矩陣.

2.3 考慮政府干預(yù)的模型

3 優(yōu)化的傳染病傳播新模型仿真計算結(jié)果

3.1 無政府干預(yù)情況

（1）城市傳染情況.運算結(jié)果輸出了4個240*10的矩陣，橫軸表示城市，縱軸表示天數(shù)，繪制出每個城市傳染病情的實時變化情況，如圖2～圖5所示（圖中為易感人群，為感染人群，為抗體人群，為死亡人群）.

圖2 第4個城市的實時情況

圖3 第3個城市的實時情況

圖4 第2個城市的實時情況

圖5 第1個城市的實時情況

從圖2～圖5可以看出，在每一個城市，四種人群的比例都趨于一定值，易感人群為0，感染人群在 0～0.1，抗體人群在 0～0.2，死亡人群在0.7～0.9，不同城市的變化趨勢基本相同，但在時間上有延遲，不同城市的病情的差異也源于這個時間上的延遲，而這個時間上的延遲和城市的空間分布有關(guān).

（2）距離因素分析.定義兩城市的距離為其最短路徑，最短路徑可以通過鄰接矩陣用dijkstra算法求出.通過這個矩陣，我們可以研究易感率、感染率、抗體率、死亡率和到發(fā)病城市的距離的關(guān)系，如圖6～圖9所示.

從圖6～圖9可以看出，距離為1的城市病情已經(jīng)趨于穩(wěn)態(tài)，感染人群很少，預(yù)計很快會降為0，易感人數(shù)已經(jīng)為0，抗體人數(shù)較少，之后會有小幅度上升，死亡率較大，但之后基本不變.距離為2的城市情況和前者類似，之后的一段時期會重復(fù)前者的情況，但病情沒有那么劇烈.距離為3的城市由于空間距離遠，時間的延遲性強，體現(xiàn)了低感染率，高易感率，低抗體率，低死亡率，預(yù)計之后的情況和前兩種的變化趨勢一致，但由于其他城市的影響較小，現(xiàn)有病情較輕，所以達到穩(wěn)態(tài)后的病情也不會比前兩者嚴重.

圖6 易感率與網(wǎng)絡(luò)上傳染范圍的柱形圖

圖7 死亡率與網(wǎng)絡(luò)上傳染范圍的柱形圖

圖8 抗體率與網(wǎng)絡(luò)上傳染范圍的柱形圖

圖9 感染率與網(wǎng)絡(luò)上傳染范圍的柱形圖

研究四種人群數(shù)量和網(wǎng)絡(luò)的關(guān)系，xy平面是每個城市的位置，連線表示兩個城市間有交流，z軸為某一個人群的人數(shù)，在各個城市中，藍色為初始的爆發(fā)城市.結(jié)束時死亡數(shù)量、結(jié)束時抗體數(shù)量、結(jié)束時感染數(shù)量、結(jié)束時易感數(shù)量與地圖形狀關(guān)系的三維點狀圖如圖10～圖13所示.

從圖10～圖13可以看出，病情最嚴重的城市主要是和爆發(fā)城市直接相連的城市和各個城市都有連接的大城市，嚴重來說后者的病情比前者更嚴重，這就可以解釋現(xiàn)實中大城市是傳染病的高發(fā)區(qū)這一現(xiàn)象了.

圖10 結(jié)束時死亡數(shù)量和地圖形狀關(guān)系的直觀表達

圖11 結(jié)束時抗體數(shù)量和地圖形狀關(guān)系的直觀表達

圖12 結(jié)束時感染數(shù)量和地圖形狀關(guān)系的直觀表達

圖13 結(jié)束時易感數(shù)量和地圖形狀關(guān)系直觀表達

3.2 有政府干預(yù)情況

每個城市傳染病情的實時變化情況，如圖14～圖17所示.易感率、感染率、抗體率、死亡率和到發(fā)病城市的距離的關(guān)系，如圖18～圖21所示.結(jié)束時抗體數(shù)量、死亡數(shù)量、感染數(shù)量及易感數(shù)量和地圖形狀關(guān)系如圖22～圖25所示.

圖14 第四城市的實時情況

圖15 第三城市的實時情況

圖16 第二城市的實時情況

圖17 第一城市的實時情況

圖18 抗體率與網(wǎng)絡(luò)上傳染范圍的柱形圖

圖19 感染率與網(wǎng)絡(luò)上傳染范圍的柱形圖

圖20 死亡率與網(wǎng)絡(luò)上傳染范圍的柱形圖

圖21 易感率與網(wǎng)絡(luò)上傳染范圍的柱形圖

圖22 結(jié)束時抗體數(shù)量和地圖形狀關(guān)系的直觀表達

圖23 結(jié)束時死亡數(shù)量和地圖形狀關(guān)系的直觀表達

圖24 結(jié)束時感染數(shù)量和地圖形狀關(guān)系的直觀表達

圖25 結(jié)束時易感數(shù)量和地圖形狀關(guān)系的直觀表達

考慮了政府的隔離措施后，很明顯，病情得到了抑制，對一些連接較多的大城市和直接與爆發(fā)城市相連的城市而言，依然迅速爆發(fā)了疫情，但對于其他城市，疫情隔了接近200天才爆發(fā)，此時疫情比起不隔離時更弱了，更關(guān)鍵的是，經(jīng)過了200多天的疫情，人們對疫情有了基本的認識，有了一定的治療方法，當200天爆發(fā)時，真實的病情將遠弱于此模擬結(jié)果.可見政府的隔離很有效.

4 結(jié)論

在隨機地圖上模擬的疾病傳播過程，利用倉室基本模型，結(jié)合SIR模型，并且利用數(shù)組迭代方法模擬城市間的傳播，獲得了詳細的數(shù)據(jù)，并且從理論上說，可以適用于多個倉室的模擬，可以大大增加SIR模型的復(fù)雜度.新引入的一個動態(tài)模型，滿足馬爾科夫鏈的條件，我們構(gòu)造兩個狀態(tài)間的概率轉(zhuǎn)換矩陣，每運行一步都乘以這個矩陣，然后通過極大似然法處理之后的狀態(tài)會如何變化.模擬表明城市采取加強隔離措施，并且對疾病的惡化傳播有明顯抑制.