基于實(shí)時(shí)再生數(shù)的COVID-19預(yù)測(cè)模型及可視化平臺(tái)開(kāi)發(fā)

何沅宸 孫傳瑞 陳胤孜 陳曉葳 王錫玲△

（1復(fù)旦大學(xué)公共衛(wèi)生學(xué)院衛(wèi)生統(tǒng)計(jì)學(xué)教研室 上海 200032； 2復(fù)旦大學(xué)公共衛(wèi)生安全教育部重點(diǎn)實(shí)驗(yàn)室 上海 200032）

COVID-19變異株Omicron毒株于2021年11月在南非等地被發(fā)現(xiàn)［1］，因傳染性和免疫逃避能力顯著超過(guò)其他變異株［2］，其迅速在全球成為主要流行毒株［3］。實(shí)時(shí)再生數(shù)（time-varying reproduction number，以下簡(jiǎn)稱(chēng)Rt）即t時(shí)刻一個(gè)病例可傳染所致續(xù)發(fā)病例的平均數(shù)量，在流行病學(xué)研究中常被用來(lái)快速評(píng)估某傳染病的實(shí)時(shí)平均傳播水平。當(dāng)Rt＞1時(shí)認(rèn)為疫情會(huì)繼續(xù)發(fā)展，Rt＜1時(shí)則認(rèn)為疫情逐漸受到控制。新冠肺炎疫情Rt的實(shí)時(shí)估計(jì)已在國(guó)際上廣泛應(yīng)用，英國(guó)、日本和香港等國(guó)家和地區(qū)均向公眾實(shí)時(shí)展示Rt的變化［4-6］以研判疫情走勢(shì)和調(diào)整干預(yù)措施。2022年3月上海市出現(xiàn)新冠本土病例并發(fā)生社區(qū)傳播，流行毒株為Omicron BA.2和BA.2.2型［7］，政府采取了一系列應(yīng)對(duì)措施。目前少有研究評(píng)價(jià)上海市此輪疫情Rt值的實(shí)時(shí)變化并基于此預(yù)測(cè)疫情后續(xù)發(fā)展［8-9］。

新冠疫情可視化平臺(tái)的開(kāi)發(fā)具有重要意義，其可將信息圖表化以提高復(fù)雜數(shù)據(jù)的可讀性，滿足大眾及科研人員的信息需求［10-11］。目前有多個(gè)網(wǎng)站進(jìn)行新冠數(shù)據(jù)可視化展示，如丁香園、騰訊新聞等［12-13］。然而，這些平臺(tái)內(nèi)容以數(shù)據(jù)描述為主，相似度高，缺乏對(duì)上海市各區(qū)的詳細(xì)描述，沒(méi)有進(jìn)行新冠病毒實(shí)時(shí)傳播能力Rt的評(píng)估和預(yù)測(cè)［14］。本研究將基于Rt構(gòu)建貝葉斯分層模型，評(píng)估和預(yù)測(cè)上海市新冠疫情發(fā)生發(fā)展并開(kāi)發(fā)新冠疫情可視化平臺(tái)。

資料和方法

數(shù)據(jù)來(lái)源和病例定義2022年3月1日—5月6日上海市COVID-19疫情每日確診病例、無(wú)癥狀感染者和新冠疫苗接種數(shù)據(jù)［15］來(lái)源于上海市衛(wèi)生健康委員會(huì)（以下簡(jiǎn)稱(chēng)“衛(wèi)健委”）官方網(wǎng)站（https：//wsjkw.sh.gov.cn/）。確診病例的定義參照國(guó)家衛(wèi)健委修訂的《新型冠狀病毒肺炎診療方案（試行第九版）》［16］。無(wú)癥狀感染者未在第九版診療方案中明確定義，一般指除確診病例以外的病原學(xué)檢測(cè)呈陽(yáng)性的感染者。

