Potential landscape-scale pollinator networks across Great Britain: structure, stability and inﬂuence of agricultural land cover

Understanding spatial variation in the structure and stability of plant – pollinator networks, and their relationship with anthropogenic drivers, is key for maintaining pollination services and mitigating declines. Constructing sufﬁcient networks to examine patterns over large spatial scales remains challenging. Using biological records (citizen science), we constructed potential plant – pollinator networks at 10 km resolution across Great Britain, comprising all potential interactions inferred from recorded ﬂoral visitation and species co-occurrence. We calculated network metrics (species richness, connectance, pollinator and plant generality) and adapted existing methods to assess robustness to sequences of simulated plant extinctions across multiple networks. We found positive relationships between agricultural land cover and both pollinator generality and robustness to extinctions under several extinction scenarios. Increased robustness was attributable to changes in plant community composition (fewer extinction-prone species) and network structure (increased pollinator generality). Thus, traits enabling persistence in highly agricultural landscapes can confer robustness to potential future perturbations on plant – pollinator networks.


INTRODUCTION
Insect pollinators face many threats that may jeopardise the crucial ecosystem service they provide to crops and wild plants (Vanbergen et al. 2013;Gill et al. 2016;Potts et al. 2016). The stability of pollinator communities and the service they deliver is mediated by the structure of ecological networks formed by interactions between pollinator and plant species (V azquez et al. 2009;Vanbergen et al. 2017). Understanding such networks is important to predict the risks associated with threats to pollinators (Gill et al. 2016). Analysis of plant-pollinator networks has provided insights into their structure and potential stability under actual or simulated environmental change, including extinctions (e.g. Memmott et al. 2004;Kaiser-Bunbury et al. 2010), climate change (e.g. Memmott et al. 2007), habitat change (e.g. Forup et al. 2008;Vanbergen et al. 2017) and restoration (Kaiser-Bunbury et al. 2017).
Studies traditionally rely on obtaining well-characterised networks from field surveys, which are time consuming and costly to construct (V azquez et al. 2009;Burkle & Alarc on 2011). Constructing networks replicated across larger spatial scales remains a daunting prospect (Burkle & Alarc on 2011), particularly at the regional and national scales relevant to land use and conservation policy-making. Although broad geographical patterns in plant-pollinator network properties have been identified across biomes (Olesen & Jordano 2002;Welti & Joern 2015) or within landscapes (Burkle & Alarc on 2011;Carstensen et al. 2014;Trøjelsgaard et al. 2015;Kaiser-Bunbury et al. 2017), these still rely on a comparatively limited number of empirical plant-pollinator networks.
Of particular interest in understanding spatial variability in plant-pollinator networks is the contrast between the benefits of insect pollinators to agricultural crops (Kremen et al. 2002;Winfree 2008;Eilers et al. 2011) and the negative impacts of intensive agriculture on pollinators (Kluser & Peduzzi 2007;Potts et al. 2010;Gill et al. 2016). However, we have very little knowledge of how plant-pollinator networks are affected by agriculture at landscape scales (e.g. > 1 km 2 ) or whether networks comprising species that pollinate agricultural crops are representative of the wider pollinator community (Kleijn et al. 2015;Gill et al. 2016).
Lack of information on ecological interactions across larger spatial, temporal and taxonomic scales, termed the 'Eltonian shortfall', represent key gaps in our large-scale knowledge of biodiversity (Hortal et al. 2015). Moreover, there are limitations on the extent to which different data sources can be combined to analyse multiple networks because data collection methods can introduce potential biases (Hortal et al. 2015). While there are exciting possibilities for molecular techniques to increase the speed and accuracy with which plant-pollinator networks can be constructed (Keller et al. 2015;Richardson et al. 2015;Bohan et al. 2017;Pornon et al. 2017), these are yet to be realised across larger spatial scales.
Biological records (i.e. records submitted to voluntary recording schemes, a form of 'citizen science') provide a valuable resource for analysing large-scale patterns in time and space (Bishop et al. 2013;Tulloch et al. 2013;Powney & Isaac 2015). Records consist of species' identification, date and location (hereafter 'occurrence' data) and provide large volumes of data over a wide spatial coverage, equivalent to innumerable hours of field survey. Methods to control for variation in recorder effort and to infer ecological signals from occurrence data are rapidly emerging (e.g. Isaac et al. 2014;Dyer et al. 2016), but hitherto their potential as a source of data on ecological networks is untapped (Gray et al. 2014).
Here, we constructed potential plant-pollinator networks for every 10 km-by-10 km grid square ('hectad') in Great Britain (GB) using interactions from a 30-year, long-term national data set of occurrence records of pollinating insects (bees, butterflies and hoverflies). Instead of inferring species interactions from spatial co-occurrence (Morales-Castilla et al. 2015; Morueta-Holme et al. 2016) we used metadata from records that detailed flower visitation as a proxy of pollination (Ballantyne et al. 2015). These networks are 'potential' in that we acknowledge their limitations in terms of assumptions that constrain their biological realism. However, while the structure of each potential network may be subject to errors, we aimed to minimise bias affecting comparisons across replicate networks. We used these potential networks to address three questions. First, does network structure and stability vary spatially across GB? Second, is network stability reduced by greater agricultural land cover, a major driver of plant and pollinator declines (Kluser & Peduzzi 2007;Potts et al. 2010;Vanbergen et al. 2013;Ollerton et al. 2014)? Finally, are the structure and stability of networks comprising crop-pollinator species consistent with those of the wider pollinator community?

