Modelling COVID-19 outbreak on the Diamond Princess ship using the public surveillance data.

The novel coronavirus disease 2019 (COVID-19) outbreak on the Diamond Princess (DP) ship has caused over 634 cases as of February 20, 2020. We model the transmission process on DP ship as a stochastic branching process, and estimate the reproduction number at the innitial phase of 2.9 (95%CrI: 1.7-7.7). The epidemic doubling time is 3.4 days, and thus timely actions on COVID-19 control were crucial. We estimate the COVID-19 transmissibility reduced 34% after the quarantine program on the DP ship which was implemented on February 5. According to the model simulation, relocating the population at risk may sustainably decrease the epidemic size, postpone the timing of epidemic peak, and thus relieve the tensive demands in the healthcare. The lesson learnt on the ship should be considered in other similar settings.

Introduction

The coronavirus disease 2019 (COVID-19), caused by the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), was first documented in Wuhan, China in the end of 2019, and spread to over 100 foreign countries in a short period of time. The World Health Organization (WHO) declared the outbreak to be a public health emergency of international concern on January 30, 2020 (2019World Health Organization). The COVID-19 outbreak on the Diamond Princess, a British-registered cruise ship that contains 3711 tourists and crew members, caused 634 confirmed cases as of February 20, 2020 (2020Ministry of Health of Japan; Mizumoto & Chowell, 2020; Mizumoto et al., 2020; Rocklov et al., 2020; Zhang et al., 2020).

Epidemiology studies had considered the situation in Diamond Princess (DP) ship as a “natural experiment” of SARS-CoV-2 exposure, and examined the epidemiological characteristics of COVID-19 (Mizumoto et al., 2020; Russell et al., 2020; Yamagishi et al., 2020). Modelling studies also attempted to uncover the transmission dynamics of COVID-19 in closed community (Azimi et al., 2021; Mizumoto & Chowell, 2020; Rockl ö v et al., 2020), as well as the effect of disease control intervention. Although recent modelling study by Emery et al. also modelled the transmission dynamics on the DP ship (Emery et al., 2020), a relatively long mean generation time of COVID-19 (at 9 days or more) was assumed, which appears larger than existing estimates from 4 to 6 days (Ali et al., 2020; Du et al., 2020; Ferretti et al., 2020; Ganyani et al., 2020; He et al., 2020; Nishiura, Linton, & Akhmetzhanov, 2020a, 2020b; Ren et al., 2021; Xu et al., 2020; Zhao, 2020; Zhao et al., 2020a, 2021a). As discussed in previous studies (Wallinga & Lipsitch, 2007; Zhao et al., 2020b), the setting of the generation time (or serial interval) is likely to affect the estimate of the disease transmissibility (Tang et al., 2021).

In this study, we model the COVID-19 outbreak on DP ship, estimate the key epidemiological parameters of this outbreak. Different types of control measures were usually implemented simultaneously, and thus their individual effects were difficult to be disentangled. We explore several hypothetical control measures, and how the effects of these measures could be reflected on epidemic curve.

Methods

We consider the population on the DP ship as a close cohort with N = 3711 individuals. Following (Althaus, 2015; Cori et al., 2013), we model the transmission of COVID-19 as a branching process, which considers the epidemic as a reproduction process across generations. Each individual case was modelled explicated in terms of his time of infection and number of offspring cases. In each generation, for one primary case with symptom onset (or reported for a case without symptoms) on the t-th day, the number of secondary cases generated by this primary case follows a negative binomial (NB) distribution with mean at Reff(t) and dispersion parameter at k (Adam et al., 2020; Riou & Althaus, 2020). The Reff(t) is the effective reproduction number, and Reff(t) = R0[N − C(t)]/N, where C(t) is the cumulative number of cases at the t-th day, and R0 is the basic reproduction number of COVID-19 to be estimated. The time delay between symptom onset dates of a primary case and its associated secondary cases is the serial interval (SI). Hence, given a primary case with symptom onset on the t-th day, the symptom onset time of its associated secondary case is expected to be the summation of t and SI of COVID-19. We model the SI to follow a Gamma distribution with mean at 4.5 days and standard deviation (SD) at 3.1 days (Du et al., 2020; Nishiura et al., 2020a, 2020b; You et al., 2020; Zhao et al., 2020c). All secondary cases in one generation will be considered as primary case in the next consecutive generation, and thus the simulation framework can be repeated iteratively.

