Assessing the potential risk of Zika virus epidemics in temperate areas with established Aedes albopictus populations

G Guzzetta 1 , P Poletti 1 2 , F Montarsi 3 , F Baldacchino 4 , G Capelli 3 , A Rizzoli 4 , R Rosà 4 , S Merler 1 1. 1. Fondazione Bruno Kessler, Trento, Italy 2. Dondena Centre for Research on Social Dynamics and Public Policy, Bocconi University, Milan, Italy 3. Istituto Zooprofilattico Sperimentale delle Venezie, Padova, Italy 4. Department of Biodiversity and Molecular Ecology, Research and Innovation Centre, Fondazione Edmund Mach, San Michele all’Adige (Trento), Italy

Based on 2015 abundance of Aedes albopictus in nine northern Italian municipalities with temperate continental/oceanic climate, we estimated the basic reproductive number R 0 for Zika virus (ZIKV) to be systematically below the epidemic threshold in most scenarios. Results were sensitive to the value of the probability of mosquito infection after biting a viraemic host. Therefore, further studies are required to improve models and predictions, namely evaluating vector competence and potential non-vector transmissions.
In 2015, the largest recorded epidemic of Zika virus (ZIKV) started in Brazil and has since then expanded progressively to most countries in Central and South America [1]. We provide estimates of the basic reproduction number (R 0 ) of ZIKV in northern Italy, based on estimates of the mosquito abundance from entomological surveillance data.

Entomological surveillance in northern Italy
Mosquito monitoring was carried out fortnightly from May to October 2015 in the provinces of Belluno and Trento, Italy ( Figure 1).
Aedes albopictus mosquitoes were collected using 54 Biogents Sentinel traps (Biogents AG, Regensburg, Germany, hereafter abbreviated as BG) baited with BG lures and CO 2 from dry ice, running for 24 hours and placed by entomologists at selected locations in nine municipalities ( Figure 1) at altitudes ranging from 74 m a.s.l. to 650 m a.s.l and geographical coordinates between 10°49'04.9"E and 12°12'54.2"E longitude and 45°53'26.9"N and 46°09'59.4"N latitude. Temperatures at trap locations were obtained from land surface temperature satellite data with a resolution of 250 m [2] (Figure 2).

Mosquito population dynamics
We developed a population model representing the developmental cycle of mosquitoes by means of temperature-dependent parameters ( Figure 3) and fitted it to capture data in order to estimate the density of female adult mosquitoes per hectare over time at each municipality. For two towns (Belluno and Feltre), human landing captures were carried out (seven and five sessions, respectively) where BG traps were positioned. Two experts performed the catches, rotating between the two sites, acting as human baits and collectors. The mosquitoes were collected by a handheld aspirator during the three hours preceding sunset. Human landing data were used for independent validation of the local mosquito abundance predicted by the model.
The four main stages of the Aedes albopictus life cycle (eggs, larvae, pupae and adults) are modelled. Biological parameters encoding mortalities, developmental rates and the length of the gonotrophic cycle depend on the average daily temperature recorded at the site of capture, according to equations provided in [6] and based on experimental data [20]. The sitespecific density-dependent factors and the capture rate (common to all sites) are free model parameters estimated by fitting model outputs to experimental capture data.
Given the model-predicted daily number of mosquitoes N V and the number of bites per mosquito per day k, the following relation should hold: where HLR is the hourly human landing rate estimated from data and T is the average duration of biting activity during a day (set to 12 hours, based on several studies on daily landing patterns, e.g. [3]).

Basic reproduction number of Zika virus
We assumed that the only route of transmission for ZIKV is via mosquito bites. R 0 can be calculated from densities of human and mosquito populations and several epidemiological parameters according to the following Formula [4]: Symbols, interpretations, values and literature references are reported in the Table. When R 0 < 1 (epidemic threshold), the probability of observing sustained transmission of ZIKV after importation of a case is negligible. When R 0 > 1, the outbreak probability is given by the following Formula [5]: where and Using baseline parameter values (Table), the expected value of R 0 stayed far below the epidemic threshold of 1 at all sites and times in our simulations ( Figure 6A), resulting in a low risk of autochthonous transmission of ZIKV.
We re-computed the values of R 0 under a range of worst-case scenarios for parameter values and model assumptions. In all scenarios, all epidemiological parameters but one were fixed at their baseline values and sensitivity was assessed against variations of the selected parameter. Firstly, we set the mosquito biting rate (k) to the largest estimate for the 2007 Italian chikungunya virus outbreak (k = 0.16 days −1 [6]). In this scenario, the peak value of R 0 never exceeded 0.8. Secondly, we assumed daily temperatures in the upcoming mosquito season to be 2 °C higher than those recorded in 2015 (an extreme scenario in climatological terms) under baseline parameter values. This resulted in an increase of the peak mosquito abundance of 17% to 95%, depending on the town; however, even in this case, R 0 remained far from the epidemic threshold (peak values below 0.4 at all sites). Thirdly, R 0 remained below 1 even when considering 100% human susceptibility to infection given a bite from an infected mosquito (p H ) [7,8].