Constructing a plant-pollinator interactions database
We constructed a national-scale (GB) plant-pollinator interaction database defining which species of pollinator visit which species of plant (Fig. 1). Data were mostly (73%) sourced from biological records. Specifically, these were species observations submitted to the Bees, Wasps and Ants Recording Society (BWARS), Butterflies for the New Millennium (BNM, Asher 1997) and Hoverfly Recording Scheme (HRS), with plant interactions recorded as incidental metadata. Interactions were inferred by algorithmically screening metadata for valid scientific or vernacular plant names (or widely used synonyms or abbreviations thereof), followed by data cleaning (see Appendix S1). Remaining interaction data were obtained from books (e.g. Morris 1998), papers (e.g. Carvell 2002) and unpublished experimental data. Where interactions were recorded only to plant genus we assumed, given the rarity of pollinators specialised to the level of individual plant species (Waser et al. 1996;Minckley & Roulston 2006), that these were indicative of interactions with all plant species within the genus that were present in the data set (full details in Appendix S1). These inferred interactions comprised 6487 unique interactions (39%) within the full data set.
Our final plant-pollinator interactions database contained 16,712 unique interactions, involving 485 pollinator species (206 bees, 56 butterflies and 223 hoverflies) and 499 plant species. This total comprises approximately 76, 92, 81 and 55% of GB bee, butterfly, hoverfly and insect-pollinated plant species respectively (Fitter & Peat 1994;Stubbs & Falk 2002;Thomas 2010;Falk 2015). We explored the completeness of our interactions database by calculating interaction accumulation curves across all records used to construct the database (i.e. pollinator occurrences where we were able to identify a valid plant interaction) and for each plant and pollinator species separately (Appendix S2). Results suggested that our database captured around 60% of estimated total interactions (mean 62% for pollinators, 57% for plants), comparable to studies which performed high-effort, multitemporal field sampling of individual networks (Chacoff et al. 2012;Falcão et al. 2016).

Modelling plant and pollinator occurrence
For all species in the interactions database, we obtained occurrence data from BWARS, HRS, BNM and, for plants, the Botanical Society of Britain and Ireland (Fig. 1). Data records were restricted to 1985 onwards, covering the vast majority of records while excluding occurrences of species that may have been more widespread prior to major changes in GB land use (Robinson & Sutherland 2002;Ollerton et al. 2014). Occurrence data were modelled to account for spatial bias in recorder effort using the FRESCALO algorithms (Hill 2012), implemented in the SPARTA (v0.1.30 August et al. 2015b) package of R (v3.4.0 R Core Team 2017). FRES-CALO weights by recorder effort to estimate trends and probability of occurrence in under-recorded areas (for validation of FRESCALO for different groups and through simulation see Hill 2012;Fox et al. 2014;Isaac et al. 2014;Dyer et al. 2016). We used the CEH Land Cover Map (LCM2007, Morton et al. 2011) as input data for FRESCALO's calculation of neighbourhoods of ecologically similar hectads (see August et al. 2015b;Dyer et al. 2016). For each species, FRESCALO produces a probability of occurrence per hectad. To transform this to presence/absence, we assigned a species as present in a hectad if its probability of occurrence was greater than a set threshold (see Appendix S3).

