文章專區

2020-03-112019冠狀病毒病何時結束?以動力學模型估算「疫情結束日」 459 期

Author 作者 徐丞志/臺灣大學化學系助理教授。
(原文發表於作者個人臉書

本文主要作教育用途,嘗試以高中與大學教科書中的觀點,讓有基礎理化背景的民眾理解疫情的可能進展。本文所提出的模型非常簡陋,許多假設可能與實際狀況差異甚遠,不能當作真正的流行病學預測,也未經同儕審查(peer review),亦不能作為正式的學術報告。

 

「老師,請問武漢肺炎疫情什麼時候才會結束?」

這是我學生在今(2020)1月29日問我的一個問題。當時全世界都還稱2019冠狀病毒病(COVID-19)作武漢肺炎。

我當下回應是,我不是公衛和流行病學專業,所以無法回答。不過學生問的這個問題,從化學課本中教過的反應動力學知識或許可以嘗試推敲出來。也因為這個問題,我隔(1月30)日在個人的Facebook上發表了一篇以科學教育為出發點的文章,討論該疫情在中國大陸境內的可能走勢,其中一個預測曲線因與後來官方公布的確診數字相當接近,引起國內外媒體的討論,甚至讓《經濟學人》(The Economist)在2月13日的報導中引用該數據(圖一)。在此,我把當日發表的原文稍作整理與改寫,並與學術界的同儕們分享討論。

 

圖一:徐丞志以簡單的動態模型,預估中國大陸在COVID-19 疫情中確診人數的走勢,並受《經濟學人》引用刊登。

流行病傳播時,隨著時間變化,被感染的人數可能是多少?關於這個問題,流行病學家從上世紀前期就已參考過往經驗,提出各 種數學模型來描述。以嚴重急性呼 吸道症候群(SARS)為例,2003年就有至少十篇的統計模型被提出來,隔年還有人寫了文獻回顧整理。這些模型所需要的數學工具較複雜,非相關專業背景的人不容易知其究竟,筆者在此取用流行病學早期最經典的易感-感染-恢復模型(susceptible-infected-recovered model, SIR)的前半部,並從化學反應的碰撞理論作為起點,以較直觀的方式描述傳染病的傳播方式以及被感染人數隨時間的變化。接著,再嘗試以過去經驗和其他專家的意見,來「預測」本次新冠肺炎在中國大陸境內可能的走向。
 

簡化動力學模型

首先可想像,在一個空間裡有兩種人:一種是健康的人(用A表示),另一種是不幸被感染的人(用B表示);A和B在這個空間中會隨機亂走,因此有一定機率會彼此相碰。當A碰A或B碰B時,並不會有反應,但是當A碰到B的時候,則被感染的B有一定機率會將健康的A也變成B(將病原傳給A)。

把上述的情況用一條不可逆的化學反應式來代表:
A(一位健康者)+B(一位被感染者)→2B(兩位被感染者)

此反應速率會跟空間中健康者A和被感染者B的數量分別呈正相關,因此當一群A裡面忽然有了B,隨著時間演進,整個群體會逐漸都變成B。而B的數量增加最快的時候,並不是最剛開始出現B時,而是A和B數量差不多的時候(因為兩者相碰的頻率最高)。因此,起初B人數不多時,人數增加得並不快,而當跨過最高峰後因為可被感染的人變少了,B所增加的速率也就會接著降低。

理論上,在不考慮健康者A自然死亡的情況下,A+B的數量會等於區域內的總人口數。但若在疫情發生後進行隔離防疫,相當於限制B可以觸及A的機會和範圍,因此A+B的總數量會隨著防疫的落實而大幅減少,形成一個初始的邊界條件,並決定最後的總感染人數。

在此極為簡化的模型中,每日新增的被感染者(delta [B])即化學的二級反應式,寫成數學式就是:
delta [B] = 「當日患者數」x 「當日健康者數」 x 常數([A]*[B]* Const) 

用這樣的式子就可以計算,隨著疫情演進在不同時間點被感染者的增加速率。SIR模型的前半段,其實就是化學的「自催化反應」,兩者間的動力學行為完全相同。但是有幾點需要特別注意:

1.此模型假設人一旦被感染就具有傳染力,即使是在潛伏期也是。
2.因為從感染後進入潛伏期到出現症狀,以及檢體送樣後到確診,都需要一點時間,故確診日期與實際被感染日中間會有誤差,以COVID-19疫情來說,約有10~15天的延遲。
3.為方便起見,假定每個人的延遲時間相同,筆者直接將官方統計的「確診日期」來套用到模型上面,用來估算每日新增「確診數」,而非無法被實際量測到的「感染人數」。
4.先前的SARS或MERS都證實存在少數超級傳染者,對疫情傳播有關鍵作用,但在這個簡易的模型中筆者假設每個人的傳播能力是一樣的。

 

過去案例探討:從SARS、登革熱案例驗證

理論建好後,就要進行驗證。

以2003年的SARS和2014~2016年的臺灣登革熱疫情來看,大致上都能有不錯的擬合(圖二)。從圖上來看,不管是登革熱或是SARS,每週或每日的「新增人數」都是兩側低、中間高的曲線。SARS疫情因為是逐日登載,所以上下跳動的幅度較大;而登革熱數據因為是每週記錄,相當於做了七點平均,故相對平滑許多(R2值甚至都可以超過0.95,意即此模型的解釋能力很高)。

 

圖二:數據來源包括WHO網站和臺灣傳染病統計資料查詢系統。圖中的校正曲線為模型預估的理論值,圓圈則為實際值。

為方便理解和比較不同傳染病間疫情的差別,先定義校正曲線上半高寬所橫跨的期間為疫情的「高峰期」,作為一個比較疫期長短的基準,如2015年臺南登革熱高峰期(約六週)就遠比2014年的高屏區登革熱疫情(約十一週)要短。而同樣是SARS,在香港疫情高峰期約40天,但在臺灣就大概只有18天(儘管當年4月24日和平醫院已封院,但從事後數據看來,真正的高峰要到5月10日之後,而此時香港的疫情已近乎結束)。

總確診數則為整個曲線的線下面積,因此曲線涵蓋的範圍愈大,感染人數愈多。特別注意的是,圖上的縱軸代表的是單位時間(天或週)內「新增」的個案數,因此這裡說的疫情結束,是指感染人數不再新增,而非感染人數歸零。即使感染人數不再新增,可能還是存在大量感染人數和新增死亡人數。

 

未來疫情走勢預判

回頭再看此次的COVID-19疫情。在用這個模型作數據推估前,可先看一下其他流行病學專家怎麼說。

各家統計的COVID-19感染人數,從幾萬人到上億都有人估計,估最高的是約翰霍金大學(JohnsHopkins University)資深學者東爾(Eric Toner)在去年10月就預測的新型冠狀病毒會大流行並殺死6500萬人,但是這個數據提出當時還尚未發現有引發COVID-19的病毒,加上這個數字是在毫無防疫措施的前提下推估的,因此先暫且擱一邊。

另一個比較可信的,是英國蘭卡斯特大學(Lancaster University)資深學者里德(Jonathan Read)等人在1月24日放上medRxiv的一篇論文。它的第一版中提到,在沒有有效的防疫措施下,到2月4日前,光是武漢地區就會有約20萬人被感染,且疫情在中國大陸各大城市會隨著疫情遍地開花。不過在他們1月27日上傳的第二版論文中,已刪去該段20萬人染病的預測,並在修正參數後更新了其中一項數字,推估1月22日前在武漢的總感染約是兩萬人。不過須注意的是,這篇並沒有討論何時會是疫情高峰。

此外,在1月底時,有兩位醫界專家不約而同分別作出預測,第一位是當年抗SARS的大將鍾南山,他認為疫情大概在2月9日達頂峰,但是沒有預估最終的確診人數。第二位是香港大學醫學院院長梁卓偉,他推估在1月27日武漢已有四萬多人感染,且因中國大陸其它大城市也會出現大規模傳染,因此整體疫情要到4月底才封頂。

綜合上述預測的可能感染人數和疫情高峰,搭配截至1月30中國官方公布的確診數據,接著嘗試以此簡化的動力學模型來推估從當天之後每天可能「新增」確診數。筆者在1月30日的文章裡給了三個預測狀況,分別以感染人數從多到少、每日確診人數高峰從遠到近來排序,如圖三。(因篇幅關係,模型參數的
設定、優化等細節就不在此贅述。)


圖三:以上述模型預測中國大陸境內武漢肺炎未來每日「新增」確診人數。

狀況一,也就是圖一上的悲觀情景(pessimistic scenario),接近里德在第一版論文中提到防疫效果低、感染人數在2月4日破20萬的情況(注意可被觀察到的「確診數」要推遲到2月18以後),另外參考當年SARS在臺灣高峰期大約18天,即可得出狀況一的曲線。在這個推估趨勢中,確診數高峰要到2月9號左右之後才開始,並且持續約三週,最高峰時每天新增的確診病例會達到一萬人,最後總人數約24萬。此外,這個預測狀況也會很類似2003年臺灣的SARS,和平醫院封院時距離高峰期還有兩個禮拜;而這次武漢封城距離2月9號也大約是兩週多一點的時間。不過這個狀況的前提是前面封城的防疫效果不佳,到2月後才停止大規模人傳人的情況,而里德那篇論文的第二版也拿掉了這個預測。(從2月中之後的數據判斷並沒有走這個曲線)

狀況二就是武漢等城市封城隔離後有達到防疫效果,因此在1月23日後人傳人的情況大幅收斂,往後推2峰期每天會有2000~4000個新增確診病例,並持續到2月中旬。這個趨勢線即為圖一上的樂觀情景(optimistic scenario),接近鍾南山的預測,並在套用模型後可推算確診人數總數,推估約六萬人左右。要區分狀況一與狀況二,則可觀察從2月中旬起是否新增數有下降趨勢,如果下降幅度明顯,到了月中後低於每日兩千人,那可能就代表高峰期已經結束,大概可以判斷這波疫情有機會在二月底之後逐漸落底。

狀況三則為立即停止(Immediate halt)曲線,狀況最樂觀,即在模型推估日1月30日達到高點的情況!但此版本並未發生,就不再敘述。

在此特別強調,圖上的預測都是在新冠肺炎疫情沒有在中國境內其它城市大規模群聚感染的情況下才成立,因為若從其它城市開始新的散播,會爆發出新的獨立疫情,並且和舊的曲線交疊,會讓數據判讀相對複雜。若要考慮這樣的情況,就會有一個圖上沒畫出來的狀況四,也就是梁卓偉講的,其它大城市遍地開花,疫情就真的要落到3、4月才能善了。(《刺胳針》(The Lancet)在2月1日發表一篇論文,根據他們的流行病學模型推估,疫情要到4月之後才會達到高峰)。

 

數據更新與討論

本文原文發表於1月30日,從圖四可以看出大致上的趨勢與模擬的狀況二有相當程度的吻合。要特別注意的是,此前的確診數都是實驗室檢測的確診人數(即黑色與灰色粗圈),然而從2月12日開始,疫情最嚴峻的湖北省加入了「臨床診斷(clinical diagnosis)」作為確診的依據,因此當天的新增確診數暴增約一萬5000例,但若單只看realtime PCR確診的人數,依然落在狀況二的曲線範圍內。而後在2月16日起,官方公布的數據就不再將「實驗室檢測」和「臨床診斷」的確診數分開報導,因此數據點包含了大量的臨床診斷確診而無法再使用預測曲線去做擬合。然而,就算直接以總確診數來看,整體的疫情趨勢依然比較接近狀況二的預測。


圖四:1月30日前的預測數據點。其中,黑色粗圈為1月30日前用於預測的數據點,灰色粗圈為1月31日起官方所公布的中國大陸境內「每日新增確診人數」。

特別注意的是,截至2月26日截稿當天,包含「臨床診斷」的累積確診數為78064人,若扣除2月12~15日等四天的18450例臨床診斷數,總確診人數大約為六萬人,此即為狀況二曲線最初所預測的六萬人總數。在沒有大量新增確診數的情況下,大致上可以判斷本次疫情的走勢應最接近狀況二的模擬。

在此特別補充,在預測趨勢公布後,收到很多來自民眾和學界的意見,認為中國官方的數據不可信,因此這篇的趨勢預測可說是「廢料進、廢品出(garbage-in, garbageout)」。無可否認地,這個爭議點確實存在。本次武漢地區因為疫情嚴重程度遠超過當地醫療體系可以承擔的量能,有不少出症狀的患者無法獲得專業的醫療照護,加上疫情初期受限於檢測量能的不足,因此可預期存在大量被感染且有症狀者未列入官方確診人數中。這部分的誤差值是用這個模型推估上的盲點,但用於整體疫情趨勢的判讀上,相信仍具有一定的參考意義。

 

結論與聲明

本文中筆者以簡化後的反應動力學模型來描述傳染病的傳播模式,試著讓具有高中以上化學基礎的人以直觀的方式理解傳染性疾病的爆發與疫情特徵。此外,截至2月26日的中國大陸境內確診數趨勢來看,若後續在武漢以外的中國大陸境內城市沒有大規模發生疫情傳播的情況下,此次中國的COVID-19疫情有機會在2 月結束後漸趨尾聲。然而,儘管新增確診數不再大幅度增加,仍然存有大量的被感染者,以及日漸增加的死亡病例數。而近日在中國大陸境外(如南韓、伊朗、義大利、以及鑽石公主號)也出現大規模的傳染,雖然不樂見到其它國家也發生類似的疫情,但在中國境內疫情走勢對其它發生疫情的國家來說,依然具有高度的參考價值。

再次強調,本篇的目的是作為教育用途,希望以相對精簡且直觀的方式去理解疫情可能的趨勢。實際上,流行病學中經典的SIR模型,就是在健康者和被感染者之後,多了一個R(代表康復者)。對此模型有瞭解的人會發現,本篇的模型只有S和I(對應到A和B),所以沒有辦法計算基本傳染數(reproductive number, R0),相當於被感染者一直具有傳染力,這點很不符合真實世界的情況。而化學動力學(chemical kinetics)過去在醫藥、生態、傳染病等領域大量地被應用,也在此希冀透過此文,讓所有的學子與學術工作者更瞭解基礎學科的重要性。

〔註〕由於WHO 所公告的數據均將臺、港、澳數據列入中國,為避免混淆,本文多折衷使用「中國大陸」作說明。

參考資料
1. Novel Coronavirus (2019-nCoV) situation reports, World Health Organization.
2. 中華人民共和國國家衛生健康委員會
3. Jonathan M. Read et al., Novel coronavirus 2019-nCoV: early estimation of epidemiological parameters and epidemic predictions, medRxiv, 2020.
4. 林則宏,〈大陸專家鍾南山預估 武漢肺炎約再十天達到高峰〉,聯合新聞網,2020年1月29日。


延伸閱讀
Compartmental models in epidemiology, Wikipedia.