Consequences of Uncertainty in Global-Scale Land Cover Maps for Mapping Ecosystem Functions : An Analysis of Pollination Efficiency

Mapping ecosystem services (ESs) is an important tool for providing the quantitative information necessary for the optimal use and protection of ecosystems and biodiversity. A common mapping approach is to apply established empirical relationships to ecosystem property maps. Often, ecosystem properties that provide services to humanity are strongly related to the land use and land cover, where the spatial allocation of the land cover in the landscape is especially important. Land use and land cover maps are, therefore, essential for ES mapping. However, insight into the uncertainties in land cover maps and how these propagate into ES maps is lacking. To analyze the effects of these uncertainties, we mapped pollination efficiency as an example of an ecosystem function, using two continental-scale land cover maps and two global-scale land cover maps. We compared the outputs with maps based on a detailed national-scale map. The ecosystem properties and functions could be mapped using the GLOBCOVER map with a reasonable to good accuracy. In homogeneous landscapes, an even coarser resolution map would suffice. For mapping ESs that depend on the spatial allocation of land cover in the landscape, a classification of satellite images using fractional land cover or mosaic classes is an asset.


Introduction
Global and continental-scale integrated assessments, like the Millennium Ecosystem Assessment and The Economics of Ecosystems and Biodiversity increasingly emphasize the importance of the provision of ecosystem services (ESs) for humanity [1,2].During the past decades, ESs have become increasingly important in policies and conservation, and sustainable use of ESs is now part of the targets set by the Convention on Biological Diversity (CBD) and the European Union [3].
Ecosystem services are commonly divided into provisioning services, regulating services, cultural services and supporting services [1].The provision of several regulating, supporting and cultural services depends strongly on the spatial allocation of land cover.For example, a network of hedgerows can regulate erosion and influences the spatial variability of soil fertility [4,5].Vegetation and relief control the water storage capacity of the landscape [6,7].Natural elements provide habitats for wild pollinators and for predators that prevent or diminish crop pests and diseases [8].Finally, the landscape structure can influence the attractiveness of the landscape for tourists and holidaymakers [9,10].However, land cover provides several services only to the direct surroundings, like erosion protection directly downslope of a hedgerow.
Mapping and modeling of ESs is an increasing field of scientific research at scales varying from landscapes [9,11] to the national scale [12,13] and globally [14,15].As a common approach to mapping ESs, relationships between ecosystem properties and ecosystem functions (EFs) or services are quantified and applied on a map of the ecosystem properties, often based on a land use or land cover map [16].For example, the distance between (semi-)natural elements and agricultural fields or the fraction of (semi-)natural elements in an agricultural area can determine the capacity of the ES provision by these natural elements.
To map ESs that depend on the spatial allocation of land cover, detailed data such as locations of hedgerows or small patches of natural growth are needed.At the landscape scale, the spatial organization of natural elements can be inventoried in detail [17] and translated into landscape measures.At the continental or global scales, direct field inventories are unfeasible, and landscape measures have to be calculated from continental or global-scale land cover maps, which are commonly derived from remote sensing (satellite imagery).Remote sensing-derived land cover maps are, however, uncertain in two ways: the shape or location of objects (geometric uncertainty) and the attribute values of objects (thematic uncertainty).
Geometric uncertainty in land cover maps arises from the spatial data resolution.Global land cover maps generally aim to map and monitor the spatial distribution of land cover [18,19].Because of the extent, the resolution is coarse.Major uncertainties caused by a coarse resolution are the poor representation of small landscape elements and minor land cover types [20].Thematic uncertainty arises because of errors and problems with translating the reflectance signature of the satellite image into a land cover classification [21][22][23].Both geometric and thematic uncertainty may propagate into erroneous outputs when remote sensing-derived land cover maps are used for ES mapping.Additionally, the resolution and extent used for mapping ESs influence the explaining factors and mapping patterns for soil organic matter contents and the level of herbivory and parasitism [16,24,25].Upon mapping ESs, the scale of input data should, therefore, match the scale of the processes studied.
In reality, especially at larger scales of analysis, the available data define the mapping scale [26]; however, the uncertainty and the effects of uncertainties are unclear.General information exists on the accuracy of the global-scale land cover maps used in global-scale assessments of ESs.Nevertheless, insight into the effects of geometric and thematic uncertainty of these input data for the mapping of ESs is lacking.In this study, we assess the utility of global and continental land cover maps by analyzing the uncertainties that arise upon mapping ecosystem properties and functions.As an example of an EF that depends on the spatial allocation of land cover, we assessed the efficiency of pollination by wild pollinators as derived from land cover patterns in the Netherlands.

