Open Science

Mathematical modelling to inform New Zealand’s COVID-19 response

Published:22 February 2021

ABSTRACT

Between February and May 2020, New Zealand recorded 1504 cases of COVID-19 before eliminating community transmission of the virus in June 2020. During this period, a series of control measures were used including population-wide interventions implemented via a four-level alert system, border restrictions, and a test, trace, and isolate system. Mathematical modelling played a key role in informing the government response and guiding policy development. In this paper, we describe the development of a stochastic mathematical model for the transmission and control of COVID-19 in New Zealand. This includes features such as superspreading, case under-ascertainment, testing and reporting delays, and population-wide and case-targeted control measures. We show how the model was calibrated to New Zealand and international data. We describe how the model was used to compare the effects of various interventions in reducing spread of the virus and to estimate the probability of elimination. We conclude with a discussion of the policy-modelling interface and preparedness for future epidemic outbreaks.

Introduction

The COVID-19 outbreak originated in Wuhan China in November 2019 (WHO 2020a) before spreading globally to become a pandemic in March 2020 (WHO 2020b). The human population lacks immunity to COVID-19, a viral zoonotic disease with reported case fatality rates that are of the order 1% (Verity et al. 2020). Many countries have experienced community transmission after undetected introductions of the disease by travellers exposed in China. This led to rapid growth of new infections in many countries, even as China, through the use of strong controls and rapid testing and tracing, managed to contain the spread of the virus.

As of early July 2020, more than 11 million confirmed cases of COVID-19 have been recorded globally and over 500,000 have died. Strict lockdowns have been implemented by many countries to varying degrees and at different times. Several countries have experienced a resurgence in case numbers after seemingly managing to contain the virus. For example, Singapore and later Australia have experienced outbreaks associated with rapid transmission within high-density housing. The USA and the UK have experienced significant national or regional outbreaks after relaxing lockdown measures.

In New Zealand, border restrictions and strong interventions to maintain physical distancing were implemented early in the outbreak. Population-wide restrictions were encoded via an alert system ranging from Alert Level 1 through to Alert Level 4, which involved strict, legally enforceable stay-at-home orders, and closure of all schools and non-essential businesses. The Alert Level was progressively raised to Level 2 on 21 March, Level 3 on 23 March, and Level 4 on 26 March, at which time there were 205 confirmed and probable cases. Together with widespread testing, case isolation and contact tracing, these measures were successful in containing the outbreak and eventually eliminating community transmission.

Restrictions were relaxed in stages starting in late April and ending on 9 June with the move to Alert Level 1 effectively removing all restrictions on domestic activities, movements and gatherings. Elimination of community cases has so far been maintained via strict border measures, namely closure of the border to non-residents, and mandatory testing and 14-day quarantine in government-managed facilities for all arrivals.

New Zealand’s first COVID-19 outbreak

Between 26 February and 22 May 2020, New Zealand recorded 1504 confirmed and probable cases of COVID-19 (). Of these, 573 had a recent history of overseas travel and were considered imported cases, 454 were import-related cases, 389 were domestically transmitted cases epidemiologically linked to another confirmed or probable case, and 88 were domestically transmitted from an unknown source. The last domestically transmitted case was reported on 22 May and the last domestically transmitted case that could not be traced to a confirmed or probable contact was reported on 30 April. Between 22 May and 11 August, there were no reported cases of COVID-19, outside of government-managed quarantine for international arrivals. On 11 August 2020, a second outbreak was detected in Auckland which is ongoing at the time of writing.

Figure 1. COVID-19 cases reported in New Zealand between 26 February and 22 May 2020, showing cases with a recent international travel history (blue) and without a recent international travel history (orange) along with changes in alert level and border restrictions.

In total, there were 22 deaths from COVID-19 in the first outbreak, giving a crude overall case fatality rate (CFR) of 1.46%. The ages of the people who died from COVID-19 ranged from 62 to 99 years. Of the 1504 cases, 95 (6.3%) required hospitalisation and 10 (0.66%) were admitted to ICU. These rates may have been affected by the relatively young, fit demographic of the internationally imported cases and their close contacts: only 7.7% of all cases were aged 70+, compared to 10.9% of the population as a whole (see ). Children are also under-represented in the case data: only 5.5% of cases were aged under 15, compared to 19.4% of the population. For the 89 hospitalised cases with admission and discharge dates recorded, the mean length of hospital stay was 7.9 days (25th percentile 1 day, 75th percentile 10 days). Out of all cases, 32 confirmed cases (2.1%) were recorded as asymptomatic. This is likely to be an underestimate because testing for COVID-19 during this outbreak was predominantly symptom-based and there was no routine testing of asymptomatic individuals.