As the quarantine program starting on February 5, 2020, we account for the effect of this event by reducing R0 to λ∙R0, where λ (0 ≤ λ ≤ 1) is the change ratio of transmissibility due to quarantine program. We estimate both R0 and λ simultaneously.

We collect the number of confirmed cases time series, including both symptomatic and asymptomatic infections, on the DP ship released in the situation reports of COVID-19 infections in Japan (2020Ministry of Health of Japan). To set up the model, we simulate the courses of the COVID-19 outbreaks stochastically starting from one infectious index case on January 20, 2020. We fix k at 0.43 as estimated in (Adam et al., 2020; Wang et al., 2021; Zhao et al., 2021b), which is considered as the main results. For sensitivity analysis, we also consider k = 1, with which the NB distribution is reduced to geometric distribution, and k = 1000000, with which the NB distribution is reduced to Poisson distribution.

We fitted model generated number of cases to the observed number of cases with Poisson-distributed likelihood framework. We adopted a Bayesian fitting procedure with Metropolis-Hastings Markov chain Monte Carlo (MCMC) algorithm with noninformative uniform distribution ranging from 1 to 9 as prior distribution for R0. The MCMC is conducted with 10 chains and 100000 iterations for each chain, including 40000 as for the burn-in period, to obtain the estimates. The median estimate and 95% credible intervals (95%CrI) of R0 are calculated.

To explore an alternative scenario with a larger SI, we repeat the fitting and estimation procedures above with mean SI at 7.5 days and SD at 3.4 days, which was estimated in (Li et al., 2020). We evaluate the modelling performance with mean SI at 4.5 or 7.5 days by the Akaike information criterion (AIC). Under the latter scenario with mean SI at 7.5 days, we find that the R0 is larger than 5, which is not in line with WHO estimates and unlikely under a ‘quarantine’ setting in DP ship. In particular, the latter scenario yields a higher AIC by 10 units than the former scenario with mean SI at 4.5 days, and thus we argue the former scenario is most likely the reality.

Furthermore, we explore several hypothetical COVID-19 control measures, and how these measures may be reflected by epidemic curve. We consider five ‘what-if’ scenarios, and they were.

scenario (0): the best-fit R0 but λ fixed at 1, i.e., baseline scenario;

scenario (1): the best-fit model, i.e., main results;

scenario (2): based on (1), 1000 susceptibles were removed on February 11, 2020;

scenario (3): based on (1), the R0 was set at 2; and

scenario (4): based on (1), combining scenarios (2) and (3).

Scenario (0) was the worst scenario in terms of disease control measures, where we assume the quarantine program was never places and thus reproduction number was not decreased. Scenario (1) was our main results, where a quarantine program was placed on February 5, 2020 and reproduction number decreased from R0 to λ∙R0. Under scenario (1), both model parameters λ and R0 were set to be the best estimates by fitting to real-world dataset, which were summarized in Table 1. Scenarios (2)–(4) were based on scenario (1), and different control measures that different from the actual situation were considered. We simulated the model under five different scenarios, and estimate the key epidemiological metric that quantify the timing and size of the outbreak. We estimate the cumulative number of cases on February 20, the outbreak final size, the epidemic doubling time, and the peaking time to capture the patterns of outbreak.

Table 1

Summary of the reproduction number (R0) and its change ratio (λ) estimates of the outbreak with dispersion parameter k fixed at 0.43, 1, and 1000000. The scenarios where k was fixed at 0.43 was considered as the main results.

dispersion	reproduction number	change ratio
0.43	2.86 (1.66, 7.68)	0.66 (0.11, 0.96)
1	2.53 (1.86, 6.60)	0.78 (0.12, 0.95)
1000000	2.99 (2.05, 3.67)	0.85 (0.37, 0.92)

The analysis in this study was carried out using R statistical software (version 3.6.3).

Results

With the dispersion parameter k fixed at 0.43, we estimate the R0 to be 2.9 (95%CrI: 1.7–7.7), which is consistent with previous estimates (Riou & Althaus, 2020; You et al., 2020; Zhao et al., 2020d). We find that the quarantine program (since February 5, 2020) in DP ship was likely to reduce [(1 − λ) × 100% = ] 34% of R0 (95%CrI: 4–89), see Fig. 1. The fitting results are shown in Fig. 2(a), which match the observed data well. For the sensitivity analysis, the estimating results of R0 and λ with different dispersion term (k) were largely consistent, see Table 1 and Fig. 1, which implies the virus is hard to go extinction that other two coronavirus. Note that low k implies high chance of superspreading events and high chance of extinction.