Case Study
A number of crops providing a significant portion of the human food supply depend on pollination performed by insects like Honeybees and their wild relatives.Recently, the European Honeybees in Europe and North America have suffered large losses caused by colony collapse disorder [27].Wild pollinators are, therefore, seen as potentially important for the production of food crops [28][29][30].For efficient wild pollination, a landscape should provide a suitable habitat for wild pollinators.To map this EF, we first quantified the relation between the land cover structure and pollination based on a literature review.Previous research indicates that the abundance and species richness of wild pollinators decreases with increasing distance from the closest edge of natural elements or decreases with a decreasing percentage of natural land cover in the landscape.This can be translated into a percentage yield loss relative to a situation with optimal pollination [31][32][33][34][35].
We mapped the distance to natural elements in agricultural land (ecosystem property) and subsequently calculated the pollination efficiency (EF) in the Netherlands based on land cover maps at different scales (Figure 1).To analyze the effects of both thematic and geometric uncertainty on EF mapping, we used both a series of "real" land cover maps (Table 1) and a series of land cover maps based on aggregating the best available land cover map of the study area.

Land Cover Data
Five land cover maps are available for the Netherlands.The land use map of The Netherlands (LGN) was considered the most accurate representation of the actual land cover and was, therefore, used as a base map for all analyses.The LGN is a 25-m resolution map, based on a stratified multi-temporal classification of satellite imagery from several sources (Table 1).Ancillary data were used to classify agricultural land into land use.The overall classification accuracy of LGN at class level is 85-90% [22] (Table 1).The other maps that were analyzed are continental or global-scale maps (Table 1).All maps attribute a dominant land cover type to each grid cell.The extent of the Netherlands was extracted from each map, and the maps were projected using the Dutch National grid coordinate system.
To analyse the effects of geometric uncertainty, we created aggregated versions of the LGN with 100, 250, 825, and 1,100-m resolutions.First, we reclassified the LGN into a map indicating nature (1) and other land (0) and a map indicating agricultural land (1) and other land (0) following the classification in Table 2.Then, the area percentage of nature and agriculture were calculated for each grid cell at the output resolution (100, 250, 825, and 1,100 m).Grid cells with over 60% of nature or agriculture were classified as such.Grid cells with 20%-60% nature and agriculture were classified as mixed land.These thresholds are consistent with the description of mixed land cover classes in GLOBCOVER and GLC2000 (Table 2).The remaining area was classified as other land (Figure 1).The effects of geometric and thematic uncertainty are analyzed using area fractions of the main land cover types (Table 2) and the patchiness of natural areas.Patches of natural land were identified and the number of nature patches and the average patch size were calculated.For the LGN map, additionally, the median patch size and the number of patches were calculated for 10 × 10 km grid cells to provide a picture of the actual pattern of the land cover structure across the Netherlands.

Mapping of Ecosystem Properties and Functions
Pollination efficiency was mapped as a function of the distance to nature in agricultural land (D nature ), using each "real" land cover map (Table 1) and each aggregated LGN map (Figure 1(b)).To map the pollination efficiency, we classified the land cover maps into agricultural land, nature, mixed land and other land (Table 2).To map the D nature , each land cover map was reclassified into a map indicating nature (1) and other land (no data) and a map indicating agriculture (1) and other land (no data).Mixed land cover classes (Table 2) are included in both the nature map and the agriculture map.Then, the Euclidean distance to the centre of the nearest nature grid cell was calculated.In mixed land cover, the D nature was assumed to be zero.
Pollination efficiency was defined as the yield loss caused by diminished pollination compared to a situation with optimal animal pollination conditions.Klein et al. [29] provide a conceptual model of the relation between the D nature and the pollination efficiency, indicating maximum pollination efficiency close to nature.The pollination efficiency decreases to an asymptote with increasing D nature .To quantify this relation, we recalculated data on the percentage of fruits set as a function of D nature [33] into an average percentage yield loss relative to the yield at a zero D nature (yield reduction fraction, YRF): We used data for crops that depend on 25% of total pollination by animal pollination [39].The data given by Steffan-Dewenter and Tscharntke [33] range to 1,200 m.Beyond 1,200 m, the YRF is kept constant at 40%, based on the conceptual model of Klein et al. [29].Equation (1) explains 51% of the variance in the data [33].