Figure 2. Age distribution of COVID-19 cases in New Zealand between 26 February and 22 May 2020 in 5-year age bands, split into cases with (blue) and without (orange) a recent international travel history. Black curve shows the number of cases there would be in each 5-year age band if the 1504 cases were distributed in proportion to the New Zealand population.

Up to 9 June 2020 (when New Zealand moved to Alert Level 1 and testing of international arrivals was significantly ramped up), New Zealand conducted 298,532 tests for SARS-CoV-2 (approximately 60 tests per 1000 people) with a positivity rate of 0.4% (Ministry of Health 2020a). The mean time between symptom onset (according to case recall) and a positive test result being returned was 7.3 days ( $N = 1108$ , $S D = 5.9$ days).

Contact tracers in regional public health units interviewed confirmed and probable cases about their contacts during their infectious period, which is usually considered to be up to 2 days prior to symptom onset. Close contacts are household members and those with at least 15 min of either face-to-face contact with a case or being within 2 m of a case in a closed environment; see Ministry of Health (2020b) for a complete definition. A total of 7683 contacts of the 1504 confirmed and probable cases up to 22 May were identified (an average of 5.1 contacts per case). Of the 7683 identified contacts, 6055 were followed up (average 4.0 per case).

Despite the small sample size, the crude observed CFR of 1.46% is consistent with international data allowing for demography and case under-ascertainment (Verity et al. 2020). Since the rates of clinically severe infection and death are higher in older age groups (CDC 2020; Williamson et al. 2020), it is possible that the hospitalisation rate and the CFR would have been higher if community transmission had become more widespread in New Zealand.

New Zealand’s largest case clusters were associated with a wedding (Bluff, 98 cases), a school (Marist College Auckland, 96 cases), a hospitality venue (Matamata, 77 cases), aged residential care facilities (Christchurch, 56 cases; Auckland, 51 cases), a private function (Auckland, 40 cases), and a conference (Queenstown, 39 cases). The settings for these clusters are typical of international patterns of COVID-19 transmission and superspreading (Adam et al. 2020; Leclerc et al. 2020). They also exhibit one or more of the ‘three Cs’ (closed spaces, crowded places and close contact) suggested to be associated with superspreading events (Nishiura et al. 2020). It is notable that in the cluster associated with a school, the majority of cases (55 out of 96) were in individuals aged over 20 and, of the 41 cases in under 20s, only 6 were in under 10s. This is consistent with international evidence on the reduced incidence and transmission of COVID-19 among young children (Boast et al. 2020; Li et al. 2020).

Mathematical model for COVID-19 transmission and control

We developed a stochastic, continuous-time branching process model for the spread of COVID-19 (). Stochastic models are a class of mathematical model that includes a random element and are defined in terms of the probabilities of certain events occurring. Stochastic models of disease spread have several advantages over their deterministic counterparts:

Deterministic models break down when the number of cases is relatively small and so cannot be used to look at questions relating to the elimination of the virus.
Stochastic models intrinsically allow for random variations in the transmission process, for example superspreading events or variations in the timing of secondary infections, symptom onset, testing and isolation. This enables uncertainty in model outputs to be quantified and the probability of elimination to be estimated.
Stochastic individual-based models track individual infected cases, so they are more compatible with data on the number of cases and individual attributes of those cases such as age, time of symptom onset, hospitalisation status.
Stochastic individual-based models allow the structure of the transmission tree () to be explicitly tracked. This is essential for models of contact tracing, which require information on who infected whom.

Figure 3. Stochastic, continuous-time branching process model of COVID-19 transmission starting from single infected seed case (top centre). Each infected individual (blue circles) infects a random number of individuals generated from the offspring distribution. The diagram shows an example of a ‘superspreading event’, where a single individual infects a large number of other individuals at one time.

We developed a stochastic model for COVID-19 in New Zealand that used data on cases with a recent history of international travel as seed cases. We did not attempt to predict cases arriving into New Zealand, only the subsequent chains of transmission originating from each imported case. Many imported cases did not infect any other people in New Zealand, or only infected one or two others, but some imported cases initiated large transmission trees. The following is a brief description of the model and its assumptions; see Appendix for mathematical details.

