article: Thermal biology predicts that vector-borne disease transmission peaks at intermediate temperatures and declines at high and low temperatures. However, thermal optima and limits remain unknown for most vector-borne pathogens. We built a mechanistic model for the thermal response of Ross River virus, an important mosquito-borne pathogen in Australia, Pacific Islands, and potentially at risk of emerging worldwide. Transmission peaks at moderate temperatures (26. 4°C) and declines to zero at thermal limits (17. 0 and 31. 5°C). The model accurately predicts that transmission is year-round endemic in the tropics but seasonal in temperate areas, resulting in the nationwide seasonal peak in human cases. Climate warming will likely increase transmission in temperate areas (where most Australians live) but decrease transmission in tropical areas where mean temperatures are already near the thermal optimum. These results illustrate the importance of nonlinear models for inferring the role of temperature in disease dynamics and predicting responses to climate change.
  The natural history of RRV is complex: transmission occurs across a range of climates (tropical, subtropical, and temperate) and habitats (urban and rural, coastal and inland) and via many vertebrate reservoir and vector species (Claflin and Webb, 2015). Marsupials are generally considered the critical reservoirs for maintaining the virus between human outbreaks, but recent work has argued that placental mammals and birds may be equally important in many locations (Stephenson et al. , 2018). The virus has been isolated from over 40 mosquito species in nature, and 10 species transmit it in laboratory studies (Harley et al. , 2001; Russell, 2002). However, four species are responsible for most transmission to humans (Culex annulirostris, Aedes [Ochlerotatus] vigilax, Ae. [O. ] notoscriptus, and Ae. [O. ] camptorhynchus), with two additional species implicated in outbreaks (Ae. [Stegomyia] polynesiensis and Ae. [O. ] normanensis). The vectors differ in climate and habitat niches, leading to geographic variation in associations with outbreaks. We assembled and mapped records of RRV outbreaks in humans attributed to different vector species (Figure 1, Figure 1—source data 1) (Rosen et al. , 1981; Campbell et al. , 1989; Russell et al. , 1991; Yang et al. , 2009; Lindsay et al. , 1993b; Lindsay et al. , 1993a; Lindsay et al. , 1996; Lindsay et al. , 2007; McManus et al. , 1992; Merianos et al. , 1992; Whelan et al. , 1992; Whelan et al. , 1995, Whelan et al. , 1997; McDonnell et al. , 1994; Russell, 1994; Russell, 2002; Dhileepan, 1996; Ritchie et al. , 1997; Brokenshire et al. , 2000; Ryan et al. , 2000; Harley et al. , 2000; Harley et al. , 2001; Kelly-Hope et al. , 2004a; Frances et al. , 2004; Biggs and Mottram, 2008; Jacups et al. , 2008b; Schmaedick et al. , 2008; Lau et al. , 2017). Ae. vigilax and Ae. notoscriptus were more commonly implicated in transmission in tropical and subtropical zones, Ae. camptorhynchus in temperate zones, and Cx. annulirostris throughout all climatic zones. Freshwater-breeding Cx. annulirostris has been implicated in transmission across both inland and coastal areas, while saltmarsh mosquitoes Ae. vigilax and Ae. camptorhynchus have been implicated only in coastal areas (Russell, 2002) and inland areas affected by salinization from agriculture (Biggs and Mottram, 2008; Carver et al. , 2009). Peri-domestic, container-breeding Ae. notoscriptus has been implicated in urban epidemics (Russell, 2002). The vectors also differ in their seasonality: Ae. camptorhynchus populations peak earlier and in cooler temperatures than Ae. vigilax, leading to seasonal succession where they overlap (Yang et al. , 2009; Russell, 1998). This latitudinal and temporal variation suggests that vector species may have different thermal optima and/or niche breadths. If so, temperature may impact disease transmission differently for each species. Transmission depends on a suite of vector, pathogen, and human traits, including mosquito density (M). Our main model (‘full R0 Model, ’ Equation 1) assumes temperature drives mosquito density and includes the relevant life history trait thermal responses (Parham and Michael, 2009; Mordecai et al. , 2013; Mordecai et al. , 2017). We initially compare this model to an alternative (‘constant M model, ’ Equation 2) where mosquito density does not depend on temperature. We make this comparison because many transmission models do not include the thermal responses for mosquito density, assuming it depends primarily on habitat availability. Here, we focus on the relative influence of temperature on transmission potential, recognizing that absolute R0 also depends on other factors. Accordingly, we scaled model output between zero and one (‘relative R0’). Relative R0 describes thermal suitability for transmission. Combined with factors like breeding habitat availability, vector control, humidity, human and reservoir host density, host immune status, and mosquito exposure, relative R0 can be used to predict disease incidence. In this approach, only the relative thermal response of each trait influences R0, which is desirable since traits can differ substantially due to other factors and in laboratory versus field settings (particularly mosquito survival: Clements and Paterson, 1981). Relative R0 does not provide a threshold for sustained disease transmission (i. e. where absolute R0 = 1), since this threshold is not controlled solely by temperature. Instead, relative R0 preserves the temperature-dependence of R0 to provide three key temperature values: upper and lower thermal limits where transmission is possible (R0 >0; a conservative threshold where transmission is not excluded by high or low temperatures) and the temperature that maximizes R0.
  Vector and pathogen traits that drive transmission consistently responded to temperature (Figure 2), though data were sparse (McDonald et al. , 1980; Mottram et al. , 1986; Russell, 1986; Rae, 1990; Kay and Jennings, 2002). Although we exhaustively searched for experiments with trait measurements at three or more constant temperatures in the Australian vector species (Cx. annulirostris, Ae. vigilax, Ae. camptorhynchus, Ae. notoscriptus, and Ae. normanensis), no species had data for all necessary traits (Figure 1). Thus, we combined traits from two species to build composite R0 models. We used mosquito life history traits measured in Cx. annulirostris: fecundity (as eggs per female per day, EFD), egg survival (as the proportion of rafts that hatch, pRH, and the number of larvae emerging per viable raft, nLR), the proportion surviving from larvae-to-adulthood (pLA), mosquito development rate (MDR), adult mosquito lifespan (lf), and biting rate (a). We used infection traits measured in Ae. vigilax: vector competence (bc) and parasite development rate (PDR). For comparison, we also fit traits for other mosquito and virus species: MDR and pLA from Ae. camptorhynchus and Ae. notoscriptus, and PDR and bc from Murray Valley encephalitis virus (another important pathogen transmitted by these mosquitoes in Australia) in Cx. annulirostris (Figure 2—figure supplements 2,3 and 4) (Kay et al. , 1989; Barton and Aberton, 2005; Williams and Rau, 2011). We used sensitivity analyses to evaluate the potential impact of this vector mismatch. However, all spatial and temporal predictions of R0 (Figures 5–7) use the full R0 model parameterized with mosquito life history traits from Cx. annulirostris and infection traits from Ae. vigilax (as shown in Figure 2). Thermal optima ranged from 23. 4°C for adult lifespan (lf) to 33. 0°C for parasite development rate (PDR; Figure 2). The data supported unimodal thermal responses for most traits, though declines at high temperatures were not directly observed for biting rate (a) and parasite development rate. Data from other mosquito species and ectotherm physiology theory imply these traits must decline at very high temperatures, so we used strong priors to make them decline near ~40°C. Because our approach is designed to identify which traits constrain transmission at thermal limits, this choice is conservative since it means R0 will be limited by other traits with better data. Accordingly, in the absence of data we preferred to overestimate upper thermal limits and underestimate lower thermal limits rather than vice versa. Transmission potential (relative R0 from the full R0 model) peaked at 26. 4°C, and was positive from 17. 0–31. 5°C (Figure 3). Removing the temperature-dependence of mosquito density [M] did not substantially affect the peak, because the optima for transmission and mosquito density were closely aligned (constant M model: 26. 6°C, M: 26. 2°C). By contrast, the range of temperatures suitable for transmission is much larger when mosquito density does not depend on temperature because M (T) constrains transmission at the thermal limits (constant M model positive from 12. 9–33. 7°C). The thermal constraints that mosquito density imposes on transmission are important because, although demographic traits are well-known to vary with temperature in the laboratory, many temperature-dependent transmission models do not assume that temperature influences mosquito density (Martens et al. , 1997; Craig et al. , 1999; Paull et al. , 2017; Caminade et al. , 2017; Hamlet et al. , 2018, but see Parham and Michael, 2009; Mordecai et al. , 2013, Mordecai et al. , 2017; Johnson et al. , 2015). The moderate optimal temperature for RRV (26–27°C) fits within the range of thermal optima found for other diseases: malaria transmission by Anopheles spp. at 25°C, and dengue and other viruses by Ae. aegypti and Ae. albopictus at 29 and 26°C, respectively (Figure 4) (Mordecai et al. , 2013; Mordecai et al. , 2017). At the upper thermal limit fecundity (EFD) and adult lifespan (lf) constrain R0, while at the lower thermal limit fecundity, larval survival (pLA), egg survival (raft viability [pRH] and survival within rafts [nLR]), and adult lifespan constrain R0 (Figure 3—figure supplement 2). All of these traits (except adult lifespan) only occur in, and adult lifespan is quantitatively more important in, the full R0 model, illustrating the importance of incorporating effects of temperature on vector life history. Correspondingly, uncertainty in these traits generated the most uncertainty in R0 at the respective thermal limits (Figure 3—figure supplement 2C). The optimal temperature for R0 was most sensitive to the thermal response of adult lifespan. Near the optimum, most uncertainty in R0 was due to uncertainty in the thermal responses of biting rate and egg raft viability. For comparison, substituting larval traits from alternative vectors or infection traits for Murray Valley Encephalitis virus did not substantially alter the R0 thermal response, since Cx. annulirostris life history traits strongly constrained transmission (Figure 3—figure supplement 1). Temperature suitability for RRV transmission varies seasonally across Australia, based on the full R0 model (Equation 1) using monthly mean temperatures from WorldClim. In subtropical and temperate locations (Brisbane and further south), low temperatures force R0 to zero for part of the year (Figures 5A and 6). Monthly mean temperatures in these areas fall along the increasing portion of the R0 curve for the entire year, so thermal suitability for transmission increases with temperature. By contrast, in tropical, northern Australia (Darwin and Cairns), the temperature remains suitable throughout the year (Figures 5 and 6). Darwin is the only major city where mean temperatures exceed the thermal optimum, and thereby depress transmission. Because most Australians live in southern, temperate areas, country-scale transmission is strongly seasonal. Using the average (1992–2013) seasonal incidence at the national scale, human cases peak two months after population-weighted R0 (T), matching our a priori hypothesized time lag between temperature suitability and human cases (based on empirical work in other mosquito-borne disease systems, see Materials and methods and Discussion; Figure 7).
  As the climate warms, it is critical to understand effects of temperature on transmission of mosquito-borne disease, particularly as new mosquito-borne pathogens emerge and spread worldwide. Identifying transmission optima and limits by characterizing nonlinear thermal responses, rather than simply assuming that transmission increases with temperature, can more accurately predict geographic, seasonal, and interannual variation in disease. Thermal responses vary substantially among diseases and vector species (Mordecai et al. , 2013; Mordecai et al. , 2017; Tesla et al. , in press), yet we lack mechanistic models based on empirical, unimodal thermal responses for many diseases and vectors. Here, we parameterized a temperature-dependent model for transmission of RRV (Figure 2) with data from two important vector species (Cx. annulirostris and Ae. vigilax; Figure 1). The optimal temperature for transmission is moderate (26–27°C; Figure 3), and largely determined by the thermal response of adult mosquito lifespan (Figure 3—figure supplement 2). Both low and high temperatures limit transmission due to low mosquito fecundity and survival at all life stages—thermal responses that are often ignored in transmission models (Figure 3—figure supplement 2). Temperature explains the geography of year-round endemic versus seasonally epidemic disease (Figures 5 and 6) and accurately predicts the seasonality of human cases at the national scale (Figure 7). Thus, the model for RRV transmission provides a mechanistic link between geographic and seasonal variation in temperature and broad-scale patterns of disease. While the thermal response of RRV transmission generally matched those of other mosquito-borne pathogens, there were some key differences. The moderate optimal temperature for RRV (26–27°C) fit within the range of thermal optima found for other diseases using the same methods: malaria transmission by Anopheles spp. at 25°C, and dengue and other viruses by Ae. aegypti and Ae. albopictus at 29 and 26°C, respectively (Figure 4) (Mordecai et al. , 2013; Mordecai et al. , 2017). For all of these diseases, the specific optimal temperature was largely determined by the thermal response of adult lifespan (Mordecai et al. , 2013, ; Mordecai et al. , 2017; Johnson et al. , 2015). However, the traits that set the thermal limits for RRV transmission differed from other systems. The lower thermal limit for RRV was constrained by fecundity and survival at all stages, while the upper thermal limit was constrained by fecundity and adult lifespan. By contrast, thermal limits for malaria transmission were set by parasite development rate at cool temperatures and egg-to-adulthood survival at high temperatures (Mordecai et al. , 2013). As with previous models, the upper and lower thermal limits of RRV transmission are more uncertain than the optimum (Figure 3) (Johnson et al. , 2015; Mordecai et al. , 2017), because trait responses are harder to measure near their thermal limits where survival is low and development is slow or incomplete. Overall, our results support a general pattern of intermediate thermal optima for transmission where the well-resolved optimal temperature is driven by adult mosquito lifespan, but upper and lower thermal limits are more uncertain and may be determined by unique traits for different vector–pathogen systems. Additionally, upper thermal limits of mosquito-borne disease transmission (a major concern for climate change) are primarily determined by vector life history traits with symmetrical thermal performance curves (like fecundity and survival at various life stages) rather than rate-based traits with asymmetrical thermal performance curves (like biting rate or pathogen development rate). The trait thermal response data were limited in two keys ways. First, two traits (fecundity and adult lifespan) had data from only three temperatures. We used priors derived from data from other mosquito species to minimize over-fitting and better represent the true uncertainty (Figure 2, versus uniform priors in Figure 2—figure supplement 1). However, data from more temperatures would increase our confidence in the fitted thermal responses. Second, no vector species had data for all traits (Figure 1), so we combined mosquito traits from Cx. annulirostris and pathogen infection traits in Ae. vigilax to build composite relative R0 models. Geographic and seasonal variation in vector populations suggests that Ae. camptorhynchus and Ae. vigilax have different thermal niches (cooler and warmer, respectively) and Cx. annulirostris has a broader thermal niche (Figure 1) (Russell, 1998). We need temperature-dependent trait data for more species to test the hypothesis that these niche differences reflect the species’ thermal responses. If true, the current model, parameterized primarily with Cx. annulirostris trait responses, may not accurately predict the thermal responses of transmission by Ae. camptorhynchus and Ae. vigilax. Hypothesized species differences in thermal niche could explain why RRV persists over a wide climatic and latitudinal gradient. Thus, thermal response experiments with other RRV vectors are a critical area for future research. The temperature-dependent R0 model provides a mechanistic explanation for independently-observed patterns of RRV transmission across Australia. As predicted by the model (Figures 5 and 6), RRV is endemic in tropical Australia, with little seasonal variation in transmission potential due to temperature, and seasonally epidemic in subtropical and temperate Australia (Weinstein, 1997). The model also accurately predicts disease seasonality at the national scale (Figure 7), reproducing the a priori predicted lag (8–10 weeks, or 2 months) for temperature to affect reported human cases (Hu et al. , 2006; Jacups et al. , 2008b; Stewart Ibarra et al. , 2013; Mordecai et al. , 2017). This lag between temperature and reported human cases arises from the time it takes for mosquito populations to increase, bite humans and reservoir hosts, acquire RRV, become infectious, and bite subsequent hosts; for pathogens to incubate with vectors; for humans to potentially develop symptoms, seek treatment, and report cases. Further, RRV transmission by Cx. annulirostris in inland areas often moves south as temperatures increase from spring into summer (Russell, 1998), matching the model prediction (Figure 6). Although temperature is often invoked as a potential driver for these types of patterns, it is difficult to establish causality from statistical inference alone, particularly if temperature and disease both exhibit strong seasonality and could both be responding to another latent driver. Thus, the mechanistic model is a critical piece of evidence linking temperature to patterns of disease. In addition to explaining broad-scale patterns, the unimodal thermal model explains previously contradictory local-scale results. Specifically, statistical evidence for temperature impacts on local time series of cases is mixed. RRV incidence is often—but not always—positively associated with warmer temperatures (Tong and Hu, 2001; Tong et al. , 2002a; Tong et al. , 2004; Hu et al. , 2004; Hu et al. , 2010; Jacups et al. , 2008b, Williams et al. , 2009; Werner et al. , 2012; Koolhof et al. , 2017). However, variation in the effects of temperature on transmission across space and time is expected from an intermediate thermal optimum, especially when observed temperatures are near or varying around the optimum. The strongest statistical signal of temperature on disease is expected in temperate regions where mean temperature varies along the rapidly rising portion of the R0 curve (~20–25°C). If mean temperatures vary both above and below the optimum (as in Darwin), important effects of temperature may be masked in time series models that fit linear responses. Additionally, if temperatures are always relatively suitable (as in tropical climates) or unsuitable (as in very cool temperate climates), variation in disease may be due primarily to other factors. A nonlinear mechanistic model is critical for estimating temperature impacts on transmission because the effect of increasing temperature by a few degrees can have a positive, negligible, or negative impact on R0 along different parts of the thermal response curve. Although field-based evidence for unimodal thermal responses in vector-borne disease is rare (but see Mordecai et al. , 2013; Perkins et al. , 2015; Peña-García et al. , 2017), there is some evidence for high temperatures constraining RRV transmission and vector populations: outbreaks were less likely with more days above 35°C in part of Queensland (Gatton et al. , 2005) and populations of Cx. annulirostris peaked at 25°C and declined above 32°C in Victoria (Dhileepan, 1996). Future statistical analyses of RRV cases may benefit from using a nonlinear function for temperature-dependent R0 as a predictor instead of raw temperature (Figure 6B versus Figure 6A). Breeding habitat availability also drives mosquito abundance and mosquito-borne disease. Local rainfall or river flow have been linked to the abundance of RRV vector species (Barton et al. , 2004; Tall et al. , 2014; Jacups et al. , 2015) and RRV disease cases (Tong and Hu, 2001; Tong et al. , 2002a; Hu et al. , 2004; Kelly-Hope et al. , 2004b; Tong et al. , 2004; Gatton et al. , 2005; Jacups et al. , 2008b; Bi et al. , 2009; Williams et al. , 2009; Werner et al. , 2012), as have high tides in coastal areas with saltmarsh mosquitoes, Ae. vigilax and Ae. camptorhynchus (Tong and Hu, 2002b; Tong et al. , 2004; Jacups et al. , 2008b; Kokkinn et al. , 2009). Overlaying models of species-specific breeding habitat with temperature-dependent models will better resolve the geographic and seasonal distribution of RRV transmission. Relative R0 peaked at similar temperatures whether or not we assumed mosquito abundance was temperature-dependent (Equation 1 versus Equation 2); however, the range of suitable temperatures was much wider for the model that assumed a temperature-independent mosquito population (Figure 3). Since breeding habitat can only impact vector populations when temperatures do not exclude them, it is critical to consider thermal constraints on mosquito abundance, even when breeding habitat is considered a stronger driver. Nonetheless, many mechanistic, temperature-dependent models of vector-borne disease transmission do not include thermal effects on vector density (Martens et al. , 1997; Craig et al. , 1999; Paull et al. , 2017; Caminade et al. , 2017; Hamlet et al. , 2018). Our results demonstrate that the decision to exclude these relationships can have a critical impact on model results, especially near thermal limits. Several important gaps remain in our understanding of RRV thermal ecology, in addition to the need for trait thermal response data for more vector species. First, the relative R0 model needs to be more rigorously validated using time series of cases to determine the importance of temperature at finer spatiotemporal scales. These analyses should incorporate daily, seasonal, and spatial temperature variation, including aquatic larval habitat and adult microhabitat temperatures (Paaijmans et al. , 2010; Cator et al. , 2013; Carrington et al. , 2013; Thomas et al. , 2018). They should also integrate species-specific drivers of breeding habitat availability, like rainfall and tidal patterns, infrastructure (e. g. drainage), and human activities (e. g. deliberate and accidental water storage). Second, translating environmental suitability for transmission into human cases also depends on disease dynamics in reservoir host populations and their impact on immunity. For instance, in Western Australia heavy summer rains can fail to initiate epidemics when low rainfall in the preceding winter depresses recruitment of susceptible juvenile kangaroos (Mackenzie et al. , 2000). By contrast, large outbreaks occur in southeastern Australia when high rainfall follows a dry year, presumably from higher transmission within relatively unexposed reservoir populations (Woodruff et al. , 2002). Third, as the climate changes, long-term predictions should consider potential thermal adaption of vectors, since transmission at upper thermal limits is currently limited by vector life history traits. To date, we know very little about standing genetic variation for thermal performance or existing local thermal adaptation in vectors for any disease system. Building vector species-specific R0 models and integrating thermal ecology with other drivers are important next steps for forecasting variation in RRV transmission. These more advanced models are necessary to translate our relative R0 results into predictions of absolute R0 (i. e. estimating the secondary cases per primary case, and where and when R0 >1 for sustained transmission). Nonlinear thermal responses are particularly important for predicting how transmission will change under future climate regimes. Climate warming will likely increase the geographic and seasonal range of transmission potential in temperate, southern Australia where most Australians live. However, climate change will likely decrease transmission potential in tropical areas like Darwin, where moderate warming (~3°C) would push temperatures above the upper thermal limit for transmission for most of the year (Figure 5). However, the extent of climate-driven declines in transmission will depend on how much Cx. annulirostris and Ae. vigilax can adapt to extend their upper thermal limits and whether warmer-adapted vector species (e. g. Ae. aegypti and potentially Ae. polynesiensis) can invade and sustain RRV transmission cycles. Thus, we can predict the response of RRV transmission by current vector species to climate change based on these trait thermal responses. However, future disease dynamics will also depend on vector adaptation, potential vector species invasions, and climate change impacts on sea level and precipitation that drive vector habitat availability.
  The ‘full R0 model’ (Equation 1) assumes temperature drives mosquito density and includes vector life history trait thermal responses (Parham and Michael, 2009; Mordecai et al. , 2013; Mordecai et al. , 2017). The ‘constant M model’ (Equation 2) assumes mosquito density (M) does not depend on temperature (Dietz, 1993). There is disagreement in the literature over whether the equation for R0 should contain the square root (Dietz, 1993; Heffernan et al. , 2005; Smith et al. , 2012). We use the version derived from the next-generation matrix method (Dietz, 1993) in order to be consistent with our previous work in other mosquito-borne disease systems (Mordecai et al. , 2013; Mordecai et al. , 2017; Johnson et al. , 2015). R0 (T) = (a (T) 2bc (T) e− μ (T) PDR (T) EFD (T) pEA (T) MDR (T) N r μ (T) 3) 1/2R0 (T) = (a (T) 2bc (T) e− μ (T) PDR (T) MN r μ (T) ) 1/2 In both equations, (T) indicates that a parameter depends on temperature, a is mosquito biting rate, bc is vector competence (proportion of mosquitoes becoming infectious post-exposure), µ is adult mosquito mortality rate (adult lifespan, lf = 1/µ), PDR is parasite development rate (PDR = 1/EIP, the extrinsic incubation period), N is human density, and r is the recovery rate at which humans become immune (all rates are measured in days−1). The latter two terms do not depend on temperature. In the full R0 model, mosquito density (M) depends on fecundity (EFD, eggs per female per day), proportion surviving from egg-to-adulthood (pEA), and mosquito development rate (MDR), divided by the square of adult mortality rate (µ) (Parham and Michael, 2009). We calculated pEA as the product of the proportion of egg rafts that hatch (pRH), the number of larvae per raft (nLR, scaled by the maximum at any temperature to calculate proportional egg survival within-rafts), and the proportion of larvae surviving to adulthood (pLA). We digitized previously published trait data (Figure 2—figure supplement 1; McDonald et al. , 1980; Mottram et al. , 1986; Russell, 1986; Rae, 1990; Kay and Jennings, 2002) using the free web-based tool Webplot Digitizer available at: https: //automeris. io/WebPlotDigitizer/. We fit thermal responses of each trait using Bayesian inference with the ‘r2jags’ package (Plummer, 2003; Su and Yajima, 2009) in R (R Core Team, 2017). Traits with asymmetrical thermal responses were fit as Brière functions: qT (T–Tmin) (Tmax–T) 1/2 (Briere et al. , 1999). Traits with symmetrical thermal responses were fit as quadratic functions: -q (T–Tmin) (T–Tmax). In both functions of temperature (T), Tmin and Tmax are the critical thermal minimum and maximum, respectively, and q is a rate parameter. For priors we used gamma distributions with hyperparameters derived from thermal responses fit to data from other mosquito species (Figure 2—source data 2) (Davis, 1932; Jalil, 1972; McLean et al. , 1974; Watts et al. , 1987; Rueda et al. , 1990; Focks et al. , 1993; Joshi, 1996; Teng and Apperson, 2000; Tun-Lin et al. , 2000; Alto and Juliano, 2001; Briegel and Timmermann, 2001; Kamimura et al. , 2002; Calado and Navarro-Silva, 2002; Focks and Barrera, 2006; Wiwatanaratanabutr and Kittayapong, 2006; Lardeux et al. , 2008; Delatte et al. , 2009; Beserra et al. , 2009; Yang et al. , 2009; Westbrook et al. , 2010; Muturi et al. , 2011; Carrington et al. , 2013; Tjaden et al. , 2013; Eisen et al. , 2014; Xiao et al. , 2014; Ezeakacha, 2015; Morin et al. , 2015). These priors allowed us to more accurately represent the fit and uncertainty. Our data did not include declining trait values at high temperatures for biting rate (a) and parasite development rate (PDR). Nonetheless, data from other mosquito species (Mordecai et al. , 2013; Mordecai et al. , 2017) and principles of thermal biology (Dell et al. , 2011) imply these traits must decline at high temperatures. Thus, for those traits we included an artificial data point where the trait value approached zero at a very high temperature (40°C), allowing us to fit the Brière function. We used strongly informative priors to limit the effect of these traits on the upper thermal limit of R0 (by constraining them to decline near 40°C). Because our approach is designed to identify which traits constrain transmission at thermal limits, this choice is conservative by allowing R0 to be limited by other traits with better data. Accordingly, in the absence of data we favored overestimating Tmax and underestimating Tmin over the alternative. For comparison, we also fit all thermal responses with uniform priors (Figure 2—figure supplement 1); these results illustrate how the priors affected the results. Bayesian inference produces estimated posterior distributions rather than a single estimated value. Because these distributions can be non-normal and asymmetric, we report and apply medians rather than means, since medians are less sensitive to outlying values in extended tails. However, we plot mean values in the figures because they show a smoother and more visually intuitive representation of where trait and R0 thermal responses go to zero at the upper thermal limit. The means and medians are not substantially different, except at this thermal limit (see example in Figure 3—figure supplement 3). We conducted sensitivity and uncertainty analyses of the full R0 model (Equation 1) to understand how trait thermal responses shape the thermal response of R0. We examined the sensitivity of R0 two ways. First, we evaluated the impact of each trait by setting it constant while allowing all other traits to vary with temperature. Second, we calculated the partial derivative of R0 with respect to each trait across temperature (∂R0/∂X · ∂X/∂T for trait X and temperature T; Equations 3-6). To understand what data would most improve the model, we also calculated the proportion of total uncertainty in R0 due to each trait across temperature. First, we propagated posterior samples from all trait thermal response distributions through to R0 (T) and calculated the width of the 95% highest posterior density interval (HPD interval; a type of credible interval) of this distribution at each temperature: the ‘full R0 (T) uncertainty’. Next, we sampled each trait from its posterior distribution while setting all other trait thermal responses to their posterior medians, and calculated the posterior distribution of R0 (T) and the width of its 95% HPD interval across temperature: the ‘single-trait R0 (T) uncertainty’. Finally, we divided each single-trait R0 (T) uncertainty by the full R0 (T) uncertainty. The partial derivatives are given below for all traits (x) that appear only once in the numerator of R0 (bc, EFD, pRH, nLR, pLA, MDR; Equation 3), biting rate (a, Equation 4), parasite development rate (PDR, Equation 5), and lifespan (lf, Equation 6). (3) ∂R0∂x= R02x (4) ∂R0∂a= R0a (5) ∂R0∂PDR= R02 lf PDR2 (6) ∂R0∂lf= R0 (1+3PDR) 2 PDR lf2 We took monthly mean temperatures from WorldClim for seven cities spanning a latitudinal and temperature gradient (from tropical North to temperate South: Darwin, Cairns, Brisbane, Perth, Sydney, Melbourne, and Hobart) and calculated the posterior median R0 (T) for each month at each location. We also compared the seasonality of a population-weighted R0 (T) and nationally aggregated RRV cases. We used 2016 estimates for the fifteen most populous urban areas, which together contain 76. 6% of Australia’s population (Australian Bureau of Statistics, 2017). We calculated R0 (T) for each location (as above) and estimated a population-weighted average. We compared this country-scale estimate of R0 (T) with data on mean monthly human cases of RRV nationwide from 1992 to 2013 obtained from the National Notifiable Diseases Surveillance System. We expected a time lag between temperature and reported human cases as mosquito populations increase, bite humans and reservoir hosts, acquire RRV, become infectious, and bite subsequent hosts; after an incubation period, hosts (potentially) become symptomatic, seek treatment, and report cases. Empirical work on dengue vectors in Ecuador identified a six-week time lag between temperature and mosquito oviposition (Stewart Ibarra et al. , 2013). Subsequent mosquito development and incubation periods in mosquitoes and humans likely add another 2–4 week lag before cases appear, resulting in an 8–10 week lag between temperature and observed cases (Hu et al. , 2006; Jacups et al. , 2008b; Mordecai et al. , 2017). With monthly case data, we hypothesize a two-month time lag between R0 (T) and RRV disease cases. To illustrate temperature suitability for RRV transmission across Australia, we mapped the number of months for which relative R0 (T) >0 and >0. 5 for the posterior median, 2. 5 and 97. 5% credibility bounds (Figure 5—figure supplement 1) for the full R0 model (Equation 1). We calculated R0 (T) at 0. 2°C increments and projected it onto the landscape for monthly mean temperatures from WorldClim data at a 5 min resolution (approximately 10 km2 at the equator). Climate data layers were extracted for the geographic area, defined using the Global Administrative Boundaries Databases (GADM, 2012). We performed map calculations and manipulations in R with packages ‘raster’ (Hijmans, 2016), ‘maptools’ (Bivand and Lewin-Koh, 2017), and ‘Rgdal’ (Bivand et al. , 2017), and rendered GeoTiffs in ArcGIS version 10. 3. 1. In 2000 there was a proposed shift in mosquito taxonomy: several subgenera within the genus Aedes were elevated to genus status (Wilkerson et al. , 2015). This change affected Aedes vigilax and Aedes camptorhynchus, which were called Ochlerotataus vigilax and Ochlerotatus camptorhynchus for a time by some researchers. More recently, there has been a consensus to return to the previous naming system, so we use Aedes here, although many of the papers we cite use Ochlerotatus instead. The fecundity and adult survival data in McDonald et al. (1980) were published as time series of one experimental population at each temperature. The resulting data needed to be transformed to fit the corresponding trait thermal responses. For survival, McDonald et al. reported the percent surviving approximately every other day (hereafter: ‘semi-daily’). We used these data—along with the number of female adults alive on the first day of oviposition at each temperature—to generate a semi-daily time series estimating the number surviving. To generate the dataset that we used to directly fit the thermal responses, we converted this time series into the number of female individuals who died on each day (i. e. lifespan data). For fecundity, McDonald et al. reported semi-daily fecundity data for entire population. Because the population was synchronized, and because mosquitoes lay discrete clutches of eggs separated by several days (the gonotrophic cycle duration), there were many data points when the populations did not produce offspring. These zero-inflated fecundity data are not ideal for fitting thermal responses. Therefore, after digitizing the semi-daily fecundity time series, we binned periods of several days (the bin size varied by temperature, since the gonotrophic cycle duration varies with temperature) and took a survival-weighted average within each bin (so days with more individual mothers contributing to offspring production counted more). To generate the dataset that we used to directly fit the thermal responses, we weighted the values within each bin by the mean number of surviving mothers in that bin. This approach allowed us to more accurately reflect daily fecundity averaged over a non-synchronized mosquito population. Note that the variation captured by these data and this approach is not variation between individual adult females, but rather variation by age for the whole population.