Fig. 1

The prior (green) and posterior (blue) distributions of reproduction number R0 (upper panels) and change ratio of reproduction number λ (lower panels) with dispersion parameter k = 0.43 (panel (a)), 1 (panel (b)), and 1000000 (panel (c)). In each panel, the vertical gold bar shows the median of posterior estimate. The top panels were the prior (green) and posterior (blue) distributions of reproduction number (R0) before the quarantine program, which was on February 5, 2020. The bottom panels were the prior (green) and posterior (blue) distributions of the change ratio of reproduction number (λ) since February 5, 2020, when reproduction number decreased from R0 to λ∙R0.

Fig. 2

The observed (dots) and fitted (curves) number of COVID-19 cases on the Diamond Princess ship. Panel (a) shows the daily number of cases, and panel (b) shows the estimated epidemic curves under different scenarios. In panel (a), the black dots are the observed number of cases time series, the blue curve is the simulation median, and the grey curves are 100 simulation samples. In panel (b), the bold curves are the simulation median, and the dashed curves are the 95% centiles. The black curves are for scenario (0), blue curves are for scenario (1), green curves are for scenario (2), gold curves are for scenario (3), and red curves are for scenario (4). The key settings and simulation outcomes of all scenarios were summarized in Table 2. In panels (a) and (b), the blue curves show the baseline scenario, and they are the same.

Table 2

Summary of the key epidemiological estimates of the outbreak under different scenarios. The epidemic curve of each scenario was visualized in Fig. 2.

scenario	description	cum. # on Feb 21	doubling time (day)	peak time	Visualization
(0)	baseline, λ = 1	1340 (190, 2493)	3.4 (1.8, 8.1)	Feb 16 (Feb 7, Mar 2)	black curve in Fig. 2
(1)	main result, R0 = 2.86, λ = 0.66	1118 (39, 2389)	3.4 (1.7, 12.6)	Feb 22 (Feb 7, Mar 18)	blue curve in Fig. 2
(2)	based on (1), remove 1000 susceptibles	1195 (26, 2701)	3.9 (2.0, 9.3)	Feb 18 (Feb 9, Mar 18)	green curve in Fig. 2
(3)	based on (1), R0 reduced by 0.5	562 (26, 2154)	4.1 (1.5, 9.1)	Mar 3 (Feb 14, Apr 11)	gold curve in Fig. 2
(4)	based on (1), implement (2) and (3)	512 (31, 1706)	3.5 (1.7, 9.3)	Mar 1 (Feb 12, Mar 30)	red curve in Fig. 2

Without the quarantine program on February 5, we estimate that the COVID-19 outbreak is likely to cause 1340 COVID-19 cases (95%CI: 190–2493) as of February 21, and the epidemic curve is likely to peak around February 16, 2020 with a doubling time at 3.4 days (95%CI: 1.8–8.1), see Fig. 2(b) and Table 1. Under scenario (1), we find the quarantine program mitigated the epidemic size at 1118 (95%CI: 39–2389) as of February 21, and postponed the timing of epidemic peak on February 22. Under scenario (2), timely reduction in the susceptible population could lower the final size and sustainably reduce the daily incidences. Under scenario (3), the lower the R0, decreasing from 2.9 to 2.0, the lower the epidemic size and the later the peaking time will be. If the susceptible pools and R0 are reduced simultaneously, the COVID-19 outbreak on the ship will be likely mitigated and postponed largely, see Fig. 2(b).

Discussion

The R0 of COVID-19 on the Diamond Princess ship is likely to be lower than that of the severe acute respiratory syndrome (SARS) ranging from 2.2 to 3.6 (Lipsitch et al., 2003), but lower than estimates in previous studies assuming relatively long SI (Emery et al., 2020; Rockl ö v et al., 2020). The situation on the Diamond Princess cruise implies that the virus could spread rapidly, most likely with a short SI than previously estimated (Li et al., 2020). Note that a recent study on the 1099 patients found that the median of incubation period is only 3 days, thus a short effective SI is also possible (Guan et al., 2020). With a shorter SI, a relatively lower R0 could also result in a rapid growth of the epidemic size (Tuite & Fisman, 2020), and with a shorter epidemic doubling time as we show in Table 1. Therefore, timely contact tracing and effectively quarantine were crucial to shutoff transmission chains.

