COVID-19 Data Analysis

Stima infetti reali

Stima infetti reali con fast backcasting.

Max Pierini


ITALIA

Possiamo definire il tasso di rilevamento $\mathrm{DR}_t$ (Detection Rate) al tempo $t$ come la percentuale di eventi $n(\hat{\alpha}_t)$ stimati sugli eventi $n(\alpha_t)$ osservati in $t$

$$ \mathrm{DR}_t = \frac{ \sum_{i=0}^{t} n(\hat{\alpha_i}) }{ \sum_{i=0}^{t} n(\alpha_i) } $$

Un $\mathrm{DR}_t$ pari al 100% significa che tutti i casi vengono correttamente rilevati. Al contrario un $\mathrm{DR}_t$ minore del 100% significa che una parte dei casi, ovvero $1-\mathrm{DR}_t$, sfuggono al rilevamento.

L'attuale $\mathrm{DR}_t$ stimato è quindi pari a 46.1% (95% CI: 25.8%, 70.3%), ovvero si stima che il 53.9% dei casi sarebbe sfuggito al rilevamento (95% CI: 29.7%, 74.2%).

Dal grafico dei nuovi casi e del tasso di rilevamento, notiamo che la maggior parte dei mancati rilevamenti è avvenuta nella prima fase della curva epidemica, indicativamente tra Febbraio e Giugno. In seguito, verosimilmente grazie alle campagne di screening e contact-tracing, la differenza tra casi reali (stimati) e osservati risulta molto meno marcata.

METODO

Concetto di backcasting

Il metodo di fast backcasting di EpiData.it si basa sul backcasting per simulazione Monte Carlo di Phipps et al. 2020 ma implementato per esecuzione veloce senza uso di simulazioni.

Siano dati due eventi $\alpha$ e $\omega$.

Sia definita Alpha-Omega-Rate $\mathrm{A}\Omega\mathrm{R}$ la probabilità condizionata non nulla del verificarsi dell'evento $\omega$ in seguito all'osservazione dell'evento $\alpha$

$$ P(\omega | \alpha) = \mathrm{A}\Omega\mathrm{R} \;\,\; >0 $$

Sia posta pari a 1 la probabilità condizionata che si sia verificato $\alpha$ data l'osservazione di un evento $\omega$

$$ P(\alpha | \omega) = 1 $$

ovvero

  • l'evento $\omega$ ha probabilità $\mathrm{A}\Omega\mathrm{R}$ di accadere in seguito ad un evento $\alpha$

  • l'osservazione di un evento $\omega$ comporta sicuramente l'accadimento di un evento $\alpha$ nel passato

Pertanto, data l'osservazione di $n(\omega)$ eventi $\omega$, si può stimare che il numero $n(\alpha)$ di eventi $\alpha$ occorsi sia

$$ n(\alpha) = \frac{n(\omega)}{\mathrm{A}\Omega\mathrm{R}} $$

Sia data una distribuzione di probabilità $p(\omega_{t_0 + \tau})$ nel tempo $t$ del verificarsi dell'evento $\omega_{t_0+\tau}$ dopo $\tau$ unità di tempo dall'osservazione dell'evento $\alpha_{t_0}$ al tempo $t_0$

$$ p(\omega_{t_0 + \tau}) = \mathbf{PDF}[\omega]_{t_0 + \tau} $$

Noto il numero $n(\omega_t)$ di eventi $\omega$ nel tempo $t$, possiamo stimare che il numero $n(\alpha_{t_0})$ di eventi $\alpha$ occorsi in $t_0<t$ siano pari alla somma degli eventi $n(\omega_{t_0 + \tau})$ moltiplicati per la densità di probabilità $\mathbf{PDF}[\omega]_{t_0 + \tau}$ al tempo $t_0 + \tau$ divisa per la probabilità condizionata $p(\omega|\alpha)$

$$ n(\alpha_{t_0}) = \sum_{\tau=1} \frac{ n(\omega_{t_0 + \tau}) }{ \mathrm{A}\Omega\mathrm{R} } \mathbf{PDF}[\omega]_{t_0 + \tau} $$

Indicando con $t=T$ l'ultimo giorno di osservazione, possiamo affermare che soprattutto negli ultimi giorni della serie temporale, il numero di eventi $n(\alpha)$ sarà ovviamente sottostimato perché i futuri eventi $\omega$, di cui potrebbero essere causa, non sono ancora accaduti.

Supponiamo, ad esempio, di avere i seguenti dati:

