The Relative Role of Climate Variation and Control Interventions on Malaria Elimination Efforts in El Oro, Ecuador: A Modeling Study

Malaria is a vector-borne disease of significant public health concern. Despite widespread success of many elimination initiatives, elimination efforts in some regions of the world have stalled. Barriers to malaria elimination include climate and land use changes, such as warming temperatures and urbanization, which can alter mosquito habitats. Socioeconomic factors, such as political instability and regional migration, also threaten elimination goals. This is particularly relevant in areas where local elimination has been achieved and consequently surveillance and control efforts are dwindling and are no longer a priority. Understanding how environmental change, impacts malaria elimination has important practical implications for vector control and disease surveillance strategies. It is important to consider climate change when monitoring the threat of malaria resurgence due to socioeconomic influences. However, there is limited assessment of how the combination of climate variation, interventions and socioeconomic pressures influence long-term trends in malaria transmission and elimination efforts. In this study, we used Bayesian hierarchical mixed models and malaria case data for a 29-year period to disentangle the impacts of climate variation and malaria control efforts on malaria risk in the Ecuadorian province of El Oro, which achieved local elimination in 2011. We found shifting patterns of malaria between rural and urban areas, with a relative increase of Plasmodium vivax in urbanized areas. Minimum temperature was an important driver of malaria seasonality and the association between warmer minimum temperatures and malaria incidence was greater for Plasmodium falciparum compared to P. vivax malaria. There was considerable heterogeneity in the impact of three chemical vector control measures on both P. falciparum and P. vivax malaria. We found statistically significant associations between two of the three measures [indoor residual spraying (IRS) and space spraying] and a reduction in malaria incidence, which varied between malaria type. We also found environmental suitability for malaria transmission is increasing in El Oro, which could limit future elimination efforts if malaria is allowed to re-establish. Our findings have important implications for understanding environmental obstacles to malaria elimination and highlights the importance of designing and sustaining elimination efforts in areas that remain vulnerable to resurgence.