summary: Mosquitoes cannot control their body temperature, so their survival and performance depend on the temperature where they live. As a result, outside temperatures can also affect the spread of diseases transmitted by mosquitoes. This has left scientists wondering how climate change may affect the spread of mosquito-borne diseases. Predicting the effects of climate change on such diseases is tricky, because many interacting factors, including temperatures and rainfall, affect mosquito populations. Also, rising temperatures do not always have a positive effect on mosquitoes – they may help mosquitoes initially, but it can get too warm even for these animals. Climate change could affect the Ross River virus, the most common mosquito-borne disease in Australia. The virus infects 2,000 to 9,000 people each year and can cause long-term joint pain and disability. Currently, the virus spreads year-round in tropical, northern Australia and seasonally in temperate, southern Australia. Large outbreaks have occurred outside of Australia, and scientists are worried it could spread worldwide. Now, Shocket et al. have built a model that predicts how the spread of Ross River virus changes with temperature. Shocket et al. used data from laboratory experiments that measured mosquito and virus performance across a broad range of temperatures. The experiments showed that ~26°C (80°F) is the optimal temperature for mosquitoes to spread the Ross River virus. Temperatures below 17°C (63°F) and above 32°C (89°F) hamper the spread of the virus. These temperature ranges match the current disease patterns in Australia where human cases peak in March. This is two months after the country’s average temperature reaches the optimal level and about how long it takes mosquito populations to grow, infect people, and for symptoms to develop. Because northern Australia is already near the optimal temperature for mosquitos to spread the Ross River virus, any climate warming should decrease transmission there. But warming temperatures could increase the disease’s transmission in the southern part of the country, where most people live. The model Shocket et al. created may help the Australian government and mosquito control agencies better plan for the future.
keywords: ecology
  epidemiology and global health
title: Temperature explains broad patterns of Ross River virus transmission