Analyses
Each D nature map and each YRF map resulting from each land cover map was aggregated to average values for a 10 × 10 km grid, only considering the grid cells within agricultural land.Then, each 10 × 10 km map was compared with the base map, i.e., the D nature and YRF maps calculated from the original LGN map.Also, for each resolution a comparison was done between the D nature or YRF map based on the real land cover map (Table 1) and the D nature or YRF map based on an aggregated LGN map.Maps were compared visually, and a fuzzy numerical similarity statistic was calculated with the Map Comparison Kit (MCK; [40]).For each grid cell, the similarity between map a and map b was calculated as: Fuzzy numerical similarity is considered the best method to distinguish the spatial patterns of similarities and errors between maps [41].The overall similarity between a pair of maps was calculated as the average similarity statistic over all grid cells.A similarity statistic of 0.7 or higher was judged as good [42].

Land Cover Characteristics
The Dutch land cover is dominated by agricultural land, although approx.15% of the country is covered by natural or semi-natural areas (Table 3).Large areas of agricultural land are grasslands used as pasture or under rotation with arable crops.The land cover pattern is strongly related to the soil and parent material.The clay areas (Figure 2(a)) are characterized by homogenous, open landscapes with few, small patches of nature.In the sand areas, traditionally, mixed agriculture occur, which has resulted in a highly heterogeneous land cover.Peatlands are mainly used as fen meadow areas and are characterized by a high density of very small patches of nature.Main natural areas are found in the sand areas.Here, patch size and density are highly variable (Figure 2(b,c)).For the aggregated LGN maps, the percentage of nature decreases upon coarsening resolution (Table 3), while the percentage mixed land increases.The number of patches decreases and patch size increases.
In the "real" land cover maps (i.e., CORINE, GLOBCOVER, GLC2000 and PELCOM), the same trend is seen, but the decrease of patchiness is stronger than for the aggregated LGN maps.Additionally, the effects of thematic differences are observed.In GLOBCOVER, the percentages of nature and mixed land are large because some pasture areas are classified as a mosaic of shrubs and sparse vegetation.Consequently, more and smaller patches are seen in GLOBCOVER compared to the other maps.

Distance to Nature in Agricultural Areas
Large values of D nature are observed in the clay and fen meadow landscapes (Figure 2(a)) and in the northeast of the Netherlands (Figure 3).For the aggregated LGN maps, both the average D nature in agricultural areas and the variance increase with coarsening resolution (Table 3).
Although the spatial pattern of D nature remains largely intact, D nature increases throughout the Netherlands upon aggregating the LGN (Figure 3), resulting in large absolute differences with the base map.The D nature map based on an LGN aggregation at 100 m shows a high similarity with the base map throughout the country (similarity statistic 0.78; Table 4).Whereas, for the map based on aggregating the LGN to 250 m, the similarity is high in the clay landscapes and in the northeast of the study area (Figure 4(a)).In the heterogeneous land cover on the sand soils, areas with a large D nature emerge at a resolution of 250 m or coarser.In the real land cover maps, D nature is larger and more varied than compared to the aggregated LGN maps, except for GLOBCOVER (Table 3).The CORINE map results in spatial patterns similar to the base map (Figure 3), but values of D nature deviate strongly, resulting in a lower overall similarity statistic (Table 4).GLOBCOVER is more similar to the base map than the LGN aggregated at 250 m (similarity statistic 0.69).Overall, at each resolution the similarity between the real map and the aggregated LGN map is around 0.50, which is reasonable (Table 5).In the sand areas, GLOBCOVER performs better than the aggregated LGN, whereas in the fen meadow areas the aggregated LGN is more similar to the base map.GLC2000 and PELCOM result in large D nature values throughout the country and low similarities to the base map.Especially in the fen meadow areas and southwestern clay areas (Figure 2(a)) the D nature pattern from the base map are poorly reproduced by GLC2000 and PELCOM (Figure 3).The differences in the D nature patterns between the real maps are influenced by thematic differences.GLOBCOVER classifies parts of the fen meadow area as sparse vegetation or shrubland, resulting in a D nature that can be somewhat lower than the base map.Many heterogeneous areas in the sand landscape are classified as mosaic land cover in CORINE and GLOBCOVER, resulting in small D nature values.

