Climate-driven spatial mismatches between British orchards and their pollinators: increased risks of pollination deficits

Understanding how climate change can affect crop-pollinator systems helps predict potential geographical mismatches between a crop and its pollinators, and therefore identify areas vulnerable to loss of pollination services. We examined the distribution of orchard species (apples, pears, plums and other top fruits) and their pollinators in Great Britain, for present and future climatic conditions projected for 2050 under the SRES A1B Emissions Scenario. We used a relative index of pollinator availability as a proxy for pollination service. At present, there is a large spatial overlap between orchards and their pollinators, but predictions for 2050 revealed that the most suitable areas for orchards corresponded to low pollinator availability. However, we found that pollinator availability may persist in areas currently used for fruit production, which are predicted to provide suboptimal environmental suitability for orchard species in the future. Our results may be used to identify mitigation options to safeguard orchard production against the risk of pollination failure in Great Britain over the next 50 years; for instance, choosing fruit tree varieties that are adapted to future climatic conditions, or boosting wild pollinators through improving landscape resources. Our approach can be readily applied to other regions and crop systems, and expanded to include different climatic scenarios.

Climate-driven spatial mismatches between British orchards and their pollinators: increased risks of pollination deficits.

Supporting Information
Material and Methods Pollinator data Table S 1 lists the wild insect pollinator species used to predict the potential pollinator availability to orchard crop flowers in Great Britain. Available records indicate the number of sites of known presence, mapped on 5 by 5 km grid cells. The honeybee (Apis mellifera) was excluded due to the difficulties of separating records from managed colonies from those of wild colonies.
Climate data for future projections Monthly averages for future projections were derived from the daily data generated by the HadRM3-PPE-UK experiment (Murphy et al., 2010), run at the UK Met Office Hadley Centre (MOHC). The dataset contained the outputs from an ensemble of eleven variants of the MOHC Regional Climate Model (HadRM3). The HadRM3 model run from 1950-2099 and was used to dynamically downscale global climate model results, as part of the climate change experiments carried out for the UK Climate Projections report (UKCP09). The data were obtained online after registration (http://badc.nerc.ac.uk/home/). For our study, we used the unperturbed model variant of the HadRM3-PPE-UK experiment (identified by the name HadRM3Q0).
The HadRM3-PPE-UK data are located on a grid characterised by a rotated-pole and a spatial resolution of approximately 25 by 25 km. We used the accompanying metadata to convert the rotated grid into a regular grid of latitude and longitude. To match the resolution of the baseline data, we then rescaled the data to a 5 by 5 km British National Grid, giving to each of the 5 by 5 km grid cell the value of the coarser 25 by 25 km grid. Correlation between selected climatic predictors Table S 2 shows the Pearson's correlation between the variables used to predict orchards distribution. Variables were selected from the original set of 20 climatic predictors, using literature (Franklin et al., 2013, Sork et al., 2010, Termansen et al., 2006, Thuiller, 2004, Warren et al., 2013. The values refer to present conditions, which for orchards corresponded to the 30-year period from 1977 to 2006. Variables are defined in Table 1 of the main text.  (Nakićenović et al., 2000). We used the 30-year period from 2040 to 2069 ("M2050"). These variables were selected on present-day conditions using literature (see Table S 2) and are defined in Table 1 of the main text.  Table S 4 shows the Pearson's correlation between the variables used to predict pollinators' distribution. We used Jolliffe's Principal Component Analysis with the rejection method "B2" and λ 0 = 0.70 (Jolliffe, 1972, Jolliffe, 1973, to select a set of variables that minimised multicollinearity, from the original 19 bio-climatic predictors. The values refer to present conditions, which for pollinators corresponded to the 10-year period from 1990 to 1999. Variables are defined in Table 1 of the main text.  Table S 5 shows the Pearson's correlation between the variables used to predict pollinators' distribution, based on the predictions for the Medium Emission Scenario (SRES A1B storyline), for the M2050 period. These variables were selected on present-day conditions using PCA (see Table S 4) and are defined in Table 1 of the main text. Pollinator distribution models Detailed settings for the Maxent pollinators' distribution models (PDM) follows Polce et al.(2013) . Models were set to run with Hinge features and with modified values of prevalence, to reflect the differences in number of available sightings at any resolution, from the original databases. The background was chosen to reflect the non-random distribution of the species' sightings in relation to the environmental range available for Great Britain. It corresponded to all localities where crop pollinators had been recorded (Phillips et al., 2009).

Contribution of different predictors
The contribution of each predictor to the final Maxent model was measured using the permutation importance. Within Maxent, permutation importance is determined for each variable by randomly permuting the values of the variable among the presence and background training points and evaluating the resulting decrease in training AUC. The drop in AUC is then normalised to percentage allowing comparison across models; a larger percentage (a larger drop) indicates that the original model depended heavily on that variable. Average and confidence interval for the importance of the different predictors were derived through 10,000 bootstrap replicates.

Results
Climatic and bioclimatic variables Each chart in Fig. S 1 shows the distribution of predictors' values for present (in black) and future M2050 (in grey) climatic conditions. The predictors are the ones used to model pollinators and orchards distribution. They are defined in Table 1 of the main text. The y-axis provides an estimate of the number of grid cells characterised by each predictor's value. Model performance Figure S 2 show the model performance for the pollinator distribution models, as measured by the the AUC of model testing. Error bars show the standard deviation of the null models (10 sets for each pollinator species, each modelled with 10-fold cross-validation). The number of available records is used to plot different species along the x-axis. Figure S 3 shows the contribution of the different predictors used for PDMs, measured by the arithmetic and bootstrap mean of each predictor's importance, pooled across pollinator species. Confidence interval shows the 95% biased-corrected accelerated percentile, based on 10000 replicates. TAR was significantly more important than the other predictors, whilst Isoth and MTDQ were equally the least important. The significance of multiple pairwise comparisons was tested using Tukey's post-hoc test (Table S 6 and Table S 7). Predictors are defined in Table 1 of the main text.

Figure S 3: Importance of predictors used for PDMs
Table S 6 shows the significance of the differences in permutation importance of the six predictors used for the pollinator distribution models. We used linear mixed effects model, with model run nested within species and group (bees or hoverflies) as a random factor, and predictor as fixed effect. Before the analysis, we applied an angular transformation to the permutation importance percent, to improve normality of residuals. We used Isoth as a baseline. Predictors are defined in Table 1 of the main text. Fixed effects only are shown here.  (Hothorn et al., 2013). These results and those in Table S 6 show that TAR was the most important predictor, significantly more important than all others (e.g. MTCQ, which followed next). On the opposite end, with equal importance, MTDQ and Isoth were the least important predictors. Predictors are defined in Table 1 of the main text. Figure S 4 shows the contribution of the different predictors used for ODMs, measured by the arithmetic and bootstrap mean of each predictor's importance, from different runs of the orchard distribution model. Confidence interval s show the 95% biased-corrected accelerated percentile, based on 10,000 replicates. TSeasSD was significantly more important than all other predictors. The significance of multiple comparisons was tested using Tukey's post-hoc test (Table S 8 and Table S 9). Predictors are defined in Table 1 of the main text. Table S 8 shows the significance of the differences in permutation importance of the five predictors used for the crop distribution models. We used linear mixed effects model with model run as a random factor and predictor as fixed effect. Before the analysis, we applied an angular transformation to the permutation importance percent, to improve normality of residuals. We used Temperature Seasonality as a baseline. Predictors are defined in Table 1 of the main text. Fixed effects only are shown here.  (Hothorn et al., 2013). These results and those obtained from the mixed effects model in Table S  8 show that TSeasSD was significantly more important than all other predictors, including mTCM which followed next. The remaining predictors followed in a significantly decreasing importance, according to the rank shown in Figure 5 of the main text. . Values of mTCM are taken from present conditions. Greatest p is predicted around 1.5 °C. The red line shows the mean from 10 Maxent models, the blue shade shows 1 standard deviation from the mean.

Future projections
The two histograms in Fig. S 6 show current and future distribution of Spearman's correlations between orchards extent and probability of occurrence (p), based on 9999 samples with replacement. The dashed lines indicate the observed correlation in each period. Probability of occurrence was estimated using climatic predictors, with Maxent. From the position of the dashed line we conclude that the correlation between crop extent and p is significantly different than what it could be expected by chance alone. In particular, there is a positive correlation between crop extent and p in the present (ρ = 0.153), and a negative correlation (ρ = -0.233) between the current crop extent and the future p.

Figure S 6: Observed and simulated correlation between orchards extent and probability of occurrence for present and future climatic conditions
Fig. S 7 shows the predicted pollinator availability (PA) currently available to orchards. PA is used as a proxy for pollination service, and measured with a relative index from 0 to 1. PA is mapped from red to green using intervals. Red is used only to indicate areas where PA is predicted to be 0 (i.e. where pollinators are predicted to be absent). The chart shows the frequency of the predicted PA, following the same intervals used for the map. Frequency is measured as number of grid cells with a certain value of PA. The intervals are based on a common scale, used for Figs S 7 to S 9.

Figure S 7: Current pollinator availability to orchards
Fig. S 8 shows the predicted future PA available to orchards, in areas identified as the most suitable to crop growth, based on future climatic conditions (i.e. M2050 scenario). Red is used to map areas where PA is predicted to be 0, due to absence of wild pollinators. The chart shows the frequency of the predicted PA, following the same intervals used for the map. Frequency is measured as number of grid cells with a certain value of PP. The intervals are based on a common scale, used for Figs S 7 and S 9. shows the predicted PA where orchards are currently grown, based on future climatic conditions (i.e. M2050 scenario). PA is mapped using interval classes, following the same scale adopted for Figs S 7 and S 8. The map suggests that, if orchards persisted where they are currently planted, they would all receive some pollination service. Green shades are used to map areas with greater PA. The chart shows the frequency of the predicted PA, following the same intervals used for the map. Frequency is measured as number of grid cells with a certain value of PA. The intervals are based on a common scale, used also for Fig. S 7 and S 8.