Vogliamo stimare il numero di eventi $n(\alpha)$ con metodo di backcasting dagli eventi $\omega$ noti e sapendo che la probabilità di osservare $\omega$ in conseguenza di $\alpha$ nel tempo $t$ è una distribuzione Gamma con i seguenti parametri di shape e scale (scale è definita come 1/rate) e il seguente $\mathrm{A} \Omega \mathrm{R}$

AΩR and Gamma parameters
AΩR shape scale
value 6.00% 6.223439 3.462716

Notiamo la sottostima nell'ultima parte della serie temporale dovuta agli eventi $\omega$ non ancora occorsi e la "caduta a zero" dei nuovi eventi $\alpha$ stimati in $t=T$: per l'ultimo giorno della serie temporale non si ha alcun dato futuro da cui stimare i probabili eventi $\alpha$.

Correzione temporale

Nota la distribuzione di probabilità del verificarsi dell'evento $\omega$ in seguito ad $\alpha$ nel tempo $t$ e indicato con $t=T$ l'ultimo giorno di osservazione, la probabilità $p(\alpha|\omega)_{t_0}$ di aver osservato tutti gli eventi $\omega$ causati dagli eventi $\alpha$ nel tempo $[t_0...T]$ è esprimibile come la densità di probabilità cumulativa di $\omega$ al tempo $t=T-t_0$

$$ p(\alpha|\omega)_{t_0} = \mathbf{CDF}[\omega]_{T-t_0} $$

possiamo perciò correggere il numero di eventi $n(\alpha_{t_0})$ stimati dividendo per la probabilità $p(\alpha|\omega)_{t_0}$

$$ n(\alpha_{t_0}) = \sum_{\tau=1}^{T} \frac{ n(\omega_{t_0+\tau}) }{ \mathrm{A}\Omega\mathrm{R} } \frac{ \mathbf{PDF}[\omega]_{t_0+\tau} }{ \mathbf{CDF}[\omega]_{T-t_0} } $$

Notiamo che nonostante la correzione abbia evitato la "caduta a zero" dei nuovi eventi $\alpha$ e la sottostima, ora siamo di fronte ad una sovrastima: la correzione effettuata non è in grado di "prevedere" i futuri eventi $\omega$ ma solamente di stimare la probabilità di aver osservato in $t_0$ tutti gli eventi $\omega$ causati da $\alpha$.

Aggiustamento

Possiamo perciò ottenere un'aggiustamento $n(\alpha_{t_0})_{\mathrm{adj}}$ della sovrastima dovuta alla correzione, prendendo la media degli stimati grezzi $n(\alpha_{t_0})_{\mathrm{raw}}$ e corretti $n(\alpha_{t_0})_{\mathrm{cor}}$

$$ n(\alpha_{t_0})_{\mathrm{raw}} = \sum_{\tau=1}^{T} \frac{ n(\omega_{t_0+\tau}) }{ \mathrm{A}\Omega\mathrm{R} } \frac{ \mathbf{PDF}[\omega]_{t_0+\tau} }{ \mathbf{CDF}[\omega]_{T-t_0} } $$$$ n(\alpha_{t_0})_{\mathrm{cor}} = \sum_{\tau=1}^{T} \frac{ n(\omega_{t_0+\tau}) }{ \mathrm{A}\Omega\mathrm{R} } \mathbf{PDF}[\omega]_{t_0+\tau} $$$$ n(\alpha_{t_0})_{\mathrm{adj}} = \frac{ n(\alpha_{t_0})_{\mathrm{raw}} + n(\alpha_{t_0})_{\mathrm{cor}} }{2} $$

COVID-19 in Italia

Se definiamo

  • $\alpha$: infezione da SARS-CoV-2
  • $\omega$: decesso per COVID-19
  • $\mathrm{A}\Omega\mathrm{R}$: Infection Fatality Rate $\mathrm{IFR}$ per COVID-19
  • $\mathbf{PDF}[\omega]$: distribuzione di probabilità del tempo dall'infezione all'eventuale decesso per COVID-19

data l'osservazione dei decessi giornalieri per COVID-19 possiamo stimare le infezioni occorse precedentemente.

I dati dei decessi e dei casi osservati dal repository GitHub COVID-19 del Dipartimento di Protezione Civile.

L'IFR stimato è tratto da Imperial College Report 34 Brazeau, Verity et al. 2020 (tabella pag 11, High Income Country).

I parametri per la distribuzione dell'intervallo dall'infezione al decesso sono tratti da Phipps et al. 2020 e Flaxman et al. 2020.

Sappiamo dunque che:

