Predicting Abundances of Aedes mcintoshi, a primary Rift Valley fever virus mosquito vector

Summary

Rift Valley fever virus (RVFV) is a mosquito-borne zoonotic arbovirus with important livestock and human health, and economic consequences across Africa and the Arabian Peninsula. Climate and vegetation monitoring guide RVFV forecasting models and early warning systems; however, these approaches make monthly predictions and a need exists to predict primary vector abundances at finer temporal scales. In Kenya, an important primary RVFV vector is the mosquito Aedes mcintoshi. We used a zero-inflated negative binomial regression and multimodel averaging approach with georeferenced Ae. mcintoshi mosquito counts and remotely sensed climate and topographic variables to predict where and when abundances would be high in Kenya and western Somalia. The data supported a positive effect on abundance of minimum wetness index values within 500 m of a sampling site, cumulative precipitation values 0 to 14 days prior to sampling, and elevated land surface temperature values ~3 weeks prior to sampling. The probability of structural zero counts of mosquitoes increased as percentage clay in the soil decreased. Weekly retrospective predictions for unsampled locations across the study area between 1 September and 25 January from 2002 to 2016 predicted high abundances prior to RVFV outbreaks in multiple foci during the 2006–2007 epizootic, except for two districts in Kenya. Additionally, model predictions supported the possibility of high Ae. mcintoshi abundances in Somalia, independent of Kenya. Model-predicted abundances were low during the 2015–2016 period when documented outbreaks did not occur, although several surveillance systems issued warnings. Model predictions prior to the 2018 RVFV outbreak indicated elevated abundances in Wajir County, Kenya, along the border with Somalia, but RVFV activity occurred west of the focus of predicted high Ae. mcintoshi abundances.

Keywords:

Kenya – Livestock – Epizootics – Mosquitoes – Surface temperature – Viral vectors – Rift Valley fever virus

Introduction

Rift Valley fever virus (RVFV) is a mosquito-borne zoonotic disease of great economic, livestock, and human health importance in Africa and the Arabian Peninsula [1, 2]. Transmission of RVFV to humans generally involves direct contact with infected tissues or body fluids of animals or bites of infected mosquitoes [3]. Human illness from RVFV often goes unnoticed, or results in flu-like symptoms, but a more severe form of the virus may present, resulting in ocular disease, meningoencephalitis, or hemorrhagic fever, the latter with a case-fatality rate of 50% [4]. Livestock infected with RVFV are less likely to be asymptomatic. High numbers of simultaneous, spontaneous abortions among ruminants (so-called “abortion storms”) and high mortality rates among young animals accompany epizootics [5]. Effects of epizootics on domestic livestock herds are devastating and result in tremendous economic losses due to imposed quarantines and embargo, and food insecurity for communities whose livelihoods depend on livestock [6]. RVFV infection based on serology has also been detected in a wide variety of wild ruminants, from African buffalo to giraffes, but without the pronounced symptoms displayed in livestock [7].

Current RVFV forecasting models use persistence of above-average rainfall, positive Normalized Difference Vegetation Index anomalies, elevated cloud coverage measurements, and the occurrence of El Niño conditions to guide early warning systems [8–13]. Although this approach has seen some success in predicting RVFV risk to general locations and time periods [8, 9, 14–18], model outputs are updated at monthly time steps and do not identify the specific timing or abundances of primary vectors, and thus, potential virus emergence. RVFV epizootics in Kenya are associated with flooding events that result in the emergence of high abundances of Aedes mosquitoes, which are considered virus reservoirs [3, 16, 19, 20]. Identifying drivers of high abundances of primary vector species is a logical step toward a better understanding of RVFV ecology, while improving current prediction methods. While high abundances of primary vectors are not the only driver of an RVFV epizootic, information gained from monitoring this component, used in conjunction with existing surveillance efforts, has the potential to contribute to more precise response times when monitoring risk.