New infections occur as a result of contact between infectious and susceptible individuals and the population is assumed to be homogeneous. We assume that 33% of infections are subclinical (Byambasuren et al. 2020; Lavezzo et al. 2020; Pollán et al. 2020) and, relative to clinical cases, these subclinicals are assumed to result in 50% as many secondary infections on average (Davies et al. 2020). To model random variation among individuals in the number of secondary infections (e.g. superspreading), each case is assigned a randomly generated individual reproduction number $R_{i}$ from a gamma distribution with mean $R_{c l i n}$ for clinical cases or $R_{s u b} = 0.5 R_{c l i n}$ for subclinical infections, and dispersion parameter $k = 0.5$ (Lloyd-Smith et al. 2005; Endo et al. 2020). In the absence of control measures, the number of secondary infections from individual $i$ is a Poisson random variable with mean $R_{i}$ .

The generation times (i.e. times between infection of an individual with COVID-19 and onward transmission to secondary cases) are independent, identically distributed random variables from a Weibull distribution with mean $μ = 5.05$ days and standard deviation $σ = 1.94$ days. This distribution is based on the estimated times between infection of 40 source-recipient case pairs (Ferretti et al. 2020). The shape of this distribution means that cases are most likely to transmit the virus to close contacts around 5 or 6 days after infection, and less likely to transmit in the first few days after infection or more than around 12 days after infection (). The incubation period (time from infection to onset of symptoms) is approximately 5 days on average (Lauer et al. 2020), so this corresponds to substantial levels of pre-symptomatic transmission (Ganyani et al. 2020; Tindale et al. 2020).

Case-targeted control was modelled by assuming that clinical cases self-isolate a random time after infection that is the sum of the time from infection to onset, modelled as a gamma distributed random variable ( $μ = 5.5$ days, $σ = 2.3$ days, Lauer et al. 2020), and the time from onset to isolation, modelled as an exponentially distributed random variable ( $μ = 2.2$ days, based on New Zealand case data). Self-isolation reduces onward transmission by 35% (, green curve). After self-isolation, there is a further delay until a positive test result is returned, which is exponentially distributed ( $μ = 3.48$ days). We assumed subclinical infections did not self-isolate and were not tested or reported. Population-wide social distancing measures were modelled by reducing the contact rate between infectious and susceptible individuals by a factor of $C (t)$ , with smaller values of $C (t)$ for higher COVID-19 Alert Levels.

Values of model parameters were estimated from published studies or from New Zealand data (see Appendix ). Using these parameter values, the model’s effective reproduction number (i.e. the average number of people infected by a single case) in the absence of population-wide control ( $C (t) = 1$ ) is $R_{e f f} = 2.4$ ; with population-wide control, $R_{e f f} = 2.4 C (t)$ .

Effective reproduction number

The effective reproduction number $R_{e f f}$ , is a key measure of a disease’s transmission potential and the most important parameter in any epidemic model. If $R_{e f f} > 1$ , the virus has the potential to cause a large outbreak, while if $R_{e f f} < 1$ , the virus will eventually die out with probability 1. The basic reproduction number $R_{0}$ for COVID-19, for a fully susceptible population in the absence of any control measures depends on both the biological characteristics of the virus and the contact rates of the population in which it is spreading, which vary between and within countries. Reported estimates range from $R_{0}$ [95% CI] = 0.48 [0.26, 0.88] to 6.91 [6.60, 7.23], with a global average of 3.17 [2.62, 3.84] (Billah et al. 2020). However, interventions including school closures, population-wide social distancing, and case-targeted controls, reduce $R_{e f f}$ to varying degrees.

Directly estimating $R_{e f f}$ empirically in real time is difficult because it requires information on the number of secondary infections from each infected case, which is usually not available. The delay between infection and symptom onset, and between onset and reporting can be over 2 weeks, meaning that the signal from changes in Alert Level in the official case data is heavily lagged. The size of the reduction in $R_{e f f}$ during New Zealand’s lockdown, in particular, was difficult to predict in advance or estimate in real-time. Indications of the reduction in activity could be obtained from real-time electronically collected data, such as telecommunications data, transit telemetry data, point-of-sale bank transaction data, Facebook co-location data, and Apple/Google mobility data. However, none of these directly measures the outcome of interest, which is close physical contact between pairs of individuals. We therefore used international data on case numbers to estimate the effect of various population-wide interventions on $R_{e f f}$ .