Yield Reduction Fraction
In the base map, the YRF is low (indicating potential yield loss caused by a lower level of pollination) in the clay and fen meadow landscape, whereas a higher YRF is expected in the sand landscape (Figure 5).For the aggregated LGN maps, the average YRF decreases with coarsening resolution (Table 3; Figure 5), following the change in D nature upon coarsening resolution (Figure 3).As the YRF remains constant beyond a D nature of 1200 m, the increasing variance upon coarsening resolution in D nature (Section 3.2) is not seen in the YRF maps.In most of the Netherlands, all YRF maps based on aggregating the LGN compare well with the base map (Similarity statistic >0.7;Table 4; Figure 4(c); Figure 5).In the east of the study area, the 825 m and 1,100 m resolution maps show a lower similarity to the base map (Figure 4(c)).
In the real land cover maps, the decrease of the YRF with coarsening resolution is strong compared to the aggregated LGN maps (Table 3).Similarities with the base map also decrease with coarsening resolution.For GLOBCOVER, the high similarity of D nature results in a highly similar YRF (Table 4).Furthermore, the average YRF is high (Table 3).At each resolution, the similarity between the real land cover map and the aggregated LGN map is high.The similarities are all within the same order of magnitude (Table 5).
The spatial patterns of the YRF base map remain intact but become more pronounced in the CORINE map (Figure 5).Moreover, GLOBCOVER shows a similar pattern, except for the high YRF in the fen meadow areas (Figure 2(a)).In the YRF map calculated using GLC2000, some areas with a high YRF can still be identified, although the absolute YRF clearly deviates from the base map.The YRF map based on PELCOM shows hardly any spatial variation.
In most of the clay and fen meadow area, all YRF maps based on real land cover maps are similar to the base map (Figure 4(d)).In the remainder of the Netherlands, two of the maps (CORINE and GLOBCOVER) are similar to the base map.These patterns are opposite to the similarity patterns of the D nature maps.In large parts of the study area, aggregating the LGN did not affect the YRF maps (Figure 4(c)), while using another real map resulted in a considerably different YRF (Figure 4(d)).Only in the eastern part of the sand area (Figure 2(a)), using another map had the same effects as aggregation.