Constructing potential networks
We used lists of modelled plant and pollinator species presence per hectad to filter the interactions database (as derived from plant associations in biological records) and create a potential plant-pollinator network for each hectad (Fig. 1). Networks were unweighted (i.e. interaction matrices consisting of ones and zeros), this being the most conservative interpretation of our interaction data because the frequency with which an interaction was recorded is unlikely to provide reliable quantitative information on abundance, due to differences in detectability, recorder bias and data sources.

Network structure metrics
From the networks constructed in each of 2823 GB hectads, we used the R package bipartite (v2.07, Dormann et al. 2009) to calculate the following metrics (Bersier et al. 2002;Dunne et al. 2002;Tylianakis et al. 2007): (1) Species richness: total number of plant and pollinator species in the network (2) Connectance: proportion of possible links which are realised (3) Pollinator generality: mean number of plants per pollinator (4) Plant generality: mean number of pollinators per plant While other, more complex metrics of network structure (e.g. nestedness, modularity) have been implicated in stability they can be comparatively insensitive to spatial or temporal change (Kaartinen & Roslin 2012;Morris et al. 2014;Kemp et al. 2017). Preliminary analyses confirmed that nestedness and modularity showed little variation even for networks greatly differing in the metrics listed above.

Robustness to simulated extinctions
To derive a metric of network stability, we assessed the impact of simulated extinctions of plants. This was a measure similar to robustness (following Memmott et al. 2004;Burgos et al. 2007), but differing in that sequential simulated extinctions were ordered according to a complete 'global' list of plants (i.e. across all hectads) and not just those in each 'local' network (i.e. individual hectad) (Fig. 2). This approach meant the same extinction scenario was universally applied across hectads and enables comparisons across networks with different plant communities. We term this approach and the resulting metric 'robustness to global simulated extinctions' (R g from here on) to avoid confusion with the usual approach. For comparison, we also calculated 'local' robustness (R l from here on) following Memmott et al. (2004) with randomised extinction of plants within each hectad. We focussed on simulating plant extinctions because many of the major impacts of agriculture indirectly affect pollinators via altered plant communities (Potts et al. 2010; Vanbergen et al. 2013Vanbergen et al. , 2017, as would restoration of agricultural plant-pollinator networks in practice (Kremen et al. 2002;Forup et al. 2008;Menz et al. 2011;Kaiser-Bunbury et al. 2017).
As well as randomised plant extinctions from the global list (R g R ), we conducted simulated extinctions by ordering the complete list of 499 plant species according potential predictors of future plant declines under three scenarios: (1) historic distribution trend   For each hectad-level network, plants were sequentially extirpated from the global list in the order determined by each scenario. After each plant extinction, any pollinator species with no remaining links were removed from the network. We assessed R g as the area under the curve (Burgos et al. 2007) of pollinators remaining in the local network against plants removed from the global sequence ( Fig. 2). This process was repeated 100 times per scenario, with random ordering of plant species with tied trends or Ellenberg values (or of the entire list for randomised plant extinctions; R g R ). Following individual plant extinctions, there is potential for pollinators to switch between plants and rewire networks (Thomsen et al. 2017). Other authors have used regional information to inform the likelihood of rewiring in local networks (Kaiser-Bunbury et al. 2010). However, our potential networks already implicitly incorporate some of this capacity for rewiring because each local network contains information from all recorded interactions across GB over three decades. Another approach is to create putative novel interactions from plant-pollinator traits but, given the assumptions and uncertainties involved, it is difficult to assess whether such rewiring scenarios are more ecologically meaningful than using only observed interactions. We conducted supplementary analyses exploring additional trait-based network rewiring scenarios but found that spatial patterns in robustness metrics were largely unaffected (Appendix S4).
Two main sources of error in our hectad-level potential networks are the methods used to model occurrence and the database used to assign interactions. Consequently, we assessed the impact of these sources of uncertainty independently and in combination by performing, for each hectad, 100 randomised resamples of species according to FRES-CALO probability of occurrence, and of interactions proportional to the number of times they were recorded (see Appendix S5).

Crop-pollinators
We repeated our analyses by calculating hectad-level network metrics, including R g and R l , for potential networks consisting solely of interactions involving known crop-pollinators. This allowed us to explore whether the structure and stability of crop-pollinator networks was similar to the wider plant-pollinator networks in which they were embedded. Bees are generally considered the most important contributors to crop pollination (Free 1993; but see Rader et al. 2016) and their predominance in the pollination of GB crops is well supported (Woodcock et al. 2013;Garratt et al. 2014). Crop-pollinators were determined from a published list of bee species with the highest contribution to crop production value (Kleijn et al. 2015) for major GB insect-pollinated crops (oilseed rape, field bean, apple and strawberry). Our interactions database included 23 such species from 5 genera (see Fig. S5 for full species list). We then compared metrics for crop-pollinator networks, overall plant-pollinator networks and bee-only networks. Because crop-pollinator networks are considerably less speciose, we resampled the bee-only network for each hectad 100 times, with a random selection of pollinators equal in number to crop-pollinator species in the hectad, and then calculated mean resampled network metrics for comparison. The number of plants in each resampled network was allowed to vary depending on interactions with the selected pollinators, as attempting to limit plants to the number in the crop-pollinator networks would severely restrict the number of resampled networks and constrain resultant network metrics.

Statistical analysis
Network metrics were modelled independently against agricultural coverage, using linear mixed-effects models in the nlme R package (v3.1 Pinheiro et al. 2015). We derived coverage of agricultural land (arable + improved grassland) per hectad from LCM2007 and used this as a fixed effect explanatory variable, along with an optional quadratic term, which was retained in models if significant. To account for the potential influence of other environmental variables on network structure and response to agricultural coverage, we assigned each hectad to an environmental zone, using a pre-existing classification (Bunce et al. 2007). Environmental zone was then included as a random factor in all models, with variable slope and intercept. Some of the network metrics we used are sensitive to the size of the network (Jordano 1987;Olesen & Jordano 2002;Forup et al. 2008;Morris et al. 2014), so models were compared with and without total species richness as a fixed covariate. Environmental zones represented by < 30 hectads were considered to have insufficient sample size for robust analysis and were excluded, as were hectads with > 50% coverage of sea, giving a final sample of 2290 hectads. All variables were standardised to mean of zero and standard deviation of one, to facilitate comparison of model coefficients. Each model was compared using a likelihood-ratio test against a model consisting only of the random effect and species richness, to determine the impact of incorporating agricultural coverage on model fit. We also applied randomisation tests (Fortin & Jacquez 2000) to account for potential complex spatial autocorrelation patterns arising from how FRES-CALO defines neighbourhoods based on spatial proximity and biological similarity (see Appendix S6).

Spatial patterns in network properties
Variation in plant and pollinator species richness conformed to known clines across GB (i.e. higher richness in the south and at lower altitudes, Fig. 3). Spatial patterns for connectance, pollinator generality and plant generality showed similar latitudinal and altitudinal trends to species richness (Fig. 3), and a significant correlation with total plant-pollinator species richness (Pearson's r = À0.95; 0.82; 0.97, respectively, n =2290, P < 0.001 in all cases).
Robustness to global simulated extinctions (R g ) also showed spatial variation across GB, with more variation than would be expected under simple conformity to species richness, latitude or altitude (Fig. 4a, c, e, g and i). The different extinction sequences gave mean R g scores across hectads of 0.84, 0.92, 0.85 and 0.93 for extinctions ordered by trend (R g Trend ), Ellenberg N (R g N ), Ellenberg F (R g F ) and randomised extinctions (R g R ), respectively (range across hectads 0.66-0.97 across all four R g measures), but with varying spatial patterns (Fig. 4). Robustness to randomised local extinctions (R l ) showed a very similar range of values and spatial patterns to R g R (Fig. 4i).

Effect of agricultural land cover on network properties and robustness
Pollinator generality and all five measures of robustness to simulated extinctions showed significant positive relationships with agricultural coverage (Table 1). All relationships apart from R g F included a significant negative quadratic term (Table 1) indicating a levelling off of the relationship as agricultural coverage approaches 100% (Fig. 5a, b, c, e, f). These results suggest that pollinator communities in more highly agricultural landscapes are more generalist and that, under all our extinction scenarios, hectads with a higher coverage of agricultural land lost pollinators less quickly than other hectads in the same environmental zone. This effect appeared most pronounced for R g Trend and R g N (Fig. 5b, c) . Neither plant generality nor connectance showed a significant relationship with proportion of agricultural coverage (Table 1).
While species richness was a significant covariate in models for all network metrics, its inclusion did not qualitatively change the relationships with agricultural coverage. For all models, likelihood ratio tests and randomisation tests generally corroborated the significance of the individual model coefficients for agricultural coverage (Table 1). There was no significant relationship between agricultural coverage and species richness of plants, pollinators or both combined, once environmental zones were accounted for (Table 1).
There were significant relationships between agricultural coverage and the mean values across the plant community per hectad of the trait and trend values used to order extinction sequences (Table 1). This indicates that agricultural coverage influences the relative position of the plant community in our global extinction sequences. Standard deviations in these values within hectads showed a significant, negative relationship only for Ellenberg N, suggesting hectads with higher agricultural cover not only host communities with higher average fertility tolerance, but also show significantly less variation in fertility tolerance between plant species.

Crop-pollinator networks
Subsets of networks consisting of only crop-pollinating bees and the plants they visit showed significant differences in their properties from complete hectad-level networks, or from randomly resampled networks of equivalent bee species richness. Crop-pollinator networks showed significantly higher plant species richness than the randomly resampled bee networks (pairwise t-test; t = 133.57, P < 0.001, d.f. = 2750), as well as higher connectance and pollinator and plant generalities (pairwise t-test; t = 159.54, 155.48 and 116.35 for connectance, pollinator generality and plant generality, respectively; d.f. = 2750 and P < 0.001 in all cases). Robustness values for crop-pollinator networks were significantly higher than for full or resampled networks (Supplementary Material, Fig. S5), for all simulated extinction scenarios (pairwise t-test; t = 152.22, 121.83, 132.14, 161.45 and 155.72 for R g Trend , R g N , R g F , R g R and R l , respectively; d.f. = 2750 and P < 0.001 in all cases). Crop-pollinator species were among the most widely occurring species in the database (median occurrence for crop-pollinators = 60% of hectads, for all bees = 30%, for all pollinators = 43%).

Spatial patterns of plant-pollinator networks and relationships with agricultural land cover
Our results revealed that national-scale spatial patterns were clearly evident in all network metrics. Those for pollinator and plant generality and for connectance largely reflected well-known latitudinal gradients in GB plant and invertebrate species richness (e.g. Woodcock et al. 2014). This is Figure 4 Spatial patterns in robustness to global simulated extinctions (R g ) as measured by area under extinction curve. Panels a, c, e, g and i show R g with extinctions ordered by historic plant occurrence trend (R g Trend ), fertility tolerance (R g N ), drought tolerance (R g F ), globally randomised extinctions (R g R ) and locally randomised extinctions (R l ) respectively. Panels b, d, f, h and j show residuals from linear mixed models of R g Trend , R g N , R g F , R g R and R l , respectively, against species richness and environmental zone. Grey-shaded cells indicate environmental zones with < 30 cells excluded from mixed models. For all panels, darker colours indicate higher values, with a linear colour stretch between maximum and minimum values. unsurprising as these metrics are a function of the connectedness between the two levels of a network and thus fundamentally affected by network species richness (Jordano 1987;Olesen & Jordano 2002;Th ebault & Fontaine 2010).
Of greater interest in terms of implications for network stability were the more complex spatial patterns of our metrics of robustness to global simulated extinctions and the positive relationship with coverage of agricultural land evident under all our extinction scenarios. This positive relationship may, at first sight, seem surprising. Highly agricultural landscapes are often considered to have depauperate plant and pollinator communities (Kremen et al. 2002;Potts et al. 2010;Ollerton et al. 2014). As noted by Kleijn et al. (2015), this is often used as justification for pollinator conservation efforts under the assumption that continued crop pollination depends upon a diverse pollinator community. While our results cannot directly shed light on the provision of pollination services, it is clear that plant-pollinator networks in landscapes with relatively high agricultural cover can exhibit higher robustness to extinction scenarios.

Explaining higher robustness of plant-pollinator networks in agricultural landscapes
Although it is somewhat counterintuitive that increased levels of anthropogenic disturbance (coverage of agriculture here) can lead to increased resilience to future perturbations (as estimated by our robustness metrics), similar relationships have been observed in other ecosystems (e.g. coral reefs, Côt e & Darling 2010). This might be due to positive correlations between traits that confer tolerance to past and future disturbance (Vinebrooke et al. 2004). Exposure to previous stressors therefore acts as a filter either extirpating vulnerable species or favouring resistant ones to produce a community more resilient to future stress. Our results suggest that the positive relationships between agricultural coverage and robustness may arise in this way from two interacting properties of the plant-pollinator networks.
Under extreme circumstances, where networks have completely extinction-prone plant communities or completely resistant ones, differences in robustness to global simulated extinctions might be driven by the second effect alone, regardless of network structure. However, this is unlikely in our data given the significant relationship with generality and robustness under randomised extinctions. Also, all hectads possessed plant communities with varying positions in the extinction sequences. For example, Figs. 2a and 2c show extinction curves for the two hectads which were, respectively, most and least robust to extinctions ordered by plant trend. From the distribution of the tick marks denoting plant extinctions from the hectad on the horizontal axes it is clear that, while these two hectads have plant communities consisting of species with differing positions in the extinction sequence, neither hectad has all its plant species at either extreme.
Of course, the observed tendency of agricultural networks to require extreme plant extinction scenarios to collapse pollinator network structure does not mean that agriculture is without detrimental effects. Simple network metrics are insufficient to capture the myriad aspects of ecological stability Figure 5 Relationships between proportion of agricultural land cover for six network properties: a) pollinator generality b) robustness to simulated global extinctions ordered by historic plant trend (R g Trend ), c) robustness to simulated global extinctions ordered by Ellenberg N (R g N ), d) robustness to simulated global extinctions ordered by Ellenberg F (R g F ), e) robustness to simulated global extinctions ordered at random (R g R ) and f) robustness to simulated local extinctions ordered at random (R l ). All relationships are statistically significant (see Table 1). Slopes were back transformed from the full model (i.e. with agricultural land cover, species richness and environmental zone) and show the effect of agricultural land coverage on response variables once the effects of species richness and environmental zone are accounted for. (Grimm & Wissel 1997). While the networks of agricultural landscapes may be more robust to the scenarios we examined, they may also have lower levels of functional diversity. Potentially, they may also have lower functional resilience due to a homogenisation of species traits in response to the selective pressures of intensive agriculture (Woodcock et al. 2014;Oliver et al. 2015;Kaiser-Bunbury et al. 2017), as seen in our results for plant Ellenberg N.
In reality, extinctions are unlikely to proceed in a rigid linear sequence according to a single predictor. Extinction cascades (Vanbergen et al. 2017), rewiring (Thierry et al. 2011;Ramos-Jiliberto et al. 2012), climate change (Chen et al. 2011), disease (Smith et al. 2006) or invasive species (Bartomeus et al. 2008) can alter the stability of networks in unpredictable ways. However, our approach for calculating robustness to global simulated extinctions is sufficiently flexible that, where information on such effects exists, these could be incorporated into the extinction sequences.