政府所采取的疫情防控措施及實(shí)施日期來(lái)源于上海市官方微信平臺(tái)“上海發(fā)布”和疫情防控工作新聞發(fā)布會(huì)。根據(jù)所影響的人群和地區(qū)，本研究將各類(lèi)措施劃定為低、中、高三級(jí)強(qiáng)度。低強(qiáng)度包括3月12日起幼兒園及中小學(xué)線下停課和3月15日起的切塊式、網(wǎng)格化篩查，中強(qiáng)度包括3月28日起以黃浦江為界分批、分區(qū)管控，高強(qiáng)度則包括4月3日起全市范圍進(jìn)行多輪抗原和核酸篩查。常住人口年齡構(gòu)成來(lái)自上海市統(tǒng)計(jì)局第七次人口普查數(shù)據(jù)［17］。

模型構(gòu)建與主要參數(shù)設(shè)定本研究使用基于實(shí)時(shí)再生數(shù)構(gòu)建的貝葉斯分層模型進(jìn)行預(yù)測(cè)，模型中假定每日?qǐng)?bào)告的感染者數(shù)It服從以It為均數(shù)和It+為方差的負(fù)二項(xiàng)式分布，考慮感染到報(bào)告延遲時(shí)間分布π和感染者診斷比例Piar來(lái)推導(dǎo)每日實(shí)際感染者數(shù)it，如公式（1）所示。公式（2）中每日實(shí)際感染者數(shù)it由該時(shí)刻的實(shí)時(shí)再生數(shù)Rt和該時(shí)刻前感染者的離散卷積和決定，離散卷積和由此前時(shí)刻的感染者數(shù)通過(guò)代際間隔分布g（即一代病例和二代病例感染的時(shí)間間隔分布）加權(quán)計(jì)算得到。實(shí)時(shí)再生數(shù)Rt則受到主要流行毒株的基本再生數(shù)R0、人群免疫水平（包括自然感染和接種新冠疫苗所獲得）和非藥物干預(yù)措施強(qiáng)度的影響，S為經(jīng)人群免疫水平校正后的易感者比例，pin和pva分別指由自然感染和接種疫苗獲得的免疫效果，見(jiàn)公式（3）。

上海市疫情發(fā)生社區(qū)傳播后病例數(shù)激增，醫(yī)療資源不足，存在病例錯(cuò)分可能，研究中以陽(yáng)性感染者而非確診病例作為擬合指標(biāo)，每日?qǐng)?bào)告的陽(yáng)性感染者數(shù)由每日?qǐng)?bào)告的確診病例數(shù)+無(wú)癥狀感染者數(shù)-由無(wú)癥狀感染者轉(zhuǎn)歸的確診病例數(shù)計(jì)算得到。為減少數(shù)據(jù)的隨機(jī)波動(dòng)性，本研究將每日?qǐng)?bào)告數(shù)據(jù)進(jìn)行移動(dòng)平均處理，取所計(jì)算日期前3天至后3天數(shù)據(jù)（n=7）的算術(shù)平均值。研究中所用的主要參數(shù)來(lái)源于國(guó)內(nèi)外相關(guān)研究和合理假定，包括流行病學(xué)參數(shù)、疫苗接種和自然感染后預(yù)防感染的免疫效果及人群年齡別疫苗覆蓋率，具體設(shè)定見(jiàn)表1?？紤]到當(dāng)病例數(shù)激增時(shí)報(bào)告的延遲性會(huì)進(jìn)一步增加，因此以每日?qǐng)?bào)告感染者數(shù)是否超過(guò)5 000例為界，使其服從均數(shù)不同的Gamma分布。

模型擬合與預(yù)測(cè)性能評(píng)估本研究以2022年3月1日—4月29日為模型擬合段，以4月30日—5月6日為預(yù)測(cè)性能評(píng)估段。模型擬合段以最大似然法［25］估計(jì)不同強(qiáng)度的防控措施對(duì)Rt的影響，并選擇最佳參數(shù)預(yù)測(cè)疫情的后續(xù)發(fā)展，擬合偏差比例的計(jì)算如公式（4）所示，其中kt指實(shí)際觀察值，Pt，median是t時(shí)刻預(yù)測(cè)概率分布的中位數(shù)。鑒于感染到報(bào)告的延遲性，未納入4月23—29日的數(shù)據(jù)用于預(yù)測(cè)，假定此時(shí)間段內(nèi)的Rt值和4月22日的數(shù)值相同。