Effects of Resolution
The consequences of the input data resolution for mapping ecosystem properties and functions are derived by comparing D nature and YRF maps calculated from the series of aggregated LGN maps.For an accurate reproduction of the patterns of D nature in the base map, in heterogeneous landscapes a 100-m resolution map was required, while in homogeneous landscapes a resolution of 250 m was sufficient (Figure 4(a)).For an accurate YRF map, a high resolution was required only in the eastern part of the Netherlands (Figure 4(c)).These different input data requirements are caused by the variation in land cover structure throughout the study area and the process range of pollination.
Generally, when aggregating high-resolution land cover data, minor land cover types will disappear, while major land cover types will expand [43].Moody and Woodcock [44] found, however, that this effect is moderated by the spatial organization of the landscape.Minor but highly clustered land cover types will disappear more slowly upon aggregation than a dispersed class of the same size.Especially linear landscape elements, like ditch banks and tree rows, are sensitive to rapid disappearance in aggregation of land cover data [17,20].
In the clay and fen meadow areas, D nature is, in reality, large because natural land is organized into small, mainly linear, patches (Figure 2(c)).Consequently, natural elements quickly disappear upon aggregating in the landscapes, resulting in an increasing D nature .In the clay areas, the patch density is low (Figure 2(b)); thus the real D nature is large.The relative change of D nature upon aggregation is therefore small, resulting in high similarities up to a resolution of 250 m.In the fen meadow areas, patch density is higher (Figure 2(b)), resulting in a lower real D nature .Because the natural elements are too small and isolated to be supported at a coarser resolution, many natural elements disappear, resulting in considerable changes of the D nature .Therefore, only the 100-m resolution map had a high similarity.Throughout the sand area, D nature is very small because of the high percentage of natural area.This results in relatively large changes of D nature upon aggregation.
The YRF is less sensitive to map aggregation, as a consequence of the process range of pollination [26].Wild pollinators have a maximum foraging range of approx.1.5 km; therefore, wild pollination is assumed to be only effective within this range [29,45].Consequently, increases of D nature beyond 1,200 m will not propagate into changes of the YRF.
In the clay and fen meadow landscape, the D nature , in reality, sometimes exceeds 1,200 m.In these areas, the YRF does not change very much when D nature increases because of the aggregation of the land cover map.Based on the land cover structure, the sand area can be subdivided in two parts; in the east natural elements are small and highly dispersed while in the north, south and central sand area nature is more clustered into large patches.In the eastern sand area, the small natural elements are only sufficiently supported for mapping the YRF up to a resolution of 250 m.In the remainder of the sand area, the larger natural elements are supported to 1,000 m resolution upon aggregation.Consequently, a coarser land cover map suffices for mapping the YRF.

Effects of Satellite Imagery Resolution and Classification
The spatial resolution of the real land cover maps used in this study is defined by the resolution of the satellite images used for producing the maps and ranges from 25 m to 1,100 m (Table 1).At a coarser resolution, there is a larger chance that within one pixel more land cover types are found in reality.In such pixels, the dominant land cover type will dominate the reflectance characteristics and, consequently, classification into dominant land cover classes will result in the underrepresentation of minor land cover types.This explains the high percentages of agricultural land in GLC2000 and PELCOM.The underrepresentation is stronger than when aggregating a map with a higher resolution (Table 3) because small natural landscape elements in landscapes dominated by agriculture are supported in the original pixels.Upon aggregating, the small natural elements are partly included in the mixed land cover.The better match between D nature and YRF maps based on aggregated LGN maps with the base map suggests that capturing small natural elements is of great importance for mapping ecosystem properties and functions and that this can easily be achieved by including such elements into mixed land cover classes.
If no single land cover type dominates the reflectance characteristics of a pixel, classification is problematic.In the Netherlands, the land cover is highly heterogeneous with cropland and pasture being strongly intertwined [22].Consequently, in satellite images with a coarse resolution, a large amount of pixels will have a high level of heterogeneity.This effect is excluded when aggregating the finer data from the LGN map [43] and is the reason that the real coarse-resolution maps have a higher percentage of agricultural land and less patches of nature than the aggregated LGN maps.These differences in classification result in an additional deviation from the base map as can be seen in the lower similarity with the base map for both D nature and the YRF for CORINE, GLC2000 and PELCOM (Table 4).The magnitude of the deviations is similar for all the maps, as can be derived from the similarity statistics that are all in the same order of magnitude (Table 5).However, the direction of the effects might differ.Classifying heterogeneous areas with a lot of small nature patches as mixed land cover, as is done in GLOBCOVER, avoids both the underrepresentation of minor land cover types and avoids incorrect classification of pixels that are difficult to classify.This is demonstrated by the higher similarity of D nature maps based on GLOBCOVER and to a lesser extent CORINE with the base map in the highly heterogeneous landscape in the east of the Netherlands.
The real land cover maps used in this study differ in definition and classification of the land cover types (Table 2).In GLC2000 for example, forest is defined as areas with a tree cover of over 15%, while CORINE uses a 30% tree cover threshold [46].The classification is related to the extent of the map.Global-scale maps have a generalized legend where details on the land cover at the landscape scale are lost, while the LGN classification is specific to the Netherlands.To calculate landscape measures and pollination efficiency, the land cover maps have been classified into the main land use classes (Table 2).This decreases the thematic differences between the maps [21].The correct reclassification of the land cover map into the main classes required knowledge of the land use and land cover to correct errors.For example, many Dutch grasslands are classified as natural land cover in GLOBCOVER and GLC2000 while they are, in reality, used as agricultural pasture.The required knowledge of the study area might be unavailable at the global scale and might be location-specific.
Furthermore, the maps have different goals and priorities that influence the classification procedure.For example, CORINE is strongly ecology and landscape oriented [18].In the interpretation procedure, first, landscape units are delineated in classification, and ecological features are given priority [47].This results in a large area of nature (Table 3) that is organized in few, but large patches.Consequently, D nature is large and in large parts of the study area a low YRF is expected.