Estimating the effective reproduction number from international data

Numerous methods have been developed for estimating the reproduction number using mathematical models (Obadia et al. 2012). We used a method by Wallinga and Lipsitch (2006) to estimate $R_{e f f}$ from data on the number of new daily cases and number of daily deaths (Dong et al. 2020) in 25 countries or states/provinces. This method calculates $R_{e f f}$ by relating it to the exponential growth rate r of daily new local cases (which we inferred using log-linear regression of daily new cases over a 10-day window) by $R_{e f f} = 1 / \int_{a = 0}^{\infty} e x p (- r a) w (a)$ where $w (a)$ is the assumed generation time distribution shown in . It assumes a randomly mixed population undergoing exponential growth of new local infections arising from domestic transmission – see Binny et al. (2020a) for details.

We broadly grouped interventions in each location into four levels of increasing stringency, approximately equivalent to New Zealand’s Alert Level framework. We then estimated the reduction in $R_{e f f}$ before and after various interventions were implemented (). In the early phase before interventions were implemented, $R_{e f f}$ was greater than 1 (indicating a growing outbreak) in all 25 countries, ranging from 1.5 to 5.4. This range corresponds to doubling times of 8.7 days to 1.8 days. In the late phase after interventions were implemented, $R_{e f f}$ decreased in 24 out of 25 locations, the exception being Singapore, where $R_{e f f}$ declined over the first 30 days of the outbreak before returning to values above 1 during a second outbreak. This was a result of the virus spreading in a subgroup of the population with a substantially higher value of $R_{e f f}$ , in this case migrant workers. This factor is likely important in the dynamics of COVID-19 spread in other countries, including New Zealand, and this illustrates one of the limitations of homogeneous population models.

Figure 5. Interventions reduced Reff in 24 of the 25 countries or states/provinces investigated. Greatest reductions are observed in locations with longer durations of sustained Alert Level 3–4 interventions. (a) Reff estimates for the early phase of an outbreak before interventions were implemented (red) and late phase after interventions were implemented (blue) (each Reff estimate is obtained by fitting to daily new case counts in a 10-day window and plotted on the last day of the 10-day window). Threshold for outbreak Reff=1 indicated by black dashed line. Note, countries (or states/provinces) experienced their first confirmed case of COVID-19 infection on different dates; units on horizontal axes are number of days since the 100th case of confirmed infection in an outbreak (by which time an outbreak is typically established and driven mainly by community transmission). (b) Reduction in Reff (difference between early-phase and late-phase Reff estimates for each location) against the length of time (days) spent under the restrictions, corresponding to Alert Level 1 (green diamonds), Level 2 (light blue triangles), Level 3 (dark blue squares) or Level 4 (purple circles), that were in place in that location during the late phase.

The effectiveness of interventions varied considerably between countries, with estimates for $R_{e f f}$ after interventions ranging from 0.3–2.1. Strong interventions (comparable to Alert Levels 3-4) reduced $R_{e f f}$ below the threshold for outbreak ( $R_{e f f} < 1$ ) in 16 locations. At the time of analysis, these countries had spent varying durations under restrictions and, in some locations, these may not have had sufficient time to achieve their full potential effect. For instance, the delays between an individual being infected with the virus to the onset of symptoms, and to being confirmed by testing, mean the full extent of the reduction in $R_{e f f}$ after interventions will not be observable in case data until around 2 weeks after interventions are begun. The largest reductions in $R_{e f f}$ were achieved in locations (e.g. New York State and Quebec) that had sustained Alert Level 3–4 restrictions for periods of more than 25 days.

Our results are in general agreement to other studies on $R_{e f f}$ for COVID-19 (e.g. Alimohamadi et al. 2020; Cowling et al. 2020; Flaxman et al. 2020; Price et al. 2020). We obtained similar late-phase $R_{e f f}$ estimates but higher early-phase estimates than a global study by Abbott et al. (2020), which used an approach based on Thompson et al. (2019). Our estimates of $R_{e f f}$ from the exponential growth rate are based on an assumed generation time distribution, which has mean 5.05 days. Other groups have estimated longer mean generation times (Flaxman et al. 2020; Li et al. 2020). This would result in our method returning higher estimates for $R_{e f f}$ in the early phase and potentially lower estimates for $R_{e f f}$ in the late phase.