Crop-pollinator network properties
Our results showed that crop-pollinator networks are significantly more robust to simulated extinction scenarios than the overall networks of which they are a subset. This is probably due to the observed ubiquity and high generality of crop-pollinator species. These characteristics might be expected, as GB crop-pollinators are by definition those species pre-adapted to exploit the resource of non-native agricultural crop species growing in highly modified landscapes. Our results support the contention of Kleijn et al. (2015) that strategies and initiatives based on conserving crop-pollinators will provide insufficient protection for wild pollinator communities overall. More generally, our results suggest caution where such functionally specific taxa are studied in isolation of the wider communities of which they are often only a small fraction. Obviously, crop-pollinators can be threatened by a wide variety of factors other than loss of nectar sources (Vanbergen et al. 2013;Gill et al. 2016;Potts et al. 2016). For example, preferential loss of crop-pollinators could be triggered if association with crops results in detrimental exposure to pesticides (Stanley et al. 2015;Woodcock et al. 2016).

Limitations of the potential network approach
Constructing potential networks from biological records has a variety of limitations and assumptions that constrain their biological realism (hence 'potential' networks) and affect the uncertainty of results. Perhaps the most obvious limiting factor in our networks are the biological records from which they are constructed. In particular, our data are affected by shortfalls in our knowledge of species occurrence and of their interactions (Hortal et al. 2015).
Regarding occurrence, although FRESCALO accounts for variation in recorder effort, there are likely to be remaining inaccuracies for rare or under-recorded species, while conversion of FRESCALO's probabilistic outputs to binary presence/absence values may also introduce errors, particularly at species' range boundaries. Regarding interactions, we know that our coverage of GB plants, pollinators and the interactions between them are incomplete (see Appendix S2). Our potential networks may exhibit either 'missing' or 'forbidden' links in some hectads (Olesen et al. 2010) as they do not account for variation in interactions due to flower phenology (Basilio et al. 2006;Rafferty & Ives 2011), pollinator life-history (Vieira & Almeida-Neto 2015;Vanbergen et al. 2017) or pollinator resource-switching (Thomsen et al. 2017).
Our exploration of some of these sources of uncertainty (see Appendices S4 and S5) suggests that uncertainties arising from occurrence and/or interaction data affect hectad-level networks in ways that are relatively consistent across space. While both sources of uncertainty affect the accuracy of individual potential networks they are far less likely to introduce a systemic bias which would affect our observed spatial patterns and relationships with agricultural land. Therefore, despite these limitations, we suggest that our potential networks properties and the spatial patterns we observe are broadly representative of real-world networks (see Appendix S5).