通過(guò)對(duì)比預(yù)測(cè)結(jié)果與官方報(bào)告數(shù)據(jù)來(lái)評(píng)估模型預(yù)測(cè)性能。預(yù)測(cè)所得的每日感染者數(shù)是由非負(fù)整數(shù)集合的概率分布，其評(píng)估首先需要驗(yàn)證預(yù)測(cè)概率分布和實(shí)際觀察值在統(tǒng)計(jì)學(xué)上具有一致性，再衡量集中度［26-28］。本研究中選擇隨機(jī)化概率分布變換（probability integral transform，PIT）和分級(jí)概率評(píng)分（ranked probability score，RPS）來(lái)分別對(duì)預(yù)測(cè)結(jié)果統(tǒng)計(jì)學(xué)的一致性和集中度進(jìn)行評(píng)估。

PIT常被用來(lái)探究實(shí)際值和連續(xù)性的預(yù)測(cè)概率分布的統(tǒng)計(jì)學(xué)一致性，而對(duì)于每日新增感染數(shù)這樣的離散性資料，隨機(jī)化PIT是一種新的替代方法［27-28］，計(jì)算方法如公式（5）所示，其中kt指實(shí)際觀察值，Pt(x)是t時(shí)刻的預(yù)測(cè)概率分布，根據(jù)定義Pt(-1)=0，v是標(biāo)準(zhǔn)均勻分布且與k無(wú)關(guān)。若Pt確為實(shí)際觀察值kt的預(yù)測(cè)概率分布，那么ut將是標(biāo)準(zhǔn)均勻分布，再用Anderson-Darling檢驗(yàn)來(lái)驗(yàn)證該假設(shè)在統(tǒng)計(jì)學(xué)上是否具有意義。當(dāng)檢驗(yàn)結(jié)果P≥0.05時(shí)認(rèn)為沒(méi)有充分證據(jù)表明實(shí)際值和預(yù)測(cè)概率分布在統(tǒng)計(jì)學(xué)上不具有一致性；當(dāng)0.01＜P＜0.05時(shí)認(rèn)為有證據(jù)表明兩者可能不具有統(tǒng)計(jì)學(xué)的一致性，而當(dāng)P≤0.01時(shí)認(rèn)為有足夠的證據(jù)表明兩者不具有一致性，即預(yù)測(cè)概率分布相對(duì)實(shí)際觀察值是失準(zhǔn)的。RPS計(jì)算如公式（6）所示，其中k，Pt(x)含義同公式（5），當(dāng)計(jì)算結(jié)果越接近于0，則認(rèn)為實(shí)際觀察值越有可能來(lái)自于預(yù)測(cè)概率分布且預(yù)測(cè)概率分布的集中度越好。研究分析均用R 4.1.3軟件完成，其中模型構(gòu)建使用RStan 2.29［29］。

可視化平臺(tái)開(kāi)發(fā)可視化平臺(tái)的開(kāi)發(fā)主要包括數(shù)據(jù)提取與整理、數(shù)據(jù)可視化、預(yù)測(cè)模型構(gòu)建和平臺(tái)搭建等4個(gè)步驟，基本架構(gòu)如圖1所示。使用Python 3.10的Request、Re和Datetime第三方庫(kù)提取數(shù)據(jù)。使用R 4.1.3軟件完成數(shù)據(jù)整理、分析和可視化。使用網(wǎng)頁(yè)開(kāi)發(fā)語(yǔ)言和R shiny 1.7.0框架［30］交互實(shí)現(xiàn)平臺(tái)搭建。R shiny框架已將網(wǎng)頁(yè)開(kāi)發(fā)中的Ajax技術(shù)、Bootstrap框架和jQuery功能框架進(jìn)行封裝以實(shí)現(xiàn)數(shù)據(jù)通信，而為滿足頁(yè)面布局的美化需求，后續(xù)再分別使用HTML進(jìn)行網(wǎng)頁(yè)前端設(shè)計(jì)，CSS進(jìn)行網(wǎng)頁(yè)布局優(yōu)化以及Java Script提高網(wǎng)頁(yè)互動(dòng)效果。平臺(tái)使用Echarts實(shí)現(xiàn)條形圖、線圖、地圖、熱圖等數(shù)據(jù)可視化，并支持鼠標(biāo)懸浮或點(diǎn)擊等交互。該平臺(tái)的系統(tǒng)架構(gòu)搭建并運(yùn)行于Ubuntu服務(wù)器。