Malaria is a vector-borne disease of significant public health concern. Despite widespread success of many elimination initiatives, elimination efforts in some regions of the world have stalled. Barriers to malaria elimination include climate and land use changes, such as warming temperatures and urbanization, which can alter mosquito habitats. Socioeconomic factors, such as political instability and regional migration, also threaten elimination goals. This is particularly relevant in areas where local elimination has been achieved and consequently surveillance and control efforts are dwindling and are no longer a priority. Understanding how environmental change, impacts malaria elimination has important practical implications for vector control and disease surveillance strategies. It is important to consider climate change when monitoring the threat of malaria resurgence due to socioeconomic influences. However, there is limited assessment of how the combination of climate variation, interventions and socioeconomic pressures influence long-term trends in malaria transmission and elimination efforts. In this study, we used Bayesian hierarchical mixed models and malaria case data for a 29-year period to disentangle the impacts of climate variation and malaria control efforts on malaria risk in the Ecuadorian province of El Oro, which achieved local elimination in 2011. We found shifting patterns of malaria between rural and urban areas, with a relative increase of Plasmodium vivax in urbanized areas. Minimum temperature was an important driver of malaria seasonality and the INTRODUCTION Malaria is the most important vector-borne disease worldwide, with 228 million cases reported in 2018 (WHO, 2019). Despite global elimination and eradication efforts, progress toward elimination is stalling in some endemic countries. In particular, increases in malaria incidence have been observed in some South American countries since 2014 including Venezuela, Peru, Colombia, and Ecuador, where movement of infected individuals between neighboring countries in the region can result in a resurgence of cases in local populations, threatening progress toward malaria reduction for other countries (Grillet et al., 2019;WHO, 2019). Malaria is highly sensitive to environmental conditions, including climate variability and land use practices (Moreno, 2006;Caminade et al., 2014), which along with socioeconomic influences can act as significant barriers to elimination. In addition, the lapse in control and surveillance efforts combined with global environmental changes also poses a threat to malaria elimination efforts in South America (Alimi et al., 2015;Conn et al., 2018).
Variation in meteorological conditions, particularly temperature, and rainfall, determines the spatiotemporal patterns of malaria through their effects on both the Plasmodium parasite and the Anopheles vector. Warmer temperatures can decrease the extrinsic incubation period (EIP), the time taken for the malaria parasite to complete its development inside the mosquito. A shorter EIP increases the transmission intensity of malaria by allowing mosquitoes to become infectious more quickly. The EIP for the Plasmodium falciparum parasite is longer (∼10 days) than that for Plasmodium vivax (∼8 days) and P. vivax can also develop at relatively cooler temperatures (minimum 15 • C), compared to P. falciparum (minimum 18 • C) (Patz and Olson, 2006). Temperature also impacts the mosquito vector, by affecting both larval and adult survival and longevity (Beck-Johnson et al., 2013), mosquito development (Bayoh and Lindsay, 2003), and vector population dynamics (Beck-Johnson et al., 2013). However, the effect of temperature on malaria transmission is non-linear, with the optimal temperature for malaria transmission estimated to be 25 • C with transmission declining above 30 • C (Paaijmans et al., 2009;Mordecai et al., 2013). Rainfall is also important for malaria transmission as it determines the availability of mosquito breeding habitat and thus mosquito abundance, although heavy and persistent rainfall can wash out larval habitats (Galardo et al., 2009;Olson et al., 2009;Prussing et al., 2019;Wolfarth-Couto et al., 2019).
Anthropogenic activities, such as urbanization and agricultural land use can also change malaria transmission patterns over longer time scales by creating more or less favorable microclimatic conditions and habitats for malaria vectors, depending on the species (Conn et al., 2002;Wayant et al., 2010;Tucker Lima et al., 2017). For example, urban development can eliminate mosquito habitat subsequently decreasing malaria incidence or in some cases increases in peri-urban malaria can occur as a result of economic migration and the creation of novel breeding habitats (Padilla et al., 2015). In addition, deforestation can facilitate malaria transmission due to the creation of new mosquito habitats on forest fringes (Vittor et al., 2006;Stefani et al., 2013). Increasing environmental suitability of malaria transmission, associated with climatic change and land use modifications, has important implications for malaria elimination efforts in South America, as these changes could compromise elimination efforts (Recht et al., 2017).
Previous studies have identified several environmental risk factors for malaria transmission that may impede elimination efforts in South America. For example, large-scale climate patterns such as the El Niño Southern Oscillation have been associated with large outbreaks of malaria in Peru, due to flooding conditions and in Colombia and Venezuela due to drought conditions (Gagnon et al., 2002). In highland areas of Colombia, warmer years were associated with an increased number of malaria cases occurring at higher altitudes (Siraj et al., 2014). Studies assessing spatial heterogeneity in malaria risk in South America have found deforestation and agricultural development have altered anopheline vector ecology in the Amazon region through creation of new Anopheles habitats and increased human-mosquito interactions (Moreno et al., 2007;Vittor et al., 2009;Hahn et al., 2014;Burkett-Cadena and Vittor, 2017).
When assessing the role of environmental factors on malaria incidence it is important to simultaneously consider the relative impact of control interventions, especially in elimination settings. Studies that evaluate the effectiveness of malaria control programs often consider environmental factors, most commonly rainfall and vegetation indices. Some studies used geostatistical modeling and after adjusting for climate conditions estimate the contribution and effectiveness of malaria control measures, including indoor residual spraying (IRS) in sub-Saharan Africa, relative to the impacts of climate variability (Graves et al., 2008;Giardina et al., 2014). Other studies often treat climate and environmental variables as confounders and do not apportion variation in malaria incidence due to climate (Aregawi et al., 2011;Bennett et al., 2014). There has been a limited assessment of how the combination of environmental factors, elimination efforts and socioeconomic pressures influences long-term trends in malaria and the route toward elimination. Understanding how malaria elimination efforts might be hindered by environmental changes has implications for targeting of vector control and disease surveillance, and is important especially given the possibility of shifting vector distributions and increasing climate suitability in South America (Pinault and Hunter, 2012).
According to the Pan-American Health Organization (PAHO), Ecuador is in the malaria pre-elimination phase, but has fallen short of its target to achieve elimination by 2020 (WHO, 2018). Since 2016 malaria has increased significantly in the country, mainly in the Amazon region, with more cases being detected than expected for a country on the verge of malaria elimination (PAHO, 2017;Sáenz et al., 2017). Recent surges of malaria in other South American countries, primarily Venezuela, are threatening current elimination efforts in Ecuador due to regional migration (Daniels, 2018;Grillet et al., 2019). This study focuses on Ecuador's southern province, El Oro, a high-risk region, which previously eliminated malaria but recently recorded indigenous cases, in addition to being vulnerable to the introduction of infections via imported cases due to high rates of human movement associated with its location on a strategic migration route (Jaramillo-Ochoa et al., 2019). Here we seek to identify threats to malaria elimination, with a view to prevent lapses in control efforts and the reestablishment of malaria in an area currently considered malaria-free but historically endemic for malaria. We advance the existing knowledge of the malaria elimination success in El Oro (Krisher et al., 2016) by quantifying the associations between malaria, climate and control interventions in a region vulnerable to resurgence.

Study Area
El Oro province (latitude: 3 • 5 45.20 S-4 • 11 3.06 S, longitude:-79 • 43 10.92 W-80 • 50 37.96 W) is located in southern Ecuador on the Pacific Coast and shares a border with the Tumbes region of northern Peru ( Figure 1A). El Oro covers 5,870 km 2 and is divided into 14 cantons, which range in size from approximately 70 km 2 up to 900 km 2 . Population densities per canton range from 7 people per km 2 in rural areas and up to 760 people per km 2 in the province capital, Machala.
Plasmodium vivax and P. falciparum are the malaria-causing parasites in the region (WHO, 2019) and the primary mosquito vectors transmitting malaria in El Oro are Anopheles albimanus and Anopheles punctimacula (Ryan et al., 2017).
During the mid-1990s public health authorities in El Oro and neighboring Tumbes unified to implement an effective binational collaboration to tackle the surge in malaria transmission in the region (Krisher et al., 2016). El Oro has been free of locally acquired malaria infections since 2011 although malaria has been increasing elsewhere in Ecuador (WHO, 2019). In 2018 seven malaria cases were recorded in El Oro (6/7 were imported cases) during a period of elevated migration of Venezuelan citizens associated with the social and economic crisis in their country (Jaramillo-Ochoa et al., 2019). El Oro's border location, along the Pan American highway, a route used by many migrants, makes it particularly vulnerable to malaria resurgence. In addition, El Oro is one of the most hazardous coastal zones in Latin America and the Caribbean (Calil et al., 2017). Vulnerable populations are susceptible to the effects of the El Niño Southern Oscillation, which intensifies annual flooding events during the rainy season (Lowe et al., 2017;Tauzer et al., 2019). Lapses in surveillance and control efforts and reductions in funding following local elimination also means malaria resurgence is likely (PAHO, 2017;Sáenz et al., 2017).