Symbol Description
Site-specific density-dependent factor The four main stages of the Aedes albopictus life cycle (eggs, larvae, pupae and adults) are modelled. Biological parameters encoding mortalities, developmental rates and the length of the gonotrophic cycle depend on the average daily temperature recorded at the site of capture, according to equations provided in [6] and based on experimental data [20]. The site-specific density-dependent factors and the capture rate (common to all sites) are free model parameters estimated by fitting model outputs to experimental capture data.

Figure 4
Comparison between observed and predicted numbers of mosquitoes captured over time at all sites, Italy, 2015

Observed Predicted
Model predictions shown as average and 95% confidence intervals over 10,000 stochastic simulations Finally, we considered the variability of R 0 with respect to the probability of a mosquito being infected upon biting of a viraemic human host, p V . The very low baseline value (6.7%) was suggested by a recent experimental study [7], but previous work had estimated a value of 100% [8]. Resulting predictions were very sensitive to the value of this parameter. When we used the value provided by the latter study (p V = 100%), the peak value of R 0 exceeded the epidemic threshold in seven of nine towns, with values as high as 3.8 in the highly mosquito-infested towns of Feltre and Riva del Garda ( Figure 6B). In Strigno and Belluno, R 0 remained systematically below the epidemic threshold because of the low ratio of mosquitoes per human (V/H in equations above). In all other towns, the minimum value of p V required to have R 0 above 1 ranged from 25% to 50%.
We call epidemic season the time of the year when the local mosquito abundance is sufficiently high for R 0 to exceed the epidemic threshold. According to model estimates, the epidemic season in the worst-case scenario of p V = 100% was predicted to last between two and three months in the seven towns at higher risk ( Figure 6C). In this scenario, for every ZIKV case imported within the epidemic season, the average probability of observing an outbreak of local transmission ranged from 18% in Tenno to 39% in Feltre and Riva del Garda ( Figure 6D).

Discussion
Although ZIKV infection in humans is generally asymptomatic or very mild, there is growing evidence of association with Guillain-Barré' syndrome [9] and congenital neuronal defects in newborns [10,11]. The flow of international travellers to and from Latin America raises potential concerns for the occurrence of outbreaks also in Europe during the summer months when the mosquito activity is higher [12]. The Latin-American epidemic is likely to be driven by Ae. aegypti, a mosquito species that is currently present in Europe only in Madeira (Portugal) and around the Black Sea [13]. However, in many European countries Ae. albopictus is now endemic [13]. This species has been demonstrated to transmit ZIKV both in the laboratory [7,8] and in the  wild [14], although its estimated transmission efficiency is much lower than measured for Ae. aegypti [7]. Based on data-driven estimates of the abundance of Ae. albopictus mosquitoes in nine municipalities of northern Italy, we expect a low risk of autochthonous mosquito-borne transmission of ZIKV.  ZIKV genotype involved in the large South American outbreak [7], and one using Ae. albopictus from a tropical rainforest climate (Singapore) infected with an Ugandan ZIKV genotype [8]. We chose the first study as a baseline, because the viral strain used in the experimental setting was more relevant for the ongoing epidemic, and the second study to set the worst-case scenario (p V = 100%).
Additional sources of uncertainty come from sexual transmission of ZIKV [16]; much higher viral loads have been found in semen than in blood [17]; however, the relative contribution of non-vector transmission is currently not quantified. To improve models and predictions, further studies are required that evaluate the vector competence and capacity of European populations of Ae. albopictus for the circulating strain of ZIKV and the potentially related contributions of non-vector transmission.