Estimating the effective reproduction number for New Zealand

The method we used to estimate $R_{e f f}$ from international data is not suitable in the very early stages of an outbreak, when daily numbers of new cases are low and exponential growth is not fully established (i.e. growth in case numbers is initially erratic). In addition, the method cannot distinguish between domestically transmitted cases and imported cases, which accounted for a significant proportion of all cases in New Zealand’s first outbreak. Applying the method to both domestic and imported cases would yield artificially inflated estimates for $R_{e f f}$ because the model assumes that all new cases on a given day arose by transmission from existing active cases. Other methods that do account for imported cases have been developed, but these also rely on sufficiently high case numbers (Thompson et al. 2019).

For these reasons, we did not apply these approaches to New Zealand case data. Instead, we fitted the stochastic model described above to the number of reported new daily cases and found the value of the relative transmission rate $C (t)$ , and hence the reproduction number $R_{e f f}$ , at different Alert Levels that produced the best match. We did this by simulating the model for a range of values of the relative transmission rates $C (t)$ : (i) for the period prior to Alert Level 4 (26 February to 25 March 2020); and (ii) for the period spent in Alert Level 4 (26 March to 27 April 2020). We compared the simulated number of reported cases per day to actual reported cases per day and calculated the best-fit values of $C (t)$ by minimising the mean-square error of square root-transformed data, over a time window from 10 March to 27 April and averaged over 1000 identically initialised realisations of the stochastic model. Bootstrap confidence intervals (95%) were obtained by re-estimating $C (t)$ as described above (and corresponding $R_{e f f}$ ) for each of 10,000 simulations from the best-fit model. Using this approach, we estimated $R_{e f f} = 1.8 [1.44, 1.94]$ prior to moving to Alert Level 4 and $R_{e f f} = 0.35 [0.28, 0.44]$ during Alert Level 4. These estimates provided a good match between model simulations and data (). These results were relatively insensitive to changes in other model parameters.

Figure 6. Simulated (red curves) and actual (black points) numbers of new daily reported cases in New Zealand, along with simulated number of new daily infections (blue curves). Note that the move to Alert Level 4 leads to an immediate drop in the number of simulated new infections, but there is a lag before the number of reported cases begins to decrease. Light curves show individual realisations of the stochastic model; bold curves show the average of all realisations. Effective reproduction number Reff=1.8 at Alert Level 2, Reff=0.95at Level 3, and Reff=0.35 at Level 4.

Our results suggest that New Zealand’s Alert Level 4 was highly effective at reducing $R_{e f f}$ below one, meaning that new infections would continue to decline as long as the Alert Level remained in place. The Alert Level 4 estimate $R_{e f f} = 0.35$ is comparable to the smallest reproduction numbers we estimated for other countries where comparable interventions were particularly effective (). For example, lockdowns in six Australian states successfully reduced $R_{e f f}$ to approximately 0.3-0.5 (Price et al. 2020). Other countries (e.g. UK, USA, Sweden) remained at $R_{e f f} \approx 1$ for prolonged periods (Binny et al. 2020a), meaning that new infections remained steady or declined only slowly. Our estimate for the reproduction number before Alert Level 4, $R_{e f f} = 1.8$ , is relatively low compared to those reported globally. This could be due, in part, to measures put in place in early- to mid-March, for example the cancellation of mass gatherings, employees being encouraged to work from home, self-isolation of international arrivals, and tracing and testing of close contacts of cases associated with international travel.

The estimated 81% reduction in $R_{e f f}$ as a result of the lockdown is in good agreement with other studies, which employed a range of modelling approaches to assess effects of physical distancing on $R_{e f f}$ . For example, Jarvis et al. (2020) reported a 74% reduction in $R_{e f f}$ (from a mean of 2.6–0.62), as a result of physical distancing during lockdown in the UK. Their approach used survey data to compare age-specific population contact rates before and during lockdown. Another study by Flaxman et al. (2020) fitted a hierarchical Bayesian model to data from 11 European countries and estimated an 81% reduction in $R_{e f f}$ due to lockdowns. A third study used a Bayesian SEIR-type model to infer a 78% reduction in contacts in British Columbia due to the introduction of physical distancing measures (Anderson et al. 2020).

Probability of elimination