Mosquito vectors in the genus Aedes are responsible for RVFV maintenance (i.e., primary vectors), while those in the Culex, Mansonia, and Anopheles genera are responsible for RVFV amplification (secondary vectors) in natural environments [21, 22]. Primary vector ecology is a key component in the RVFV disease system, warranting special consideration when inferring risk. Evidence suggests that adult Aedes mosquitoes in the subgenera Neomelaniconion and Aedimorphus transmit the virus transovarially to their eggs, strongly implicating these species as primary disease vectors [3, 19]. Adult Aedes primary vector mosquitoes can emerge already infected with the virus and become infectious before feeding on wild or domestic ungulates, thereby establishing low levels of virus activity within a geographic area. If suitable environmental conditions persist, Culex, Mansonia, and Anopheles secondary vector mosquitoes emerge and amplify the virus broadly across vulnerable populations [21]. In Kenya, where 11 epizootics occurred between 1951 and 2007 [23], the mosquito species Ae. mcintoshi (included within the subgenus Neomelaniconion) has been implicated as a major RVFV primary vector, and it is the only species to demonstrate the capacity for transovarial transmission [19, 24, 25]. These factors, combined with high RVFV prevalence in the species during the 2006–2007 epizootic warrant special focus on Ae. mcintoshi to investigate primary vector abundances and RVFV ecology in this region [3].

Early field studies provided fundamental information regarding Ae. mcintoshi population biology and ecology. Ae. mcintoshi prefer dambo habitats: shallow depressions in the landscape that become flooded following heavy rainfall [26, 27]. Linthicum et al. [28] flooded a dambo artificially for 18 continuous days and found that female emergence occurred at ~14 days, and blood feeding at ~18 days; mosquito life expectancy was <45 days, female dispersal was low (~0.15 km), and direction of dispersal did not always correspond to surface winds. Logan et al. [29] conducted a sequence of artificial flooding events and found that ~90% of Ae. mcintoshi eggs hatched during the initial flooding of a dambo habitat, with fewer eggs hatching in subsequent flooding periods. Our study uses the biological knowledge gained from these initial field studies to predict Ae. mcintoshi adult abundances across broader landscapes and at different time periods, using georeferenced mosquito abundance data, remotely sensed environmental variables, and a predictive modelling approach.

Specifically, we investigated effects from land surface temperature, cumulative precipitation, compound topographic wetness index values, and percent clay in the soil on adult Ae. mcintoshi abundances, using zero-inflated negative binomial regression and a multimodel averaging approach. We hypothesized that wetness index values, absolute values of land surface temperatures, and cumulative precipitation would all have positive effects on subsequent adult abundance. High wetness index values indicate greater potential for pooling of water during a precipitation event, a factor important to the initiation of Ae. mcintoshi hatching and development. Positive correlations between elevated temperatures and mosquito development warranted the inclusion of temperature as a variable in our analyses, and greater cumulative precipitation promotes flooding and standing water. Additionally, we hypothesized that because soils with more clay retain water better, higher percent clay would be associated with lower probabilities of so-called structural zeros in the zero-inflated model framework (see Methods). Having parameterized models with available Ae. mcintoshi survey data, we used the model to predict Ae. mcintoshi abundances retrospectively across Kenya and western Somalia. We then observed visual correlations between the timing and locations of high predicted Ae. mcintoshi abundances and RVFV activity reported in the literature and through public and veterinary health reporting systems, and we compared retrospective predictions across a 16-year time period (2002 to 2018) to observe general differences between predicted values during epizootic and inter-epizootic time periods.

Methods

Rainfall varies greatly across Kenya, with more moisture in the West. Kenya experiences two rainy seasons referred to as the short rains (October to December) and the long rains (March through May), that coincide with movement of the intertropical convergence zone (ITCZ) [30]. More rainfall occurs from October to December during warm El Niño Southern Oscillation (ENSO) anomalies [31].

Mosquito abundance data for Ae. mcintoshi at 23 locations across Kenya were acquired by RS and JL as part of broader ongoing studies in the United States Army Medical Research Detachment in Kenya (Fig 1). Access to locations did not require official permits but was obtained through informal permission from nearby community leaders. Locations selection included areas where arbovirus disease incidence, outbreaks or serosurveillance were reported through previous reports or publications. Insets of sampling locations at a finer scale are available in supplementary materials (S1 Fig). The data consisted of repeat daily sampling at each location for up to 10 sampling days usually twice yearly during the short rains and long rains from 2007 to 2012, except for 2011. One CO₂-baited CDC light trap was placed overnight at each sample site, and trapped mosquitoes were identified to species by Kenya Medical Research Institute entomology personnel using taxonomic identification keys [32–35].