IFR values (Imperial College)
IFR
2.5% 0.78%
mean 1.29%
97.5% 1.79%
infection-to-onset & onset-to-death values
infection-to-onset onset-to-death
2.5% 4.10 12.80
mean 5.55 16.00
97.5% 7.00 19.20

Da questi parametri e conoscendo il coefficiente di variazione cov delle distribuzioni Gamma degli intervalli del tempo infection-to-onset $io$ (periodo d'incubazione) e del tempo onset-to-death $od$ (dai sintomi al decesso), possiamo facilmente calcolare i parametri di shape e scale della distribuzione Gamma dell'intervallo infection-to-death $id$ (dall'infezione al decesso).

Sappiamo infatti che

$$ \mathbf{cov} = \frac{\sigma}{\mu} $$$$ shape = \frac{\mu^2}{\sigma^2} = \frac{1}{\mathbf{cov}^2} $$$$ scale = \mu \cdot \mathbf{cov}^2 $$

Da Phipps et al. 2020 e Flaxman et al. 2020 sappiamo che i coefficienti di variazione stimati per COVID-19 dell'intervallo infection-to-onset (periodo d'incubazione) $\mathbf{cov}_{io} = 0.86$ e dell'intervallo onset-to-death (dai sintomi al decesso) $\mathbf{cov}_{od} = 0.45$ da cui, nota la media $\mu$, possiamo ricavare la deviazione standard $\sigma$ dei due intervalli

$$ \sigma_{io} = \mu_{io} \cdot \mathbf{cov}_{io} $$$$ \sigma_{od} = \mu_{od} \cdot \mathbf{cov}_{od} $$

e le distribuzioni di probabilità dei due intervalli.

Da Stewart et al. 2006 sappiamo che media e deviazione standard della convoluzione di due distribuzioni Gamma possono essere stimate come

$$ \mu_{1,2} = \mu_1 + \mu_2 $$$$ \sigma_{1,2} = \sqrt{\sigma_1^2 + \sigma_2^2} $$

Pertanto la distribuzione $\mathbf{Gamma}_{id}$ dell'intervallo infection-to-death data dalla convoluzione delle distribuzioni $\mathbf{Gamma}_{io}$ infection-to-onset (periodo d'incubazione) e $\mathbf{Gamma}_{od}$ onset-to-death, sarà

$$ \mu_{id} = \mu_{io} + \mu_{od} $$$$ \sigma_{id} = \sqrt{\sigma_{io}^2 + \sigma_{od}^2} $$

da cui possiamo ricavare il coefficiente di variazione $\mathbf{cov}_{id}$ dell'intervallo infection-to-death

$$ \mathbf{cov}_{id} = \frac{\sigma_{id}}{\mu_{id}} $$

e infine i parametri shape e scale dell'intervallo infection-to-death

$$ shape_{id} = \frac{1}{\mathbf{cov}_{id}^2} $$$$ scale_{id} = \mu_{id} \cdot \mathbf{cov}_{id}^2 $$

Pertanto, la distribuzione $\mathbf{PDF}[\omega]$ per l'intervallo infection-to-death risulta

$$ \mathbf{Gamma} \left( \alpha = \frac{1}{\mathbf{cov}_{id}^2} \;,\; \beta = \mu_{id} \cdot \mathbf{cov}_{id}^2 \right) $$

Con i valori sopra definiti, otteniamo i seguenti parametri per la distribuzione Gamma dell'intervallo infection-to-death

infection-to-death
Gamma parameters
shape scale
2.5% 6.261966 2.698833
mean 6.223439 3.462716
97.5% 6.190279 4.232443

Calcolate le tre distribuzioni per l'intervallo infection-to-death

  • CI 2.5%
  • media
  • CI 97.5%

le combineremo con i tre valori di IFR

  • CI 2.5% di IFR
  • media di IFR
  • CI 97.5% di IFR

in modo da ottenere una matrice di nove possibili stime

Backcasting parameters matrix
id 2.5% id mean id 97.5%
IFR 2.5% 1 4 7
IFR mean 2 5 8
IFR 97.% 3 6 9

da cui stimeremo mediana e intervallo di confidenza al 95%.

Dal risultato del backcasting notiamo la differenza tra la sottostima dei valori grezzi, la sovrastima della correzione temporale e l'aggiustamento con metodo EpiData.it.

Confronto con Monte Carlo

Confrontiamo dunque i risultati dell'algoritmo di backcasting veloce implementato da EpiData.it con le simulazioni Monte Carlo dell'implementazione di Phipps et al. 2020.