11  Imported cases

Code 
# Setting up
using Distributions, Plots, Measures
include("../src/loadData.jl")

We suggested in section ?sec-smc-data that separating imported from local cases may be critical when modelling the example data (the first 100 days of the COVID-19 pandemic in New Zealand). Separating the graph of reported cases into local and imported cases demonstrates why:

Code 
Y = loadData("NZCOVID")
bar(Y.date, Y.border + Y.local, label="Local", color=:darkorange)
bar!(Y.date, Y.border, label="Imported", color=:darkblue)
plot!(xlabel="Date", ylabel="Reported cases", size=(800,300), margins=3mm)
Figure 11.1: Reported cases of COVID-19 in the first 100 days of the pandemic in New Zealand, separated by imported cases (dark blue) and local cases (orange).

The problem of modelling imported cases has previously been covered by Thompson et al. (2019), …, and …

This chapter demonstrates one way in which local and imported cases can be distinguished in the sequential hidden-state framework.

11.1 Without quarantine

We retain the hidden-state model from Chapter 8:

\[ \log R_t \sim \text{Normal}(\log R_{t-1}, \sigma) \tag{11.1}\]

and now assume that only local cases \(L_t\) are infected by past local and imported \(M_t\) cases:

\[ L_t | R_t \sim \text{Poisson}\left(R_t \Lambda_t^{(m)}\right) \tag{11.2}\]

where

\[ \Lambda_t^{(m)} = \sum_{u=1}^{u_{max}} \omega_u \left(L_{t-u} + M_{t-u}\right) \tag{11.3}\]

The bootstrap filter for this model is nearly identical to the simple model, we simply change one line (#TODO: install highlight extension and highlight):

Code 
function importedModel(σ, Y::DataFrame, opts::Dict)

    # Extract frequently used options
    T = opts["T"]
    N = opts["N"]
    L = opts["L"]

    # Initialise output matrices
    R = zeros(N, T) # Using R instead of X to highlight we're estimating Rt
    W = zeros(N, T)

    # Sample from initial distribution
    R[:,1] = rand.(opts["pR0"], N)

    # Run the filter
    for tt = 2:T

        # Project according to the state-space model
        R[:,tt] = exp.(rand.(Normal.(log.(R[:,tt-1]), σ)))

        # Weight according to the observation model
        Λ = sum(Y.Ct[tt-1:-1:1] .* ω[1:tt-1])
        W[:,tt] = pdf.(Poisson.(R[:,tt] .* Λ), Y.local[tt]) # <- This line is the only line that has changed!

        # Resample
        inds = wsample(1:N, W[:,tt], N; replace=true)
        R[:, max(tt - L, 1):tt] = R[inds, max(tt - L, 1):tt]

    end

    return(R, W)

end

Fitting the model and plotting \(R_t\) against our original estmiates reveals substantial differences:

[Single figure to go here]

11.2 With quarantine

The model above assumed that imported cases are just as infectious as local cases. In

11.3 Concluding remarks