CONCLUSIONS
Our results demonstrate the ability of potential networks constructed from biological records to provide new insights into spatial patterns of ecological networks across national scales that would be impossible to monitor using conventional direct observation approaches. The positive relationship between agricultural cover and robustness to a range of extinction scenarios supports previous observations that anthropogenic disturbance can result in ecological networks which are more robust to further perturbation. Furthermore, from our results, crop-pollinator networks are not representative of wider plant-pollinator networks, such that targeting landscape management for the retention of crop pollination may be entirely insufficient to conserve wider biodiversity (Kleijn et al. 2015).
Our findings suggest potentially productive fields of further investigation, including further investigation of the mechanisms underpinning spatial patterns in network properties, validation of potential networks against those constructed from large-scale molecular data and exploration of more complex scenarios of extinction, invasion or restoration. In the future, the production of potential networks from biological records is likely to become easier and more accurate, as new technology and methods increase the quality and quantity of biological records (Tulloch et al. 2013;Gray et al. 2014;August et al. 2015a;Powney & Isaac 2015) and novel molecular techniques increase the potential for wide-scale validation (Keller et al. 2015;Richardson et al. 2015;Bohan et al. 2017;Pornon et al. 2017). assistance with compilation of the plant-pollinator interactions database; Oli Pescott and Tom Humphries for assistance with BSBI data; Gary Powney and Tom August for guidance on the use of FRESCALO in R; and Mike Edwards (BWARS) for advice on plant-pollinator interactions. We also extend our thanks to four reviewers who provided insightful and constructive comments on previous versions of this paper. This work was supported by NERC national capability funding to CEH and the Natural Environment Research Council (NERC) and the Biotechnology and Biological Sciences Research Council (BBSRC) under research programme NE/ N018125/1 LTS-M ASSIST -Achieving Sustainable Agricultural Systems (www.assist.ceh.ac.uk). AUTHORSHIP JWR performed the data handling, conducted spatial and statistical analyses, wrote the first draft of the manuscript and undertook revisions. THO and AJV coordinated the project under which this work was performed. BAW and MJOP supported ecological network analyses. RFP conceived and coordinated the compilation of plant-pollinator interaction records. All authors contributed to writing of the manuscript.

DATA ACCESSIBILITY STATEMENT
Plant-pollinator network data and metadata are available from the NERC Environmental Information Data Centre (EIDC) at doi.org/10.5285/6d8d5cb5-bd54-4da7-903a-15bd4bb d531b. Biological records are available on request from the relevant recording schemes.