Predicting malaria vector distribution under climate change scenarios in China: Challenges for malaria elimination

Projecting the distribution of malaria vectors under climate change is essential for planning integrated vector control activities for sustaining elimination and preventing reintroduction of malaria. In China, however, little knowledge exists on the possible effects of climate change on malaria vectors. Here we assess the potential impact of climate change on four dominant malaria vectors (An. dirus, An. minimus, An. lesteri and An. sinensis) using species distribution models for two future decades: the 2030 s and the 2050 s. Simulation-based estimates suggest that the environmentally suitable area (ESA) for An. dirus and An. minimus would increase by an average of 49% and 16%, respectively, under all three scenarios for the 2030 s, but decrease by 11% and 16%, respectively in the 2050 s. By contrast, an increase of 36% and 11%, respectively, in ESA of An. lesteri and An. sinensis, was estimated under medium stabilizing (RCP4.5) and very heavy (RCP8.5) emission scenarios. in the 2050 s. In total, we predict a substantial net increase in the population exposed to the four dominant malaria vectors in the decades of the 2030 s and 2050 s, considering land use changes and urbanization simultaneously. Strategies to achieve and sustain malaria elimination in China will need to account for these potential changes in vector distributions and receptivity.


Mosquito presence data
The presence records for malaria vectors in this study came from an exhaustive and systematic search of formal as well as informal publications. Records were available from six years (2005)(2006)(2007)(2008)(2009)(2010)  sinensis  , An. lesteri 20,24, , were searched for and confirmed by a technical advisory group of Anopheles experts including malaria epidemiologists, entomologists and ecologists. We also searched formally published literatures from 2000 to 2010. From these searches, 247 articles were identified and the full articles were downloaded. We removed the articles that did not contain information relating to these four malaria vectors occurrence. Records of the presence of the dominant malaria vectors at a particular site and survey date were entered into the database so that information collected at different times from a locality was documented. In order to ensure the quality of malaria vectors presence data, we only kept the administrative unit indicating confirmed more than three times occurrences of malaria vectors in a given searching period. Finally, data from a total of 120 published articles from 2000 to 2010 were compiled (Supplementary Table S1). We recorded the county names and reported Anopheles species. These data were then matched with county level administrative maps in order to assign a location to each presence observation.
A total of 27 (An. dirus), 33 (An. minimus), 59 (An. lesteri) and 95 (An. sinensis) county level presence records from 2000 to 2010 were identified. For visualization purposes, the centroid of each county was used to map current presence data ( Supplementary Fig. S2).

Supplementary Notes Environmental variables contributions
The relative contributions of the environmental variables are shown in Supplementary Fig. S2. Based on the Maxent model's internal jacknife test of variables importance, the current distributions of ESA for An. dirus were largely affected by the annual temperature range (bio7) and precipitation of wettest quarter (bio16). The mean temperature of coldest quarter (bio11) and annual temperature range (bio7) largely controlled the spatial distribution of ESA for An. minimus. The current distribution of ESA for An. lesteri and An. sinensis were both largely affected by the precipitation of driest quarter (bio17) and fraction of urban area within grid cell (gurbn).

Relationship between vector occurrence and environmental variables
A negative near-linear relationship is observed between the mean temperature annual range (bio7) and the probability of An. dirus presence, peaking at 20 °C. In contrast, the precipitation of the wettest quarter (bio16) seems to be restrictive and a relative high probability of An. dirus presence is noted where precipitation is between 900 mm to 1,000 mm ( Supplementary Fig. S3-1). Above 5mm, the precipitation of the driest quarter was negatively related to predicted probability of occurrence of An. minimus, while the mean temperature of the coldest quarter was positively related to predicted probability of occurrence (Supplementary Fig. S3-2). The predicted probability of occurrence for An. lesteri was positively related to precipitation in the wettest quarter (bio16) from 0mm to 50mm, while there was a negative relationship when greater than 50mm ( Supplementary Fig. S3-3). The response curve of An. sinensis indicated that both the precipitation of the driest quarter (bio17) and the minimum temperature of the coldest month (bio6) were positively related to probability of occurrence ( Supplementary Fig. S3-4). As our model predicted, the probability of An. dirus occurrence was non-linear negatively related to gurban. While predicted probability of occurrence of An. minimus was negatively related to gurban where gurban is below 0.03. However, both the probability of An. lesteri and An. sinensis was non-linear positively related to gurban.  Fig. S4-S12). In general, the ESA for An. sinensis estimated by CCCma_CanESM2 were larger than the other two GCMs.

Assumptions, limitations, and evaluation of point sampling approach
The point sampling approach is a useful method to predict finer-resolution species distributions when only coarse presence data (e.g. county level vector presence data) are available 121 . In our study, point sampling involved assigning the presence data to a random location (random point) within each county. The set of random points were then associated with fine resolution environmental variables via an ecological niche model (e.g., Maxent) 122 . However, the approach inherits several inevitable assumptions and limitations. This approach assumes all areas are suitable for the modelled species within the county boundary 123 . This assumption may not be true due to fine scale heterogeneities in suitable habitats. As a degree of uncertainty arises in predictions derived from different random sampling iterations 124