德國的 COVID,時間序列的 LOO-CV
Science [1] 中的論文推斷了德國 COVID 傳播的變化點。假設一個(紅色)、兩個(橙色)和三個(綠色)變化點,作者擬合了每日病例數。每個更改點都會向模型添加兩個參數。
很難相信三個變化點模型捕捉到了一個變化點模型中缺少的一些基本物理現實。“檢測到三個對應的變化點”的結論是基於留一法交叉驗證(LOO-CV)分數的比較:
[loo log-score] [standard error] [effective number of parameters] three points 787 15 13 two points 796 17 13 one point 819 17 13
pymc3.compare(..., ic='LOO', scale='deviance')
返回(d_loo
是相對差異,dse
是每個模型與排名靠前的模型之間得分差異的標準誤差):loo p_loo d_loo weight se dse three points 786.543 13.3241 0 0.933612 15.2098 0 two points 795.797 12.5467 9.25366 0.0662461 16.6689 4.88424 one point 819.280 13.3403 32.737 0.000141764 17.106 8.25306
pymc3.plot_elpd
顯示這個情節:LOO-CV 使用是否正確?
有與論文 [4] 相關的電子信件交流和作者的技術報告 [5]。
- Dehning, J., Zierenberg, J., Spitzner, FP, Wibral, M., Neto, JP, Wilczek, M., & Priesemann, V. (2020)。推斷 COVID-19 傳播的變化點揭示了乾預措施的有效性。科學。 http://dx.doi.org/10.1126/science.abb9789
代碼和數據:https ://zenodo.org/record/3780722 2. Vehtari, A.、Gelman, A. 和 Gabry, J. (2017)。使用留一法交叉驗證和 WAIC 的實用貝葉斯模型評估。統計和計算,27(5),1413-1432。
https://doi.org/10.1007/s11222-016-9696-4(PyMC3 中使用相同的引用) 3.
pymc3.loo
以及pymc3.compare
文檔和代碼https://docs.pymc.io/api/stats.html
https://github.com/arviz-devs/arviz/blob/18797b81/arviz/stats/stats.py 4. https://science.sciencemag.org/content/early/2020/05/14/science.abb9789/tab-e-letters 5. https://github.com/Priesemann-Group/covid19_inference_forecast/blob/aa2516680c1f3800225e5e7abce07607ad03a795/technical_notes_dehning_etal_2020.pdf
概覽 速記
- 具有三點的模型確實更適合。
- 與三點的配合只是稍微好一點。
- 只有一分的模型還不錯。locv 分數的差異可能表明具有更多點的模型是顯著/可能/可能的改進,但影響大小只是很小。
- 即使三點模型很合適,它也可能不需要是物理現實。
- 更好的擬合應該被解釋為確認具有一個轉折點的零假設 SIR 可能不正確(在“不完全正確”的意義上,它可能仍然是一個相當好的描述)。它不能確認具有三個點的替代模型是正確的(在物理意義上)。正確的模型(真實模型)實際上可能是不同的模型(例如,平滑過渡而不是變化點)。它只是確認替代模型表現更好。
很難相信三個變化點模型捕捉到了一個變化點模型中缺少的一些基本物理現實。
三個變化點的擬合確實更準確
不難相信具有三個變化點的模型會做得更好。一個簡單的 SIR 模型(假設所有人的均勻混合)並不完全符合現實。這些變化點將有助於彌補這一缺點(使其更靈活,能夠適應更廣泛的不同曲線)。
但它可能無法捕捉到物理現實
然而,你懷疑它是否捕捉到了物理現實是正確的。SIR 模型被設計為機械模型。然而,當它不夠準確時,它實際上就變成了一個經驗模型。
基本參數可能不一定代表某些物理現實。(如果你願意,你可以擬合一個顯然沒有任何物理現實的機械模型)
在不改變流行病學參數的情況下,有很多方法可以降低增長率。在空間和網絡 SIR 模型中,這可能是由於局部飽和(例如,參見此處的示例)。
因此
- 與 SIR 模型的擬合將低估 $ R_0 $ 值(因為較低 $ R_0 $ 值往往更適合曲線中的偏轉)。
- 當 SIR 模型通過變化點變得更加靈活時, $ R_0 $ 最初可能會更高,但擬合將表明生長參數下降 $ \beta $ 這實際上可能不存在。
一變點
那麼,這些變化點是虛構的嗎?我想不是。的價值 $ \beta $ 在那個模型中確實改變了很多。
我不認為這種增長率的下降不會發生,而是由於對 SIR 模型的奇怪調整使其自動下降。
雖然當 $ N $ 較低,我認為它不包括在模型參數之一中並且似乎是固定的,那麼在不改變流行病學參數的情況下,增長率可能會急劇下降。
$$ \frac{dI}{dt} = \overbrace{\frac{S}{N}}^{\substack{ \llap{\text{If N or}}\rlap{ \text{ S = N-I}} \ \llap{\text{are over/un}}\rlap{ \text{der estimated} }\ \llap{\text{then the dro}}\rlap{ \text{p in this term}} \ \llap{\text{becomes un}}\rlap{ \text{derestimated}} \ }} \underbrace{\beta}_{\substack{ \llap{\text{In that case}}\rlap{ \text{ $\beta$ will get}} \ \llap{\text{underestimate}}\rlap{ \text{d in order to}} \ \llap{\text{correct for the w}}\rlap{ \text{rong S/N term} }\ }} I - \mu I $$
如果錯了 $ N $ 使用,則模型將被推送以糾正此問題。當我們錯誤地假設所有案例都在測量時也是如此(因此低估了案例的數量,因為我們沒有包括漏報)。
但無論如何,我想可以合理地說,有轉折點/下降 $ \beta $ 有許多流行病學曲線顯示增長率迅速下降。我認為,這不是由於飽和(增強免疫力)等自然過程,而是主要由於參數的變化。
兩三點
這些模型的效果實際上只是非常微妙的。這些額外的變化點所做的是使從增長到下降的變化更加平滑,而且這只發生在很短的時間內。因此,在 3 月 8 日至 22 日之間,您將獲得三個小步驟,而不是一大步。
不難相信你會得到一個平穩的下降 $ \beta $ (許多機制可能會產生這種變化)。更難的是解釋。變化點與特定事件相關。
例如,請參閱摘要中的此引用
“專注於 COVID-19 在德國的傳播,我們發現了與公開宣布的干預措施的時間密切相關的有效增長率變化點”
或者在文中
第三個變化點……在 3 月 24 日推斷 $ (CI [21, 26])] $ ; 這個推斷的日期與第三次政府乾預的時間相匹配
但這只是猜測,可能只是虛構。尤其如此,因為他們將先驗準確地放置在這些日期上(標準偏差或多或少與可信區間的大小相匹配,我們有“後驗分佈 $ \approx $ 先驗分佈”,這意味著數據沒有添加太多有關日期的信息):
所以這不像他們做了一個三變化點模型,結果巧合地匹配了特定乾預的日期(這是我快速瀏覽文章後的第一個解釋)。他們沒有檢測到變化點,更像是該模型具有與特定乾預措施密切相關的傾向,並將“檢測到”點放置在干預措施的日期附近。(此外還有報告延遲的免費參數,它允許在曲線變化日期和乾預變化日期之間有幾天的靈活性,因此無法精確定位/檢測/推斷變化點的日期非常準確和總體上它更模糊)
留一出交叉驗證。
LOO-CV 使用是否正確?
我相信 LOO-CV 已正確應用。(但解釋很棘手)
我必須深入研究代碼才能確切知道,但我沒有理由懷疑它。這些分數意味著具有三個變化點的函數沒有過度擬合,並且能夠更好地捕捉模型的確定性部分(但並不是說具有三個點的模型比具有一個點的模型好太多,只是一個小的改進)。
- 函數沒有過度擬合併不奇怪。有相當多的數據點可以消除噪聲並防止擬合函數捕獲太多噪聲而不是潛在的確定性趨勢。
- 這三個變化點能夠更好地捕捉確定性模型並不奇怪。開箱即用的標準 SIR 模型並不是很合適。您可以使用高階多項式擬合或樣條曲線獲得類似的改進,而不是更改點。更改點改進模型可能不需要是因為機械的潛在原因。
您可能會想,嘿,但是紅色、橙色、綠色三條曲線之間的細微差別呢?
是的,確實差異很小。變化點僅在一小段時間內出現。雖然 LOO-CV 分數的差異(從 819 到 796 到 787)可能表明具有某種意義,但這可能不需要與“大”效應相關,替代模型的效應也不需要與某些現實機制。例如,請參見下圖中的示例,其中還有一個 $ x^2 $ 項能夠顯著提高擬合度,但效果差異很小,“真實”效果是 $ x^3 $ 術語而不是 $ x^2 $ 學期。但是對於那個例子,對數似然分數有很大不同:
> lmtest::lrtest(mod1,mod2) Likelihood ratio test Model 1: y ~ x Model 2: y ~ x + I(x^2) #Df LogLik Df Chisq Pr(>Chisq) 1 3 15.345 2 4 19.634 1 8.5773 0.003404 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
此外,微小的差異也可能是有問題的。它可能不是很重要,尤其是當您認為噪聲可能相關時。因此,某種程度的過度擬合可能不會在留一式 CV 中受到懲罰。
示例圖像和代碼:
set.seed(1) x <- seq(0,1,0.02) ydeterministic <- x + 0.5*x^3 y <- ydeterministic + rnorm(length(x),0,0.2) mod1 <- lm(y~x) mod2 <- lm(y~x+I(x^2)) plot(x,y, main="small but significant effect", cex.main = 1, pch = 21, col =1, bg = "white", cex = 0.7, ylim = c(-0.2,1.7)) lines(x,mod1$fitted.values,col="red", lty = 2) lines(x,mod2$fitted.values,col="blue", lty =2) lines(x,ydeterministic, lty = 1 ) lmtest::lrtest(mod1,mod2) legend(0,1.7,c("true model: y = x + x³", "fit 1: y = x", "fit 2: y = x + x²"), col = c("black","red","blue"), lty = c(1,2,2), cex = 0.6)
此示例適用於線性模型,而不是貝葉斯設置,但它可能有助於直觀地了解“顯著但很小的影響”的情況,以及如何根據對數似然值而不是影響大小進行比較,與此相切。