Methodological Issues
We assumed the LGN land cover map to be an adequate representation of the Dutch landscape.Although the LGN has an accuracy of 92% [36], small landscape elements are not correctly mapped at a 25-m resolution.In the fen meadow area, structures like ditches are underrepresented [17].Similar behavior is expected for hedgerows and tree lines.The "real" D nature might be smaller than calculated from the LGN.
All land cover maps have accuracy issues related to the satellite imagery equipment used for data acquisition.Additionally, there are differences in the time of data acquisition between the maps (Table 1) that will have caused deviations because of land use changes between the mappings.The global area-weighted accuracy of GLC2000 is 68.6%, and the accuracy is low especially in shrubland and herbaceous cover [21].Consequently, especially isolated heath areas have been wrongly classified in GLC2000.Because such areas are important habitats for wild pollinators, good classification is important for accurate mapping of pollination efficiency.This error has probably caused the underestimation of the D nature and YRF in the centre of the Netherlands.Clevers et al. [22] evaluated MERIS satellite images for land cover mapping for the Netherlands.Particularly arable land was very difficult to classify from MERIS data and was confused with all other land cover types, particularly grassland and built-up area.This might be a result of the lack of a short-wave infrared spectral band in the MERIS satellite; A SWIR band could improve the distinction of barren and sparsely vegetated areas.MERIS data have been used for the GLOBCOVER map.We suspect that parts of the areas wrongly classified in the study by Clevers et al. [21] are included in the mosaic land cover in the final GLOBCOVER map.Schmit et al. [26] assessed the accuracy of CORINE and PELCOM in Belgium, by comparing the maps with detailed land cover inventory data.CORINE had an agreement of almost 80%, whereas the accuracy of PELCOM was 50-60%.
Finally, a challenge in mapping EFs is distinguishing both land cover and land use.Satellite-derived land cover maps provide information about the actual cover of the land and not about the way the land cover is used [48].Although there is often a strong overlap between land cover and land use, it is not possible to make a one-to-one translation of land cover into land use.The use of the land, e.g., the intensity of management or the choice for a certain crop type, can influence the capacity of a land cover to provide an EF.For a more detailed mapping of EFs, thus, spatial information on both the land cover and the land use is needed.In our study, the lack of distinction between land cover and land use posed particular problems in classifying the grassland areas.An area that is classified as grassland based on satellite imagery can be natural grassland that provides pollinator habitat, while many Dutch grasslands are intensively managed and grazed.