Data Sources
Monthly counts of blood smear confirmed cases of P. falciparum and P. vivax malaria for each canton in El Oro, January 1990-December 2018, were provided by the Ecuadorian Ministry of Health, where malaria is a mandatory notifiable disease. Cases were recorded at local clinics across El Oro and collated by the Ministry of Health. Population data, available for each canton, were sourced from the national census in Ecuador, from the Instituto Nacional de Estadística y Censos (INEC 1 ) for 1990, 2001. Data for 2011 were provided by INEC as annual population projections. Population values between census years 1990, 2001, and 2010 were estimated by interpolating assuming linear growth, to provide annual population estimates for each canton. To estimate the proportion of the population in poverty per canton in El Oro, we used Unmet Basic Needs (UBN) sourced from the 2010 census, an indicator based on measures including housing quality, education and access to water and sanitation. We used UBN data from 2010 as a socioeconomic covariate in all our models.
Between 2001 and 2015 an intensive period of vector control was carried out across El Oro following a resurgence of malaria cases and monthly canton-level data for three control measures were available from the Ministry of Health. Monthly estimates per canton of the number of households treated by IRS with insecticides (deltamethrin 5% concentrated suspension, deltamethrin 2.5%, malathion 50%, alphacypermethrin 10% concentrated suspension, and betacypermethrin 2.5%) were available for January 2001-September 2013. The number of neighborhoods that were treated with insecticide via ultra-low volume (ULV) fumigation, which is performed by spraying entire neighborhoods with 96% malathion from trucks were available between January 2004 and May 2015. Finally, the number of households space-sprayed with 2.5% deltamethrin concentrated emulsion, using a backpack fogger that creates a fog insecticide to treat both inside and outside the home, were available between January 2004 and May 2015 (Supplementary Figure S1). For each control measure we tested for lagged relationships, up to 3 months, in order to account for delays between implementation and impact on measurable malaria cases (Supplementary Table S1).
Other malaria control interventions including fumigation with DDT (until banned in 1996), elimination of larval habitats, and bed net provision occurred during the study period, but no detailed data for control measures were available prior to 2001.
Monthly climate data, for temperature and precipitation, between January 1990 and December 2018 were sourced from TerraClimate, a high-spatial resolution (∼4 km) dataset that uses climatically aided interpolation to produce monthly estimates (Abatzoglou et al., 2018). To obtain climate estimates for each canton in El Oro, TerraClimate variables were aggregated by taking mean grid cell values using the raster package in R.
Annual land cover maps at a 300 m spatial resolution generated by the European Space Agency (ESA) Climate Change Initiative (CCI) were compiled for 1990-2018 2 (European Space Agency, 2019). For each annual map, the number of grid cells in each canton that were classified in the urban category were extracted, and from this the proportion of urban grid cells per canton were calculated to give a continuous variable of urbanization (Supplementary Figure S3). All continuous variables were scaled by subtracting the variable mean from each value and dividing this by the standard deviation. This transformation allowed for the direct comparison of continuous variables and enabled their relative importance to be determined in the model.

Statistical Analysis
Bayesian hierarchical mixed models are powerful statistical models that can be used to estimate the marginal posterior distributions of multiple covariates whilst simultaneously accounting for multiple sources of uncertainty that can arise from disease count data. Here, spatiotemporal models were constructed, for each malaria parasite, P. falciparum and P. vivax, to explore the contributions of climate and control interventions to variations in monthly malaria incidence in El Oro between January 1990 and December 2018 (348 months).
We fitted zero-inflated negative binomial models that take into account excess zeros that cannot be explained by a standard negative binomial distribution (Supplementary Appendix S1.1). Where µ st is the mean number of malaria cases in each canton (s = 1,. . .,14) for each timestep (t = 1,. . .,348), the hierarchical spatiotemporal model is defined as: Annual population estimates per 1,000 inhabitants for each canton in El Oro, log(P st ), were included as an offset in the model to account for canton-level differences in the annual parasite incidence (API). The estimated API, log(ρ st ) is then made up of a combination of climate covariates, including mean monthly minimum temperature and precipitation for each canton. Other explanatory variables considered included socioeconomic factors (poverty rates), level of urbanization and vector control interventions (IRS, ULV fumigation and space-spraying).
To account for the seasonality in malaria incidence, a monthly random effect was introduced using a first order autoregressive prior (m t ), which allows malaria in 1 month to depend on incidence in the previous month. Independent random effects for each year (y t ), 1990-2018, were included to allow for additional sources of variation due to unobserved confounding factors such as variations in healthcare access, case reporting and diagnostics.
Spatially structured random effects, υ s , were introduced into the model to allow for correlated heterogeneity in malaria incidence across cantons in El Oro. This spatial dependency structure was accounted for by assuming a conditional intrinsic Gaussian autoregressive (CAR) model prior distribution for the spatial effects, which takes into account the neighborhood structure of the area (Besag et al., 1991). We also included spatially unstructured random effects in the model to allow for additional uncorrelated spatial variation across cantons, ν s , which were assigned independent diffuse Gaussian exchangeable prior distributions (Lowe et al., 2016). Separate models were constructed for both P. falciparum and P. vivax malaria due to their intrinsic differences, which include diagnostic potential, development time inside the mosquito vector (extrinsic incubation period, EIP), the differential impact of vector control on transmission, infection reservoirs of P. vivax hypnozoites and the presence of asymptomatic cases.
Detailed data on interventions was only available for the period January 2001-December 2015. Therefore, we formulated separate sub-models for each malaria type to investigate the relative impact of the control measures implemented between 2001 and 2015 (see Supplementary Appendix S2), herein referred to as "intervention models." The intervention models, fitted to data for the period 2001-2015 included the same explanatory variables and random effects as the "full" models (fitted to data for the period 1990-2018) with the addition of the three control measures. We examined differences between the random effects structures for the full models and the intervention models, to see how much variation in malaria due to interventions could be accounted for by random effects in the absence of intervention data for the entire time period.
Posterior distributions of model parameters were estimated using Integrated Nested Laplace approximations (INLA) (Rue et al., 2009). INLA provides a computationally more efficient alternative to Markov Chain Monte Carlo (MCMC) methods, by using numerical approximations of model parameters (see Supplementary Appendix S1.2). Covariate time lags and the most parsimonious models were selected using the deviance information criterion (DIC) (Spiegelhalter et al., 2002) and the Watanabe-Akaike information criterion (WAIC) (Watanabe, 2010), which are Bayesian methods of model comparison that trade off model adequacy against model complexity. Models with lower DIC and WAIC indicate a more parsimonious model. We also used the logarithmic score to assess model fit, which is based on the conditional predictive ordinate (CPO) leave-one-out cross-validation score, where a smaller value indicates a greater predictive power of the model (Gneiting and Raftery, 2007).
Covariates were added iteratively to the models starting with a baseline model, which included spatial and temporal structured and unstructured random effects and retained if model fit was improved, assessed through a decrease in DIC and WAIC. The most appropriate temperature variable, maximum or minimum temperature was selected using DIC and WAIC (Supplementary Table S2). We tested different time lags (from 0 to 3 months) to take into account the lagged effect of climate factors on malaria transmission, including development time of mosquitoes and parasites as well as the time between malaria diagnosis and case recording (Supplementary Table S2). We also tested non-linear functions of temperature and precipitation (Supplementary Table S3). We used the root mean square error (RMSE) to assess the extent to which models of malaria incidence for 2001-2015 were improved by the inclusion of each control measure in turn. RMSE is a measure of the standard deviation of the model residuals, with smaller values indicating better model fit to the observed data.