Various definitions of elimination of an infectious disease have been suggested (e.g. Dowdle 1998; Heymann 2006; Wilson, Parry, et al. 2020). A reasonable definition is that there are no active cases that could contribute to future community transmission (i.e. excluding cases that are no longer infectious or are in strictly managed isolation, for example border quarantine facilities). In a global context, this is actually ‘local’ elimination because COVID-19 transmission is still occurring in many other countries. This is impossible to guarantee empirically, because there is always a risk of active cases that are asymptomatic or otherwise undetected. However, it is straightforward to test in a model and stochastic models are ideally suited to this purpose. The probability of ultimate extinction of a branching process can be inferred directly from the distribution of the number of secondary cases when the latter is fixed and does not change over time (e.g. Lloyd-Smith et al. 2005). However, we wanted to investigate scenarios where population-wide control measures change the relative transmission rates over time, and in which estimates of the probability of elimination could be updated using real-time observations. We therefore took a simulation-based approach as follows.

We estimated the probability that COVID-19 would eventually be eliminated in New Zealand, conditional on observing a specified number of consecutive days with no new reported cases (Binny et al. 2020b). Such estimates are important to inform decisions on timings for the easing of restrictions. We simulated the branching process model up to 31 July 2020, using estimated values for $R_{e f f}$ of 0.35 at Alert Level 4 (26 March to 27 April 2020), 0.95 at Alert Level 3 (28 April to 13 May 2020), and 1.78 at Alert Level 2 (14 May 2020 onwards). We ran 1000 realisations of the model and defined elimination as having no remaining cases within 30 days of infection at the end of the simulation (i.e. no new infections after 1 July). Results are not highly sensitive to the choice of end date which was chosen for convenience. To obtain the probability of elimination conditional on having $n$ consecutive days with no new reported cases, we restricted the sample to realisations with $n$ consecutive days with no new reported cases and calculated the proportion of these realisations that resulted in elimination. We used a longer delay from isolation to reporting (mean 6 days) and, under an optimistic scenario with high detection and reporting of clinical cases (75%), we found that reaching 95% probability of elimination required 10 consecutive days with no reported cases (). For a pessimistic scenario with low detection of clinical cases (20%) and lower effectiveness of Alert Levels 2 and 3 ( $R_{e f f}$ of 2.3 and 2.2 respectively), 95% probability of elimination required 22 days with no reported cases (). These results are similar to those of Wilson, Parry, et al. (2020) who used a different model.

Figure 7. Probability of eventual elimination as a function of the observed number of consecutive days with no new reported domestic cases. This is defined as the proportion of model realisations with n consecutive days with no new reported cases that result in elimination. Optimistic scenario with high detection and reporting of clinical cases (75%) and moderate effectiveness of Alert Level 2 (Reff=1.8) and Level 3 (Reff=0.95). Pessimistic scenario with low detection and reporting of clinical cases (20%) and low effectiveness of Alert Level 2 (Reff=2.3) and Level 3 (Reff=2.2). Alert Level 4 Reff=0.35.

The last reported case of COVID-19 in New Zealand in the first outbreak was 22 May. This means that New Zealand reached a 95% probability of elimination on 1 June 2020 under the optimistic scenario, or 13 June 2020 under the pessimistic scenario. The relaxation to Alert Level 1 occurred on 9 June 2020.

Estimated inequities by ethnicity and healthcare access

Previous pandemics have shown that Māori and Pacific people are at greater risk of negative outcomes. During the 1918 influenza pandemic, Māori experienced death rates seven times higher than New Zealand European/Pākehā (Pool 1973; Summers et al. 2018). In the H1N1 influenza pandemic in 2009, Māori were infected at twice the rate of Pākehā, and with increased severity (Wilson et al. 2012). International literature suggests that COVID-19 follows similar patterns. Black, Asian and Minority Ethnic (BAME) groups in the United Kingdom have experienced higher mortality rates, even after controlling for comorbidities (Williamson et al. 2020). In the United States, Pacific populations are experiencing greater rates of infection and hospitalisation (Jackson 2020).

In a New Zealand context, Te Tiriti o Waitangi requires the Crown to ensure that Māori have equitable access to and outcomes from healthcare (Waitangi Tribunal 2019). However, as Talamaivao et al. (2020) have shown in their systematic review of quantitative studies ‘experience of racial discrimination is an important determinant of health in New Zealand, as it is internationally’.

Risk of outbreaks in disadvantaged communities