圖1 新冠肺炎疫情預(yù)測(cè)預(yù)警平臺(tái)基本架構(gòu)Fig 1 Architecture of the COVID-19 prediction and early warning platform

結(jié) 果

疫情概況2022年3月1日，上海市衛(wèi)健委報(bào)告首2例本土新冠感染者。截至4月29日，累計(jì)報(bào)告確診病例52 576例（其中18 686例由無(wú)癥狀感染者轉(zhuǎn)歸而來(lái)），無(wú)癥狀感染者530 642例（圖2A）。上海市本土疫情發(fā)展呈現(xiàn)先緩慢增加后迅速上升的態(tài)勢(shì)，在4月13日達(dá)到疫情高峰，當(dāng)日新增27 605例感染者。此后疫情放緩，每日?qǐng)?bào)告的感染者數(shù)波動(dòng)下降，4月27日起每日新增感染者數(shù)已小于10 000例。

圖2 上海市2022年COVID-19疫情發(fā)展趨勢(shì)、實(shí)時(shí)再生數(shù)及預(yù)測(cè)結(jié)果曲線Fig 2 COVID-19 epidemic curve， Rt values and forecast result of Shanghai in 2022

根據(jù)疫情發(fā)展，上海市的非藥物干預(yù)措施逐漸從“精準(zhǔn)防控”策略過(guò)渡到3月15日開(kāi)始的切塊式、網(wǎng)格化防控篩查，又調(diào)整為3月28日起以黃浦江為界分批、分區(qū)實(shí)行封控，再到4月初對(duì)全市范圍進(jìn)行核酸篩查和防控管理。隨著防控措施升級(jí)，Rt曲線總體呈現(xiàn)波動(dòng)下降，在4月5日其均值首次降至1以下（0.94，95%CI：0.87～1.03），結(jié)合延遲時(shí)間分布，與4月13日疫情達(dá)到高峰相符（圖2B）。在4月13—14日Rt的95%CI上限再次超過(guò)1，達(dá)到1.01，與4月22日感染者新增數(shù)出現(xiàn)新的小高峰的趨勢(shì)貼合?？紤]感染到報(bào)告的延遲時(shí)間，4月23—29日（即圖中灰色區(qū)域）的Rt估計(jì)值仍存在較大不確定性。

模型擬合結(jié)果擬合時(shí)間段中，經(jīng)移動(dòng)平均處理后的每日新增感染者數(shù)與預(yù)測(cè)概率分布中位數(shù)擬合良好，擬合偏差比例為14.41%±15.20%。在4月20—26日模型的預(yù)測(cè)中位數(shù)明顯高于對(duì)應(yīng)的移動(dòng)平均值（如圖2C藍(lán)線和黃線所示），其原因可能是4月22日新增感染者數(shù)又出現(xiàn)一個(gè)小高峰而模型擬合未能體現(xiàn)。各級(jí)防控措施的調(diào)整對(duì)Rt值的相對(duì)效應(yīng)如表2所示，上海市所有防控措施可將Rt值降低至1以下（0.30，95%CI：0.28～0.32），而在放松部分區(qū)域的防控后，Rt值增加至0.43（95%CI：0.41～0.46），仍未超過(guò)關(guān)鍵閾值水平，但使得本輪疫情中每日新增感染者數(shù)的下降速度在一定程度上放緩。