Consequences for Mapping Ecosystem Services at a Global Scale
We mapped one ecosystem property and one EF for a small case study that does not contain the full variety of land cover as displayed on global land cover maps.We do, however, consider the case study representative of a larger area and for other EFs and ESs that depend on the spatial organization of the land cover within a landscape for several reasons.First, the Dutch land cover is highly heterogeneous [22], and, although not all land cover types occur in the Netherlands, the variation in D nature might be representative for larger areas.Second, in our model of pollination, we only distinguished between natural land cover and agricultural land cover and used a simple relation with distance to nature for mapping pollination.Although this is a strong simplification of the processes controlling pollination, it does capture the main effects of land cover on pollination efficiency correctly [28][29][30][31][32][33][34][35].A similar model applies for other EFs and ESs that depend on the spatial organization of the land cover, including pest control, mitigation of erosion and floods and air quality regulation.
Heterogeneous landscapes are hotspots for the provision of pollination.For mapping EFs that depend on the land cover structure, this heterogeneity should be captured in the land cover maps used for mapping EFs.Currently, the heterogeneity of a land cover map is strongly related to the accuracy and accurate mapping of heterogeneous areas is seen as a key challenge in land cover mapping [21].GLOBCOVER might be a map that gives an accurate picture of the heterogeneity; the GLOBCOVER map in Africa shows high levels of patchiness and heterogeneity [49].In this study, GLOBCOVER performed well in areas that are, in reality, heterogeneous.
To translate maps of EFs into maps of ESs, a link to human use of the EFs needs to be made [16].Thus far, we only included this by mapping the pollination efficiency in agricultural land.For a complete model the dependency of the crops on pollination by animals needs to be included, along with other factors that influence pollination like pesticide use and weather conditions.

Conclusions and Recommendations
In this case study, GLOBCOVER was the best land cover database for mapping ESs that depends on the land cover structure.GLOBCOVER has a patchier land cover in areas that are, in reality, highly heterogeneous, resulting in a large variety in D nature that matches the spatial pattern of D nature and YRF that is seen in the base map in most of the Netherlands.
For mapping the D nature , in a heterogeneous landscape, high-resolution maps were needed.In a more homogeneous landscape, the requirements for the resolution were somewhat lower.To distinguish between heterogeneous and homogenous landscapes, however, a high-resolution map was needed.
To identify areas that risk insufficient pollination, resolution was less important.In highly heterogeneous landscapes, a resolution of 250 m was sufficient for accurate mapping of both D nature and YRF.In homogenous landscapes for the YRF, a coarser resolution sufficed.Additionally, differences in input satellite data and goals of the map propagate into the utility of the land cover maps for mapping ecosystem properties and functions.
In our study, heterogeneous areas are hotspots for the provision of pollination.In heterogeneous areas, the natural land cover that provides EFs and the agricultural land use that uses the services are strongly intertwined.For accurate mapping of ecosystem services, the land cover map should accurately distinguish the different land cover types in heterogeneous areas or should, alternatively, provide a measure of the heterogeneity.Therefore, a classification with fractional land cover or including mosaic land cover classes or measures for patchiness would be an asset for mapping ecosystem functions and services at a global scale.

Figure 1 .
Figure 1.(a) Flowchart of the analyses; and (b) Procedure for aggregating the 25-m resolution Land Use Map of the Netherlands (LGN).

Figure 2 .
Figure 2. Landscape characteristics of the Netherlands.

Figure 3 .
Figure 3. Distance to nature in agricultural areas.(Top left) Base map; (Top) Maps based on aggregated LGN; (Bottom) Maps based on real land cover.

Figure 4 .
Figure 4. Number of maps that show a high similarity to the basemap.(a) D nature , aggregated LGN maps; (b) D nature , real maps; (c) Yield reduction fraction, aggregated LGN maps; (d) Yield reduction fraction, real maps.

Figure 5 .
Figure 5. Yield reduction fractions.(Top left) Base map; (Top) Maps based on aggregated LGN; (Bottom) Maps based on real land cover.

Table 1 .
Characteristics of the land cover maps.
* Only the number of land cover categories that occur in the Netherlands is given.

Table 2 .
Classifications of land cover maps into main land cover types.For the mixed land cover types, the percentage nature is given between brackets.
inside/outside urban area with high housing density; Marshes; Open/closed dune vegetation; Open drift sands; Heathland in dunes with no/some/strong invasion of grasses;Peat moors; Forest on peat moors; Other swamp vegetation; Reedlands; Forest in swamps; Other open natural vegetation; Bare soil in nature areas.

Table 3 .
Characteristics of land cover and ecosystem properties and functions for the study area in all land cover databases.

Table 4 .
Fuzzy numerical similarity statistics for comparison with the base map.

Table 5 .
Fuzzy numerical similarity statistics for comparison between real maps with aggregated LGN at each resolution.