Several countries initially managed to control the spread of COVID-19 successfully, but experienced major outbreaks when the virus reached disadvantaged groups, including ethnic minorities. For example, Singapore appeared to have contained the spread of COVID-19 until a rapidly growing outbreak was detected in migrant worker populations in predominantly substandard housing. An antibody prevalence study in England found infection rates were higher in BAME people, groups with high deprivation score, and for those living in large households (Ward et al. 2020). Lack of access to primary care and testing facilities, inequities and racism within the healthcare system, and crowded living conditions are likely factors contributing to the rapid spread within disadvantaged communities.

We developed a structured model of COVID-19 transmission to describe groups with different levels of access to healthcare and testing facilities (James, Plank, Binny, et al. 2020). We modelled a reduction in transmission corresponding to Alert Levels 3 and 4, but assumed that transmission rates dropped less in the group with low access to healthcare. This reflects factors such as substandard housing, financial stress, lack of paid sick leave, insecure and high-contact employment, and inability to work from home. For example, Pacific people in Auckland experience a high level of segregation by ethnicity (Manley et al. 2015; Salesa 2017), higher rates of overcrowding (Schluter et al. 2007) and higher rates of unmet healthcare needs (Ryan et al. 2019). These factors mean that there is a higher risk of a major outbreak in this group. The modelling results showed that if contact rates between groups with different levels of access to healthcare are relatively low, the outbreak could grow very large before detection (James, Plank, Binny, et al. 2020). This implies that these inequities in access to healthcare and structural factors affecting ability to undertake physical distancing need to be urgently addressed (McLeod et al. 2020).

Relative fatality rates by ethnicity

Early in the pandemic, we sought to quantify the potential risk of fatality from COVID-19 infection faced by Māori and Pacific people in New Zealand (Steyn, Binny, Hannah, et al. 2020). We used international data on COVID-19 infection fatality rates by age (Verity et al. 2020) alongside life expectancy statistics (StatsNZ 2020) to estimate age-specific infection fatality rates for each ethnicity. This reflects the significantly increased infection fatality rates experienced by older populations. We rescaled these age-structured fatality rates to reflect differing levels of unmet healthcare need and the prevalence of various comorbidities, weighted by their prevalence in each age-ethnicity group.

We estimated that the risk of fatality for Māori is around 50% greater than the risk for the New Zealand European and other group. Sensitivity analysis shows this may be as high as 300% greater if the impact of greater unmet healthcare needs and lower life expectancy have a larger effect than is reflected in officially recorded data. These findings are consistent with evidence from other countries. A large study in the United Kingdom (Williamson et al. 2020) found that Black people, even after controlling for age, comorbidities, healthcare and deprivation, were 48% more likely to die from COVID-19 infection.

These results only consider the risk of fatality following infection. As described above, other factors may increase the rate of spread in disadvantaged communities, where Māori and Pacific people are typically overrepresented. These communities therefore face elevated levels of risk on two fronts: higher risk of infection and higher risk of clinically severe outcomes and fatality following infection. Appropriately engaging with affected communities to understand and respond to the inequities in transmission and clinical outcomes of COVID-19 is critical for any government strategy or plan to be successful.

Policy-modelling interface

In New Zealand, there were relatively few modelling tools available early in the outbreak to support government decision-making. The University of Otago’s School of Public Health made use of a deterministic compartment model, developed with collaborators in Germany (Wilson, Telfar Barnard, et al. 2020) to investigate long-term control strategies. The Ministry of Health commissioned modelling from this group to support the Cabinet decisions that led the country into Alert Level 4 in late March 2020. However, this type of deterministic model is not well suited to supporting the types of decisions the government faced in April and May 2020, as case numbers came under control and the country progressed towards elimination. In mid-March, our team began work on the stochastic model described above, and from early April, this model was used extensively to advise government, both on operational requirements and on policy settings.

On 27 March 2020, a ‘COVID-19 Modelling Workstream’ was established within the National Crisis Management Centre (NCMC), with one of us (Shaun Hendy) reporting directly to this workstream. It was subsequently agreed to divide the workstream into two parts. The first would focus on operational advice, supported by a team at Orion Healthcare, including forecasting potential case numbers and patient loads on hospitals and intensive care units. The second would investigate possible scenarios that might result from future policy decisions.