表2 上海市疫情防控措施對(duì)Rt值的效應(yīng)估計(jì)Tab 2 The estimated relative effectiveness of nonpharmaceutical interventions on Rt in Shanghai

疫情走勢(shì)預(yù)測(cè)模型預(yù)測(cè)結(jié)果如圖2C所示，本輪疫情預(yù)計(jì)累積感染者數(shù)將達(dá)到625 040例（95%CI：600 465～670 211），預(yù)計(jì)將于5月17日（5月12—22日）單日新增感染者數(shù)＜1 000例，將于6月5日（5月31日—6月11日）單日新增感染者數(shù)＜100例。對(duì)4月30日—5月6日的預(yù)測(cè)結(jié)果進(jìn)行評(píng)估，隨機(jī)化PIT的結(jié)果是P=0.12，尚無(wú)足夠證據(jù)表明結(jié)果概率分布和實(shí)際觀察值無(wú)統(tǒng)計(jì)學(xué)一致性，分級(jí)概率評(píng)分指標(biāo)為382.53±105.29。

可視化平臺(tái)可視化平臺(tái)實(shí)現(xiàn)了上海市疫情變化趨勢(shì)描述和疫情預(yù)測(cè)兩大功能。疫情變化趨勢(shì)描述包括上海市及上海市各區(qū)新增、累積陽(yáng)性感染者和死亡人數(shù)變化趨勢(shì)、各區(qū)累積感染率和感染者來(lái)源（隔離管控或風(fēng)險(xiǎn)人群排查）。疫情預(yù)測(cè)包括實(shí)時(shí)再生數(shù)、新增陽(yáng)性感染者預(yù)測(cè)曲線和預(yù)測(cè)感染人數(shù)關(guān)鍵節(jié)點(diǎn)。為滿足用戶端的不同需求，平臺(tái)內(nèi)置圖表適應(yīng)電腦和手機(jī)端多種主流瀏覽器。該平臺(tái)于2022年4月19日起開(kāi)放訪問(wèn)，訪問(wèn)地址：http：∥124. 222.178.35：3838/COVID19-tracker/。

討 論

本研究使用基于實(shí)時(shí)再生數(shù)構(gòu)建的模型，評(píng)估和預(yù)測(cè)了2022年3月以來(lái)上海市本土疫情的發(fā)展態(tài)勢(shì)，并完成了新冠肺炎預(yù)測(cè)預(yù)警平臺(tái)的搭建。研究中考慮上海市現(xiàn)有免疫水平并設(shè)置隨時(shí)間變化的模型參數(shù)，包括感染者從感染到報(bào)告的延遲時(shí)間分布和各級(jí)非藥物干預(yù)措施的相對(duì)效應(yīng)，預(yù)測(cè)結(jié)果更貼合疫情的實(shí)際發(fā)展。而在預(yù)測(cè)段中雖假定Rt值不變，模型評(píng)估指標(biāo)和后續(xù)實(shí)際報(bào)告數(shù)據(jù)依然證實(shí)了預(yù)測(cè)結(jié)果的準(zhǔn)確性。2022年底，我國(guó)新冠疫情防控策略已有所調(diào)整，但調(diào)整相關(guān)參數(shù)后，此模型也可用于預(yù)測(cè)其他新發(fā)呼吸道傳染病以評(píng)估公共衛(wèi)生風(fēng)險(xiǎn)和干預(yù)措施有效性。

