First outbreak of Zika virus in the continental United States: a modelling analysis.

Since 2015, Zika virus (ZIKV) has spread throughout Latin and Central America. This emerging infectious disease has been causing considerable public health concern because of severe neurological complications, especially in newborns after congenital infections. In July 2016, the first outbreak in the continental United States was identified in the Wynwood neighbourhood of Miami-Dade County, Florida. In this work, we investigated transmission dynamics using a mathematical model calibrated to observed data on mosquito abundance and symptomatic human infections. We found that, although ZIKV transmission was detected in July 2016, the first importation may have occurred between March and mid-April. The estimated highest value for R0 was 2.73 (95% confidence interval (CI): 1.65-4.17); the attack rate was 14% (95% CI: 5.6-27.4%), with 15 (95% CI: 6-29) pregnant women involved and a 12% probability of infected blood donations. Vector control avoided 60% of potential infections. According to our results, it is likely that further ZIKV outbreaks identified in other areas of Miami-Dade County were seeded by commuters to Wynwood rather than by additional importation from international travellers. Our study can help prepare future outbreak-related interventions in European areas where competent mosquitoes for ZIKV transmission are already established. This article is copyright of The Authors, 2017.

Introduction

Zika virus (ZIKV) has recently emerged as a significant threat to public health worldwide. Transmission of ZIKV to humans is thought to occur mainly via bites of Aedes aegypti mosquitoes; however, other mosquito species, such as Ae. albopictus, have been demonstrated as potential vectors [1,2], and other routes of infection are possible, including sexual transmission [3], blood transfusion [4] and vertical transmission in both humans [5] and mosquitoes [6]. Although infection is usually asymptomatic or mild [7], a causal link with congenital birth defects has been established, and a strong association exists with Guillain-Barré syndrome [8], a severe neurological condition. The low frequency of severe complications has been counterbalanced by a rapid geographic spread of the virus, which led the World Health Organization to declare a public health emergency of international concern (PHEIC) from 1 February to 18 November 2016. ZIKV was originally identified in Africa in 1947 [9], but the first large outbreak was not reported before 2007 in Micronesia [9], followed by further outbreaks in 2013 and 2014 in French Polynesia [10] and other Pacific Islands [7]. After the first cases were notified in Brazil in March 2015, ZIKV spread throughout South and Central America, with the notable exception of Uruguay and continental Chile, within 18 months [11,43].

In Europe, locally transmitted cases of ZIKV have not occurred to date. Ae. aegypti is currently present only on Madeira (Portugal) and around the Black Sea [12]; the much more widespread Ae. albopictus has a similar competence for ZIKV at 27°C, but none at 18 °C [2]. For these reasons, the risk of local transmission in temperate climate regions has been estimated to be minor [13,14], but the possibility of outbreaks transmitted by Ae. albopictus cannot be ruled out in warmer areas, allowing both greater vector densities and higher competence. Indeed, 11% of the European population is estimated to live in areas suitable for a ZIKV epidemic, mostly concentrated in the Mediterranean basin [14].

In the continental United States (US), the first outbreak of ZIKV was recorded at the end of July 2016 in the Wynwood neighbourhood of Miami-Dade County, Florida [15]. In mid-August, mid-September and mid-October, three additional outbreaks in other areas of the same county (South Miami Beach, North Miami Beach and Little River) were identified. On 9 December 2016, the Florida Department of Health declared to have cleared active ZIKV transmission from all identified areas, which counted ca 250 locally acquired infections overall [16], and the area has been Zika-free since then.

In this work, we analysed the transmission dynamics of the Wynwood outbreak using a mathematical model calibrated to outbreak data, and we assessed the efficacy of the implemented vector control measures in containing viral transmission. Results from this analysis provide useful insights for prevention and control of possible future outbreaks in European areas.

Methods

The outbreak under study involved an area of approximately 2.6 km2 [15], with a population of 7,725 inhabitants [17]. Health authorities identified 21 locally transmitted symptomatic cases, with onset of illness ranging from 26 June to 5 August 2016.

We modelled the mosquito abundance over time M(t) by considering a fixed mortality rate μ and a time-dependent recruitment rate ψ(t) (defined as the average number of adult female mosquitoes produced by a single female adult) previously estimated on Ae. aegypti populations from Miami [18]. A density D of ca 1,560 female mosquitoes per hectare was estimated on  = 26 July, before vector control treatments [15], in the following equation, using parameter values reported in Table 1: D = K/aπr.