The operational workstream produced daily reports, distributed to the NCMC, the Department of Prime Minister and Cabinet (DPMC), the Ministry of Health, and District Health Board (DHB) planners based on model simulations run using the latest case numbers, clinical data, and parameter estimates. These operational reports had to make assumptions concerning future policy settings (e.g. the dates on which lockdowns might be eased), known as ‘working scenarios’ that were developed in discussions with NCMC and Ministry of Health. The New Zealand Treasury also implemented a version of the model to use in its forecasting of healthcare demand.

The policy workstream provided advice on policy settings through a process we called the ‘scenario sandpit’. As decision points approached or particular issues arose, policy-makers would request simulations to explore the outcome of their decisions. For instance, in the lead up to the 20 April 2020 Cabinet decision on the date for exiting Alert Level 4, DPMC requested simulations of ‘optimistic’ and ‘pessimistic’ scenarios for exits on 23, 27 April, and 7 May. The 23 April and 7 May scenarios were both explicitly included in the 20 April Cabinet Paper (DPMC 2020a). In the end, Cabinet split the difference by deciding to set the exit date as 27 April. The three scenarios provided to DPMC were developed via an iterative process that involved intense engagement and feedback between the modelling team, departmental science advisors, and policy analysts. Similar scenarios were also included in the 11 May and 8 June Cabinet Papers (DPMC 2020b; DPMC 2020c) reviewing Alert Levels 3 and 2 respectively.

The branching process model that formed the basis of policy advice was subject to several rounds of rapid peer review. During April and May, the model was updated weekly after review by several panels and advisory groups. Te Pūnaha Matatini established a technical peer review panel, chaired by Dr Matthew Parry, and also made use of reviews by the Ministry of Health’s Epidemiology Technical Advisory Group, chaired by Associate Professor Patricia Priest. This provided a degree of quality assurance (although not to the standard that would normally be expect from blind peer review), as well as important feedback and guidance on further model developments. Once reviewed and revised, descriptions of the models were released as working papers via Te Pūnaha Matatini’s website. Considerable efforts were also made to communicate results publicly via mainstream media and social media. Several data visualisations were produced in partnership with the New Zealand Herald and an interactive version of the model with a graphical user interface was made available to the public in late April (‘COVID-19 Take Control Simulator’). This layer of quality assurance and transparency allowed for informed public debate and enabled policy-makers to make decisions based on independent, publicly available science.

Discussion

Mathematical modelling has played a key role in informing New Zealand’s response to the COVID-19 pandemic. This included the decision to ‘go hard, go early’ and issue stay-at-home orders from 25 March 2020, as well as the timing of relaxation of alert levels (DPMC 2020a, 2020b, 2020c). Public communication of the science behind the response, in the mainstream media, social media and via interactive web apps, helped build public understanding and trust in government decision-making.

In this paper, we have focused on the first COVID-19 outbreak in New Zealand between February and May 2020. Mathematical modelling has continued to play a role in preparedness and decision-making subsequently. This includes modelling the effectiveness of manual and digital contact tracing systems (Plank et al. 2020; James, Plank, Hendy, et al. 2020), and the risks from international arrivals and government-managed quarantine facilities (Steyn, Binny, Hendy et al. 2020). The models were also used to provide advice about the likely number of infected cases and the probability of regional spread at the time of detection of the second outbreak on 11 August 2020.

The capacity of modelling tools and the data pipelines on which they rely has been built up over a very short space of time. Ongoing investment in mathematical modelling expertise and the structures, systems, and relationships needed to provide model-based policy and operational advice will help increase New Zealand’s capacity to manage to infectious diseases and increase preparedness for future pandemics.

Acknowledgements

The authors acknowledge the support of StatsNZ, ESR, and the New Zealand Ministry of Health in supplying data in support of this work. The authors are grateful to Samik Datta, Nigel French, Markus Luczak-Roesch, Anja Mizdrak, Fraser Morgan, Matt Parry, Tom Cunningham, Ann Brower, and Mick Roberts for discussions on model development and informal peer review comments on draft manuscripts. The authors also thank the Ministry of Health’s Technical Advisory Group Epidemiology Subcommittee chaired by Patricia Priest for their comments and suggestions. The authors are grateful to Pieta Brown, Kevin Ross and the team at Orion Health for their role in establishing the COVID-19 Modelling Update and running model simulations for the daily report.

Disclosure statement

No potential conflict of interest was reported by the author(s).

This work was supported by Ministry of Business, Innovation and Employment [grant number UOAX1941]; Te Pūnaha Matatini, New Zealand’s Centre of Research Excellence for complex systems.