Urban and Rural Malaria Incidence in El Oro
Between 1990 and 2018, 62,120 cases of malaria were recorded in El Oro, with 54% of infections resulting from P. falciparum malaria and 46% from P. vivax. Malaria transmission in El Oro was historically high, reaching up to 2,000 total annual cases in 1999 and an API of 61 cases of P. falciparum and 35 P. vivax cases per 1,000 people recorded in Huaquillas canton along the Peruvian border. Between 1998 and 2002 large outbreaks, particularly of P. falciparum malaria occurred across the province before declining to low levels of incidence after the period of more intensive vector control ( Figure 1B). A small outbreak of P. vivax malaria occurred between 2007 and 2009 before malaria incidence was dramatically reduced and remained at a low and stable level ( Figure 1B). The highest malaria incidence was concentrated in the western part of the province, along the Peru-Ecuador border (Supplementary Figure S4).
Between 1990 and 2018, outbreaks of both P. falciparum and P. vivax malaria occurred roughly at the same time. Prior to the decline in cases in 2001-2015 there was a greater difference between P. falciparum and P. vivax incidence in rural areas, compared to urbanized areas (Supplementary Figure S5). During the elimination phase (2001-2010) the incidence of P. falciparum and P. vivax malaria in rural areas was reduced much more than malaria in urbanized areas, where outbreaks of P. vivax were still observed. There was also a greater overall reduction of P. falciparum malaria in rural and urbanized areas.
In the models of malaria incidence in El Oro between 1990 and 2018, we found that areas that were more urbanized were associated with greater incidences of P. falciparum malaria (Figure 2A). To examine how the epidemiology of malaria in El Oro had changed with the three vector control interventions, we tested whether there was a difference in the relationship between malaria and level of urbanization in El Oro, before and after the interventions were implemented. We found a statistically significant interaction between level of urbanization in El Oro and the period the control measures were implemented (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015) for P. vivax malaria and no significant interaction for P. falciparum malaria (Supplementary Figure S6). Prior to 2001, more urbanized areas of El Oro were associated with decreased P. vivax malaria and after 2001 this relationship reversed with urbanized areas associated with more P. vivax malaria.
During the period of malaria decline (post 2001) there was also a reduction in the seasonal pattern of both P. falciparum and P. vivax malaria in El Oro especially for P. vivax malaria incidence, which ceased to peak between June and July (Supplementary Figure S7). For P. falciparum a seasonal peak in incidence during this period was evident between June and July, but less distinct than before 2001.

