Effects of changes in temperature on Zika dynamics and control

When a rare pathogen emerges to cause a pandemic, it is critical to understand its dynamics and the impact of mitigation measures. We use experimental data to parametrize a temperature-dependent model of Zika virus (ZIKV) transmission dynamics and analyse the effects of temperature variability and control-related parameters on the basic reproduction number (R0) and the final epidemic size of ZIKV. Sensitivity analyses show that these two metrics are largely driven by different parameters, with the exception of temperature, which is the dominant driver of epidemic dynamics in the models. Our R0 estimate has a single optimum temperature (≈30°C), comparable to other published results (≈29°C). However, the final epidemic size is maximized across a wider temperature range, from 24 to 36°C. The models indicate that ZIKV is highly sensitive to seasonal temperature variation. For example, although the model predicts that ZIKV transmission cannot occur at a constant temperature below 23°C (≈ average annual temperature of Rio de Janeiro, Brazil), the model predicts substantial epidemics for areas with a mean temperature of 20°C if there is seasonal variation of 10°C (≈ average annual temperature of Tampa, Florida). This suggests that the geographical range of ZIKV is wider than indicated from static R0 models, underscoring the importance of climate dynamics and variation in the context of broader climate change on emerging infectious diseases.


Introduction
Vector-borne viruses (arboviruses) are emerging threats to both human and animal health. The global expansions of dengue virus (DENV), West Nile virus (WNV), chikungunya (CHIKV) and most recently Zika virus (ZIKV) are prominent examples of how quickly mosquito-transmitted viruses can emerge and spread through naive host populations. Currently, 3.9 billion people living within 120 countries are at risk of mosquito-borne arboviral diseases [1] with effects on human well-being that can be devastating (e.g. death, illness, as well as social and human ramifications of Zika-induced microcephaly and other congenital disorders) [2]. Anticipating and preventing outbreaks of emerging mosquito-borne viruses across these host populations is a major challenge.
Despite growing research to develop new therapeutics and vaccines, mitigating arbovirus disease spread still depends on conventional mosquito control methods, often with mixed success. Developing tools that allow us to successfully predict outbreaks of these viruses and efficiently target current and future interventions to specific times and locations can aid effective mosquito and disease control. Such efforts are often limited by gaps in knowledge on the relationships among mosquito vectors, pathogens and the environment, especially for emerging arboviruses such as CHIKV and ZIKV. Even in well-researched disease systems (e.g. malaria and DENV), key transmission parameters are only estimated from a few studies [3][4][5].
Several modelling frameworks have been used to predict environmental suitability for vector-borne disease transmission, including, most recently, temperature-dependent R 0 models [3][4][5]18,19] and compartmental models of vectorborne disease dynamics [8,34,35]. The parameter R 0 is broadly considered to be the most important summary statistic in epidemiology and disease ecology. It is defined as the expected number of new human (respectively, mosquito) cases generated by a single infectious human (respectively, mosquito) introduced into a fully susceptible human (respectively, mosquito) population throughout the period within which that human (respectively, mosquito) is infectious [36]. As a simple metric, it can easily incorporate the nonlinear influence of multiple temperature-dependent mosquito and pathogen traits, and has been applied to define the thermal optimum and limits for malaria [4,5,37], DENV, CHIKV [3,38,39], ZIKV [3,19] and Ross River virus [18]. However, temperature-dependent R 0 formulations only define the relative risk of disease emergence and do not predict the final epidemic size (or incidence). The derivation, interpretation, and validation of R 0 models are thus problematic in highly variable systems [40]. Dynamical models of transmission that track densities of infectious individuals over time, on the other hand, can more readily capture the impact of varying environmental conditions.
To better understand potential climate effects on control strategies for ZIKV, we developed a temperature-dependent dynamical model based on recent experimental work characterizing temperature-trait relationships between ZIKV vector competence, extrinsic incubation rate, and the per capita daily mosquito mortality rate [41]. Unlike other published results in the literature (e.g. R 0 model in Tesla et al. [19]), we model ZIKV transmission dynamics between humans and vectors and the flow of humans and vectors between various classes explicitly through a compartmental SEIR model for the human population and SEI model for the vector population. Because the model is dynamic, we are also able to account for seasonal temperature variation. The model and analysis differ from Huber et al. [8] in that, through numerical and sensitivity analyses, we explicitly analyse the simultaneous effects of parameters that are influenced by control measures (including vaccination) on both the basic reproduction number R 0 and the final epidemic size. Our analysis thus addresses the following questions: (1) How do the thermal optima and ranges for R 0 compare to those for the human final epidemic size? (2) How does seasonal temperature variation affect the final epidemic size relative to a constant temperature environment? (3) Which parameters have the greatest impact on R 0 and the final epidemic size that can inform control efforts? (4) Are different thermal environments more or less suitable for specific control strategies?
Our results show that R 0 and the final epidemic size were largely driven by different sets of parameters, with the exception of temperature being the dominant driver of both. Furthermore, the human final epidemic size was maximized across a wider range of temperatures than what would have been predicted from the temperature-dependent R 0 model. The human final epidemic size was highly sensitive to seasonal temperature variation, suggesting the potential invasion map of ZIKV may be wider than previously reported. Furthermore, the effectiveness of potential control strategies (e.g. vaccines, drug treatment, and insecticides, assessed through model parameters that are influenced by those strategies) is predicted to be sensitive to such differences in seasonal temperature variation.