模型擬合結(jié)果顯示，上海市所采取的非藥物干預(yù)措施將Rt值降低并保持在1以下，且本輪疫情預(yù)測(cè)將于6月中下旬實(shí)現(xiàn)“動(dòng)態(tài)清零”目標(biāo)。另有研究通過(guò)上海市每日?qǐng)?bào)告的新增感染者數(shù)據(jù)估計(jì)Rt值，同樣證實(shí)了上海市防控措施的效果［31］。劉可伋等［32］使用時(shí)滯動(dòng)力學(xué)模型預(yù)測(cè)疫情發(fā)展并估計(jì)實(shí)時(shí)再生數(shù)，但與本研究不同的是未充分考慮感染到報(bào)告的延遲時(shí)間，預(yù)測(cè)結(jié)果具有滯后性。疫情后期，連續(xù)14天無(wú)新增感染者區(qū)域的管控放松，而存在傳播風(fēng)險(xiǎn)的地區(qū)仍在嚴(yán)格管控中，防疫政策的變化對(duì)結(jié)果影響不大。但當(dāng)跨區(qū)及跨省市人員流動(dòng)增加時(shí)，新冠傳播風(fēng)險(xiǎn)增加，感染規(guī)模可能進(jìn)一步擴(kuò)大。

研究中假定上海市自2022年3月1日?qǐng)?bào)告首2例感染者后已暫停全市的新冠疫苗接種，而部分區(qū)域在4月下旬已開(kāi)始恢復(fù)對(duì)老年人等重點(diǎn)人群的疫苗接種［33］，這部分人群的免疫水平未在模型中體現(xiàn)，考慮到日均接種千余劑次和疫苗產(chǎn)生有效抗體水平需要2～3周的時(shí)間，可認(rèn)為其對(duì)預(yù)測(cè)結(jié)果影響有限。但若老年人群的接種率進(jìn)一步增高時(shí)，其所提供的保護(hù)效果將使得疫情控制速度快于預(yù)測(cè)結(jié)果。鑒于上海市內(nèi)各區(qū)域人員流動(dòng)頻繁，且非藥物干預(yù)措施的調(diào)整多是基于全市范圍，研究使用全市數(shù)據(jù)進(jìn)行預(yù)測(cè)，區(qū)級(jí)疫情預(yù)測(cè)需要更多數(shù)據(jù)來(lái)支撐分析。

預(yù)測(cè)預(yù)警平臺(tái)使用Python進(jìn)行公開(kāi)數(shù)據(jù)實(shí)時(shí)爬取，基于實(shí)時(shí)再生數(shù)進(jìn)行疫情走勢(shì)預(yù)測(cè)，使用Echarts和R shiny實(shí)現(xiàn)數(shù)據(jù)可視化和交互，具有實(shí)時(shí)性、科學(xué)性、頁(yè)面友好等優(yōu)點(diǎn)。通過(guò)數(shù)據(jù)可視化，高效展示了上海市及各區(qū)新冠疫情變化的時(shí)空差異，降低公眾捕獲并理解公開(kāi)疫情數(shù)據(jù)的時(shí)間成本。目前已有預(yù)測(cè)預(yù)警平臺(tái)無(wú)法提供Rt值以估算疫情實(shí)時(shí)傳播速度［34-35］。本研究基于實(shí)時(shí)再生數(shù)的預(yù)測(cè)模型可視化，填補(bǔ)了已有平臺(tái)對(duì)疫情數(shù)據(jù)深層次挖掘的空白，提供更加科學(xué)的防控措施效果評(píng)價(jià)和疫情預(yù)警預(yù)判。該可視化平臺(tái)可快速擴(kuò)展到其他城市的疫情數(shù)據(jù)解讀及預(yù)測(cè)預(yù)警，乃至其他可能發(fā)生的具有類(lèi)似流行病學(xué)特征的感染性疾病，具有較好的應(yīng)用前景。

作者貢獻(xiàn)聲明何沅宸 研究設(shè)計(jì)，數(shù)據(jù)分析，論文撰寫(xiě)和修改。孫傳瑞 數(shù)據(jù)分析，論文撰寫(xiě)和修改。陳胤孜 數(shù)據(jù)整理和分析。陳曉葳 數(shù)據(jù)采集。王錫玲 研究設(shè)計(jì)和指導(dǎo)，論文修改，經(jīng)費(fèi)支持。

利益沖突聲明所有作者均聲明不存在利益沖突。