The simulation results under scenarios (1)–(4) indicate that the public health control measures could effectively mitigate the COVID-19 outbreak on the ship in terms of the epidemic size, see Table 2. The estimated 34% reduction in R0 was associated with the quarantine program on February 5, and similar reduction in transmissibility was also reported in (Emery et al., 2020; Rockl ö v et al., 2020). Decreasing the disease transmissibility in terms of R0 could postpone the peak, which may gain valuable time to prepare and allocate the resources in response to incoming patients. Relocating the population at risk (if possible) could sustainably decrease the daily incidences, and thus relieve the tensive demands in the healthcare, and improve in the treatment outcome.

From Japan official reports (https://www.mhlw.go.jp/stf/newpage_09542.html), we obtain the following facts:

“2/15/2020 陽性が確認されたのは、延べ930名の検査中285名(うち無症状病原体保有者延べ73名)となりました。2/20/2020陽性が確認されたのは、延べ3063名の検査中634名(うち無症状病原体保有者延べ328名)となりました。”

The English translation from Google translation:

“February 15, 2020: Positive cases were confirmed in 285 of the 930 people tested (including 73 asymptomatic pathogen carriers). February 20, 2020: A total of 634 people (including 328 asymptomatic pathogen carriers) tested positive for a total of 3063 people.”

From the above official announcement of Japan, we calculate the ratio of asymptomatic cases out of all confirmed cases increased from (73/285 = ) 25.6% to (328/634 = ) 51.7% between February 15 and 20, 2020. It was unlikely that such a large increase in asymptomatic cases were infected before the quarantine starting on February 5, 2020, but unascertained in February 15. As such, the quarantine was less effective to prevent the emergence of asymptomatic cases, while it successfully reduced the occurrence of symptomatic cases, which is also suggested by studies only including symptomatic cases (Mizumoto & Chowell, 2020). Recent modelling study also inferred that there were 74% of infections proceeded asymptomatically, and the prevalence of asymptomatic infections appears higher after February 5 than before (Emery et al., 2020). They also estimate that asymptomatic individuals were the source for 69% of all infections. In other words, 74% of all infections are asymptomatic and may generate 69% offsprings, and thus asymptomatic or symptomatic infections made roughly equal contributions to the COVID-19 transmission.

In this study, we avoid separating the asymptomatic and symptomatic cases, and we model all laboratory confirmed infections regardless their symptoms. We capture the overall transmissibility of COVID-19 on DP ship, and evaluate the effect of the quarantine program in terms of the reduction in R0. One of the contributions is to explore several ‘what-if’ scenarios, under which the epidemic outcomes were assessed and projected, see Table 2. We suggest maintaining the enhancement in both populational level public health control as well as the individual level self-protection actions in combating the COVID-19 outbreak. With detailed information on public intervention, the analytical framework in this study can be extended to a complex context and used for evaluating the effects of certain control measures.

Ethics approval and consent to participate

The data were collected via public domain (2020Ministry of Health of Japan), and thus neither ethical approval nor individual consent was not applicable.

Availability of materials

All data used in this work were publicly available via (2020Ministry of Health of Japan).

Consent for publication

Not applicable.

Funding

DH was supported by General Research Fund (Grant Number 15205119) of the Research Grants Council (RGC) of Hong Kong, China and an Alibaba-Hong Kong Polytechnic University Collaborative Research project. WW was supported by (Grant Number 61672013) and Huaian Key Laboratory for Infectious Diseases Control and Prevention (Grant Number HAP201704), Huaian, Jiangsu, China.

Disclaimer

The funding agencies had no role in the design and conduct of the study; collection, management, analysis, and interpretation of the data; preparation, review, or approval of the manuscript; or decision to submit the manuscript for publication.

Authors’ contributions

SZ and DH conceived the study, carried out the analysis, and drafted the first manuscript. All authors discussed the results, critically read and revised the manuscript, and gave final approval for publication.

Declaration of competing interest

This study was partially supported by a grant from the Research Grants Council of the Hong Kong Special Administrative Region, China (HKU C7123-20G).