Table 1

Model parameters for Zika virus outbreak, Wynwood, 2016

Parameter	Interpretation	Value	Source
K	Female mosquitoes per trap per day before interventions	30	[15]
a	Capture rate (%/day)	2.46	[40]
r	Flight range for Aedes aegypti (in m)	50	[41]
ρ	Reduction in mosquito abundance following treatments (%)	75	[15]
pMH	Probability of transmission from mosquito to human per bite	0.214	[1]
pHM	Probability of transmission from human to mosquito per bite	0.767	[1]
θM	Mosquito incubation period (days)	Gamma distribution(μ = 10.5, σ = 0.5)	[19]
LH	Latency before symptom development (days)	Gamma distribution(μ = 5.9, σ = 0.5)	[19]
αH	Symptomatic period (days)	Normal distribution(μ = 4.5, σ = 0.75)	[19]
ti	Duration of asymptomatic infection (days)	Uniform(range 7–13)	[42]
q	Probability of being asymptomatic	0.8	[4,7]
μM	Ae. aegypti death rate (1/days)	0.1	[18]
ψ(t)	Ae. aegypti recruitment rate (1/days)	μM (1+0.25 cos(2π/365 (t-90.89))	[18]
d	Reporting probability for symptomatic individuals	0.1	[19]
T0	Date of index case importation	18 March(95% CI: 1 March–15 April)	Calibrated
b	Mosquito biting rate (1/days)	0.058(95% CI: 0.055–0.061)	Calibrated

CI: confidence interval.

We modelled the mosquito abundance over time M(t) by considering a fixed mortality rate μM, a time-dependent recruitment rate ψ(t) previously estimated on Ae. aegypti populations from Miami [18] and their density D. Specifically, we applied the following system to compute the mosquito abundance on each day t:

We modelled the interventions by imposing a sudden decrease of the mosquito population size by ρ = 75% after 6 August, as suggested by the entomological data presented in [15] (Figure 1, panel A).

Figure 1

Modelled mosquito abundance and transmission model for Zika virus outbreak, Wynwood, 2016

A. Predicted mosquito abundance over time in the study area. Black: 75% reduction in mosquito abundance following 6 August due to vector control interventions; orange: no vector control implemented.

B. Model flow chart for ZIKV transmission in humans (blue) and mosquitoes (yellow). Compartments:
A: permanently asymptomatic humans; I: asymptomatic humans who will develop symptoms; C: symptomatic humans; Me: exposed mosquitoes; Mi: infectious mosquitoes; Ms: susceptible mosquitoes; R: humans who recovered after symptoms; S: susceptible humans. Parameters: λH and λM are the force of infection for humans and mosquitoes respectively and are computed as λH = b pMH Mi/H and λM = b pHM Hi/H, with H being the total human population in Wynwood; Ψ(t) is the mosquito recruitment rate, given by Ψ(t) = ψ(t) (Ms + Me + Mi). See Table 1 for parameter values and references.

Transmission dynamics of ZIKV in humans and mosquitoes were modelled according to the compartmental scheme reported in Figure 1, panel B. Susceptible humans (S) contract ZIKV from bites of infectious mosquitoes. A fraction q of infected individuals will remain asymptomatic (A) for the entire duration of viraemic infection t, after which they will recover (R); the remaining fraction (1 − q) will remain temporarily asymptomatic (I) before developing symptoms (C) and will then recover (R). Recovered individuals are no longer infectious and become immune to reinfections. We assumed that infected humans (H) are infectious regardless of apparent clinical symptoms. Symptomatic individuals have an overall probability d of being detected at some point. Susceptible mosquitoes (M) can become infected (M) after biting infectious humans; in such cases, they will become infectious to humans (M) after an extrinsic incubation period and for the rest of their life. Given the very small rates of vertical transmission [6], we assumed that all mosquitoes newly emerged from their breeding site are susceptible.

To take into account the uncertainty surrounding the epidemiological parameters, we ran a stochastic version of the model proposed in Figure 1, panel B. We implemented an agent-based model representing the 7,725 inhabitants of the Wynwood neighbourhood, following an approach similar to the one presented in [19]. On each day, susceptible humans move to the infectious compartments based on a Poisson sample with a rate equal to the force of infection λH; infected humans are subject to a binomial probability q of remaining asymptomatic. Asymptomatic individuals will remain infectious for a time t sampled from a uniform distribution with a range of 7–13 days; symptomatic individuals will show their symptoms after a time LH drawn from a gamma distribution with a mean of 5.9 days and shape parameter 0.5. Symptoms will last until recovery for a time αH, sampled from a normal distribution with a mean of 4.5 days and standard deviation 0.75. Vector dynamics were implemented as a stochastic stage-structured model representing the three possible epidemiological stages of the vector: susceptible, exposed and infected. See Table 1 for details on model parameters and the corresponding references for the adopted probability distributions.

We assumed that the Wynwood outbreak was initiated by a single index case imported in the area via international travel at time T. Two free model parameters (T and the mosquito biting rate b) were estimated by fitting the model-predicted weekly number of reported symptomatic cases to the observed cases. We used a Markov chain Monte Carlo (MCMC) approach with a standard Metropolis–Hastings algorithm with 100,000 iterations, uniform priors (range for T0: 1 January–20 June; range for b: 0–1) and a Poisson likelihood. We used the parameter values accepted by the algorithm as the posterior distributions, which were used to compute model predictions.

We computed the initial reproduction number R0(t) from the following formula [20]:

,

using the posterior distribution of the mosquito biting rate b, the modelled mosquito abundance M(t) and the epidemiological parameter values for ZIKV (Table 1). represents a measure of the invasion potential of an infectious host introduced in a fully susceptible host population on day t. We also computed the effective reproduction number R(t) as the product of R(t) and the fraction of available susceptible individuals in the population [21]. R(t)represents the residual transmission potential of an epidemic as the reservoir of susceptible individuals is depleted by the process of infection and recovery.

We also estimated the number of infected pregnant women and blood donors based on the model-predicted prevalence of human infections over time. Data on fertility [17] and blood donation [22] rates by age and ethnicity were adapted to the demographic structure of Miami-Dade County [17]. We considered all blood donations occurring between T0 and 27 July, when health authorities suspended blood collection from the area [23].

We evaluated the robustness of the proposed model (‘baseline model’) by fitting alternative model structures after varying some of the key assumptions: (i) We considered a broader mosquito flight range r of 82.5 m [24] in the estimation of the mosquito abundance (model M1); (ii) we assumed a relative reduction in mosquito abundance due to the vector control treatments (ρ) as a free model parameter sampled with a uniform prior with a range of 0–100% (model M2); (iii) we let the mosquito biting rate b vary over time, according to the same temporal dependence assumed for the recruitment rate [18] (model M3); (iv) we evaluated a model where the reporting rate d was a free model parameter with a uniform prior of 0–100% (model M4). Models were ranked according to the deviance information criterion (DIC) [25].

Results

The baseline model was able to reproduce the number of cases by date of symptom onset (Figure 2, panel A), as reported by health authorities [15]. The arrival of the index case was estimated between early March and mid-April (average: 18 March; 95% confidence interval (CI): 1 March–15 April), approximately corresponding to the dates when R0(t) approached the epidemic threshold (average: 28 March; 95% CI: 8 March–22 April). We estimated R0 to vary from a minimum of 0.4 (95% CI: 0.24–0.61) during winter to a maximum of 2.73 (95% CI: 1.65–4.17) in mid-June (Figure 2, panel B).

Figure 2

Reported cases, R0, and viral prevalence in humans as predicted by the model for Zika virus outbreak, Wynwood, 2016

A. Predicted weekly number of reported symptomatic cases by date of symptom onset. Red dots represent observations [15] whereas boxplots indicate quantiles of the corresponding model-predicted distribution at 2.5%, 25%, 50% 75% and 97.5%.

B. Model estimates of the initial reproductive number for ZIKV in Wynwood over time.

C. Model-predicted prevalence of ZIKV in mosquitoes.

D. Model-predicted prevalence of ZIKV in human hosts.

Solid lines: mean values; dashed lines: 95% confidence interval; black: 75% reduction in mosquito abundance following 6 August due to vector control interventions; orange: no vector control implemented.

The total estimated number of human infections was 1,112 (95% CI: 436–2,120), corresponding to an average attack rate of 14.4% (95% CI: 5.6–27.4%). We estimated 15 (95% CI: 6–29) pregnant women to have been infected during the outbreak and a probability of 12% that at least one infectious individual donated blood before 27 July. The low occurrence of viraemic blood donors is due to smaller donation rates among the Hispanic population [22], the main ethnicity in Miami-Dade [17], and to the low (< 4%) predicted ZIKV prevalence until July.

According to data from the American Community Survey [26], incoming commuters represent ca 44% of the Miami-Dade population; of this proportion, the large majority (87%) are within-county commuters, while 13% come from outside the county (mostly from the neighbouring Broward County). This would correspond to ca 5.7% of infections (13% of 44%), i.e. ca 1 in 18, in commuters from outside the county and indeed, one of the 21 symptomatic patients reported in the outbreak was a commuter from Broward County identified via workplace investigation [15]. About 453 infections occurred in commuters resident in other parts of Miami-Dade; for comparison, only 208 travel-related ZIKV cases were recorded in Miami-Dade as of 19 September [27].

Vector control strategies were crucial for containing virus circulation by rapidly abating the average value of Re from 1.49 to 0.37, with a sharp reduction of infection prevalence in both vectors and humans (Figure 2, panels C and D). In the absence of interventions, the model suggested that Re would have remained above the critical threshold until the end of August, resulting in 55 (95% CI: 29.5–81.7) reported symptomatic infections and a total attack rate of 33% (95% CI: 16.7–49.3). Thus, according to our model, the implemented interventions have prevented ca 60% of all potential cases. Table 2 shows a sensitivity analysis of the expected total number of reported symptomatic cases, attack rate, and fraction of prevented infections for alternative scenarios on the assumed reduction in mosquito abundance allowed by the interventions. We estimated that even a treatment with low effectiveness, corresponding to a reduction of 25% of the vector population, would significantly reduce the expected number of infections by almost 40%.

Table 2

Effectiveness of alternative scenarios of vector control interventions on Zika virus infections, Wynwood, 2016

Reduction of mosquito abundance (6 August)	Total number of reported symptomatic cases		Total attack rate		Fraction of prevented infections compared with no intervention
	n	95% CI	%	95%CI	%
0% (no intervention)	55	29.5–81.7	33.0	16.7–49.3	0
25%	33.7	17.1–57.8	21.4	10.9–37.4	35.1
50%	26.5	11.4–48.8	17.0	6.2–31.7	51.5
75%	22.6	9.2–42.2	14.4	5.6–27.4	56.4
90%	21.3	9.7–35.8	11.7	5.4–20.3	64.5
75%a	18.1	9.2–33.2	11.5	5.9–20.8	65.1

CI: confidence interval.

The assumed effectiveness of implemented measures is shown in bold.

a Reduction of mosquito abundance occurring on the same day of outbreak detection (23 July).

Table 3 shows a summary of results obtained with alternative assumptions on the model structure and parameter values. While the baseline model had the best performance in terms of DIC score, the overall qualitative conclusions were robust with respect to these variations, with alternative models suggesting an even earlier introduction of the virus and a slightly higher attack rate. The assumptions of ρ = 75% and d = 10% for the baseline model M0 were compliant with the estimates computed for models M2 and M4, respectively, where these parameters were estimated through the MCMC procedure.

Table 3

Comparison between alternative models, Zika virus outbreak, Wynwood, 2016

Model	DIC	Average T0		Average b		Average    attack rate		Average d		Average ρ
		Date	95% CI	Days− 1	95% CI	%	95% CI	%	95% CI	%	95% CI
M0: baseline model	36.55	18 Mar	1 Mar–15 Apr	0.05	0.055–0.061	14.4	5.6–27.4	NA	NA	NA	NA
M1: model with increased flight range	36.84	20 Feb	23 Jan–24 Mar	0.093	0.088–0.100	17.4	4.1–36.6	NA	NA	NA	NA
M2: model with free vector control efficacy	36.92	16 Feb	21 Jan–21 Mar	0.058	0.054–0.061	15.3	4.2–28.0	NA	NA	85.4	60.1–99.9
M3: model with time-dependent biting rate	37.20	13 Feb	1 Feb–21 Mar	0.052	0.049–0.055	16.9	4.2–35.7	NA	NA	NA	NA
M4: model with free reporting rate	39.15	10 Mar	10 Feb–3 Apr	0.057	0.054–0.060	15.0	4.2–27.3	12.8	10.0–19.7	NA	NA

CI: confidence interval; DIC: deviance information criterion; NA: not applicable.

Discussion

In this work, we estimated that ZIKV was introduced to Wynwood between 3 and 5 months before the recognition of the outbreak by health authorities at the end of July. This long delay is justified because the infection prevalence in humans was small until June and because ZIKV symptomaticity [28] and reporting rates were low [19]. Our conclusion is supported by molecular-clock analyses of 32 ZIKV isolates from Florida, reported by the collaborative project nextstrain [29,30], which suggest that importation occurred between the end of February and the end of March. Model simulations indicate a likely extinction of transmission during winter even in the absence of interventions; however, we propose that undetected infections have occurred even after authorities declared the end of the outbreak on 19 September [27], given the sustained presence of infected mosquitoes. This prediction is consistent with the identification of locally transmitted cases in Miami-Dade as late as 28 December, almost three weeks after all active outbreaks in the area had been declared over [31]. Finally, our results suggest that the other three outbreaks in neighbouring areas of Miami-Dade County were probably seeded from Wynwood commuters rather than initiated from travel-related cases; a phylogenetic analysis of ZIKV isolates [29] supports this finding and shows that at most two viral clades were circulating in the four Florida outbreaks in 2016.

Our estimate of the mosquito biting rate is in line with a recent study that found parity and blood meal frequency to increase with temperature [32]. In particular, 60% of female mosquitoes were parous at 24 °C and took blood meals at an average interval of 11.7 days, whereas a parity of 86% and an interval between blood meals of 9.8 days were found at 27 °C. These figures correspond to a rate between 0.051 and 0.088 bites per mosquito per day. We estimated an average biting rate of 0.058 per mosquito per day throughout the period March to September, when the average temperature in Miami is 26.1 °C [33]. Because data on the biting rate at lower temperatures was not available, we did not include its temperature dependence as a possible driver of temporal variations in the transmission risk; however, the model assuming a seasonal dependency in the biting rate resulted in a similar quality of fit, an earlier estimated date of importation and a slightly higher attack rate. Other model parameters may be influenced by seasonal variations in temperature, with possibly shorter incubation periods (as in the case of dengue [34]) and increased transmission rates [2] during warmer months. Data for natural history parameters of ZIKV are insufficient to factor the temperature dependence in our model. Nonetheless, by analogy to what we found for the time-dependent biting rate, we expect that a further concentration of transmission in summer months would increase the estimated attack rate and push further back the date of introduction.

We estimated a peak value for R0 of ca 2.7, well within the range of other estimates for outbreaks in the Pacific and the Americas [14,19,35]. Previously published predictions for European areas with endemic Ae. albopictus populations also suggest potential values for R0 below 3 in the large majority of sites [14]. Therefore, the Wynwood outbreak is a relevant case study for potential future ZIKV transmission in Europe and can provide useful insights for prevention and control.

In Wynwood, about two weeks elapsed between recognition of local transmission on 23 July and a massive reduction in captured mosquitoes on 6 August [15]. These timely measures were successful in containing the attack rate below 15%. Even if effective interventions were immediately deployed on the day when local transmission was detected, the attack rate would still be above 11%. Because ZIKV is transmitted in a population for a long time before detection, reactive control measures are insufficient to prevent a large number of infections. A more effective approach would require the application of preventive measures with the aim of keeping R0 below the epidemic threshold. A recent study conducted on European Ae. albopictus showed that integrated vector control strategies can halve the mosquito abundance (and therefore R0) compared with sites where no intervention is implemented [36]. This reduction would not be sufficient to completely eliminate the risk of local transmission in the most exposed areas, but it would greatly limit the time window over which an epidemic is possible and, most importantly, its potential size, thereby reducing the risk of seeding further outbreaks in neighbouring areas. Furthermore, preventive control of mosquito populations would simultaneously reduce the risks of other mosquito-borne infections such as chikungunya and dengue [37].

Our conclusions are subject to some unknowns on ZIKV epidemiology, such as the role of asymptomatic infections, sexual transmission and spatial dynamics, along with uncertainties in parameter values and possible drifts in vector competence following adaptations of the viral genome [38]. The implications of this study for control are expected to be robust with respect to these uncertainties, since they depend on the silent transmission of ZIKV in the early months after importation, which is now a well-established trait of this emerging infection [39].