Climate Variability and Malaria Incidence in El Oro 1990-2018
Minimum temperature was the best temperature variable in the full models (Supplementary Table S2) and in the 3 months prior to case reporting was an important predictor of the spatiotemporal distribution of P. falciparum and P. vivax malaria in El Oro (Figure 2A). Warmer temperatures in El Oro between 1990 and 2018 were associated with increases in malaria incidence and greater increases in P. falciparum malaria than P. vivax malaria (Figure 2A). Precipitation lagged by 3 months for P. falciparum and 1 month for P. vivax were selected as the best time lags although when included in the full model, precipitation was not significant in driving malaria incidence in El Oro (Figure 2A). Introducing a non-linear relationship between malaria and climatic variables improved model fit by decreasing model DIC and RMSE for P. vivax but not P. falciparum malaria incidence ( Table 1 and Supplementary Figure S8). We also found that higher temperatures were associated with larger increases in malaria incidence ( Figure 2B).
We also investigated the influence of minimum temperature on the unexplained variation in P. falciparum and P. vivax incidence in El Oro between 1990 and 2018. We compared the monthly and interannual random effects of models with and without temperature to determine the extent to which temperature accounted for the seasonal pattern of malaria incidence. The monthly random effects of the model for P. falciparum incidence that included minimum temperature were near zero, showing that temperature accounted for all of the seasonal variation in P. falciparum malaria in El Oro 1990-2018 ( Figure 3A). In contrast there was minimal reduction in the random effects for the P. vivax model with minimum temperature, indicating the seasonal pattern was driven by other factors (Figure 3B). The proportion of the variation accounted for by the yearly random effects in a model without minimum temperature was greater than the proportion of the variation accounted for by the monthly random effects in a model without minimum temperature. This suggests that temperature is the main driver of malaria seasonality and other interannual signatures as well as temperature contribute to the interannual variability of malaria in El Oro. For example, for some years the proportion of variation accounted for by the yearly random effects decreased when minimum temperature was included into models of malaria incidence 1990-2018. A reduction in the value of the random effects can be seen for 1998, when a strong El Niño event occurred (Supplementary Figure S9). Reductions were also particularly apparent for P. falciparum malaria between 2001-2004 and 2009-2011, suggesting that some additional variation was explained by minimum temperature. In contrast for some years there was no reduction for example in 2008, when an outbreak of P. vivax malaria occurred suggesting the outbreak was not driven by a climate event.
We explored how climate suitability for malaria transmission had changed over the time period 1990-2018 in El Oro, Posterior mean and 95% credible intervals for minimum temperature and precipitation, included as linear terms, level of urbanization and poverty covariates for P. falciparum (pink) and P. vivax (blue) malaria. (B) Relationships between P. falciparum and P. vivax relative risk, on the log scale defined as the annual parasite incidence (API), log(ρ st ), and minimum temperature, lagged by 3 months, included in the model as a function to allow for non-linearities.  Covariates were added iteratively, starting with a baseline model that included structured and unstructured spatial random effects, and a seasonal term. The level of urbanization was interacted with a categorical variable indicating a decline in malaria incidence post 2001. The most parsimonious models included climate (temperature and precipitation) as linear terms for P. falciparum and as non-linear functions for P. vivax malaria.
specifically for minimum temperature since it was an important predictor. To do this, we calculated the number of months that had suitable temperature conditions for transmission of each malaria parasite. We used the lower temperature limits of the range for which the development of the malaria parasites inside the mosquito vector is considered most suitable for transmission, as an indicator of this suitability. For P. falciparum a lower temperature limit of 18 • C was used and for P. vivax 15 • C (Gilles, 1999;Watts et al., 2019). Between 1990 and 2018, the number of months with FIGURE 3 | Effect of minimum temperature (T min ) on the annual cycle of malaria in El Oro 1990-2018. Difference in the monthly random effect marginal posterior distributions for models of (A) P. falciparum and (B) P. vivax malaria that include minimum temperature (orange), lagged by 3 months and exclude minimum temperature (gray). Relative risk, on the log scale, is defined as the annual parasite incidence (API), log(ρ st ).

FIGURE 4 | Suitable temperature conditions for malaria transmission in El
Oro 1990-2018. Mean (solid curve) number of months per year where minimum temperature exceeds 18 • C (i.e., considered suitable for P. falciparum malaria transmission; pink) and where minimum temperature exceeds 15 • C (i.e., considered suitable for P. vivax malaria transmission; blue), logistic regression line (dashed curve) and 95% confidence intervals (gray shading).
suitable temperatures for P. falciparum transmission showed an increasing trend (Figure 4), with noticeable peaks of around 8 months of suitable temperature conditions occurring between 1997 and 1998 and between 2014 and 2016. This trend of increasing temperatures for P. falciparum transmission coincides with warmer temperatures throughout the year in El Oro and increasing temperatures between 1990 and 2018 (Supplementary Figure S10). In contrast, temperature conditions for P. vivax malaria transmission showed no increasing trend in El Oro between 1990 and 2018 (Figure 4), although similar peaks in suitability were observed in 1997-1998 and in 2015-2016.
The most parsimonious models of malaria in El Oro 1990-2018 included minimum temperature and precipitation as FIGURE 5 | Model posterior distributions with and without climate information for P. falciparum and P. vivax malaria in El Oro 1990-2018. Observed (gray solid line), posterior mean (blue dashed line) and 95% credible intervals (blue shading) for annual parasite incidence (API) for (A) models that include minimum temperature, lagged by 3 months and precipitation, lagged by 3 months for P. falciparum models and 1 month for P. vivax models and (B) without climate information. For models of P. falciparum malaria minimum temperature and precipitation were included as linear terms. For P. vivax malaria models, model fit improved when including non-linear functions of minimum temperature and precipitation.
linear terms for P. falciparum and as non-linear functions for P. vivax malaria, poverty rates, level of urbanization including an interaction with the period when the vector control measures were implemented, and spatial and temporal random effects ( Table 1). Although in these models the credible intervals for precipitation and poverty rates contained zero, the addition of these covariates to the models decreased DIC values and increased model fit for P. falciparum and P. vivax malaria incidence in El Oro between 1990 and 2018. Accounting for poverty rates (socioeconomic effects) across El Oro also decreased DIC and WAIC for models of P. falciparum and P. vivax malaria. Including meteorological variables (as linear terms for P. falciparum and as non-linear functions for P. vivax malaria) decreased the uncertainty of the model posterior distributions, in comparison to distributions from models that excluded meteorological variables (Figure 5). In particular, model uncertainty was reduced more for P. falciparum malaria, especially during the large outbreaks that occurred in 1992-1994 and 1998-2002. Interestingly during the later years of the study period, model posterior distributions show an upsurge in malaria incidence, particularly for P. vivax incidence between 2015 and 2018 although this upsurge was not observed in case reports.

