國立臺灣大學公共衛生學院 流行病學與預防醫學研究所
博士論文
Institute of Epidemiology and Preventive Medicine College of Public Health
National Taiwan University Doctoral Dissertation
貝氏非循環隨機模型於傳染病監測應用
Stochastic Processes Applications to Surveillance of Infectious Disease with Bayesian Directed Acyclic
Graphic (DAG) Approach
廖翎均 Ling-Chun Liao
指導教授：陳秀熙 博士 Advisor: Hsiu-Hsi Chen Ph.D.
中華民國 112 年 2 月 Feb 2023
誌謝
這本論文得以完成，首先要非常感謝我的指導教授陳秀熙教授，對於學生的 疼惜、百忙中指導訓練，也總在我感覺快要墜落的時候可以接住我。再來要謝謝 533 研究室的嚴明芳老師、陳立昇老師、許辰陽醫師、賴昭智醫師、林庭瑀博士 對於研究的幫助，以及其他學長學姊學弟妹們的陪伴。感謝流預所所辦小姐霽修 行政上的協助。謝謝口試委員在口試時給予的珍貴的建議，讓本論文能夠更加完 整。
在博士班一年級上學期，新冠肺炎才剛開始爆發，之後研究團隊也開始了疫 情的監測追蹤與報導，持續發表以及每週為全球大眾分析新冠肺炎相關研究，除 了無法預料疫情對於全球造成的改變外，教授帶領研究團隊的熱情與嚴謹效率的 研究態度以及學者風範在我心中留下深刻的印象，並且我也非常開心可以在這個 過程中，學習將研究成果學以致用。
回想過去碩班與博班總共六年多的時光，佩服自己的勇敢堅持，懵懵懂懂的 選擇唸生物統計，讓所有的一切從頭開始，並且於讀書會、報告修課、論文修改 期間、口試當下老師們的諄諄教誨都歷歷在目，讓曾經的挫折成為日後成長的養 分、並且除了知識之外還提點許多待人處事的道理。
另外要感謝新竹臺大分院新竹醫院同仁、藥劑部主任、各單位組長，所有藥 師夥伴、謝謝各位協助排班與換假支援，讓我得以順利安心地完成學業。謝謝小 兒部新生兒/兒童加護病房跨領域夥伴工作上給予的鼓勵支持。也謝謝陪我走過 一小段路途的好同學們、Jolie 賢伉儷、包容我一切狀態的書庭。最後謝謝我親愛 的父母、弟弟，總是提醒我要記得早點休息，懂得生活!希望未來可以繼續不忘
初衷、保持興趣、並且持續在領域中創造更多的對話。
2023/02/14
中文摘要 研究背景
新興傳染病與再現傳染病之監測，首重於在族群層次，評估疾病大規模傳播 特性與流行爆發風險。另外感染病毒於個人層次所造成之動態變化，由於病毒特 性與人體免疫機轉，亦對於受感染宿主在不同病毒量之下所造成的疾病傳播風 險、宿主在感染後由於不同病毒量負擔所可能造成之疾病傳播與感染擴散，以及 感染後產生之臨床照護需求息息相關。有效精準防治策略的發展，有賴於對前述 在個人與族群層次之疾病傳播特性，以及各種防治策略，包含公衛防疫措施、族 群疫苗施打策略，以及抗病毒藥物治療策略等，對於預防感染與疾病進展可達到 效益之評估結果而規劃。然而傳統對於傳染病之監測評估，多著重於族群感染層 次運用，如時間序列方法與傳染病傳播隔間模型，配合傳播基礎再生數指標，評 估新興及再現傳染病對族群可能造成之流行風險。然而考慮前述傳染病之病毒與 宿主交互影響所造成之動態變化，以及感染後所造成之疾病進展風險，傳統以族 群為層次之傳染病評估架構有其限制。
研究目標
(1)在疾病傳播動態模型下拓展貝氏網絡 DAG 方法，建立新興與再現傳染病爆發 流行之風險族群監測實證基礎。
(2)建構疾病進展多階段隨機過程貝氏網絡 DAG 模型運用於評估病毒以及宿主特 性在個人層次對於感染後病毒量變化之影響監測。
(3)運用貝氏網絡 DAG 方法結合廣義線性階層模式建立大流行與區域流行評估架 構，並且拓展為含括傳染病疾病進展狀態之介入效益評估方法。
材料與方法
本研究運用新冠肺炎以及猴痘全球流行為新興傳染病與再現傳染病實證資料，
結合所發展之貝氏網絡模型進行實證評估。對於新冠肺炎全球流行，本研究摘取 包含全球各國通報個案、住院、重症，以及死亡等公開資料。對於台灣之流行資 料，本研究摘取疫苗施打以及通報個案與其發展為中症、重症，以及死亡及年齡
公開資料。對於病毒量對於疾病進展之影響，本研究運用台灣某縣市之疫情調查
與檢測資料。猴痘全球流行本研究運用全球各國於2022 全球各國通報個案進行
評估。
本研究運用貝氏網絡結合有向無環圖模型在考量傳染病傳播特性下，分別運 用傳染病動態傳播模型、多階段馬可夫隨機過程，以及廣義線性迴歸模型在前述 實證資料支持下，以貝氏馬可夫蒙地卡羅演算法產生前述模型中之參數事後分 佈，並據以對新興傳染病以及再現傳染病之疾病傳播特性、不同防治策略，包含 疫苗施打、公衛防疫措施，以及抗病毒藥物可達到之保護效果，以及感染後病毒 量動態階段變化進行評估。
結果
運用貝氏網絡傳染病動態模型，結合鑽石公主郵輪侷限空間情境之新冠肺炎 評估結果，顯示新冠肺炎之基礎再生數(R0)為 5.70 (95% CI: 4.23-7.79)。施行鑽石 公主郵輪之防疫措施對於疾病傳播可達到37% (95% CI: 33-40%)之效益。在此侷 限空間情境下，若可提前5 日施行防疫措施則可達到 53% (95% CI: 44-62%)之疾 病傳播防治效益，使新冠肺炎個案數由761 人減少為 403 人。運用此貝氏網絡動
態模型於再現傳染病猴痘2022 年全球流行顯示，全球之疾病傳播再生數由初期
之1.001 (95% CI: 0.986-1.150)增加至 1.459 (95% CI: 1.370-1.507)，隨著介入措施 包含公衛防疫措施以及天花疫苗施打達到31% (95% CI: 27.0-33.6%)之介入效益，
此再現傳染病之再生數下降至1.027 (95% CI: 1.026-1.026)，傳播逐漸受到控制。
本研究亦對各國之猴痘傳播流行再生數變化以及對應之疾病傳播防治策略進行評 估。
運用貝氏網絡多階段隨機過程，本研究評估Alpha 與 Omicron 變異株之病毒 量，對於新冠肺炎感染個案由無症狀進展為臨床狀態、持續無症狀個案之比例，
以及進展發生臨床症狀之時間中位數。評估結果顯示於Alpha 時期多數感染者由
症狀前期個案期進展成為症狀之路徑，僅有少數依循持續無症狀路徑(發生率 24.9，95% CI: 15.6-35.1)。而 Omicron 時期依循無症狀路徑者則增加達 271.4 (95% CI: 240.4-303.7)。Alpha 與 Omicron 變異株感染發生臨床症狀之中位數時間
分別為4.07 天(95% CI: 3.33-4.84)以及 1.22 天(95% CI: 1.12-1.33)。在考慮宿主特
性與接觸模式後，病毒量對於受到Alpha 變異株感染個案是否產生臨床症狀以及
其病程進展皆具有顯著影響，然而對於Omicron 感染個案病毒量高低所造成的影
響則相對減低。病毒量對於Alpha 與 Omicron 感染後之疾病進展模式亦有所差 異。病毒量於Alpha 感染個案之進展呈現劑量效應，但對於 Omicron 感染個案病
毒量影響較小，多數個案在1-2 天之內由症狀前期個案期進展之發生臨床症狀。
疫苗施打與否對於疾病進展亦呈現不同的模式。完整施打疫苗者受病毒量影響較 小，但未施打疫苗者其感染病毒量對於其臨床動態進展則會有所影響。
運用貝氏網絡階層廣義線性迴歸DAG 模型，分析新冠肺炎實證資料結果顯
示，新冠肺炎之流行無法達到地方流行之平衡態，而將以四週遞延之模式受到變 異株病毒傳播特性為主要影響下爆發流行。就不同防疫措施之評估結果，疫苗施 打對於症狀感染可達到55% (95% CI: 38-67%)之保護效果，其中 mRNA 疫苗施打 可提供較佳之保護效果(63%，95% CI: 40-77%)，此保護效果在年輕族群(未滿 50 歲)優於年長族群(大於 50 歲)。疫苗對於死亡、重症，以及中症之保護方面，完 成追加劑施打分別可提供73.9% (95% CI: 72.7-75.1%)、73.9% (95% CI: 73.0- 74.8%)，以及 72.5 (95% CI: 71.8-73.2%)之效益，完成基礎劑之保護效益則分別為 52.7% (95% CI: 49.6-55.6%)、54.5 (95% CI: 50.0-54.5%)，以及 51.9% (95% CI:
50.3-53.6%)。公衛防疫措施則可提供約 12%-20%之保護效益。口服抗病毒藥物約 可降低中症風險約12%。
結論
本研究運用隨機過程與DAG 模型建立創新方法，結合多階段馬可夫模式運
用以全球以及臺灣地區實證資料，發展新興與再現傳染病由族群至個人之監測架
構與介入效益評估貝氏DAG 階層模型。所提出多層次之傳染病監測方法有助於
傳染病傳播風險完整評估以及發展個人化防治策略。
關鍵字: 隨機過程、貝氏有向無環圖模型、新興與再現傳染病、新冠肺炎、猴 痘、病毒量、臨床疾病進展
Abstract Background
The evaluation for the risk of outbreak and transmission at large scale has been the first and foremost goal in the surveillance of emerging and reemerging infectious diseases at population level. In addition to this purpose, the risk of disease evolution in terms of viral shedding and clinical severity resulting from the characteristics of pathogen, viral load, and host immunity has also been an important aspect of
surveillance at individual level. The development of effective and individual-tailored containment strategies including non-pharmaceutical interventions (NPIs), mass
vaccination, and antiviral therapies is highly dependent on the empirical evidence taking into account these factors with multilevel characteristics. Conventional approaches for surveillance of infection disease including the use of time-series models and basic reproductive number (R0) derived from compartment models, both focusing on the risk at population level. However, the characteristics of multilevel and multiple outcomes regarding the surveillance of infectious disease render the use of these conventional approaches intractable.
Aim
The aims of this study are
(1) to develop a Bayesian network (BN) analysis with DAG model supported by the compartment model to characterize the dynamics of disease evolution and to assess the risk of disease outbreak at population level;
(2) to apply the BN DAG model with stochastic process underpinning for the
surveillance of disease evolution associated with viral load and viral characteristics at individual level;
(3) to develop a BN hierarchical DAG in conjunction with generalized linear regression model to assess the risk of pandemic and endemic and the effectiveness of
containment strategies at different level.
Materials and Methods
A series of BN DAG models were developed with the support from the information contained in the empirical data of emerging (COVID-19) and reemerging (monkeypox) infectious diseases. Regarding the data on COVID-19 pandemic, global open data with the information on country and region, reported cases, hospitalized cases, cases
admitted for intensive care were collected from open data. Taiwan outbreak data on the
clinical severities including moderate, severe, and death, age, and vaccination history were abstracted from surveillance report. The empirical data on contact tracing and viral load measured by Ct level collected from a county in Taiwan were used for assessing disease evolution at individual level. The open data on global monkeypox outbreak surveillance in 2022 were used for assessing the risk of outbreak of re-emerging infectious disease.
A series of BN DAG model in conjunction with dynamic compartment model, multistate Markov process, and generalized linear regression model taking into account disease characteristics were developed. By using the Bayesian Markov Chain Monte Carlo (MCMC) algorithm, information from the empirical data mentioned above were used for the derivation of posterior distributions of the parameters of the proposed BN DAG models.
Results
For the risk of COVID-19 outbreak in the confined space of cruise ship, the R0
derived from BN DAG susceptible, exposed, infected, and recovered (SEIR) were estimated as 5.70 (95% CI: 4.23-7.79). The effectiveness of containment measures implemented on board was estimated as 37% (95% CI: 33-40%). On the basis of the posterior distributions of the BN DAG SEIR model, the implementation of containment measured 5-days earlier can enhance the effectiveness to 53% (95% CI: 44-62%) and reduce the COVID-19 cases from 761 to 403 in such a confined space of cruise ship.
Regarding the re-emerging disease of monkeypox, the R0 increased from 1.001 (95%
CI: 0.986-1.150) in the early stage to 1.459 (95% CI: 1.370-1.507) in the latter stage, suggesting the risk of transmission at large scale and outbreak. The effectiveness of containment strategies including NPIs and vaccination was estimated as 31% (95% CI:
27.0-33.6%), which bring the R0 down to 1.027 (95% CI: 1.026-1.026), indicating the containment of outbreak.
By using the BN DAG with stochastic process, a four-state disease progression model was constructed for the surveillance of COVID-19 evolution at individual level regarding Alpha and Omicron variants of concern (VOC) infection. The results
demonstrated the difference between subjects infected by these two VOCs. While the small fraction of subjects infected by Alpha VOC turned into asymptomatic (incidence:
24.9, 95% CI: 15.6-35.1), a high incidence for asymptomatic infection (271.4, 95% CI:
240.4-303.7) was estimated for Omicron infection. The median time from pre-
symptomatic to symptom phase for Alpha and Omicron VOC infection was estimated as 4.07 (95% CI: 3.33-4.84) and 1.22(95% CI: 1.12-1.33) days, a significant short period for Omicron infection. After adjusting for host characteristics and contact pattern, there was different roles of viral load for clinical evolution for Alpha and Omicron infection.
A signification impact of viral load on the clinical progression for subject infection by Alpha VOC with a dose-response pattern was observed. For those infected by Omicron VOC, although viral load remains a significant effect on clinical evolution, especially for the unvaccinated population, a lower extent was observed compared with those infected by Alpha VOC
The results derived by applying GLM with BN DAG show the outbreak of
COVID-19 will reach equilibrium in the long run. Recurrent outbreaks affected mainly by the characteristics of dominant VOC are expected with the optimal lag function spanned over four weeks. By using the BN DAG with GLM, the effectiveness of vaccination against symptomatic infection of Omicron was estimated as 55% (95% CI:
38-67%), with a higher protective effectiveness conveyed by mRNA -based vaccine (63%，95% CI: 40-77%). Regarding the effectiveness of vaccination against death, severe, and moderate disease, the estimated results were 73.9% (95% CI: 72.7- 75.1%)、73.9% (95% CI: 73.0-74.8%), and 72.5 (95% CI: 71.8-73.2%) for booster vaccination. The corresponding figures for primary series were estimated as 52.7%
(95% CI: 49.6-55.6%)、54.5 (95% CI: 50.0-54.5%), and 51.9% (95% CI: 50.3-53.6%).
The protective effectiveness resulting from NPIs and antiviral therapy were estimated as 12-20% and 12%, respectively.
Conclusion
By using Bayesian DAG approach with stochastic process underpinning, a series of novel applications to surveillance of infectious disease were developed. The proposed Bayesian DAG framework considering the effect at population level and individual level can facilitate the surveillance of emergence and re-emergence infectious disease in terms of outbreak risk assessment and the development of individual-tailored
containment strategies.
Keywords: stochastic process, Bayesian Directed Acyclic Graphs (DAGs), emerging and reemerging infectious diseases, COVID-19, monkeypox, viral load, clinical disease progression
目錄
口試委員會審定書 ... I 誌謝 ... II 中文摘要 ... III Abstract ... VI 圖目錄 ... XII 表目錄 ... XIII
第一章 研究背景 ... 1
1.1 傳統傳染流行病學模式限制 ... 1
1.2 新冠肺炎及猴痘疫情 ... 1
1.3 研究目的 ... 4
第二章 文獻回顧 ... 5
2.1 傳染病疾病特性 ... 5
2.1.2 猴痘疾病傳播特性 ... 5
2.1.1 新冠肺炎疾病傳播特性 ... 5
2.2 新冠肺炎疫苗保護效益評估 ... 6
2.2.1 新冠肺炎疫苗效益評估研究 ... 6
2.2.2 新冠肺炎疫苗效益研究設計 ... 7
2.2.3 時間對疫苗效益之影響 ... 9
2.2.4 其他的介入措施影響效益評估 ... 10
2.3 傳染病傳播與防治評估方法 ... 10
2.3.1 廣義線性模型 ... 10
2.3.2 廣義自我回歸模型 ... 11
2.3.3 兩階段馬可夫過程與傳染病流行平衡態 ... 13
2.3.4 多階段隨機過程用於疾病狀態預測 ... 17
2.3.5 貝氏有向無環圖(Directed Acyclic Graphic，DAG)模型與多階段馬可 夫過程 ... 20
第三章 研究方法 ... 21
3.1 傳染傳播隨機過程模式貝氏 DAG 模型 ... 21
3.1.1 貝氏微分方程 SEIR 模型 ... 21
3.1.2 傳染傳播隨機過程貝氏 DAG 模型與馬可夫鏈蒙地卡羅演算法 .... 23
3.2 傳染病疾病進展隨機過程貝氏 DAG 模型 ... 29
3.2.1 疾病進展多階段隨機過程貝氏網絡(BN)模型建構 ... 29
3.2.2 傳染病多階段 DAG 模型馬可夫蒙地卡羅演算法 ... 34
3.3 平穩過程之時間序列資料(Stationary Process of Time-series Data) ... 36
3.3.1 自我回歸過程 Autoregressive process ... 36
3.3.2 移動平均過程 Moving Average process (MA) ... 39
3.3.3 廣義自我回歸模型 Generalized Autoregressive Model ... 40
3.4 貝氏網絡(BN)架構模型 ... 40
3.4.1 新冠肺炎臨床疾病狀態貝氏網絡架構建立 ... 40
3.4.2 全球新冠肺炎傳播模型 ... 42
3.4.3 台灣新冠感染保護效益評估模型 ... 44
3.4.4 階層貝氏廣義線性 DAG 模型馬可夫蒙地卡羅學習演算法 ... 45
第四章 新興與再現傳染病傳播與進展資料 ... 51
4.1 猴痘國際傳播資料 ... 51
4.2 新冠肺炎疾病傳播資料 ... 51
4.2.1 國際新冠疫情評估資料 ... 51
4.2.2 台灣新冠肺炎疫苗、NPIs 與抗病毒藥保護效益實證資料 ... 51
4.2.3 新冠肺炎於局限空間傳播實證資料 ... 52
4.2.4 多階段新冠肺炎進展實證資料 ... 53
第五章 結果 ... 53
5.1 猴痘傳播評估 ... 53
5.2 新冠肺炎基礎傳播再生數 ... 84
5.3 多階段新冠感個人化動態病毒量曲線 ... 90
5.3.1 不同變異株的爆發個案情形 ... 90
5.3.2 病毒載量相依與共變數特定的狀態轉換 ... 90
5.3.3 病毒量和相關共變數的多階段感染動態過程 ... 91
5.4 全球新冠疫情貝氏網絡傳播評估 ... 104
5.6 疫苗對新冠肺炎臨床嚴重程度保護台灣效益評估 ... 111
第六章 討論 ... 116
6.1 猴痘病毒傳播 ... 116
6.2 COVID-19 病毒郵輪傳播 ... 116
6.3 多階段馬可夫感染過程 ... 118
6.3.1 以實證資料為基礎的多階段馬可夫感染過程 ... 118
6.3.2 個人化多階段感染過程 ... 118
6.4 資料模擬自我相關函數對於傳染病觀察資料之影響 ... 119
6.5. 採用貝氏分層模型評估疫苗有效性與 NPIs 與抗病毒藥物影響 ... 120
6.5.1 大規模施打疫苗對於預防 Omicron 感染重症與死亡之有效性 ... 121
6.5.2 新冠疫苗引發的 T-cell 調節免疫 ... 121
6.6 結論 ... 122
參考文獻 ... 126
圖目錄
圖一、估計傳染疾病傳播之SEIR 模型之有向無環圖模型 ... 25
圖二、估計傳染疾病傳播之SIR 模型之有向無環圖模型 ... 26
圖三、估計有無介入措施防止新冠肺炎傳播之SEIR 模型之有向無環圖模型 ... 27
圖四、估計猴痘傳播之有向無環圖 ... 28
圖五、納入調整變數之新冠肺炎疾病臨床進程 ... 32
圖六、新冠肺炎疾病進展貝氏有向無環圖模型 ... 33
圖七、傳染病貝氏網絡(BN)廣義線性階層 DAG 模型 ... 46
圖八、全球新冠肺炎傳播貝氏網絡(BN)廣義線性階層 DAG 模型 ... 47
圖九、台灣新冠疫苗感染保護效益貝氏網絡(BN)廣義線性階層 DAG 模型 ... 48
圖十、台灣新冠疫苗臨床保護效益貝氏網絡(BN)廣義線性階層 DAG 模型 ... 49
圖十一、傳染病群體與個人化之監測與介入措施評估模型架構 ... 50
圖十二、全球和不同國家的累積猴痘病例、貝氏 SIR、傳統 SIR 預測的累積病例 ... 81
圖十三、鑽石公主遊輪上新冠肺炎觀察病例與貝氏SEIR 模型預測個案數 ... 87
圖十四、鑽石公主遊輪預期自然發展病例與控制措施觀察病例 ... 88
圖十五、鑽石公主遊輪觀察累積個案與控制措施前後預測病例數 ... 89
圖十六、多階段模式之新冠肺炎個案密切接觸者主動追蹤資料建立流程 ... 98
圖十七、Alpha 變異株的多階段感染動態過程，整體以及和依病毒量(Ct 值)和年齡組分層 ... 99
圖十八、Omicron 變異株的多階段感染動態過程，整體以及和依病毒量(Ct 值)和年齡組分層 ... 102
圖十九、新冠肺炎感染率對感染以及臨床嚴重度遞延影響趨勢圖 ... 105
圖二十、資料模擬不同自我相關函數 ... 123
表目錄
表一、猴痘的累積病例數、貝氏 SIR 模型與傳統SIR 預測的累積病例數 ... 55
表二、使用貝氏 SIR模型估計猴痘爆發相關參數的結果 ... 79
表三、使用 SIR 模型估計猴痘爆發相關參數的結果 ... 80
表四、觀察累計病例數與貝氏 SEIR 模型模擬與預測病例數 ... 86
表五、兩種變異株的接觸者追踪確診個案依症狀分類的基礎特徵 ... 93
表六、某縣市ALPHA和OMICRON 變異株新冠肺炎症狀進展情形 ... 95
表七、單變項分析某縣市兩種變異株之病毒量對疾病惡化之淨效應 ... 96
表八、新冠肺炎感染率對感染遞延作用週數評估 ... 106
表九、全球新冠傳播遞延效應貝氏網絡(BN)階層模式評估結果 ... 107
表十、台灣新冠疫苗感染保護效益貝氏網絡(BN)階層模式評估結果 ... 110
表十一、以貝氏網絡(BN)階層模式評估台灣 OMICRON BA.2 流行期間疫苗效益 ... 113
第一章 研究背景
1.1 傳統傳染流行病學模式限制
由於全球化與氣候變遷趨勢之影響，新興傳染病或是再現傳染病再起，在過 去十年造成重大影響 (McMichael et al., 2006)，傳染病防治因此重新成為全球公 共衛生的重要議題，流感以及新冠肺炎即是最佳見證。過去傳統流行病學模式如 Susceptible-Exposed-Infected-Recovered (SEIR)模式應用於流行推估及進行狀況探 究，以評估各種介入措施如:非藥物公衛防疫措施(non-pharmaceutical interventions (NPIs))、接種疫苗(vaccination)、醫療照護等，雖然相當有用，但也常遭遇許多問 題，造成應用上受到許多限制如下:
(1) 充分統計量資訊不足，例如 SEIR 模式中 S 與 E 資訊經常無法區別。
(2) 很難產生以統計分佈以包含參數的不確定性，疫情評估僅有再生數 (reproductive number，R0)之點估計以及區間估計。
(3) 對於無症狀個案及症狀前期個案期個案在原有 SEIR 模型中非常困難區別，必 須改變SEIR 模式型態。
(4) 評估各種介入防疫措施如 NPIs、疫苗、抗病毒藥物的因果關係相當困難。
(5) 難以評估多種形態變種病毒流行在不同時間及不同國家資料如何交互影響也 是傳統SEIR 模式之限制。
(6) 與時間相關之傳染及疾病進展模式也是傳統 SEIR 模式無法處理
(7) 不同傳染病有不同特性。如:COVID-19 之免疫逃脫特性，是否可以使用相同 傳統傳染病模式來評估全球疫情。
為克服這些限制，本論文提出貝氏有向無環圖(DAG)隨機過程來克服前述傳統 流行病學模式之不足。
1.2 新冠肺炎及猴痘疫情
本研究第一個運用之新興傳染病為新冠肺炎全球大流行。從西元 2000 年 起，冠狀病毒已經引發三次歷史上的傳染病爆發，包含2012 至 2013 年 SARS- CoV-2 於 29 國的傳播，2012 年的 MERS-CoV，以及目前新冠肺炎( COVID-19)
的全球大流行。且因其RNA 病毒的善變特性，新冠肺炎大流行至今演變出多種
SARS-CoV-2 變異株，三年期間，全球共享流感數據倡議組織(GISAID)將新冠肺 炎關注變異株(SARS-CoV-2 VOC)區分為幾個主要的時期包含 D614G、Alpha、
Delta、Omicron 等不同變異株病毒流行時期。其傳播力、族群嚴重度以及致死率
等皆有不同，在有無疫苗情境下的NPI、隔離天數、快篩檢驗、傳播場域因應措
施等問題上需要被解決。
衛生政策領導者須因應以上情況，對於各種防疫政策與措施做出決策，包括 NPIs、疫苗、以及抗病毒藥物等。於 2020 年 D614G 流行時期之疫情初期，NPIs 配合核酸檢測(RT-PCR)是當時阻隔新冠肺炎傳染僅有的選擇。在具有高度感染與 臨床嚴重度保護效益之疫苗問世後，得以藉由疫苗基礎劑施打以建立免疫保護防
火牆，成為解除NPI 與恢復正常生活的策略。不幸的是大部份的疫苗都是在關注
變異株演變之前研發，這些mRNA 疫苗以及以載體為基礎的 RNA 疫苗由於新病 毒基因變異產生的免疫逃脫，使得無論是從Delta 或到目前的 Omicron 變異株，
保護效益皆受到挑戰。目前研究已經發現疫苗對於有症狀的感染之保護力在 Omicron 變異株流行期已減弱，因此疫苗追加劑施打，結合在特定情境下進行快
篩檢測、個人的NPI 包含口罩與社交距離，成為後疫情時代流感化防疫策略主
軸。
在此一情境下，對於疫苗有效性的評估，包含降低中症、重症以及死亡的風 險的效益，在不同的國家以及不同變異株流行時期，推展精準防疫策略扮演重要 的角色，有助於了解何時第一代疫苗的保護力會因為免疫減弱而消逝，對於目前
的Omicron 變異株來說特別重要。由於大部分的國家已有多次的感染爆發，包含
D614G 、Delta 等，透過 T-細胞調節的記憶免疫反應可由前期感染誘發，所以人 群中真正透過全劑量以及半劑量疫苗所誘發產生的免疫反應的效果至今仍難以釐 清。
台灣於Omicron BA.2 感染爆發之前未曾有大規模的社區感染，因此台灣的
感染族群屬於首次感染Omicron，提供了不受到過去病毒在全球大規模感染的影
響評估疫苗的真實效果的難得機會，無論是疫苗基礎劑或是追加劑對於預防重症 以及死亡之效益。
再來針對SARS-CoV-2 症狀分期有所了解，闡明區分症狀前期個案（症狀前 期個案）模式和無症狀感染在不同變異株的流行病學軌跡，並且了解與 SARS- CoV-2 感染病毒有關的表現與症狀表現的進展機會，對發展個人化精準防疫模式 有所助益。
伴隨新冠肺炎後疫情時代也隨之出現了猴痘的非洲外多國傳播，其疫情也在 近五十年發生變化。猴痘的介入與預防措施包含保持良好衛生的接觸隔離
(Parker et al., 2014)、出現前驅症狀早期到醫療院所確定診斷 (Suarez Rodriguez et al., 2022)，以及高風險族群透過暴露前預防或是暴露後預防天花疫苗接種，降低 其病毒量與可能之傳播。隨著後疫情時代對於新冠肺炎的防疫政策逐漸鬆綁，可 能造成大型集會後透過性行為的猴痘疫情的蔓延 (Thornhill et al., 2022)，而在這 期間猴痘病毒的傳播模式是否已經發生改變，本來的預防措施是否有效，需要進 一步評估。
1.3 研究目的
新興與再現傳染病監測包含由族群到個人層次之疾病傳播，以及病程進展過 程，為達到個人化防疫措施亦涉及不同層次之介入措施，如公衛防疫措施、大規 模疫苗施打，以及抗病毒藥等防治策略。前述傳染病源與宿主特性造成的疾病傳 播與進展，以及包含由感染、中症、重症，以及死亡的臨床進展與防治策略以及 治療皆息息相關。這些多層次且多階段的疾病特性使傳統時間序列方法有其限 制。綜合以上，本研究之主要目的在於以新興傳染病(新冠肺炎)以及再現傳染病 (猴痘)實證資料應用隨機過程貝氏 DAG 模型達到以下目標:
(1) 在疾病傳播動態模型下拓展貝氏網絡 DAG 方法建立新興與再現傳染病爆發 流行之風險族群監測實證基礎
(2) 建構疾病進展多階段隨機過程貝氏網絡 DAG 模型運用於評估病毒以及宿主 特性在個人層次對於感染後病毒量變化之影響監測
(3) 運用貝氏網絡 DAG 方法結合廣義線性階層模式建立大流行與區域流行評估 架構，並且拓展為含括傳染病疾病進展狀態之介入效益評估方法。
第二章 文獻回顧
本章節將分為三個部分，第一部分為傳染病包含本研究主要運用實證資料，
猴痘以及新冠肺炎之疾病特性，第二部分為疫苗保護效益之評估方法，第三部分 為因應不同的監測與臨床問題運用於傳染病傳播與防治效益評估之方法文獻回 顧。
2.1 傳染病疾病特性 2.1.2 猴痘疾病傳播特性
猴痘和天花都屬於正痘病毒科(orthopoxvirus)，兩者均可通過飛沫接觸、近 距離接觸受感染的皮膚或受污染的物質而傳播。在先前的研究中，天花的 R0 估 計在 3.5 和 6.0 之間 (Gani & Leach, 2001)然而，猴痘的傳播主要是經由與感染 者或動物的密切或親密接觸而傳播。目前猴痘估計的 R0 低於天花，因此，猴痘 的傳播力低於天花 (Du et al., 2022)。
雖然猴痘具有較長的潛伏期（長達 21 天）和較長的傳染期，但是依據先 前的模型估計，其隔離前的可傳染時間小於 3 天。因而可以推測經由早期診斷 和隔離，使得猴痘的傳播得到了部分的遏止。因為猴痘的潛伏期一般為 6 至 13 天，但也可能長達 5 至 21 天 (Farahat et al., 2022; Miura et al., 2022; Thornhill et al., 2022)，所以高風險人群應該接種疫苗以防止猴痘的進一步傳播。
在 2022 年 4 月之後，截至 2022 年 7 月 25 日，在六大洲超過 75 個國 家/地區報告了 16,000 多例的猴痘病例。2022 年 6 月 23 日猴痘疫情被宣布為
「國際關注緊急公共衛生事件（PHEIC）」 (WHO, 2022)。超過 90% 的猴痘病 例，主要是經由男性間的性行為、天花疫苗免疫力減弱、對於新冠肺炎介入措施 放鬆，以及大型集會後的危險性行為活動等原因，造成猴痘病毒感染的全球爆發 (Suarez Rodriguez et al., 2022)。本論文也將以猴痘疫情作為研究主題，以貝氏 SIR 模式探討其在族群下的傳播。
2.1.1 新冠肺炎疾病傳播特性
新冠肺炎病毒病例起源於2019 年末，2020 年 1 月完成病毒分型檢測確認為 新型病毒 (2019-nCoV)，國際衛生組織(WHO)於 1 月 30 日宣布為國際關注緊急公 共衛生事件(PHEIC) (WHO, 2020)，病毒於 2 月更名為 Coronavirus disease 2019 (COVID-19)。感染途徑主要為呼吸道傳播，如近距離或是密閉空間的空氣溶膠傳
播，以及接觸呼吸道飛沫，或接觸過病毒污染之表面後觸碰眼、口、鼻而感染 (Meyerowitz et al.)。
以新冠肺炎基本再生數估計高達 5.70（95% CI：4.23-7.79）與 SARS 估計的 再生數範圍為2 到 3 (Riley et al., 2003)相比，新冠肺炎的傳播力則比較高。在這 種高的 R0 值應該可以解釋韓國、伊朗和義大利等多次的新冠肺炎流行開始時，
以及隨後在歐洲和美國的最初流行的浪潮中，由於公眾群聚導致的加速傳播，造 成短期確診病例激增的原因。由於新冠肺炎具有如此高的傳播能力，該病毒導致 大量病例，進一步給醫療保健系統以及社會和經濟基礎設施帶來巨大壓力。
由於病毒的不斷變異，不同關注變異病毒株所造成的症狀與潛伏期也可能有 所不同，過去也曾有人以系統性回顧研究不同變異株之潛伏期，收錄從 2019 年 12 月至 2022 年 2 月共 142 篇研究，獲得 Alpha、Beta、Delta 與 Omicron 變異株 之潛伏期(Incubation Period)分別為 5 天、4.5 天, 4.41 天與 3.42 天(Wu et al., 2022)。然而，至今仍然須面對不斷變異之病毒株傳播。
2.2 新冠肺炎疫苗保護效益評估 2.2.1 新冠肺炎疫苗效益評估研究
以往疫苗的功效(efficacy)、有效性(effectiveness)以及與保護(protect)需要經 過許多個階段的研究來證明 (WHO TEAM, 2021)。類似於藥物的開發，疫苗在通 過動物研究確認有效引發免疫反應的抗原後，在正常的情況下，同樣會進行三期 臨床試驗，以開發有效且安全的疫苗。而在全球衛生的緊急情況之下，採取「緊 急使用清單」 (Emergency use listing, EUL)的作法，目標是加速評估可能使用的 疫苗、檢測或是診斷工具能夠提供給遭遇緊急衛生事件的人群。
作為疫苗的抗原主要有兩類，一類是採用完整的病毒或是細菌或其替代品，
如減毒疫苗、病毒載體疫苗。新冠肺炎的腺病毒載體(ChAdOx1、Ad26.COV2.)疫 苗也屬於此類，另一類則是使用病毒或是細菌，能夠產生免疫識別的結構部分，
例如表面結構的蛋白質或醣類。兒童常規疫苗多半屬於此類，例如百日咳、破傷 風、白喉和腦膜炎球菌性腦膜炎等疫苗，新冠肺炎疫苗中的高端疫苗、
NOVAVAX 疫苗都屬於此類重組蛋白疫苗。而利用核酸(DNA、mRNA)方法所開 發的疫苗屬於嶄新技術，其原理是以DNA 或是 mRNA 方法向人體傳遞特定訊 息，使身體製造出免疫系統可以識別與產生免疫反應的蛋白質。除了本身抗原的
技術不同之外，疫苗中的賦形劑也可能是產生接種後過敏反應的原因。基於快速 開發過程中未觀察到的不良反應情形屬於未知數，在人群廣泛接種疫苗之情形須 持續受到安全監測，以確保可能的不良反應與風險適時被評估，而疫苗效益、副 作用以及產地因素也與人群是否能夠接受接種疫苗的疫苗猶豫有關 (Schwarzinger et al., 2021)。
疫苗的功效或是效力(efficacy)在三期臨床試驗，通常是第三期隨機分派試 驗之研究(phase 3 randomized controlled trial)確認安全以及有效性之後，接下來便 迫切的想了解其有效性(效益)，例如是否可以在脆弱族群發揮作用。多個第一代 新冠肺炎疫苗於2020 年 4 月後開始進行臨床研究，當時的主流病毒株 D614G，
而在之後陸續出現新的病毒株如Alpha、Delta，以及 2022 年起成為主流的 Omicron VOC。雖然在疫苗技術上的發展邁進了一大步，但持續演變的新冠肺炎 變異株，卻持續的帶來新的問題: 疫苗是否對於新的變異株仍然有效? 疫苗產生 的免疫效果可以保持多久?以及所發展的疫苗是否可以有效地阻止病毒的傳播?
(Koelle et al.) 在新冠肺炎開始於各國爆發，疫苗的效益成為了一個具有天時、地 利、人和的重要問題，常無法透過大型臨床試驗來評估，觀察性研究便能夠在接 下來衡量不同人群、區域以及變異株的變化。 (Lipsitch et al., 2020)。
2.2.2 新冠肺炎疫苗效益研究設計
目前針對新冠肺炎疫苗有效性之研究，常以疫苗之接種後能夠預防感染、臨 床症狀嚴重度、死亡等做為評估效益之主要追蹤結果，透過觀察性研究設計，如 個案對照研究(case-control study) (McConeghy et al., 2022; Tenforde et al., 2022; Yan et al., 2022)、檢測陰性研究(test-negative case-control study) (Bruxvoort et al., 2021;
Chemaitelly et al., 2021; Heidarzadeh et al., 2022; Lopez Bernal et al., 2021; Rearte et al., 2022)、世代研究(cohort study) (Chanda et al., 2022; Haas et al., 2021; Intawong et al., 2023; Kiss et al., 2022; Surie et al., 2022)等，多半以勝算比(Odds ratio)或是相對 風險(Relative risk)作為疫苗效益(Vaccine Effectiveness,VE)之評估。若以檢測陰性 研究設計，可以獲得勝算比(odds ratio,OR)可以利用( 1 – OR )×100% 得到 VE。
大型資料庫之回溯性個案對照研究則可以計算 (incidence rate ratios,IRR) 並估算 ( 1 – IRR ) ×100% 作為 VE。若實證資料型態為長期追蹤之時間橫斷資料，則有 等比風險假設模型(Cox regression model)等 (Arbel et al., 2022)，也有少數研究採
用數理模式之角度評估疫苗效益，如利用傳染病隔間模型(SEIR)等相關之研究 (Pageaud et al., 2022; Watson et al., 2022)或是其他的模型 (Sinclair et al., 2022;
Sonabend et al., 2021)。
檢測陰性研究(test-negative case-control study)通常以健康數據資料庫為基礎，
評估因為感染疫苗預防之疾病後尋求醫療行為的特定人群，其研究設計假設接種 疫苗的人以及未接種疫苗的人，沒有其他因素會影響兩組人群中新冠肺炎的感染 風險，但這些因素可能未被研究所測量，有測量者也需要利用模型加權調整，如 Thompson 等人調整了個案之時間、臨床數據以及人口學變數 (Thompson et al., 2021)，以獲得疫苗對於新冠感染後進入急診以及住院之預防效益。而此類個案對 照研究設計必須公平地做到研究組以及對照組，兩組人群進入研究之機會必須相 同，透過以上之研究設計可獲得疫苗效益，但在這些回溯性之設計也存在研究限 制，主要如下:
(1) 選擇偏誤(selection-bias): 積極接種新冠肺炎疫苗者可能也會積極檢測新冠病 毒，因而獲得醫療服務的機會也比未接種疫苗者多，造成疫苗的效益被低 估，另外就是未涵蓋在醫療保健服務中的族群無法被評估到，例如失業、非 法移工等無法納保者，由於它們獲得醫療保健服務的門檻相對更高 (Dean et al., 2021)，而他們在感染疾病的風險也需要被考量。而尋求醫療行為對於重症 以及死亡可能對於效益的影響較不敏感 (Dean et al., 2021)。
(2) 分類偏誤(classcification-bias): 在確定是否感染疾病，檢測的結果以及入院時
間性非常重要，例如是否在接種疫苗前14 天直到接種疫苗後 3 天內的新冠肺
炎感染入院，患者被分類在哪一組，在分類定義上會影響最後效益結果。在 Thompson 等人的研究中也透過敏感度分析衡量因未分類錯誤之影響程度 (Dean et al., 2021)。
2.2.3 時間對疫苗效益之影響
疫苗產品間產生有效性的時間具有部分差異。以病毒載體疫苗為例，英國疫
苗接種政策建議70 歲以上之老年人接種疫苗，但研究觀察到老年人的效益具有
時效以及有效性差異: 接種 ChAdOx1-S 一劑後 14-20 天(約 2-3 週)便具有效益，
28-34 天(4-5 週)到達 60% [95%信賴區間 41%-73%]之效益，35 天後更增加到 73% [95%信賴區間:27% -90%]之效益。而接種 mRNA 疫苗後 BNT162b2 一劑 後，約在28-34 天到達效益的 61% [95%信賴區間 51%- 69%] (Lopez Bernal et al., 2021)。
以 mRNA 疫苗為例，其在推出的前幾個月，對於降低感染以及住院率上對 於公共衛生具有許多貢獻，並且也在alpha(B.1.1.7)主流時期，於以色列 (Haas et al., 2021)、卡達 (Chemaitelly et al., 2021)、加拿大 (Nasreen et al., 2022)、美國 (Dickerman et al., 2022; Lauring et al., 2022)、英國 (Lopez Bernal et al., 2021)發揮 作用。以色列於2020 年 12 月初期全面接種 BNT 疫苗，並在 12 月至 2021 年 2
月之後實施兩次的封城與兩階段的解封，於3 月解封之後逐漸到達疫苗基礎劑接
種率56%，63%至少接種過一劑，對於 alpha 感染分別在整體、無症狀、有症 狀、入院、進入加護病房、死亡之有效性分別達95%、92%、98%、97%、97%、
98% (Haas et al., 2021)。
雖然研究發現最早的武漢株病毒的感染康復者，其血清中和抗體保護力可持 續至少約9 個月，包含記憶性 B 細胞、抗體、CD4+T 細胞、CD8+T 細胞免疫等 (Dan et al., 2021)。由於變異株的免疫逃脫，使得預防不同變異株感染之效益減低 (Lopez Bernal et al.)。而疫苗的保護也在 3-6 個月後逐漸遞減 (Chemaitelly et al., 2021; Feikin et al.; Gupta & Topol, 2021)。
由前述所知，目前疫苗保護力效果有其時效性，對於效益之評估須納入此因 素綜合評估。
2.2.4 其他的介入措施影響效益評估
全球疫苗資源有限也並非平均分配，除了疫苗的接種計次與接種人次因為 區域發展與政策執行時間而異之外 (Hannah Ritchie, 2020)，在未有疫苗的情況之 下各國衛生領導者也會採取可能的公衛應對措施，疫苗效益研究常僅針對有打疫 苗以及未打疫苗者情況主體作為評估，缺少對於其他介入措施的測量以及描述，
並且若想要描述以及評估他們也不是件容易的事，這些條件也可能因為時期不同 而改變。在疫苗施打區域是否有其他的介入措施以及治療政策，此也可能造成實 際效益被高估。
Hsu 等(2022)應為目前在評估疫苗效益的觀察性研究當中，首次共同評估 NPIs 以及抗病毒藥物治療使用時期前後之下，台灣的疫苗接種情況對於新冠肺炎 感染之臨床結果，包含降低中症、重症、以及死亡之效益，其除了回溯性之研究
設計之下，也納入貝氏廣義分層模型以將NPIs 以及疫苗之效益以分層結構區分
開來 (Hsu et al., 2023)。而我的論文也將利用此類模型以評估疫苗效益。
2.3 傳染病傳播與防治評估方法 2.3.1 廣義線性模型
傳染病監測時間序列資料較多屬於計數(count)資料，其特徵是可能具有 0 且非負整數(≧0)組成。不同於過去時間序列的方法學文獻主要以平穩高斯過程一 階或是二階的自我回歸(AR)或移動平均(MA)或是混合的 ARMA 的建模為主，計 數(count)資料在近年的發展並非以單一系統的模型主宰 (Davis et al., 2021)。
廣義線性模型 (Generalized linear model (GLM)) 由 Nelder 與 Wedderburn 在1972 所發表。經歷數年研究發展，GLM 從原本的基礎模型發展出更進階且能 處理非線性關係之模型。
假設兩個條件獨立之 變項𝑥𝑖與變項𝑦𝑖，變項𝑦𝑖之條件分佈屬於指數家族。
∫ 𝑓(𝑦𝑖; 𝛽𝑖|𝑥𝑖) 𝑑𝑦 = 1 ， (2-1)
且 𝐸(𝑦𝑖|𝑥𝑖) = 𝜇𝑖 (2-2)
在基礎模型中包含三個部分，指數家族(exponential family)、連結函數 link function(response)以及設計向量(design vector)，構成如下；
𝑦𝑖 = 𝒛𝑖𝛽 + 𝜀𝑖， i=1，….，n (2-3)
𝒛𝑖 = 𝒛(𝑥𝑖) 為設計向量(design vector)或稱共變異數向量(covariate vector) 用以定
義所有共變數之關係。式中β 為模型預估計之未知參數
𝑔(𝜇𝑖) = 𝜂𝑖 = 𝒛𝑖𝛽 ， (2-4)
g(˙) 為連結函數，作為連結模型中的預測變項(或反應變項) 𝑦𝑖。預測變項可假
設為連續型分佈，如常態分佈與Gamma 分佈，若為二元之預測變項則可運用
logistic 分佈或二項分佈。預測變項若為計數型態，則可以使用卜瓦松(Poisson)分 佈。
2.3.2 廣義自我回歸模型
廣義自我回歸模型(Generalized autoregressive (GAR))由 GLM 延伸以適應具 有自我回歸關係之時間序列資料。若yt 代表於 t 時間之觀察資料，令
𝐻𝑡 = {𝑦𝑡−1，𝑦𝑡−2， … . ，𝑦1，𝑥𝑡−1，𝑥𝑡−2， … … 𝑥1} 代表過去的觀察y 與有關變項 x，如 t-1, .t-2,…….等，
而連結函數h 為 𝜇𝑡= ℎ(𝒛𝑡𝛽)，𝒛𝑡= 𝒛(𝐻𝑡)為 𝐻𝑡 之函數。
假設y 在 t 時間下之條件期望值為
𝜇𝑡= 𝐸(𝑦𝑡|𝐻𝑡) (2-5) 條件變異數
𝜎𝑡2 = 𝑉𝑎𝑟(𝑦𝑡|𝐻𝑡) (2-6) 廣義自我回歸模型可以表示為以下:
𝑦𝑡= 𝜙𝑡−1𝑦𝑡−1+ 𝜙𝑡−2𝑦𝑡−2+ ⋯ . . +𝜙𝑡−𝑝𝑦𝑡−𝑝+ 𝑒𝑡 (2-7)
𝑦𝑡 = 𝒛𝑡𝛽 + 𝑒𝑡， t=1，….，n (2-8)
𝒛𝑡 = (1，𝑦𝑡−1，𝑦𝑡−2， … . . ，𝑦𝑡−𝑝，
𝛽 = (𝜙0，𝜙1， … … ，𝜙𝑝) (2-9)
在一階廣義自我回歸模型(GAR(1) model)共變異數結構以設計向量表示如下；
∑ 𝜎2 [
1 𝜌 𝜌2 … 𝜌𝑘−1 𝜌
𝜌2
⋮
1 ⋱ ⋱
⋱ ⋱ ⋱
⋱ ⋱ 1
⋮ 𝜌2
𝜌 𝜌𝑘−1 … 𝜌2 𝜌 1 ]
(2-10)
Hsu 等人應用卜瓦松廣義自我回歸模型描述氣候因素對於日本腦炎以月份為 單位的滯後效應，其模型如下
𝑙𝑜𝑔[𝐸(𝑦𝑡|𝑦𝑡−1， … . ，𝑦𝑡−𝑙; 𝑥1. … 𝑥𝑝，𝑡，𝑚1… 𝑚11)]
= 𝑙𝑛( PY) + 𝑎 + ∑ 𝑣𝑘𝑇𝑀𝑡−𝑘 + ∑ 𝜋𝑖𝑅𝑎𝑡−𝑖+ 𝛽1𝑥1+ ⋯ . +𝛽𝑝𝑥𝑝+ 𝛿1
𝑚
𝑖=1
𝑡
𝑛
𝑘=1
+𝜔1(𝑚1) + ⋯ + 𝜔1(𝑚11) + ∑𝑙𝑗=1𝜙𝑗𝑦𝑡−𝑗， (2-11)
在本模型中等號的左邊 𝑦𝑡 代表每月的日本腦炎個案數，當時間(月)表示為 t 𝑦𝑡−𝑙 為 t 時間前 l 個月的日本腦炎個案數，
𝑥1,….𝑝 為與日本腦驗個案數有關的相關變項，如豬隻養殖密度、變異程度與地
理區域，t 為線性趨勢，𝑚1,….,11為季節效應表示為11 個虛擬變數。
等號的右邊𝑙𝑛( PY) 代表發生率(人-年) (右移為 Y 分母 offset)，
∑𝑛𝑘=1𝑣𝑘𝑇𝑀𝑡−𝑘與∑𝑚𝑖=1𝜋𝑖𝑅𝑎𝑡−𝑖 代表隨時間(t-k)與時間(t-i)變化之氣溫、濕度以及 其參數𝑣𝑘 and 𝜋𝑖 的積，其他的共變異係數則以𝛽1，….，𝑝，𝛿，𝜔，𝛾𝑗表示。
Lopmen 等人將每日實驗室檢測諾羅病毒的計數資料以時間序列適應卜瓦松 迴歸建模如下:
𝐸(𝑌𝑖) = 𝑒𝑥𝑝{𝛼 + 𝛽1𝑖 ⋅ 新變異株+𝛽2免疫因子𝑟𝑖}
+𝛽3((∑ 溫度
𝑙=0，48 𝑥−1)
𝑙 ) + 𝛽4((∑ 濕度
𝑙=0，35 𝑥−1)
𝑙 )
+[干擾因子] + [背景季節性] + [自我回歸係數] (2-12)
透過檢視局部自我回歸圖形(partial autocorrelation)的離差量殘數差(deviance residuals)後，模型中包含 21 個 自我回歸係數 (Lopman et al., 2009) 且表示為
[自我回歸係數] = ∑ 𝛿𝑖離差量
𝑖=1，21 𝑥−1 (2-13) Zeger 等人發展時間序列之線性模型，使用自我回歸(AR)於美國 1970-1983 年小兒麻痺疫情監測資料，並且假設每一筆的計數資料都具有潛在的相關性，來 自隱藏過程，且此過程不可觀察 (ZEGER, 1988)。𝑦𝑡為條件獨立計數序列且性質 如下:
𝑢𝑡=E(𝑦𝑡|𝜀𝑡)=exp(𝑥′𝑡𝛽)𝜀𝑡 ， (2-14)
𝑤𝑡 = 𝑣𝑎𝑟( 𝑦𝑡|𝜀𝑡) = 𝑢𝑡， (2-15)
在2-14 式中 𝜀𝑡 為隱藏關聯過程，導入過度擬合與𝑦𝑡的自我回歸關係，最終希 望找到𝑦𝑡的邊際分佈
𝜇𝑡 = 𝐸(𝑦𝑡) = 𝑒𝑥𝑝( 𝑥′𝑡𝛽) (2-16)
𝑣𝑡 = 𝑣𝑎𝑟( 𝑦𝑡) = 𝜇𝑡+ 𝜎2𝜇𝑡2 (2-17)
2.3.3 兩階段馬可夫過程與傳染病流行平衡態
由 Cox 與 Miller(1965)著作中所提及之兩階段馬可夫過程(Two state Markov process) (D.R.Cox, 1965)，可應用於描述傳染病狀態改變。以各國新冠肺炎開放資 料為例，以其隨時間變化的有效再生數𝑅𝑒作為傳播控制狀態之分界，以𝑅𝑒<1 為
非流行期，𝑅𝑒 ≥1 為流行期。透過真實世界數據，可以估計出在兩種狀態間變化
的重要參數，且包含不同時期流行的變異株，如D614G、 Alpha、 Delta、 以及 Omicron 之參數的變化，估計兩種狀態的平衡分佈(equilibraium distribution)。
透過有效再生數將狀態空間為非流行期與流行期。透過流行病學隔間模型 如SIR、SEIR (Kermack & McKendrick, 1927)獲得之有效再生數(𝑅𝑒 或稱 𝑅𝑡)可得 知疫情狀態:
𝑅𝑒 ≥1 ， 傳染性疾病疫情處於流行期 state (1)，
𝑅𝑒 <1 ， 傳染性疾病疫情處於非流行期 state (0)
𝑅𝑒 代表著在介入措施下，一個新發個案可以傳染個個案數。當有措施介入時，
𝑅𝑒 將下降。
在短時間內(𝛥𝑡)的兩個狀態的轉移機率將由以下表示
𝛼𝛥𝑡 + 𝑜(𝛥𝑡) 代表從非流行期 state (0) 到流行期 state (1)，
𝛽𝛥𝑡 + 𝑜(𝛥𝑡) 代表從流行期 state (1) 到非流行期 state (0)，
在 Markov property 條件下， 兩狀態的機率分佈在 t 時間為
𝑝0(𝑡) =𝛼+𝛽𝛽 + {𝑝0(0) −𝛼+𝛽𝛽 } 𝑒−(𝛼+𝛽)𝑡 𝑝1(𝑡) = 𝛼
𝛼+𝛽+ {𝑝1(0) − 𝛼
𝛼+𝛽} 𝑒−(𝛼+𝛽)𝑡 (2-18) 令 𝑝1(0) = 1 − 𝑝0(0) 當 t=0
而在t 時間夠長時，一定的 t 時間後，到達之平衡分佈為 𝑝0 = 𝛽
𝛼+𝛽， 𝑝1 = 𝛼
𝛼+𝛽 (2-19) 由上式中可見到平衡分佈並未與起始狀態有關，當𝑋(𝑡)的邊際分佈在所有 t>0
皆成立下，平衡分佈(equilibrium distribution)也可以稱為過程的穩定分佈 (stationary distribution)。
令隨機變數𝑋(𝑡) 表示 t 時間下屬於兩個狀態的其中一個狀態(state1 或 state0)，h 為一小段時間。 在時間 t 到 t+h(h>0)，維持在 State 1 轉移機率可以表 示為
𝑝𝑟𝑜𝑏{𝑋(𝑡 + ℎ) = 𝑋(𝑡) = 1}
= 𝑝𝑟𝑜𝑏{𝑋(𝑡) = 1} × 𝑝𝑟𝑜𝑝{𝑋(𝑡 + ℎ) = 1|𝑋(𝑡) = 1}
= 𝛼
𝛼+𝛽𝑝1(ℎ) (2-20) 其中𝑝1(ℎ)為在時間 h 內是 State 1 的機率
其他三種轉移機率分別依前述推導如下，
𝑝𝑟𝑜𝑏{𝑋(𝑡 + ℎ) = 𝑋(𝑡) = 0} = 𝛽
𝛼+𝛽𝑝0(ℎ)，
𝑝𝑟𝑜𝑏{𝑋(𝑡 + ℎ) = 0，𝑋(𝑡) = 1} = 𝛽
𝛼+𝛽𝑝1(ℎ)，
𝑝𝑟𝑜𝑏{𝑋(𝑡 + ℎ) = 1，𝑋(𝑡) = 0} = 𝛼
𝛼+𝛽𝑝0(ℎ). (2-21) 由於馬可夫過程的同質性(homogeneous property)，最終求得的聯合分佈(joint distribution)只會與時間間隔(Duration)的長短有關、與起始的時間點無關。
舉例如
𝑋𝑛1， … … 𝑋𝑛𝑙 (2-22) 以及
𝑋𝑛1+𝑘， … … 𝑋𝑛𝑙+𝑘 (2-23)
在時間間隔為n1,n2,……,nl 與每一段時間延後 k，具有相同的聯立分佈 兩階段馬可夫過程在連續時間中假設嚴格穩定，實際執行上可能需要足夠 的疾病傳播觀察資料進行評估。
在兩階段馬可夫過程中
𝐸{𝑋(𝑡)} = 𝑝𝑟𝑜𝑏{𝑋(𝑡) = 1} = 𝛼 𝛼 + 𝛽， 𝑉𝑎𝑟{𝑋(𝑡)} = 𝛼𝛽
(𝛼+𝛽)2. (2-24)
運用以下兩階段遞移矩陣以馬可夫特性之條件機率可估算𝜌(ℎ) 𝐸{𝑋(𝑡)𝑋(𝑡 + ℎ)}
= 𝑝𝑟𝑜𝑏{𝑋(𝑡) = 𝑋(𝑡 + ℎ) = 1}
= 𝑝𝑟𝑜𝑏{𝑋(𝑡) = 1} × 𝑝𝑟𝑜𝑏{𝑋(𝑡 + ℎ) = 1|𝑋(𝑡) = 1}
= 𝛼
𝛼+𝛽𝑝1(ℎ) (2-25)
透過式 (3-6)對於 h≥0 ，可得到自我回歸相關係數𝜌 𝛾(ℎ) = 𝐶{𝑋(𝑡 + ℎ)，𝑋(𝑡)}
= 𝛼𝛽
(𝛼+𝛽)2𝑒−(𝛼+𝛽)ℎ (2-26) 𝜌(0) =𝛾(ℎ)
𝛾(0)= 𝑒−(𝛼+𝛽)ℎ (2-27)
對於常態分佈隨機變數𝑋𝑡1， … … 𝑋𝑡𝑙而言，其聯合分佈為多元常態分佈 (Multivariate normal)，𝐸(𝑋𝑡)為 t 之函數，而自我回歸相關係數 𝐶(𝑋𝑡1，𝑋𝑡2) 為 𝑡1, 𝑡2 之函數且共變異數矩陣(covariance matrix)為 𝛾(𝑡1，𝑡2) = 𝐶(𝑋𝑡1，𝑋𝑡2).
對於可達到平衡分佈之隨機過程，
𝐸(𝑋𝑡) = 𝜇
𝐶(𝑋𝑡+ℎ，𝑋𝑡) = 𝛾(ℎ) (2-28)
上述 𝐸(𝑋𝑡) 在為恆定且對應於自我回歸相關係數的 C(Xt+h, Xt)在𝑡1− 𝑡2間存 在。
兩階段馬可夫過程在應用甚廣於探討醫學領域資料的穩定分佈(stationary distribution)或者又稱為平衡分佈(equilibrium)。無論資料是連續(continuous)或者是 離散(discrete)資料，皆可以透過馬可夫鏈(Markov chain)方法獲得平衡分佈。
Biritwum 與 Odoom (1995) 應用兩階段馬可夫過程研究嬰兒之健康(well)與 生病(sick)的行為模型轉換。在其研究中，以嬰兒 0 月到 18 個月的資料，不只找 到了兩個階段隨時間改變的轉移機率，且也能夠得知轉移機率隨時間改變之情 況，如5 個月以前與嬰兒月齡有關，6 個月後則月齡因素逐漸消失，在 12 個月大 以後，健康與生病的狀態將會維持穩定。 (Biritwum & Odoom, 1995)
Zhang(2022)等應用於研究人類乳突病毒(HPV)在男男間高風險性行為傳播的
感染發生率與病毒清除率。他們得知不同的HPV 病毒株有不同的持續時間，如
HPV6(12.6 個月)、HPV 16(11.9 個月)、HPV42 與 HPV56 則最短，綜合與階段間 的轉移速率與其共變數，他們得已成功的考量特定期間的動態感染過程，包含感
染、病毒清除以及再感染，並辨別HPV 感染發生率與清除率之風險因子。
(Zhang et al., 2022). Dona (2016, 2022)等也用同樣的方法研究男男性行為者中之 HPV 傳播發生率。 (Dona et al., 2022; Dona et al., 2016)
其他二階段馬可夫方法的臨床應用包含阿茲海默症之輕度認知困難，
(Ficiara et al., 2021)，血糖測量 (Woldaregay et al., 2020)，以及人類後天免疫不 全症候群傳染窩等 (Hawkins, 2016)。
以上回顧之方法可有效地建立對於族群層次下之傳染病模型，但對於個人 層次討論傳染病，個人的傳染狀態以及疾病狀態在群體下並無法被觀察，如新冠 肺炎的病毒量無法被觀察，但其對於傳染狀態以及疾病致死是重要因素，個人層 次的模型也是本論文之討論方向。接下來將由族群層次至多階段馬可夫模型再來 以貝氏網路探討個體層次下的傳染病疾病狀態與傳染階段之變化。
2.3.4 多階段隨機過程用於疾病狀態預測
由於資料的特性也可以利用其他的不同於SEIR 模型的觀點來了解一個傳 染病的傳播參數。Chen 等(2022)利用 Bayesian Markov Chain Monte Carlo (MCMC) method 應用於隨機過程四階段模型的方法，估計台灣入境旅客在不同新冠肺炎病 毒株流行期間以及疫苗接種、NPIs 改變下的中位數未有症狀傳播期(Median of pre-symptomatic transmission time)、以及無症狀的發生率。
參考陳等(2022)的模型，令 𝑋(𝑡)是在 t 時間下的隨機變數，其中可能為四個 狀態中的其中一種: 未感染期(state 1)、感染症狀前期個案期(state 2)、感染症狀期 (state 3)以及無症狀感染(state 4)。其中有兩個重要的假設:
(1)所有的狀態都是由未感染(state 1)出發。
(2)無症狀感染(state 4)永遠都不會有症狀，僅以檢測(如 RT-PCR)發現有感染。
他們的方法是透過競爭風險的觀點將症狀期(state 3)以及無症狀感染(state 4) 分別為兩種感染後的最終結果，無症狀感染(state 4)之個案永遠不可能會出現症 狀。依照以上定義計算感染後未出現症狀直到出現症狀的發生率以及瞬間速度(危 險函數ℎ(𝑡))，轉換成中位數無症狀傳染時間(Median of pre-symptomatic to
symptomatic time, MPTT)。
假設個體 i，在𝑡1出發，於旅遊飛行過程與感染者接觸(包含有症狀、未出 現症狀與無症狀)並且於 𝑡2時間抵達到站檢測，並且經過檢測之後可能出現以上 四種狀態，最後經過隔離至𝑡3採檢陰性之後解除隔離。
從出發到檢測這段期間(𝑡1, 𝑡2) 由未感染(state 1)到感染症狀前期個案期 (state 2)的轉移機率，在時間 s 下，
𝑃12(𝑡1, 𝑡2) = ∫ 𝑖(𝑠 − 𝑡𝑡𝑡2 1) × [1 − 𝐹(𝑡2− 𝑠)] × [1 − 𝐴(𝑠 − 𝑡1)]𝑑𝑠
1 (2-29)
𝑖(𝑡) 為由 state 1 到 state 2 的機率密度函數(pdf)
由 感染症狀前期個案期(state 2)演變到感染症狀期(state3)的累積分佈以 F(t)表 示。
由 未感染 (state 1)演變到無症狀感染(state 4)的累積分佈以 A(t)表示。
當個案由感染症狀前期個案期(state 2)到出現症狀(state 3)的轉移機率為
𝑃13(𝑡1, 𝑡2) = ∫ 𝑖(𝑠 − 𝑡𝑡𝑡2 1) × [1 − 𝐴(𝑠 − 𝑡1)] × 𝐹(𝑡2− 𝑠)𝑑𝑠
1 (2-30)
且由𝐹(𝑡2− 𝑠)可以得知此個案絕對不會變成無症狀。
若個案由未感染(state 1)之後變成感染後無症狀(state 4)，則其轉移機率為
𝑃14(𝑡1, 𝑡2) = ∫ 𝑎(𝑠 − 𝑡𝑡𝑡2 1) × [1 − 𝐼(𝑠 − 𝑡1)]𝑑𝑠
1 (2-31) 在 (𝑡1, 𝑡2)期間，維持一直在未感染狀態(state 1)的轉移機率為
𝑃11(𝑡1, 𝑡2) = 1 − 𝑃12(𝑡1, 𝑡2) − 𝑃13(𝑡1, 𝑡2) − 𝑃14(𝑡1, 𝑡2) (2-32)
(𝑡2,𝑡3)期間，感染症狀前期個案期(state 2)到感染症狀期(state 3)的轉移機率為
𝑃23(𝑡2, 𝑡3) = ∫ 𝑓(𝑟)𝑑𝑟𝑡𝑡3
2 = 𝐹(𝑡3− 𝑡2) (2-33) r 為瞬時時間。
𝑖(𝑡)與接觸感染的危險程度有關。
𝑓(𝑡) 是感染還沒出現症狀到有症狀的時間函數，與疾病的進展速度有關。
兩者有許多方式可以描述疾病進展隨時間的變化，透過指定不同的分佈，像是 單調函數或是非單調函數。Chen 等(2022)認為依據 Sethuraman(2020)等的研究，
以log-logistic 分佈比起 Weibull 分佈更為適合去描述疾病從沒有出現症狀到有症
狀的進展時間，因其描述病毒量分佈。
由未出現症狀到出現症狀的風險函數以 log-logistic 形式表示如下
ℎ(𝑡) = 𝑒𝜃𝜅𝑡𝜅−1
1+𝑒𝜃𝑡𝑘 0 ≤ 𝑡, 𝜅 > 0 (2-33)
在 (𝑡1,𝑡2)期間發展至有症狀的累積風險則由以下表示 𝐹(𝑡1, 𝑡2) = 1 − 1
(1+𝑒𝜃⋅(𝑡2−𝑡1)𝜅) (2-34) Chen (2022)等的研究，透四階段隨機模式建模，應用於入境檢疫或隔離政策 的調整可更佳精準 (S. L. Chen et al., 2022)。在他們的四階段模式當中，有幾項特 點值得注意:
(1) Omicron 相較於其他變異株的高傳播特性(從未感染(state 1)到感染未出現症狀 (state2))結果與 SEIR 模型結果接近。
(2) 將無症狀感染與感染未出現症狀分成兩個不同的路徑，避免混淆。同樣的 SEIR 模式也能夠做到。
(3) 在出現症狀的風險隨時間函數中，以 log-logistic 分佈可接近早期風險上升，
但在某時間點之後下降之變化趨勢，具有生物醫學上的合理性。
(4) 了解 SARS-CoV-2 自然感染史以及因應來自全球變異株的傳染病重要參數變 化，包含可以獲得傳染病模式常見到的有效或是基本再生數(R1/R0)。中位數 症狀前期個案傳染但時間(MPTT)能提供病毒生成期以及潛伏期的資訊。
而無症狀感染發生率也隨時間變化增加，除了與病毒造成的疾病嚴重度減輕與 傳播增加之關係之外，也可能與人群開始接種疫苗有關。
(5) 透過模擬可以知道不同的樣本數在不同感染階段的預測的檢定力情況，也反 映在目前資料中，在2020 年 3-6 月，因為入境旅客的檢測出 D614G 的樣本數 不足，在模擬研究中無法達到足夠的檢定力，因此此模型仍然需要更多資料 以支持模型參數學習過程。
本研究將運用既有之時間序列方法以前述之新興(新冠肺炎)與再現(猴痘)傳染