Methods
We construct a temperature-dependent compartmentalized model of ZIKV dynamics with and without seasonal temperature change. Where possible, model parameters are estimated from the most recent laboratory experiments on temperature effects on the life cycle of the virus [3,19]. We first compare how temperature dependence affects R 0 , the human 'final epidemic size' (total number of infected individuals over the course of the epidemic), and key ecological characteristics of the system, such as extrinsic incubation period, the probability of transmission from the mosquito to the human, the probability of transmission from the human to the mosquito, and daily rates of mosquito and egg to adult survival. To understand how temperature change can influence the effects of control measures, we then analyse the combined effects of temperature and parameters that correspond to disease control measures on R 0 and the final epidemic size. These 'control parameters' include vaccination, recovery, vector biting rates, vector-to-human transmission, vector carrying capacity, egg survival, and adult mosquito survival. Through a Latin hypercube sampling-based sensitivity analysis [42], we identify key parameters that most drive the epidemiological outcomes (R 0 and the final epidemic size). All simulations were carried out using MATLAB R2019b.

The basic dynamic model
The model apportions humans into four groups based on ZIKV infection status that changes over time, t: susceptible S h (not infected), latent E h (contracted the virus, but not yet infectious), infectious I h (contracted the virus and can transmit it), and recovered R h with lifelong immunity. The mosquito population is divided into similar classes, where the state variables have subscript v, but without an immune class since it is assumed that infectious mosquitoes do not clear the virus once it is in the salivary glands. The total human and mosquito populations are The model assumes a constant human population during the epidemic. Susceptible humans acquire the virus at rate (force of is the number of human bites per mosquito per unit time, β vh is the probability that an infectious mosquito successfully transmits the virus while taking a blood meal from a susceptible human (i.e. the transmission rate), and β hv is the probability that an infectious human successfully transmits the virus to a biting, susceptible mosquito (i.e. the infection rate). The respective average residence times of infected humans and mosquitoes in the latent classes are 1/σ h and 1/σ v , while the respective rates at which humans and mosquitoes become infectious are σ h and σ v . Humans are infectious for approximately 1/γ h days before recovering with permanent immunity (γ h is the per capita human recovery rate), while infectious mosquitoes remain infectious until they die. Mosquito recruitment occurs at a per capita rate f( where κ v is the carrying capacity (maximum number of mosquitoes a breeding site can support). Furthermore, α v = θ v ν v ϕ v /μ v consists of θ v , or the number of eggs a female mosquito produces per day; ν v , the probability of surviving from egg to adult; and ϕ v , the rate at which an egg develops into an adult mosquito. Mosquitoes die naturally at per capita rate μ v , where 1/μ v is the average lifespan of mosquitoes. See figure 1 for a schematic of the model and table 1 for descriptions and baseline values of the parameters. The dynamic model for the Zika virus is described by the equations (2:1) Dots denote differentiation with respect to time, t (in days). The dynamics of the total human population and mosquito populations are described, respectively, by the equations Without Zika virus, the mosquito population grows according to equation (2.2) . 0 for α v > μ v is the positive equilibrium obtained by setting the right-hand side of the equation to zero. Observe that N v (0) ¼ N 0 v , and that when α v > μ v , the total mosquito population relaxes on the equilibrium population (N Ã v ) in the long run. Therefore, the equilibrium point N Ã v is stable when α v > μ v and vanishes when α v < μ v . The case for which α v < μ v results in a trivial mosquito equilibrium represents a situation in which the mosquito population becomes extinct. Since To illustrate the dynamics of the system, we set N h = 1000. The next generation operator approach [36,51,52] is used to compute the basic reproduction number of the model (2.1). This involves re-writing equations (2.1) as two sub-systems-a disease-free sub-system (consisting of the equations for the susceptible human and vector classes and the recovered human class) and a disease sub-system (consisting of equations for the exposed and infectious human and mosquito classes-i.e. the second, third, sixth and seventh equations of the model (2.1)).

The disease system is expressed as the difference of two vectors-a vector of new infections
. The corresponding matrices of new infections F and transitions V (given by the Jacobians of the vectors F and V), as well as the inverse of the matrix V are royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20210165 The spectrum (set of eigenvalues) of the next generation matrix The basic reproduction number of the model (2.1) is the spectral radius, i.e. the largest eigenvalue of the matrix FV À1 . Hence, in the presence of the Zika virus, the basic reproduction number of system (2.1) is The main difference between this R 0 calculation and that from the Ross-MacDonald model [53,54] is in the probability that the mosquito survives the latent period. The Zika virus can spread when R 0 > 1 and can be contained when R 0 < 1.
Variants of the model (2.1) have been used to assess the impact of various control measures including insecticide-treated bed nets on vector-borne diseases such as malaria [55][56][57][58][59][60]. For the purposes of exploring control strategies, we also consider a variant of this model that includes vaccination, where vaccinated susceptible humans are assumed to enter the immune class directly. For this special case, the first and fourth equations of (2.1) are replaced by _ respectively, where δ h represents the per capita human vaccination rate.

Introducing temperature
The majority of the parameters associated with the mosquito vector (θ ev , ν v , ϕ v , b v , μ v ), as well as ZIKV transmission (β hv ) and replication (σ v ), are known to be influenced by environmental and climate conditions [3,19]. We investigate the effects of temperature variation on the dynamics of the mosquito population and ZIKV transmission over time. We follow the approach in [8] and model temperature-dependent parameters with the functional forms presented in equation (2.4). We rely on values and ranges of temperature-dependent parameters from recent laboratory-generated analyses for Zika virus [19] and Ae. aegypti lifehistory parameters (i.e. the biting rate of mosquitoes, the number of eggs a female mosquito lays per day, the probability of an egg surviving to an adult mosquito, and the rate at which an egg develops into an adult mosquito) [3] as specified in table 2. As in [8], the temperature dependent parameters are based on the quadratic and Briere [61] functional forms, c(T − T m )(T − T 0 ) and cT(T À T 0 )(T m À T) 1 2 , respectively, where T, c, T 0 and T m are the temperature, rate (or scaling) constant, critical thermal minimum temperature and critical thermal maximum temperature, respectively: Eggs per female mosquito per day: Vector mortality rate: The subscripts θ v , ν v , ϕ v , b v , β hv , σ v and ls on c, T 0 , and T m indicate that the temperature-related parameter (c, T 0 or T m ) is associated with the corresponding original model parameter . These relationships between some model parameters, model outcomes, and temperature are illustrated in figure 2.
We further introduce seasonal variation in the system by modelling temperature through the sinusoidal functional form: where T m is the mean annual temperature and T a is the amplitude (divergence from mean temperature or mid-point between Here, we are assuming that at the onset of an epidemic outbreak, the phase is set to zero.

Control strategies
To analyse the relationship between temperature and Zika control, we identified the following parameters of the system that correspond to potential control strategies: vaccination (δ h ) decreases susceptibility and is directly incorporated into the models as described above; recovery rates (γ h ) can, for example, be increased through treatment with antiviral medication; vector biting rates (b v ) can be reduced through decreasing exposure to mosquitoes with personal protection or household improvements; vector-to-human transmission probability (β vh ) can decrease with transmission-blocking Wolbachia; the vector carrying capacity (κ v ) can be reduced by eliminating vector breeding grounds near human habitats; egg-adult survival probability (ν v ) can be reduced through larvicides; and adult mosquito survival rate (μ v ) can be decreased through indoor spraying and the use of adulticides. We investigate the effects of the interaction between temperature and these parameters that are influenced by common disease control methods on R 0 and the final epidemic size (total I h ).

Sensitivity analysis
Two types of sensitivity analyses-local and global-were used to explore the impact of temperature and selected parameters that are affected by disease control measures on the basic reproduction number (R 0 ) and the human final epidemic size (total I h ). The local sensitivity analysis was conducted by varying only one parameter while holding all other parameters fixed, or varying both temperature and a control-sensitive parameter while holding the other parameters fixed. Each parameter that was varied was divided into 50, 100 and 250 equally spaced points within biologically feasible bounds. See figures 2-6 for results. As the human vaccination rate (δ h ) does not appear explicitly in the expression of R 0 , we could not assess how temperature modified the effect of this parameter on the basic reproduction number. However, we explore the impact of temperature on the human vaccination rate and associated implications for I h , the human final epidemic size ( figure 3). Global sensitivity analysis is presented in figure 7. The analysis is carried out using the Latin hypercube sampling and partial rank correlation coefficient (PRCC) technique [42]. The process involves identifying a biologically feasible mean, minimum and maximum value for each of the parameters (e.g. [3,19]) and subdividing the range of each parameter into 1000 equal sub-intervals, assuming a uniform distribution between the minimum and maximum values of each parameter. We then sample at random and without replacement from the parameter distributions to generate an m × n Latin hypercube sampling matrix, whose m rows (i.e. 1000 rows) consist of different values for each of the model parameters and the n columns (corresponding to the number of parameters in the system) consist of different values for the same parameter. Thus, each row of the Latin hypercube sampling matrix provides a parameter regime that is used for computing the basic reproduction number, solving the dynamic system, and computing the human final epidemic size. The parameters, basic reproduction number, and the human final epidemic size are then ranked with partial correlation coefficients estimated for each parameter along with corresponding p-values. PRCCs range from −1 to 1 and are used to examine the correlation between model parameters and model outputs (R 0 and the final epidemic size). This method thus identifies parameters with the most significant influence on model outputs; it does not quantify the effect of a change in a parameter on the output.

Mapping seasonal control
We mapped the R 0 as a function of monthly mean temperature (figure 5). Globally gridded monthly mean current temperatures were downloaded from WorldClim.org [62], at a 5 arc-minute resolution (approx. 10 km 2 at the equator), and predicted rates as a function of temperature at 0.2°C were mapped to the global grids. All raster calculations and graphics were conducted in R, using package raster [63].

Impact of temperature on model parameters and key outputs
The models show unimodal relationships between temperature and the temperature-dependent parameters, resulting in an optimal temperature that maximizes parameter values and a critical minimum and maximum temperature at which parameter values go to zero. Figure 2 presents the effects of temperature on mosquito and pathogen parameters, the final epidemic size in humans (total number of infected Table 2. Parameter values for temperature-dependent functional forms. The temperature in degrees Celsius is denoted by T, the rate scaling factor is denoted by c (main Column 2), the minimum critical thermal temperature is denoted by T 0 (main Column 3), and the maximum critical thermal temperature is denoted by T m (main Column 4). The subscripts denote the corresponding traits, e.g. c u v , T 0 u v and T m u v are for the number of eggs laid by a female mosquito per day, θ v . See equation ( royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20210165 individuals over the course of the epidemic) and mosquitoes, and the basic reproduction number, R 0 (via temperature effects on mosquito and ZIKV parameters). The response of R 0 to temperature is strongly peaked (with an optimum around 30°C) as has been demonstrated for Zika and other systems (e.g. dengue, malaria, Ross River virus) [3,5,18,19]. By contrast, the relationship between the final epidemic size and temperature is flat; i.e. within the thermal range of disease transmission, the final epidemic size does not change (figure 2i versus 2j ). At temperatures associated with lower epidemic peaks, there are longer epidemic periods, resulting in the same total number of infected individuals over the course of the epidemic (figure 2l).

Zika virus control
The effects of static temperature (i.e. not including seasonality) and control-related parameters on final epidemic size are presented in figure 3. The most striking feature of these plots is that the difference between very small and very large epidemics (represented by blue and red areas, respectively) is discrete for most parameters. Crossing these thresholds, the most significant changes occur not from incremental changes in control parameters at a given temperature, but when temperatures move into the suitable thermal range, in which case the models predict widespread transmission (this concept is more rigorously explored in the sensitivity analysis below). However, vaccination (δ h ) and recovery (γ h )   have more incremental effects (within their feasible ranges) on final epidemic size even within the suitable thermal range; i.e. proportion of the population that needs to be vaccinated is higher at optimal temperatures than at sub-optimal temperatures to achieve a given reduction in the overall final epidemic size (figure 3d). Thus, warming temperatures (for most countries) will require greater vaccination coverage and treatment rates in order to maintain control of Zika.
The effects on the basic reproduction number (R 0 ) of most parameters that are sensitive to control measures were more dependent on temperature (figure 4) than was the case for the final epidemic size (total I h ), with the greatest effect involving the clearance rate of infection (γ h ), the probability of transmission from an infectious human to a susceptible mosquito (β hv ), mosquito carrying capacity (κ v ) and the mosquito mortality rate (μ v ).  Figure 4. Effect on the basic reproduction number, R 0 , of temperature and parameters that are sensitive to control measures for the model without seasonal temperature variability: (a) the mosquito biting rate, b v ; (b) the human recovery rate, γ h ; (c) the probability of transmission from the mosquito to the human, β vh ; (d ) the probability of transmission from the human to the mosquito, β hv ; (e) the mosquito carrying capacity, κ v ; and ( f ) the mosquito mortality rate, μ v . (Vaccination δ h and the egg-to-adult mosquito survival probability, ν v , do not appear explicitly in the R 0 model.) The effect of these parameters on R 0 is more dependent on temperature than their effects on the final epidemic size (e.g. the colour bands are less vertical and more diagonal in parts of the parameter space).  Figure 5. Estimated basic reproduction number R 0 , for different mean monthly temperatures, globally. The colour bar represents increasing R 0 values with the lowest values denoted by dark blue and the highest values denoted by light yellow. A value of R 0 below 1 means the disease cannot take off. Note the month-to-month variation in locations such as Australia and the southern USA, demonstrating the impact of seasonal temperature variation on R 0 , and therefore on the control level needed to reduce or eradicate disease.

Seasonal variation
Seasonal temperature variation affects outcomes by providing transient temperatures (variation from the annual mean; equation (2.5)) where the basic reproduction number can rise above (or fall below) 1, allowing for transmission to occur (or cease). At constant temperatures, epidemics only occur in humans between ≈23−37°C (figure 3). Such average annual temperatures are only found in tropical countries. However, the model shows the potential for epidemics in areas with mean temperatures below 23°C if there is adequate seasonal variation. This would be the case, for example, for a subtropical area, such as Tampa, FL, with a mean annual temperature of ≈20°C and an amplitude of ≈10°C (figure 6a). In contrast to the models without seasonal variation (figure 3), the models with seasonal variation ( figure 6) indicate that the effectiveness of parameters on the final epidemic size is generally sensitive to changes in temperature (e.g. the colour bands in the subplots of figure 6 are diagonal in more of the parameter space than they are in figure 3). Figure 5 shows how the thermal conditions that are suitable for Zika (where R 0 > 1) change with seasonal temperature variation across the globe.

Global uncertainty and sensitivity analysis
A global sensitivity analysis using Latin hypercube sampling showed that R 0 and the final epidemic size are largely sensitive to different parameters (as indicated by differences in figures 3 and 4). However, temperature is a dominant driver of variation in both the basic reproduction number (R 0 ) and the final epidemic size (total I h ) when it is included in the model (figure 7c-e). The human recovery rate, γ h , was a consistently influential driver of the final epidemic size. By contrast, the basic reproduction number was not sensitive to recovery in the models with and without variable temperature. While R 0 was also sensitive to vector competence (β vh and β hv ), biting rate (b v ), and mosquito lifespan (1/μ), total infection burden was far less sensitive to these parameters and was mainly sensitive to human recovery rate (γ h ) (figure 7a,b).

Discussion
We are interested in what drives arbovirus epidemics, with Zika as a model, and how to reduce the burden of these diseases, focusing on temperature and key parameters that correspond to existing or potential control methods (e.g. pesticides, reduced breeding habitats, vaccines or treatment). We investigated temperature-dependent dynamic transmission models that incorporated recent empirical estimates of the relationships between temperature and Zika infection, transmission, and mosquito lifespan [19]. These dynamical models that can measure final epidemic size and account for temperature variation generate qualitatively different results from static R 0 models. Temperature had an overwhelmingly strong impact on both R 0 and the final epidemic size (total infectious individuals, equivalent to area under the I h epidemic curve), but the response was much more gradual and had a clear optimum for R 0 , while the final epidemic size responded as a threshold function (figure 2i,j ). This is because, while epidemics have a higher peak at the maximum R 0 (at optimal temperatures), the epidemics are longer at sub-optimal temperatures (lower R 0 ). Thus, Zika   Figure 6. The human final epidemic size is sensitive to the annual mean (oscillation) temperature, T m , and several important control-related parameters: (a) the seasonal divergence of the annual temperature from the mean, T a ; (b) the human recovery rate, γ h ; (c) the human vaccination rate, δ h ; (d) the scaling factor of vector biting rate, c b v ; (e) the probability of transmission from the mosquito to the human, β vh ; (f ) the scaling factor of the probability of transmission from the human to the mosquito, c b hv ; (g) the scaling factor of the vector mortality rate, c ls ; (h) the vector carrying capacity, κ v ; and (i) the egg survival probability scaling factor, c n v . The annual mean (oscillation) temperature varies along the x-axis, the seasonal divergence of the annual temperature from the mean and the other parameters that are sensitive to control measures vary along the y-axes, and the colour scale indicates the total infectious humans. Apart from (a), where the temperature amplitude (T a ) is varying, the amplitude is set at 10°C for the other plots. Plot (a) shows that there can be large epidemics even when mean temperatures are low if the seasonal variation (the amplitude) is high enough, as would be found in subtropical and temperate regions.
royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20210165 virus is capable of spreading efficiently through the host population (high I h ) across a broad range of temperatures for which R 0 > 1, spanning from 17 to 37°C in constant environments (figure 6a) [8]. This is broadly consistent with the high seroprevalence of Zika found in a number of countries [64,65]. This suitable temperature region expanded and shifted toward cooler mean temperatures under seasonally varying environments (figure 6a).
These results have two key implications. First, large epidemics can occur under realistic, seasonally varying, temperature environments even in regions where the mean temperature alone would be expected to suppress transmission, for example in a location with a mean of 20°C and a seasonal amplitude of 10°C (e.g. Tampa, FL). Second, temperature determines both upper and lower thresholds for whether or not epidemics are possible [8]. However, within the predicted suitable temperature range defined by R 0 , the final epidemic size is largely limited by the density of susceptible hosts (figures 2 and 6a) [8]. More broadly, the results highlight the important principle that metrics of transmission (e.g. R 0 ) have a nonlinear relationship with the human final epidemic size (total I h ) and contribute distinct implications for our understanding of the transmission process.
Whether or not temperature affects the potential for disease control is an important applied question for designing public health campaigns (either via vector control, reduction in host biting rate, vaccination, or drug administration). Temperature did not strongly affect the impact of most control-related parameters on the final epidemic size when the models did not include seasonal variation. In the more realistic situation where models included seasonal variation, the effectiveness of most parameters that are sensitive to disease control measures depended on temperature. In all models, human vaccination rate required to control epidemics varied strongly with mean temperature (figures 3-6). Achieving herd immunity and thereby suppressing transmission via vaccination is more difficult when temperatures are highly suitable (20−35°C under constant temperatures or 15−32°C under varying temperatures; figures 3-6). By contrast, the effects of the human recovery rate (γ h ) and the vector mortality rate (μ v ) on R 0 were sensitive to temperature, but their effects on the final epidemic size were not sensitive to temperature.
Similar to previous work on dengue [8], our results show that Zika can invade and cause large outbreaks during the summer in seasonally varying environments with lower average temperatures, such as temperate regions of the USA, Europe, and Asia. This implies that differences in the size of epidemics in tropical versus temperate locations occur not just because of differences in temperature (and its impacts on R 0 ) but also because of differences in vector breeding habitat availability, humidity, human mosquito exposure, and other socioecological factors. Much of the globe-including regions in  Figure 7. Global sensitivity analysis indicates the sensitivity of the basic reproduction number, R 0 (a and c), and the final epidemic size, total I h (b, d and e), to all model parameters. Bars indicate partial rank correlation coefficients (PRCC), illustrating the contribution of parameters to variability or uncertainty in the model outputs (R 0 and I h ). (a-b) No temperature dependence; (b-c) temperature dependence but no temperature variation; (e) temperature dependence and temperature variation. Without accounting for temperature, models show different sets of drivers for R 0 (a) than burden (b). When temperature is included, it is the dominant contributor to model sensitivity for all models (c-e).
royalsocietypublishing.org/journal/rsif J. R. Soc. Interface 18: 20210165 temperate, subtropical, and tropical climates-is already suitable for Zika transmission for all or part of the year, and climate change is likely to expand this suitability geographically and seasonally [66]. However, processes that increase the density of susceptible human populations and their exposure to mosquitoes, including urbanization and urban poverty, human population growth, and the growth and geographical expansion of vector populations, are likely to expand the burden of Zika even more dramatically in the future.
There are several limitations of such a modelling study. First, the parameters are determined by a combination of laboratory-based estimates as well as from the literature on dengue, instead of being fitted to empirical data on the spatio-temporal dynamics of Zika from the field. Such divergent approaches can generate different parameter estimates. Furthermore, the projections of the model on potential geographical distribution of Zika are based on the average of temperatures by country and season with constant parameters. In reality, there is substantial heterogeneity of temperature and parameters over time and space, which have important implications for disease dynamics. For this reason, further investigation of Zika models that are calibrated from field based dynamics will be valuable for a fuller understanding of effects of temperature variation on Zika control.

Conclusion
The unexpected emergence and global expansion of Zika in 2015-2017 and its association with Zika congenital syndrome and Guillain-Barre syndrome revealed once again how poorly prepared the global community is for the looming and expanding threat of vector-borne diseases. Given the recent history of the rapid global expansion of Aedes aegypti-transmitted viruses (e.g. DENV, CHIKV and ZIKV), as well as the challenges of controlling these epidemics, understanding the ecological drivers of transmission and their effects on potential disease control tools is crucial for improving preparedness for future vectorborne disease emergence. If a Zika vaccine becomes available, then the precisely defined temperature thresholds for large epidemics predicted in our model imply that vaccination targets should be set based on climate. By contrast, because other potential interventions that would reduce vector population sizes, biting rates, and human recovery rates act more independently of temperature, targets could be set based on other socio-ecological factors in a given epidemic setting. This dynamic temperature-dependent modelling framework, which depends most strongly on vector and host parameters that are virus-independent, may be a useful first step for responding to future Aedes-borne disease epidemics.
Data accessibility. This article has no additional data.