A total of 158 sampling days were accumulated across all locations during the study period, with abundance values ranging from 0 to 4,426 individuals; dates on which sampling took place but the species was not recorded were assigned an abundance value of zero; such zero counts represented 41% of the total counts (Fig 2).

We obtained Land Surface Temperature/Emissivity data from the Moderate Resolution Imaging Spectroradiometry (MODIS) sensor onboard the Aqua satellite at an 8-day temporal resolution and 1 km spatial resolution [36] through the Reverb ECHO NASA data portal (http://reverb.echo.nasa.gov/reverb). Climate Hazards Group Infrared Precipitation with Stations (CHIRPS) data were obtained through the University of California at Santa Barbara data portal at a daily temporal resolution and a 5 km spatial resolution (http://chg.ucsb.edu/index.html). Compound topographic wetness index values [37] were derived from Shuttle Radar Topography Mission (SRTM) version 4.0 data at a 90 m spatial resolution accessed through the Consultative Group on International Agricultural Research Consortium for Spatial Information (http://srtm.csi.cgiar.org/). Higher wetness index values occur in areas with low slopes and high flow accumulation, which are indicative of locations where water could pool during a precipitation event. Minimum and maximum wetness index values within a 500 m radius of each sampling location were calculated. Percentage of soil clay content was obtained at a 1 km spatial resolution from the Global Soil Dataset for Earth System Modeling Soils [38].

Environmental values at the geographic location of each sample site were extracted from associated raster layers using the raster package in R (Table 1) [39]. Daily cumulative precipitation data corresponding to sampling dates were constructed in three time windows: sampling date to 14 days prior, 14 to 18 days before the sampling date, and 14 to 28 days prior. Land surface temperature data were acquired in 8-day mean composites. We identified the date of sampling and then subtracted 8, 16, and 24 days from the sampling date. The composite data with dates closest to this subtracted value were used in the analysis (Table 1).

Our count data included a greater number of zeros than may be expected under a Poisson or negative binomial distribution (Fig 2). Ignoring this phenomenon would have led to large biases in estimated parameters and their standard errors, and zeros can contribute to overdispersion [40, 41], so we used a zero-inflated statistical modelling approach. Zeros may result from inevitable ecological factors or human error, including sampling error, observer error, or situations in which suitable habitat is present, but is not occupied due to essentially random events [41]. Following standard zero-inflated modelling approaches, we refer to zeros recorded owing to inevitable circumstances as structural zeros, and to those recorded by chance due to sampling variation as sampling zeros [40]. The zero-inflated regression models that we used are mixture models that fit processes for both structural and sampling zeros [42]. Yeşilova et al. [43] described zero-inflated negative binomial regression using

We evaluated the importance of relationships between mosquito abundances and our various environmental variables using an Akaike’s Information Criterion (AIC) and a Bayesian Information Criterion (BIC) approach [46–48]. Lower AIC or BIC scores indicate better-supported models. AIC model weights (denoted AIC_w) were also calculated; these values sum to 1 across all models and indicate the weight of evidence supporting a model. Likewise, BIC weights (BIC_w) were computed. Best-performing models were considered the so-called 99% confidence set of models [47], i.e., for AIC, the models with highest AIC_w values and with a cumulative sum of these weights just exceeding 0.99. Importance of predictors was assessed using a standard sum-of-weights approach [47]: the sum of AIC_w (respectively, BIC_w) values for all models in the 99% confidence set that contained a given predictor was computed, this sum being an index of importance of that predictor. Model-averaged coefficients were computed using a weighted average of the coefficients of the models in the 99% confidence set, with weights the AIC_w (respectively, BIC_w) values. We fit zero-inflated negative binomial models to a random sample of 80% of our data set (n = 127), assigning Ae. mcintoshi abundance values as the response variable and the assembled environmental data as predictor variables (with different subsets of variables used in different models). All models were fitted using the zeroinfl function specifying a negative binomial distribution with a logit link in the pscl package in R 3.13 [49].

We calculated a Pearson’s correlation matrix to assess the potential for multicollinearity between environmental variables, and found high correlation values within, but not between, sets of variables (S1 Table). Only one land surface temperature variable, precipitation variable, or wetness index variable was included in a candidate model at the same time, and we included percent clay in the zero-inflated portion of the model to investigate our hypothesis that zeros recorded in areas where water is likely to pool have a higher probability of being “sampling” zeros, while zeros recorded in areas where water is not likely to pool will have a greater probability of being a “structural” zero.

We investigated model residuals from the lowest-AIC model for evidence of spatial autocorrelation using a spline correlogram in the ncf package in R [50, 51], and found no evidence of it (S2 Fig). Although sampling dates were inconsistent across study site locations and time periods, we investigated the potential for residual temporal autocorrelation, but did not find obvious temporal patterns in residuals.

Estimated parameters from best performing models were used to predict mosquito abundances using a random sample of 20% of the abundance data (n = 31) withheld from the regression models, and a root mean square predictive error was calculated to evaluate accuracy. We used the predict function in the pscl package, with type = “response” to incorporate both the structural zero portion and the sampling portion of the zero-inflated negative binomial model when predicting values [49, 52].

For retrospective predictions, unsampled locations were generated across the study region at 10 km intervals, predictions were generated from each model in the 99% confidence set, and a weighted average of predictions was computed using the AIC_w values as weights. Predicted values were rasterized into a 10 x 10 km grid for visualization purposes. Values greater than 50,000 were considered extrapolative and excluded from visualizations, based on the fact that the greatest number of mosquitoes trapped in a single day across all species in our study locations was 47,694. Values > 50,000 constituted a very small percentage (~0.01%) of predicted values in each raster and occurred only at the fringes and across portions of large water bodies, presumably, where extreme compound topographic wetness index values are present. Data and sample R code are available in S1 Data Code.

Results and discussion

Model results corroborated several factors known from local-scale studies to be important for Ae. mcintoshi development and population ecology. The 99% confidence set of models with respect to AIC consisted of 11 models (Table 2; full model results in S2 Table), and, with respect to BIC, of 19 models (S2 Table); these are the models best supported by data according to each information criterion, with the support for each model quantified by its AIC (or BIC) weight (AIC_w or BIC_w). The same variables tended to be included as predictors in the best-supported models according to both AIC and BIC, and model-averaged coefficients were similar (S3 Table), so the choice of information criterion did not substantially affect results. We henceforth use AIC because our objectives include prediction, for which AIC is considered more suitable [53, 54].

The sum of AIC_w across the 11 best models including a given predictor indicated the importance of that predictor; these sums showed that minimum wetness index values within 500 m of a sampling site, cumulative precipitation from the sampling date to 14 days prior, land surface temperature values from 8-day mean composites 24 days prior to sampling, and percentage of clay in the soil (in the structural-zero portion of the model—see Methods) had the greatest importance for predicting abundance (Table 2).

The data strongly supported a positive relationship between minimum wetness index values within 500 m of a sampling site and mosquito abundance (Table 2). This result corroborated existing information regarding Ae. mcintoshi ecology: locations with high minimum wetness index values have low slopes and high flow accumulation, indicative of dambo habitats or landscapes likely to collect water during a precipitation event. We also found a positive effect on mosquito abundance of cumulative precipitation from the sampling date to 14 days prior, and no meaningful effect of cumulative precipitation 14 to 18 or 14 to 28 days prior to sampling. These results indicated that precipitation within a short time period prior to emergence has a greater impact on Ae. mcintoshi abundances, suggesting that forecasting may be enhanced with the addition of shorter time intervals. Model results indicated a positive effect on mosquito abundances of 8-day mean land surface temperature composites at 8, 16, and 24 days prior to sampling, with 24 days prior having the greatest effect among these variables. Temperature is known to play an important role in mosquito life stage development and may impact RVFV vector competency in some species [55, 56]. Our study results suggest that warmer land surface temperatures a few weeks prior may produce conditions that subsequently contribute to higher adult abundances, if precipitation events initiate mosquito development. Interestingly, Anyamba et al. and Linthicum et al. [57, 58] have found negative associations between anomalously lower temperatures and RVFV activity in eastern and southern Africa, suggesting a more complex association between this variable and events leading to epizootic events. We found a negative relationship between percent clay in the soil and the probability of a structural zero within the zero-inflated portion of the model, indicating that areas with low clay content had a higher probability of a structural zero, presumably because standing-water pools are less likely to form over such soil types.

Model-averaged out-of-sample predictions demonstrated the capacity of our models to predict elevated abundances (Fig 3). Our models’ predictive accuracy was low for low observed abundances: the relationship between observed and predicted values was weak on the left-most side of Fig 4. However, as the intended use of the models is to predict very high abundances, this point is of secondary importance. Prediction of high abundances was effective, although the model underestimated the observations. Model performance statistics, including root mean square error values from out-of-sample predictions for the 11 best models are provided in more detail in S2 Table.

Retrospective model-averaged vector abundance predictions were produced at 8-day intervals from 30 September 2002 through 25 January 2018, and model predictions were generated from 25 January 2018 to August 2018 because of anomalous precipitation events, resulting in a RVFV epizootic during this time period. Monitoring and human health agencies such as the World Health Organization reported [59–61] high abundances of RVFV cases at five times and locations during the 2006–2007 epizootic period (Fig 4, numbered circles). Retrospective vector abundance predictions were high pre-dating three of these, while predictions were low prior to two RVFV foci (Fig 4, heat maps). Elevated model predictions prior to the 2018 RVFV epizootic suggested two time periods in the month of April in which Ae. mcintoshi abundances may have been elevated, but the locations of theses predicted abundances were approximately 125 km from the first RVFV focus. Additionally, very high abundances were predicted along the southeastern coast of Kenya, where outbreaks were not reported. These results highlight the complexity of RVFV transmission in space and time and warrant additional investigation into the role that primary vector activity as a precursor to epizootic events, especially when immune status of animals and /or their vaccination history preclude viremia and a disease outbreak in these animals. Notable model outputs during epizootic and inter-epizootics are described in more detail, and a timeline between predicted abundances and epizootic activity is provided in supplementary materials (S4 Table).

Elevated abundances were predicted for 1 November 2006 in Bay, Bakool, Jubbuda Dhexe, and Gedo Provinces, Somalia (Fig 4, circle 1), ~7 weeks prior to reports of human cases [62]. Our results support virus circulation within Somalia independent of Kenya, with high vector abundance predictions in isolated areas prior to elevated predicted abundances in Kenya (Fig 4). Although, the ecological habitats and some of the drainage systems in Kenya and Somalia were essentially identical it is possible that flooding of mosquito habitat in Somalia occurred before habitats in Kenya due to rainfall events. High abundances were predicted for 17 and 25 November in Garissa, Ijara, and Kitui Districts, Kenya; human cases of RVFV were first reported in mid-December in Garissa, with the potential index case presenting on 30 November (Fig 4, circle 2) [63]. Additionally, several of the Kenyan districts located within high predicted abundance areas for 25 November later reported RVFV activity for the first time: Kitui, Tharaka, Mwingi, Embu, Kirinyaga, Meru South, Meru Central, and Malindi District (Fig 4, circle 3) [60]. The World Health Organization recognized cases in these districts as a new RVFV focus during the outbreak, although little information is available regarding the exact date of onset of symptoms in humans or animals [61]. It is important to note that the Food and Agricultural Organization of the United Nations published a RVFV outbreak alert in October 2016 (see [64] and www.geis.fhp.osd.mil/GEIS/SurveillanceActivities/RVFWeb/indexRVF.asp) for East Africa based upon the NASA/USDA/DoD RVFV risk prediction model (https://ars.usda.gov/saa/cmave/rvf) indicating elevated risk in September of 2016 [15].

Human RVFV cases were first reported retrospectively on 6 December 2006 in Kilifi District, near the eastern coast of Kenya (Fig 4, circle 4) [63]. Our models did not predict high Ae. mcintoshi abundances at any time in this area. Nguku et al. (2010) [63] found that illness in Kilifi District coincided with heavy rainfall, rather than emerging approximately one month after heavy rainfall as was the case in other regions, and that movement of infected livestock from the outbreak area in Northeastern Province, Kenya, into Kilifi District may have been the catalyst for the outbreak. These observations suggest that primary vector emergence in Kilifi District may not have been responsible for virus circulation in the area, but further investigation is needed.

Illness in Baringo District was estimated to begin in late December, 2006, with the first human case reported on 25 January 2007 [60, 63], but our model results did not predict elevated Ae. mcintoshi abundances in Baringo District at any time during the study period, suggesting that other species may have played a role in driving the RVF transmission (Fig 4, circle 5). Although it is possible that our model predictions did not reflect the potential for Ae. mcintoshi abundance in this region, more likely the model was accurate and suitable habitat for this species was deficient in the area. Low numbers of Ae. mcintoshi were collected in Baringo District during the 2006–2007 epizootic, and follow-up sampling also retrieved low densities of this species [3, 27]. Our data set had only two records for Ae. mcintoshi in Baringo District with 2 and 6 females recorded, respectively. Despite the low numbers of Ae. mcintoshi collected in our sample data set, this result likely had to do with the timing and sampling effort in this area, as Ae. mcintoshi is known to be present and abundant in this region at various time periods (KJL personal communication). Additional sampling is required to capture environmental conditions that produce high numbers of Ae. mcintoshi in Baringo District in order to predict abundances accurately.

In February 2018, unusually heavy rainfall activity occurred in Northeastern Kenya; County Directors were put on alert for RVF activity, and a national alert was issued 30 May 2018 [65, 66]. The first human RVF case entered a medical facility in Wajir County on 2 June, and an RVF outbreak was declared on 8 June [65]. Subsequent cases followed in Garissa, Kadjiado, Kitui, Marsabit, and Tana River Counties. Model predictions from February to June 2018 indicated high predicted Ae. mcintoshi abundances on 15 April in Wajir County on the border with Somalia, approximately 6–7 weeks prior to the RVF outbreak, but the area with highest predicted abundances was ~125 km southeast of the outbreak area (Fig 5). Slightly elevated abundances were predicted for the following week in Garissa and Tana River Counties. Very high predicted abundances on 9 May 2018 spanned southern Garissa and Lamu Counties on the coast, and east into Somalia.

Comparative predictions during the short rain season between 2002 and 25 January 2018 indicated lower predicted abundances across the study period (S1 Movie File), with the exception of one date in 2002 (Fig 6) and two dates in 2016 (Fig 7), and a few small- and medium-scale predictions of elevated abundances which apparently did not result in reported RVFV cases (Fig 8). High predicted abundances from 9 November 2002 were located to the northwest of the 2006–2007 focus in Isiolo and Laikipia Districts, Kenya (Fig 6, circle 1); virus activity was not reported there during this time period, nor was it reported in the regions shown in Figs 7 and 8. Even though conditions appeared to be suitable for high primary vector emergence, the possibility exists that ova deposited in these areas were not infected with RVFV transovarially and did not emerge with the capacity to transmit the disease; conditions were not suitable for large numbers of secondary vectors to amplify the virus; or susceptible livestock were not present in the area. In fact, some evidence indicates possible livestock movement restrictions in this area and during this time period because of rinderpest virus detected in cattle on 23 October 2002 [67], but further investigation is needed before determining whether this factor could have impacted potential RVFV circulation.

Importantly, abundances predicted by our models were low during two time periods in which RVFV warnings were issued in Kenya but no virus activity was reported, demonstrating the contribution of investigating RVFV at multiple spatial and temporal scales. Animal health authorities issued a warning on 13 November 2008 in the far northeastern region near the borders of Ethiopia and Somalia because of persistent rainfall. Our model did not predict high vector abundances in this region prior to this date. The highest abundances predicted around this date were only moderate and were south of this region (Fig 6, circle 1). Several national and international agencies issued warnings for potential RVFV activity in the 2015–2016 period because of an El Niño event in the study area [68], but we did not predict high abundances during this time period (Fig 5, right panel) and no outbreaks were reported. Interestingly, very high abundances were predicted for 26 November 2016 and 4 December 2016 (Fig 7), but again no outbreak was reported. One likely possibility is that implementation of widespread livestock vaccination program occurred in the lead up to the 2015 El Niño event provided a substantial level of herd immunity throughout 2016 [18]. Integrating model predictions with livestock restocking rates and immunity status into a broader model of RVF potential will provide a more grounded assessment of virus risk than estimates of primary vector abundances alone. Weekly predictions during the short rain season from 1 September to 25 January between 2002 and 2018 are available in the supplementary material (S1 Movie File).

Although our model can predict high mosquito abundances, our data set does not appear to capture the breadth of environmental conditions suitable for Ae. mcintoshi in central or western regions of Kenya, highlighted in the absence of high predicted values during the 2006–2007 epizootic in Baringo District and the collection of high abundances of this species in the past (Linthicum and Anyamba, personal communication). The very low number of Ae. mcintoshi collected from this region was likely due to the sampling effort available in the data set used for model calibration, demonstrating the need for additional and regular sampling to capture conditions that contribute to large emergences of Ae. mcintoshi.

In addition to this limitation, model predictions do not incorporate mosquito density dependence, which may affect the accuracy of predicted values, particularly at higher values. Additionally, we have lower confidence in predicted values greater than ~7,000 because of the potential for model extrapolation. The choice of 7,000 individuals is a user defined threshold and should be interpreted as such. Predicted values greater than the maximum 4,426 individuals in our data set could be interpreted as extrapolation, but we assume that the possibility exists for a greater number of individuals to be collected in a single trap night. This threshold allows for predicted values to exceed the maximum individuals in our data set, while providing a value at which to indicate visually when and where much higher predicted values occur in model outputs. Input data consist of irregular sampling in space and time. This sampling scheme proved challenging to assess the potential for spatial and temporal autocorrelation that may affect parameter estimates, but we did not find evidence of this phenomenon. We predict mosquito abundances across a relatively broad spatial region, which also increases the potential for extrapolation; although, the combination of environments across the prediction area are not unlike those encountered at the sampling sites, and we assume that the effect of these environments on mosquito abundances are similar to the effect of these environments on the sample data, except in the central and western portion of the study area, where we determined insufficient data for accurate model predictions. We also predict across a broad temporal period. The irregular sampling scheme, across a five-year time period helps mitigate temporal extrapolation across the predictive time period.

Predicting abundances from mosquito trap data presents challenges because of the potential for bias between individual traps, trap types, and trap bait [69]. Additionally, mosquito traps collect a fraction of the diversity and relative abundances of a mosquito community on a given night. In the case of this study, consistency in trap type and bait used across sampling sites helped to reduce bias from multiple trap types, but uncertainty remains in exact diversity and abundances of mosquitoes in collected areas. Satellite measurements of environmental variables aggregate values across a landscape to a single spatial resolution, even though the potential exists for local variations. Our data consisted of repeat sampling at locations across Kenya, but sampling was sparse in some areas. Our model predicts abundances for one primary RVFV vector, but additional primary vector species likely exist: several species were identified for further investigation during the 2006–2007 epizootic [3]. RVFV surveillance may benefit from applying our modelling framework to those species.

Our model is the first to use predicted vector abundances at an 8-day temporal resolution as a means of predicting primary RVFV mosquito vector emergence in this region. This approach provides a powerful framework to help inform current RVFV monitoring systems. Additionally, our framework has the potential to reveal new information about primary vector activity during inter-epizootic years that may contribute to a better understanding of RVFV ecology. The potential exists to extend our framework to multiple vectors and disease systems. Our approach also allows exploration of the effects of climate change on Ae. mcintoshi or other medically important vector species. For improved model performance, we recommend conducting mosquito collections at regular and more frequent temporal intervals and that these collections capture greater spatial and ecological diversity throughout RVFV potential epizootic areas.

Supporting information

S1 Fig [tif]
Insets of sampling locations at a finer scale.

S2 Fig [tif]
Spline correlogram of distance between locations and model residuals for the lowest-AIC model.

S1 Data Code [zip]
Data used in analysis and R code to execute model runs.

S1 Table [xlsx]
Pearson’s correlation matrix.

S2 Table [xlsx]
Candidate models: AIC and BIC rankings, scores, differences, and weights for all candidate models run in this analysis.

S3 Table [xlsx]
Model averaged coefficient values for AIC and BIC models corresponding to a cumulative weighted value of 0.994 (11 models for AIC and 19 models for BIC).

S4 Table [xlsx]
Timeline between predicted abundances and epizootic activity.

S1 Movie File [mp4]
Animated plots between 1 September and 25 January from 2002 to 2018.