Control Interventions and Malaria Incidence in El Oro 2001-2015
We fitted separate sub-models for incidence of P. falciparum and P. vivax malaria in El Oro between 2001 and 2015, to evaluate the impact of the three vector control measures (IRS, space spraying and ULV fumigation) that were implemented during this period. These intervention models included the same covariates as the full models (linear and non-linear climate information and socioeconomic effects) and random effects but with the addition of the intervention data that were available for 2001-2015.
By comparing the differences in the covariate parameter estimates from the full models for the entire time series 1990-2018 and the intervention models for the period 2001-2015, we found that posterior mean estimates were roughly similar across the different models (Figure 6 and Supplementary Tables S4, S5). However, there was more uncertainty in the FIGURE 6 | Posterior mean estimates for environmental, socioeconomic and malaria vector control covariates for models of (A) P. falciparum and (B) P. vivax malaria from full models (dark blue) and intervention models (green). All interventions were lagged by 3 months, apart from space spraying in the P. falciparum model, which was lagged by 2 months.
estimates for the environmental covariates in the intervention models compared to the 1990-2018 models, particularly for minimum temperature. We also found that between 2001 and 2015 areas in El Oro with greater levels of poverty were associated with higher incidences of P. vivax malaria. In the intervention models, IRS implemented in El Oro 2001-2015, 3 months prior to case detection was associated with decreased incidence of P. falciparum malaria but was not associated with significant decreases in incidence of P. vivax malaria. Space spraying was associated with decreased P. vivax malaria incidence in the following 3 months but not P. falciparum malaria in the following 2 months. ULV fumigation, also lagged by 3 months was not statistically significant (i.e., the credible intervals contained zero).
We also assessed which vector control measure implemented in El Oro between 2001 and 2015 provided more valuable information to the intervention models in explaining malaria incidence across the province. We found there was considerable heterogeneity in the model improvement for each control measure, as measured by RMSE difference (Figure 7). IRS improved the model for P. falciparum malaria in the coastal northwest of El Oro and in cantons along the Ecuador-Peru border, improving the model by up to 14%. In contrast, there was little improvement in model fit for P. vivax. Space spraying improved models of both P. falciparum and P. vivax in the province capital, Machala, with a greater improvement for P. vivax. ULV fumigation provided minimal or no model improvement for either parasite.
Finally, we explored if the variation in malaria incidence in El Oro attributed to the vector control measures in the intervention model could be captured in the random effects structure of the full models. We compared the yearly random effects of the full models (fitted to data from 1990 to 2018), which did not include intervention data, to the yearly random effects of the intervention models (fitted to data from 2001 to 2015. We found that for some years (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008) there was a reduction in the magnitude of the random effects for the 2001-2015 models, suggesting the control measures accounted for some of the unexplained variation in malaria incidence in the model for 1990-2018 during those years. However, a considerable portion of the random variation was likely due to other unmeasured factors. Between 2012 and 2015 the random effects of the P. falciparum models are close to zero, indicating that variation in malaria incidence for those years was well captured by factors included in the model (Figure 8). We found a similar pattern for P. vivax malaria, with a reduction in the value of the random effects for the intervention model, compared to the full model during the early years of the interventions 2001-2009. This reduction indicates that the control measures influenced the interannual variation of malaria risk during these years.

DISCUSSION
For regions approaching malaria elimination and where funding and disease surveillance are limited, it is important to understand the drivers of malaria incidence in order to prevent disease reestablishment and to ensure that elimination efforts are sustained. In this study, we developed a statistical modeling framework to disentangle the relative role of multiple explanatory variables in driving variation in malaria incidence, while accounting for unobserved heterogeneity. We found that minimum temperature  Distributions from the full model are shown in gray, which does not include intervention data. Relative risk, on the log scale, is defined as the annual parasite incidence (API), log(ρ st ).
was an important driver of malaria incidence in El Oro. We also found that control measures had a differential impact on the two types of malaria. After the intensive period of vector control in El Oro between 2001 and 2015 malaria incidence in rural areas declined more than malaria in urbanized areas and P. vivax malaria became more dominant. Incidence of P. falciparum malaria was also reduced more than P. vivax malaria and relatively more P. vivax cases were reported in urbanized areas during this time period. The greater reduction in P. falciparum malaria is expected given that malaria control measures deployed during this time, are primarily designed for use in rural settings and target mosquito vectors that are exophagic such as An. albimanus, a dominant malaria vector in El Oro.
Indoor residual spraying has largely been shown to be effective in reducing malaria transmission in settings with high prevalence and when applied on a large scale, especially in African countries (Pluess et al., 2010). In this study we found that IRS was effective at reducing the incidence of P. falciparum malaria in El Oro between 2001 and 2015, but not P. vivax malaria. Conventional control methods are less effective in reducing transmission of P. vivax compared to P. falciparum malaria because of the ability of P. vivax to cause multiple relapsing malaria episodes after the initial infection, from the activation of dormant liver stages (hypnozoites). In addition the P. vivax parasite develops more rapidly in the mosquito than P. falciparum and so control methods such as IRS deployed in El Oro, which aim to shorten the mosquito lifespan are less effective (Mendis et al., 2001;Bassat et al., 2016). The rapid development of P. vivax and hypnozoite stage also means that drug treatments, which are often of a longer duration are less effective in reducing transmission as mosquitoes can become infected during a pre-symptomatic period, allowing for onward transmission before drug treatment is initiated (White, 2008;McCarthy et al., 2013).
We found that IRS only reduced the incidence of P. falciparum malaria and in the west of the province. IRS is most effective against endophilic mosquitoes, which rest indoors following a blood meal. The two main malaria in vectors in El Oro, A. albimanus and A. punctimacula, are most commonly observed biting outdoors and mainly rest outdoors (Ryan et al., 2017). As a result, IRS and other conventional control methods may not be suitable for reducing malaria transmitted by these mosquitoes.
In contrast we found that space spraying was effective in reducing incidence of P. vivax in El Oro, especially in the capital, Machala and by a larger amount than IRS. Evidence for the success of space spraying is limited and it is currently only recommended during outbreaks. In our study, spacespraying is possibly targeting outdoor resting mosquitoes, such as A. albimanus, especially in peri-urban areas where housing conditions are poor and outdoor exposure is higher. There is evidence of a positive effect of space-spraying on malaria incidence in India, but the current evidence base is limited (Pryce et al., 2018).
The P. vivax malaria parasite is able to develop at lower temperatures than P. falciparum (Nikolaev, 1935;Moshkovsky, 1946;Olliaro et al., 2016;Ohm et al., 2018), enabling its persistence in less favorable environmental conditions (Mendis et al., 2001). We found a stronger association between minimum temperature and P. falciparum malaria than P. vivax, and all of the seasonality in P. falciparum malaria was explained by temperature. A greater association between minimum temperature and P. falciparum malaria compared to P. vivax has been previously identified in China (Bi et al., 2013), but to our knowledge this difference has not been previously quantified in South America. The stronger effect of temperature on P. falciparum malaria is due to the greater dependence of P. falciparum transmission on the mosquito vector, which is sensitive to climate conditions such as temperature. In contrast P. vivax transmission is less influenced by the mosquito vector and can be sustained by unpredictable relapsing infections, which have been suggested to be explained by other factors such as systemic illness in humans (White, 2011). Similarly, between 2001 and 2015 we found that areas in El Oro with higher levels of poverty were associated with increased incidences of P. vivax malaria, which suggests that socioeconomic conditions, such as limited access to healthcare, may have influenced relapsing P. vivax infections during this time.
Warmer minimum temperatures in El Oro were associated with increases in P. falciparum and P. vivax malaria incidence 3 months later. This finding is in agreement with other studies in South America that have found a similar relationship between higher temperatures and malaria incidence (Poveda et al., 2000;Basurko et al., 2011;Laneri et al., 2019). This relationship can be attributed to the physiological effects of temperature on both the mosquito vector and parasite that increases transmission. For example warmer temperatures shorten the development time of the parasite inside the mosquito, the EIP and increase mosquito larval reproduction thus increasing disease transmission (Bayoh and Lindsay, 2003;Blanford et al., 2013;Mordecai et al., 2013). We also found a trend of rising minimum temperatures in El Oro between 1990 and 2018, with an increase in the number of months that have suitable minimum temperature conditions for P. falciparum malaria transmission, which suggests a lengthening of the malaria transmission season. Large peaks in the transmission suitability of both P. falciparum and P. vivax malaria are evident between 1997-1998 and 2014-2016, when major El Niño events occurred. El Niño conditions bring warmer and wetter conditions to southern Ecuador, conditions that are favorable for malaria transmission and have previously been suggested to have caused a peak in cases observed during this time (Krisher et al., 2016). Increasingly suitable temperature conditions and large-scale climate events such as El Niño pose a real threat for the re-establishment of malaria in El Oro if cases are allowed to return, surveillance is not maintained, and interventions are not deployed appropriately. Our model posterior distributions for 2012-2018 show increases in malaria incidence despite zero or few cases being reported, which is possibly due to warmer temperatures and the large El Niño event that occurred during this time. This highlights how incidence could rise if risks, such as relaxing control efforts are not mitigated.
Rainfall is essential for providing suitable habitats for mosquito breeding (Thomson et al., 2005;Parham and Michael, 2010) and is thus considered to be a dominant factor in driving malaria transmission. Studies in the Amazon and Argentina have found that rainfall was important in determining malaria incidence (Dantur Juri et al., 2009;Olson et al., 2009). In this study, rainfall was not a statistically significant explanatory variable for malaria incidence in El Oro, possibly because ample larval habitat was already available. Two to three months prior to the peak malaria season (March-June), monthly rainfall accumulation reaches 80 mm across El Oro, which is considered to be suitable conditions for malaria transmission (Watts et al., 2019), and up to 530 mm rainfall in the southeast of the province. In addition, much of El Oro is rural with extensive mangroves in the northwest coastal regions, which provide habitat suited to A. albimanus mosquitoes (Pinault and Hunter, 2012). Therefore, our results suggest that rainfall is not a limiting factor for malaria transmission in El Oro.
The addition of temporal random effects for each year in the full models for 1990-2018 accounted for some of the variation in malaria incidence due to the control measures that were explicitly accounted for in the intervention model for 2001-2015, particularly for P. falciparum malaria, which appeared more sensitive to interventions. The use of temporal random effects also allowed us to pick up variation in malaria in El Oro due to other unobserved factors. For example, the large El Niño events mentioned previously likely caused an increased number of cases and a reduction in the unexplained yearly variation during 1997-1998 and between 2014 and 2016. However, we expect that the temperature and precipitation variables included in the models to pick up some variation to due climatic anomalies from El Niño events. The Cenepa War, a period of political instability across the Ecuador-Peru border (Krisher et al., 2016) likely hindered malaria control efforts and case reporting during 1995, and there is evidence of reduced malaria risk during this time in the interannual random effects, which is likely explained by reduced reporting. In addition we found that the outbreak of P. vivax in 2008 was not driven by a climatic event and may instead be driven by a political event or lapse in control efforts.
Despite the considerable length of the dataset, there are some limitations to consider. In Ecuador a high number of malaria cases are asymptomatic. Low levels of parasitemia, which act as a reservoir of transmission, are difficult to diagnose (Sáenz et al., 2017). It is likely that many cases of P. vivax malaria were missed and not reported during the study period. In addition, many P. vivax cases reported between 1990 and 2018 are likely to be relapses from the same initial infection, which could mask the true climatemalaria relationship. We also found that urban malaria was approximately 70% higher than incidence of malaria in rural areas. However, this may be a product of migration from neighboring countries and cases being imported from rural areas and subsequently recorded as originating in urban areas, as was found in Colombia (Molina Gómez et al., 2017). We could not quantify the amount of malaria variation explained by migration owing to the lack of detailed migration data for the whole time period. Due to El Oro's location on a key migratory route, including for Venezuelan refugees, it is likely that human movements are influential in driving case reporting and malaria incidence. In addition, outbreaks in northern Peru can increase the risk of malaria transmission in El Oro due to cross border human movement. Other factors that likely contributed to the spatiotemporal variation in malaria incidence in El Oro that were accounted for through the use of random effects include variations in case reporting and surveillance, as well as health seeking behavior in El Oro.
The incompleteness of the intervention data for the study period is also a limiting factor in this study. We were only able to evaluate the role of the three control measures implemented between 2001 and 2015 although there were other important elimination efforts that were carried out in El Oro between 1990 and 2018. For example, the use of DDT was widely used until 1996 and a campaign to distribute insecticide-treated bed nets (ITNs) was initiated in 2004, although no detailed data exist for this in El Oro, or for other forms of personal protection that would have contributed to local elimination. However, by examining the interannual random effects of our models we observe a reduction in P. falciparum malaria risk between 2004 and 2006, which could be attributed to the use of ITNs. In 2005, the Ministry of Health suspended the use of chloroquine for treatment of P. falciparum malaria, adopting the recommended change to artemisinin-based combination therapy (ACT) amid reports of drug resistance, which local partners reported led to a major decline in transmission (Krisher et al., 2016). The introduction of ACT in Colombia and Peru at similar times contributed to significant reductions in malaria case numbers (Rodríguez et al., 2011;Quispe et al., 2016). ACT introduction is also likely to have contributed to local malaria elimination El Oro as we observed declining malaria risk from this time through the use of interannual random effects in our models. At the same time in El Oro, a new shorter treatment regime for P. vivax malaria was adopted to ensure patients completed their treatment course, which was also reported to have contributed to declining transmission (Krisher et al., 2016).
In summary, we have developed a modeling framework for exploring the barriers to malaria elimination in a lowtransmission setting. Accounting for unobserved interannual variation via yearly random effects is useful for investigating malaria variation in many different contexts where detailed data, such as intervention information, may not be available. Here, we used available intervention data to explore the relative impact of three vector control measures on P. falciparum and P. vivax malaria in El Oro. We also provide an assessment of where control measures improved models of P. falciparum and P. vivax incidence, which will prove useful for targeting future malaria control efforts in the region. We have shown that warmer temperatures increase incidence of P. falciparum malaria, which is important to consider in light of global environmental change and increasing climate suitability for malaria transmission (Laporta et al., 2015;Watts et al., 2019). In addition, with a global warming of 1.5 • C above pre-industrial levels it is predicted that the western coast of South America, including southern Ecuador, will experience temperatures between 2 and 3 • C warmer and precipitation could increase by up to 20% (IPCC, 2018), providing climate conditions that are highly suitable for malaria transmission. In contrast, we show that P. vivax malaria is less sensitive to temperature variations and shows more complex transmission patterns, making it particularly challenging to eliminate (Feachem et al., 2010). An assessment of the environmental obstacles to elimination efforts in Ecuador is timely and important. Lapses in control efforts should be avoided, especially as it becomes harder to prevent and detect cases due to limited funding for surveillance and control. Recently, indigenous cases of malaria have been detected in El Oro and re-introduction of malaria parasites along with warming temperatures threaten current elimination progress.

DATA AVAILABILITY STATEMENT
The data and code used to produce the analysis are available on GitHub via the following link: https://github.com/isabelfletcher/ malaria_eloro.

AUTHOR CONTRIBUTIONS
IF and RL conceived the study. IF collated the data, fitted the models, analyzed the results, and wrote the manuscript. AS-I helped conceive the study, provided information on the vector control measures, and contributed to discussion of the study results. RS provided insight on vector-borne diseases and surveillance in El Oro. GC-E helped with the environmental data sources and extraction. MS, EB-A, TO, and JA provided the malaria case data and intervention data for analysis. FS provided insight on the malaria situation in Ecuador, including treatments, and diagnostics. KJ provided input on model formulation. KJ and CD provided input on interpretation of model results. RL and AS-I edited the manuscript. All authors critically reviewed and approved the final version of the manuscript. Special recognition of MS, posthumously (1959MS, posthumously ( -2020, a leader in the field of medical entomology and tireless advocate for community public health in El Oro. All authors contributed to the article and approved the submitted version.

FUNDING
IF was supported by the Biotechnology and Biological Sciences Research Council (grant BB/M009513/1). RL was supported by a Royal Society Dorothy Hodgkin Fellowship. RL and GC-E were supported by the National Socio-Environmental Synthesis Center (SESYNC) under funding from the National Science Foundation (DBI-1639145).