Rebecca Spake
email: r.spake@soton.ac.uk

Rémy Lasseur
email: remy.lasseur@gmail.com

Emilie Crouzat
email: emiliecrouzat@gmail.com

James M Bullock

Sandra Lavorel
email: sandra.lavorel@ujf-grenoble.fr

Katherine E Parks
email: k.e.parks@soton.ac.uk

Marije Schaafsma
email: m.schaafsma@soton.ac.uk

Elena M Bennett
email: elena.bennett@mcgill.ca

Joachim Maes
email: joachim.maes@jrc.ec.europa.eu

Mark Mulligan
email: mark.mulligan@kcl.ac.uk

Maud Mouchet
email: maud.mouchet@mnhn.fr

Garry D Peterson
email: garry.peterson@su.se

Catharina J E Schulp

Wilfried Thuiller
email: wilfried.thuiller@ujf-grenoble.fr

Monica G Turner
email: turnermg@wisc.edu

Peter H Verburg
email: peter.verburg@vu.nl

Felix Eigenbrod
email: f.eigenbrod@soton.ac.uk

Unpacking ecosystem service bundles: Towards predictive mapping of synergies and trade-offs between ecosystem services

Keywords: cross-study comparison, ecosystem services, French Alps, land use, social-ecological systems, trade-off, natural capital, biodiversity

published or not. The documents may come

Introduction

Current understanding of how multiple ecosystems services (ES) are associated across heterogeneous landscapes remains limited [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF]Qui & Turner et al. 2013;[START_REF] Bennett | Linking biodiversity, ecosystem services, and human well-being: three challenges for designing research for sustainability[END_REF]. This understanding is essential for informing better management of multi-functional landscapes. Although the idea that the spatial distribution of ES and their associations are driven by the interplay between social and ecological variables is well-established [START_REF] Reyers | Getting the measure of ecosystem services: a social-ecological approach[END_REF], the relative importance of different social and biophysical drivers of sets of ES and how these change across different socio-ecological systems remains unclear [START_REF] Bennett | Linking biodiversity, ecosystem services, and human well-being: three challenges for designing research for sustainability[END_REF]. Consequently, there have been calls to achieve a greater understanding of the drivers of ES distributions and associations [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF][START_REF] Howe | Creating win-wins from trade-offs? Ecosystem services for human well-being: A meta-analysis of ecosystem service trade-offs and synergies in the real world[END_REF][START_REF] Bennett | Linking biodiversity, ecosystem services, and human well-being: three challenges for designing research for sustainability[END_REF].

Associations among ES are understood to occur when multiple services respond to the same driver of change or ecological process, or when interactions among the services themselves cause changes in one service to alter the provision of another [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF]. Such associations are commonly referred to as ES interactions [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF], with synergies and trade-offs being routinely explored in multi-ES assessments [START_REF] Howe | Creating win-wins from trade-offs? Ecosystem services for human well-being: A meta-analysis of ecosystem service trade-offs and synergies in the real world[END_REF]. Synergies arise when multiple services are enhanced simultaneously, while trade-offs occur when the provision of one service is reduced due to increased use of another. While ES associations can be highly context-specific [START_REF] Howe | Creating win-wins from trade-offs? Ecosystem services for human well-being: A meta-analysis of ecosystem service trade-offs and synergies in the real world[END_REF], there have been calls for the development of general rules about the relationships among ES [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF][START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF]. In attempting to distinguish ES associations that are context-specific from those that are universal, several authors have emphasised the need for cross-study comparisons (e.g. [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF]Raudsepp-Hearne et al. 2014[START_REF] Queiroz | Mapping bundles of ecosystem services reveals distinct types of multifunctionality within a Swedish landscape[END_REF]. However, cross-study comparisons are hampered by differences in approaches, the services covered, spatial scale, how ES are modelled and what drivers are used [START_REF] Grêt-Regamey | On the effects of scale for ecosystem services mapping[END_REF][START_REF] Queiroz | Mapping bundles of ecosystem services reveals distinct types of multifunctionality within a Swedish landscape[END_REF].

The concept of 'ecosystem service bundles' has been operationalised to help in the search for general rules determining ES associations [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF][START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF]. While rather confusingly the use of the term varies in the literature, with bundles and synergies used interchangeably [START_REF] Berry | Ecosystem services bundles[END_REF] see Box 1 for definitions used here), the term has been widely used in conjunction with the application of a spatially explicit framework developed by [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF] for identifying and mapping ES associations based on cluster analysis. [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF] defined ES bundles as coherent sets of ES repeatable in space or time. This clustering approach has been applied across the world to facilitate cross-study comparisons of ES associations and their drivers (Table 1;Fig 1

). Maps of ES bundles delineated with this approach can indicate what

services can be expected to associate based on where we find services repeatedly occurring together across a landscape [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF]). Their distributions have been typically interpreted with regards to known distributions of principal human activities or land use within the region (Table 1), and are therefore considered useful for communicating the potential impact of management decisions to policy-makers [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF]. This qualitative interpretation of ES bundle distribution provides some information about the drivers of ES associations and whether different social-ecological systems have particular sets of ES associated with them [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF].

In addition to qualitative interpretation of ES bundles, recent studies have attempted a more mechanistic approach to understanding ES bundle distribution, based on the relative roles of different social-ecological drivers, with multi-variate approaches being increasingly used [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF] reviewed the quantitative methods that are available for such analyses. [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF] suggested that spatially explicit analyses of the social-ecological variables driving ES bundles could ultimately allow for the prediction and modelling of ES bundles and thus, critical trade-offs and synergies across regions [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF]. Studies that aim to achieve such an understanding typically infer ES associations from the analysis of spatial trends in the distribution of two or more ES, and relate these to underlying social-ecological determinants [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF].

Further, if widely accessible data on social-ecological drivers (such as land use and population density) can predict ES bundles, this could potentially overcome problems associated with complex and data-intensive models that are required to produce ES maps [START_REF] Queiroz | Mapping bundles of ecosystem services reveals distinct types of multifunctionality within a Swedish landscape[END_REF]. Indeed, an ability to use limited variables to inform about the ES context is particularly important in data scarce regions [START_REF] Meacham | Social-ecological drivers of multiple ecosystem services: what variables explain patterns of ecosystem services across the Norrström drainage basin[END_REF]).

Here, we critically assess the strengths and limitations of current approaches for explaining and/or predicting the distribution of spatial associations between multiple ES. Most studies of this type to date follow the spatially explicit ES bundle approach first outlined by [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF] (Table 1). We first review studies that have applied this approach (Table 1; Fig. 1) and synthesise the application of it into four steps (Fig. 2), that capture the plurality of methods currently used, and illustrate them with a case studya cross-study comparison of the North and South regions of the French Alps. We then use the outcomes of this case study to assess the strengths and limitations of current approaches for linking social ecological drivers to ES bundles. . Finally, we outline a roadmap for research required to enable a general understanding of ES associations.

Box 1. Definitions of key concepts surrounding ecosystem services (ES) used in this article

ES associations

Arise when two or more services respond to the same driver of change or ecological process or when true interactions among the services themselves cause changes in one service to alter the provision of another [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF]. Commonly referred to as ES interactions [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF] and are inferred from spatial overlaps or lack thereof.

ES bundle

''Sets of ES that appear together repeatedly across space or time'' [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF]. Have been delineated and mapped using cluster analysis following Raudsepp-Hearne et al. 2010 (Table 1). In a bundle, ES can be positively (synergy) or negatively (trade-off) associated [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF].

ES demand "the amount of a service required or desired by society'' (Villamagna et al., 2013).

Different sectors of society can have different, and even conflicting demands.

ES flow "the service actually received by people, which can be measured directly as the amount of a service delivered, or indirectly as the number of beneficiaries served'' (Villamagna et al., 2013).

ES supply

The capacity of the structures and processes of a particular ecosystem to provide ES within a given time period (modified from [START_REF] Burkhard | Mapping ecosystem service supply, demand and budgets[END_REF].

ES use

Refers to an ecosystem being accessed/altered/managed/protected due to ES demand (Turkelboom et al. 2015).

ES indicator Proxy measures derived from empirical data or modelled estimates of ES.

Realised ES By definition, an ES is only realised if there is a human benefit. Without human beneficiaries and demand for an ES, ecosystem functions and processes are not services [START_REF] Fisher | Defining and classifying ecosystem services for decision making[END_REF].

Socialecological system

A set of social and ecological components that interact in a constantly evolving and interdependent manner [START_REF] Berkes | Linking social and ecological systems: management practices and social mechanisms for building resilience[END_REF].

Synergy

Arises when multiple services are enhanced simultaneously by the use of an ES.

Typically inferred from positive spatial overlaps.

Trade-off

When the provision of one service is reduced as a consequence of increased use of another, such as the case of crop production diminishing water quality. Inferred from negative spatial overlaps.

Win-win A situation (or area) where a synergy occurs. Step 1: Aggregation and harmonisation of ES data governance (in most areas of Europe) where many decisions regarding planning and landscape management are taken [START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF][START_REF] Queiroz | Mapping bundles of ecosystem services reveals distinct types of multifunctionality within a Swedish landscape[END_REF]. The selected grain for multi-ES research is also likely to have been driven by data availability; municipalities often are the finest scale at which some ES (such as provisioning ES) and potential social data are available (e.g. census data).

We consider the potential limitations of municipality-level analyses in the discussion.

Following collation and aggregation of multi-ES datasets, data are usually harmonised to a common range and unit to allow for comparison prior to data analysis. The methods used such as standardisation (transformation to z-scores by centring and scaling), serve to adjust the magnitude and variability of the variables to make them compatible for analysis (Legendre & Legendre 2012).

Application of step 1 to French Alps case study

The French Alps represent a relatively large, highly socially and ecologically diverse region characterized by excellent ES data over this large extent (e.g. [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF]. Within the region, elevation, climate and vegetation gradients have had historical influenced social dynamics and economic activities, resulting in the conventional separation into the North and the South Alps [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF]; a detailed description of study system is given in SI). This social-ecological divide is also recognised by an administrative boundary at the NUTS II level (Nomenclature of Territorial Units For Statistics by Eurostat [http://ec.europa.eu/Eurostat], basic regions for the application of regional policies).

We selected nine ES that have been quantified and mapped in the French Alps previously by [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF]. These services were deemed socially, ecologically, and economically relevant to the region following consultation with scientists and local collaborators [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF] (Box & Cox, 1964), centring and scaling was applied.

Step 2: Assessment of ecosystem service associations and delineation of ES bundles

ES associations have typically been assessed by mapping multiple ES across broad regions, and any spatial overlaps (or absence of overlaps) are assumed to signify a particular type of ES association;

positively correlated ES are assumed to be synergistic, while negative correlations infer trade-offs [START_REF] Tomscha | Ecosystem service trade-offs and synergies misunderstood without landscape history[END_REF]. Spatial overlaps between multiple ES have been most commonly quantified through assessments of pairwise correlations or PCA [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF]); a correlation biplot from a PCA (scaling type 2; [START_REF] Borcard | Numerical ecology with R[END_REF]) is considered a useful way to visualise the strength of correlations between multiple ES indicators (e.g. [START_REF] Maes | Synergies and trade-offs between ecosystem service supply, biodiversity, and habitat conservation status in Europe[END_REF][START_REF] Turner | Bundling ecosystem services in Denmark: Trade-offs and synergies in a cultural landscape[END_REF].

Raudsepp-Hearne et al. ( 2010) developed an approach for identifying ES bundles based on cluster analysis, which has since been widely applied to social-ecological systems across the world (Table 1; Fig 1). In this approach, clustering algorithms (e.g. k-means, self-organizing maps) have been applied

to define groups of ES that are associated in space by delineating spatial units supplying the same magnitude and types of ES [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF][START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF]). As such, ES bundles as defined by cluster analysis are emergent properties of the maps of different ES that are used in the cluster analysis and will often result from the distribution of underlying driver variables that drive more than one ES. Following clustering, ES associations have frequently been visualized using star diagrams [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF], showing the relative delivery of different ES within each bundle.

Clustering approaches also underpin many current methodologies for mapping social-ecological systems [START_REF] Ellis | Putting people in the map: Anthropogenic biomes of the world[END_REF]Asselan and Verburg 2012;[START_REF] Levers | Archetypical patterns and trajectories of land systems in Europe[END_REF], by identifying localities that have similar sets of multiple social-ecological variables.

Application of step 2 to French Alps case study

Following the spatially explicit ES bundle approach of Raudsepp-Hearne et al. ( 2010) we used kmeans cluster analysis to delineate ES bundles across the N and S French Alps separately (Full Methods in Appendix S2). Briefly, for both the North and South regions, a two step clustering approach was adopted [START_REF] Turner | Bundling ecosystem services in Denmark: Trade-offs and synergies in a cultural landscape[END_REF]. A PCA was firstly used to quantify the main multivariate relationships between the ES variables to assess whether ES co-occur in spatial bundles. As a precursor to cluster analysis, PCA can serve to separate signal from noise and lead to a more stable clustering solution (Husson et al. 2010). We applied k-means clustering to the relevant PCA axes (selected according to the Kaiser-Guttman criterion; [START_REF] Legendre | Numerical ecology, 3rd edition[END_REF][START_REF] Turner | Bundling ecosystem services in Denmark: Trade-offs and synergies in a cultural landscape[END_REF], to delineate ES bundles with 1000 random starts and 10,000 iterations to find a solution with the lowest within-cluster sum of squares according to the relevant PCA axes. K-means clusters municipalities so that the composition of ES values are more alike within than between clusters.

Following Renard et al. (2015), we quantified the effective number of ES provided in each bundle using a transformation (H) of the Gini-Simpson's index (S): H = 1/(1 -S), [START_REF] Jost | Entropy and diversity[END_REF]; Appendix S2).

In both the North and South Alps, three ecosystem service bundles (ESBs) were identified. In both regions, bundles were identified that were characterized by high crop production and far below average levels of most other services (ESB1(N) and ESB1(S)). Crop production was negatively correlated with most services across both study regions, except for water quantity regulation in the south (Appendix S2). In both the north and south, these crop-dominated bundles had the lowest diversity (H=2.8 for the north Alps and 1.8 for the south Alps). In the north and south regions, bundles were identified that were characterised by a high delivery of forest ecosystem services (carbon storage, wood production), and relatively high provision of other services but a complete lack of crop production (ESB2(N) and ESB2(S)). These forest ES-dominated bundles had the highest diversity in both the North and South regions. In the North Alps, multifunctionality was higher (H=9.0) than in the South Alps (H=6.0). A third ESB had a more intermediate mix of ecosystem services in the north and South Alps. In the north, ESB3(N) exhibited intermediate levels of crop production while remaining relatively multi-functional, delivering other services including high levels of tourism and intermediate hunting and recreation (Fig. 3; H = 6.9). In the South, ESB3(S) was dominated by delivery of hunting, erosion mitigation, and carbon storage (H=5.9; Fig. 3).

Step 3: Identification of social-ecological determinants of ES bundles

Understanding the spatial distribution of ES associations means identifying key drivers and their interactions that produce coherent sets of ES across landscapes [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF][START_REF] Meacham | Social-ecological drivers of multiple ecosystem services: what variables explain patterns of ecosystem services across the Norrström drainage basin[END_REF]. Several studies have mapped ES associations to allow for their qualitative interpretation by association with broad social-ecological systems (Table 1). The results of cluster analysis are made spatially explicit when the spatial units (typically administrative units or grid cells, In addition to qualitative interpretation, several quantitative methods are available for analysing ES bundles in relation to potential social-ecological determinants or predictors [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF].

Widely used methods include those frequently used in community ecology to study the relationships between ecological communities and the environment, through the coupling of two data tables, a site × environmental variable table and a site × species table [START_REF] Doledec | Co-Inertia Analysis -an alternative method for studying species environment relationships[END_REF] 1 for a review of quantitative methods for identifying drivers of ES associations.

Application of step 3 to the French Alps case study

In To identify candidate variables significantly affecting the co-variation of multiple ES, we performed a preliminary redundancy analysis (RDA) with all potential social-ecological driver variables followed by forward stepwise selection to select the model with the combination of variables with the highest R² and p-value [START_REF] Legendre | Numerical ecology, 3rd edition[END_REF]. This stepwise procedure defined which variables are relevant in exploring relationships among ES. RDA and the stepwise selection of variables were performed using the "vegan" and "packfor" R packages (Oksanen et al., 2013;[START_REF] Dray | packfor: Forward Selection with permutation (Canoco p.46)[END_REF].

RDA revealed that the combinations of the following variables significantly explained the co- (2016) also assessed the spatial overlap with 'well-being bundles', as identified using cluster analysis of social and demographic factors such as income and education.We include this last step, as it represents a logical progression from testing the relative predictive power of individual socialecological variables.

variation

Application of step 4 to the French Alps case study

We followed the approach of [START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF]Hamann et al. ( , 2016) ) to ascertain whether ES bundles were congruent with social-ecological systems. Having identified the most important social-ecological determinants of ES bundles in step 3 using RDA, we used the k-means algorithm to cluster these variables into social-ecological bundles (SEBs). SEBs delineate spatial units supplying the same magnitude and types of social-ecological variables. [START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF] found that anthromes offered little predictive power for provisioning service bundles in South Africa. We therefore used the variables deemed important from the RDA to delineate SEBs, as opposed to those used in the original construction of anthromes (Ellis & Ramunkutty, 2008).

To assess whether particular ES bundles are associated with SEBs, or whether SEBs can act as proxies for ES bundles, the spatial congruence between SEBs and ES bundles was assessed using overlap analysis, a simple and intuitive way to run a spatially explicit detection of possible associations [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF]). We calculated overlap as the percentage of municipalities of a particular bundle category that overlapped with each SEB category. any SEB (Fig. 5), except in the north-east of the region (Fig. 4), dominated by high elevation grasslands and semi-natural areas with high levels of protected area coverage. However, in the South, the forest bundle (high wood production and carbon storage) (ESB2(S)) does not overlap with forest cover, but with high elevation areas with grassland and semi-natural coverage).

Discussion

A multitude of methods are available to analyse and explore ES associations relative to possible social-ecological predictors [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF]. Here, we have reviewed the application of a widely (Fig 1 ) and increasingly used (Table 1) method that analyses the spatial distribution of ES bundles, delineated by cluster analysis, in relation to possible socio-ecological predictors. A common theme across all such studies is the reliance on the spatial coincidence of ES and driver variables [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF], assuming that consistency in the spatial congruency between ES likely emerges from common social-ecological drivers. While comparison among multiple studies, such as cross-site comparisons, could help disentangle the effect of context-dependent drivers from interactions between services within bundles [START_REF] Queiroz | Mapping bundles of ecosystem services reveals distinct types of multifunctionality within a Swedish landscape[END_REF] Thereofore, studies that have bundled different ES, or measured or modelled ES in different ways, are not straightforward to compare, or necessarily generalisable to other regions. We attempted to overcome both issues in our French Alps case study by comparing two regions using the same ES and social-ecological datasets, and do so using a widely used method to analyseES bundles. However, we show that even within the French Alps, there is enormous variation in the degree to which different social-ecological variables can explain the distributions of ES bundles (See Appendix S4 for more discussion on the findings from the case study).

Importantly, our case studywhich is based on the current state of the sciencedoes not enable us to identify why the explanatory power of different social and ecological variables considered here differs so much between our two regions. As such, current approaches based on readily available data that may have little relationship to underpinning mechanisms may not provide an effective basis for predicting ES bundles across space or time, as is required for effective sustainable management of ES.

Here we discuss why current approaches for analysing ES bundles are poorly suited to enabling sound understanding and prediction of ES bundles and propose a roadmap to guide future studies aimed at understanding, mapping or predicting ES associations.

Issues of scale in understanding determinants of ES associations

Here we detail issues of scale related to the ES bundle approach. We address two key components of scale: i) grain, the size of the spatial unit of analysis; and ii) extent, the size of the study area.

Spatial unit and grain

ES associations are often analysed using municipalities or similar administrative spatial units (e.g. 1), justified by the fact that municipalities are expected to be a grain at which synergies and trade-offs between ES are observed (Rodríguez-Loinaz et al. 2012), and as while ES synergies and trade-offs can be causally linked, they do not necessarily occur in close proximity [START_REF] Berry | Ecosystem services bundles[END_REF]. However, municipality boundaries could be relevant for some ES, such as cultural ES, but arbitrary for others in management terms, such as for managing water quality.

Raudsepp-Hearne et al. 2010; Table

Boundaries may often dissect ecologically meaningful units, such as watersheds, that could be appropriate for measuring and managing some ES.

The choice of municipality-level analysis is also often driven by data availability; municipalities often are the finest scale at which some ES (namely provisioning ES) and social variables are available (census data). Despite some good reasons for municipality-scale analyses, several considerations must inform their interpretation. At such coarse scales, the identification of ES bundles relies on spatial coincidence [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF], and cannot show direct causal relationships between ES and socialecological variables. This is a key assumption with the approach; that consistency in the spatial congruency between ES likely emerges from common social-ecological drivers. In actuality, the finescale processes that some ES respond to might not be represented at this scale.

As one moves across different grain sizes, different processes are responsible for apparent synergies and trade-offs between ES and relationships to social-ecological drivers. At coarse grains such as municipalities, spatial units are highly heterogeneous, encompassing multiple LULC types. ES relationships are likely to be largely driven by fractional land cover of the large spatial units, due to its representation of i) natural conditions; e.g. natural land cover and soil conditions as well as ii) human impacts; mainly via land use [START_REF] Burkhard | Mapping ecosystem service supply, demand and budgets[END_REF]. ES relationships will, therefore, principally reflect land use distribution. For example, ES may trade-off against each other simply because they compete for space (e.g. a negative relationship between timber and crop production; [START_REF] Lautenbach | Quantifying ecosystem service trade-offs[END_REF]. At smaller grain sizes, where individual spatial units are less heterogeneous and likely to comprise a principal land cover type, the main drivers of ES variation are still likely to be land use. If

ES within a single land cover type are analysed at small grains, however, such as individual forest plots or stands, then it is possible that a more useful understanding might be obtained. By analysing a single land cover type, one can understand drivers of ES variation in relation to land use activities that result in 'land modifications', changes that occur within the same LULC type (e.g. [START_REF] Lavorel | Using plant functional traits to understand the landscape distribution of multiple ecosystem services[END_REF]). These remain much less studied than multi-ES relationships to LULC, but data are becoming increasingly available (Erb et al. 2016).

Another well-documented scale effect related to spatial unit is the modifiable area unit problem, in which statistical results can depend on the size and shape of spatial units in which a variable is aggregated [START_REF] Openshaw | A million or so correlation coefficients: three experiments on the modifiable areal unit problem[END_REF]. Grain size-dependence in the direction of correlations of ES has been demonstrated in several studies (e.g. Naidoo et al., 2008;[START_REF] Anderson | Spatial covariance between ecosystem services & biodiversity priorities[END_REF]. Various processes can cause this phenomenon. Aggregation obscures ES trade-offs particularly when ES compete for space. For example, different crop competing for productive floodplain soils could be seen as spatially concurrent in aggregated datasets, thereby suggesting a synergistic relationship [START_REF] Tomscha | Ecosystem service trade-offs and synergies misunderstood without landscape history[END_REF].

When administrative units are used, the degree of variation in the grain size among units is likely to be an issue for the interpretation of relationships, as the mechanisms essential to an ES at one grain can be less important or absent at another. Significant variation in areal size could then reduce the specificity of the measured associations, and also decrease their strength [START_REF] Arsenault | How to choose geographical units in ecological studies: Proposal and application to campylobacteriosis[END_REF]. Such a phenomenon could affect the apparent relationships between ES or social-ecological variables, e.g. population density could appear to be inversely related to landscape multi-functionality, but in actuality, this could be a function of municipality size, as densely populated areas often divided into smaller administrative units for health care and mail delivery [START_REF] Arsenault | How to choose geographical units in ecological studies: Proposal and application to campylobacteriosis[END_REF] Peterson (2016).

Study spatial extent and context-dependency

The spatial extent of the study region can impact ES relationships. At present, most studies have delineated ES bundles at regional scale (Table 1), likely due to data availability, but also due to the relevance to management of considering variation in ES bundles across municipalities within a region. However, regions will differ in the variability of both the ES and social-ecological drivers that may underpin these ES, as seen, for example in our case study, confusing our results. The relative importance of social-ecological variables in driving ES variation can change across regions, and therefore study extent. For example, [START_REF] Holland | Spatial covariation between freshwater and terrestrial ecosystem services[END_REF] found a negative relationship between agricultural production and river habitat quality at the extent of Britain, due to the negative effects of agriculture on aquatic ecosystems. However, within some heavily urbanized sub-regions of Britain, a positive relationship was observed; this was attributed to urban land cover having a larger negative effect on aquatic ecosystems than agricultural land. Variability of predictor and response variables also affects the degree of statistical power that is available to detect relationships between spatial variables [START_REF] Eigenbrod | Sub-optimal study design has major impacts on landscape-scale inference[END_REF]. Moreover, the types of social-ecological driver variables considered will likely vary with spatial extent. For example, over larger study regions, it is possible to analyse the effect of slow variables, that exhibit variation at larger extents, but remain homogeneous across spatial units at small extents. Given these issues, cross-study comparisons will not necessarily enable meaningful comparisons of the relative explanatory power of different drivers between regions, even when the same ES and the same explanatory variables are considered (as in this study).

Careful selection of ES indicators in multi-ES analyses is critical for interpretation

The studies that have delineated ES bundles based on spatial associations in indicators make attribution and prediction difficult. Firstly, because the ES indicators may be anywhere along a spectrum from ecological stocks to flows to benefits in support of human wellbeing, some ES indicators may not respond to the influence of social factors [START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF].

Indeed, supply and demand bundles are likely to exhibit very different dynamics and respond to different drivers, potentially making mixed-indicator bundles more difficult to interpret or predict, as in this and previous studies [START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF][START_REF] Meacham | Social-ecological drivers of multiple ecosystem services: what variables explain patterns of ecosystem services across the Norrström drainage basin[END_REF]

Careful selection of social-ecological variables in multi-ES analyses is critical for attribution

There have been several calls for ES analysts to improve understanding of ES associations, to allow for knowledge of how to minimize trade-offs and enhance synergies [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF][START_REF] Bennett | Linking biodiversity, ecosystem services, and human well-being: three challenges for designing research for sustainability[END_REF]. This understanding requires identifying key social-ecological variables that determine the co-variation in ES. Other authors have suggested the potential benefit of predicting ES associations from widely available social-ecological datasets, that are not necessarily causal [START_REF] Queiroz | Mapping bundles of ecosystem services reveals distinct types of multifunctionality within a Swedish landscape[END_REF]. If widely accessible data on social-ecological drivers (such as land use and population density) can predict ES associations, this may overcome problems associated with complex and data-intensive models that are required to produce ES maps in data scarce regions [START_REF] Meacham | Social-ecological drivers of multiple ecosystem services: what variables explain patterns of ecosystem services across the Norrström drainage basin[END_REF]. While causal relationships are predictive (within similar contexts), prediction of ES associations does not necessarily require causative links. We emphasise that causal social-ecological predictors for multi-ES analysis are likely to be more robust and less-context dependent (see also [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF].

Land-use change is a management intervention that can drive demand and supply in one or more ES [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF], and therefore land use/land cover (LULC) has been considered as a determinant of individual ES or ES bundles in this study and many others (e.g. [START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF][START_REF] Meacham | Social-ecological drivers of multiple ecosystem services: what variables explain patterns of ecosystem services across the Norrström drainage basin[END_REF][START_REF] Schulze | Assessing Regional-Scale Impacts of Short Rotation Coppices on Ecosystem Services by Modeling Land-Use Decisions[END_REF]. There are several issues with using LULC as a determinant in multi-ES analyses. In this study and others, land cover categories were treated as homogeneous across study regions, ignoring significant variations due to management and biophysical gradients (e.g. tree species and age structure in forests). In our study, forest cover was correlated with forest services (wood production, carbon storage) in the North (Figs 3 and S6), but not in the South (Figs 4 and S7). This is because the French South Alps experienced extensive afforestation during the last century due to both natural regeneration and deliberate planting on abandoned agricultural land. The secondary forests are not widely harvested because their uniform and dense structure makes cutting expensive, and because local populations are concerned for their conservation [START_REF] Douguédroit | Reafforestation in the French South Alps[END_REF]. By using forest cover as a driver, we gained no fine understanding of ecological processes and interactions.

Issues relevant to using cluster analysis for modelling ES associations

Cluster analysis is considered a useful first step when no prior knowledge about existing relationships in a multivariate dataset exists [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF][START_REF] Dheng | A review on trade-off analysis of ecosystem services for sustainable land-use management[END_REF]). However, its exploratory nature makes it unsuitable for understanding causality in ES assocations. Cluster analysis requires somewhat subjective decisions including the clustering algorithm and the number of clusters, which is not straightforward (Legendre & Legendre 2012). The clustering solution is also entirely dependent on the input variables, rendering the results ungeneralizable to other regions. In summary, the subjectivity of cluster analysis makes it poorly suited to cross-study comparisons that are required for understanding general socio-ecological causes of ES associations. This will likely have led to the poor congruence between ES-bundles and social-ecological bundles as found in this study (Fig. 5). Maps produced in this way should therefore be used with caution when presented to stakeholders. The 'air of authority' [START_REF] Hauck | Maps have an air of authority: potential benefits and challenges of ecosystem service maps at different levels of decision making[END_REF] imparted by these maps and their associated star diagrams completely mask any uncertainty and could lead to erroneous management decisions.

Summary: ES bundles display pattern-based multifunctionality, but not process-based multifunctionality

The visualisation of relationships among multiple ES is considered a challenge to ES analysts (Birkhofer et al., 2015) and for effectively communicating with policy makers [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF].

Maps of ES bundles are therefore useful for visualising the joint spatial distributions of multiple ES.

They can be used to identify 'pattern-based multifunctionality', the joint supply of multiple ES in space, without regard for the ecological processes underlying the pattern [START_REF] Mastrangelo | Concepts and methods for landscape multifunctionality and a unifying framework based on ecosystem services[END_REF], and help guide land management decisions, such as where to allocate urban development or prioritise conservation efforts. This is possible when the scale of analysis (spatial unit type, grain and extent) are close to the desired scale required by key stakeholders [START_REF] Scholes | Multi-scale and cross-scale assessments of social-ecological systems and their ecosystem services[END_REF]. We suggest that analyses that wish to map ES bundles compare multiple scales corresponding to a portfolio of management policies [START_REF] Qiu | Spatial fit between water quality policies and hydrologic ecosystem services in an urbanizing agricultural landscape[END_REF], focussing for example on biophysically bounded spatial units such as watersheds of different size (e.g., [START_REF] Qiu | Spatial interactions among ecosystem services in an urbanizing agricultural watershed[END_REF].

However, while such correlational analysis is a logical first step in assessing ES associations, it cannot allow for a mechanistic understanding [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF]. When ES bundles are delineated using correlation at coarse resolutions, with spatial units exhibiting high within-unit heterogeneity in land cover and thus ES, and with each ES mapped at the same resolution and extent, the approach cannot help ES analysts understand general rules of mechanistic relationships between key drivers and ES.

They therefore cannot provide 'process-based multifunctionality', the joint supply of ES in space caused by well-understood relationships [START_REF] Mastrangelo | Concepts and methods for landscape multifunctionality and a unifying framework based on ecosystem services[END_REF]. Such a mechanistic understanding of relationships between ES and management will allow the transferral of management recommendations outside the context where data were collected (Birkhofer et al., 2015).

A roadmap for predictive mapping of bundles of ecosystem services

Determining the cause of a relationship among ES based on studies that track only their spatial concordance is difficult [START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF]). Here, we outline three key requirements for improvements to current approaches to understanding and predicting ES associations. The theme that underlies all these requirements is that studies that aim to explain or predict associations between ES must be designed to have a clear mechanistic basis in order to be confident about any relationships found. recent framework and set of specific predictions about how habitat fragmentation will affect ES provides an excellent example of the types of clearly defined questions that are required for a predictive science for ES. The need for formulating specific questions and hypotheses in ES research is also relevant to the generation of policy-relevant knowledge. Indeed, designing problem-oriented ES assessments, which focus on the information demands of decision-makers, can help make ES studies more decision relevant [START_REF] Förster | Assessing ecosystem services for informing land-use decisions: a problem-oriented approach[END_REF][START_REF] Willcock | Do ecosystem service maps and models meet stakeholders' needs? A preliminary survey across sub-Saharan Africa[END_REF].

4.1

Requirement 2: The testing of specific research questions requires bespoke study designs

Observational studies of the relationships between ES and their drivers are unlike experimental studies in that the identity, crossing, replication and interspersion of variables are, by definition, outside the control of the observer [START_REF] Smart | Clarity or confusion? - Problems in attributing large-scale ecological changes to anthropogenic drivers[END_REF]). Careful study designs can help to deal with these challenges and generate meaningful tests of very specific and focused predictions about relationships between ES. Here, ES science should build on the large literature examining the effects of habitat loss and fragmentation on biodiversity (Fahrig, 2003;[START_REF] Mcgarigal | Comparative Evaluation of Experimental Approaches To the Study of Habitat Fragmentation Effects[END_REF]. Of key importance is the need to account for habitat amount before considering effects of habitat configuration when attributing effects. For example, Qiu & Turner (2015) examined whether adding configuration variables could significantly improve the explanatory power of models explaining water quality after accounting for the effect of composition. Using this two-step procedure, they found forests to be more effective at retaining nutrients when more dispersed across subwatersheds.

One major consideration in designing studies to test predictors of relationships between ES is the issue of scale (section 3.1). Multi-scale assessments of social-ecological relationships with individual ES are vital to understanding scale-dependent social and ecological processes and causality [START_REF] Scholes | Multi-scale and cross-scale assessments of social-ecological systems and their ecosystem services[END_REF][START_REF] Eigenbrod | Redefining Landscape Structure for Ecosystem Services[END_REF]. Multi-scale assessments may not be possible, for example when the highest spatial resolution of the data is the municipality as with census-derived socioeconomic variables [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF][START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF][START_REF] Queiroz | Mapping bundles of ecosystem services reveals distinct types of multifunctionality within a Swedish landscape[END_REF]. Recent developments in downscaling or disaggregating datasets hold promise for higher resolution analyses with available datasets (e.g. [START_REF] Keil | Downscaling the environmental associations and spatial patterns of species richness[END_REF][START_REF] Lamboni | Multi-scale land-use disaggregation modelling: Concept and application to EU countries[END_REF].

Requirement 3: Utilize a wider range of statistical and modelling approaches

While statistical techniques cannot compensate for poor study design (e.g. [START_REF] Hurlbert | Pseudoreplication and the design of ecological field experiments[END_REF], taking advantage of the best statistical approaches will maximize the inferential strength of a given study design. As such, a predictive science for ES should take advantage of recent advances from ecological modelling including models that take account of biases in data, confounding variables, and mechanistic relationships (e.g. Sugihara et al. 2012;[START_REF] Warton | So Many Variables: Joint Modeling in Community Ecology[END_REF]).

One approach with potential to provide major insights in refining our hypotheses about how different predictor variables may affect relationships between ES is simulation modelling. For example, the creation of artificial landscapes could enable researchers to control and tease apart variables that are inherently confounded in real landscapes. Such studies have led to major insights in landscape ecology (e.g. With and King 1997;Gardner et al. 1989), macroecology (e.g. Lennon, 2000), but also in our understanding of how landscape structure might affect ES at different spatial scales [START_REF] Mitchell | Reframing landscape fragmentation's effects on ecosystem services[END_REF]. Simulation models can also be linked with future scenarios in which effects of changing drivers, such as land-use patterns and climate, on spatial dynamics of ecosystem services are explored (e.g., [START_REF] Carpenter | Plausible futures of a social-ecological system: Yahara Watershed[END_REF].

The use of primary data or process models rather than land cover based proxies

A major issue for understanding causal drivers of relationships between ES is that most available maps of ES are themselves modelled rather than measured. For example, regulating services such as pollination and erosion mitigation are typically and necessarily quantified using models that incorporate causal relationships between social-ecological variables (Martínez-Harms & Balvanera, 2012). An element of circularity therefore exists in ours and most other studies from having assessed the relationship between social-ecological variables and modelled surfaces of ES derived from exactly such variables. As such, a true understanding of determinant predictors of ES will only come through increased availability of primary data on actual services rather than LULC surrogates, including from remote sensing [START_REF] Ayanu | Quantifying and mapping ecosystem services supplies and demands: a review of remote sensing applications[END_REF]) and field studies that measure ES indicators such as water quality and carbon storage. That said, understanding the degree to which widely accessible social-ecological data can be used to predict ES associations, composed of ES that are either data-intensive or complex to model is still useful [START_REF] Meacham | Social-ecological drivers of multiple ecosystem services: what variables explain patterns of ecosystem services across the Norrström drainage basin[END_REF], as it facilitates modelling of such ES associations in data-poor regions.

The consideration of temporal changes in ES and drivers

Inferring interactions from spatial co-incidence is loosely analogous to a space-for-time substitution in that spatial relationships are used to infer dynamics over time [START_REF] Tomscha | Ecosystem service trade-offs and synergies misunderstood without landscape history[END_REF]. A major limitation of this approach is that most spatial studies use ES snapshot data to assess ES associations and relationships with drivers. Mismatches in the timing between change in a driver (including demand) and the supply of an ES may cause relationships to be misinterpreted or overlooked, particularly in transitioning landscapes. This can also be due to mismatches in the time series of available datasets. ES are not static but spatially and temporally dynamic in terms of their delivery and associations with other services.; Municipalities have been found to change in the bundles of services they provide over time raising concerns about using snapshots of ES provision to build understanding of ES relationships in complex and dynamic social-ecological systems (Renard et al. 2015),. Long-term monitoring studies could potentially capture complex long-term ES interactions and help us avoid or minimize trade-offs and adequately track synergies that simultaneously support multiple ES [START_REF] Tomscha | Ecosystem service trade-offs and synergies misunderstood without landscape history[END_REF]. spatial context that is relevant for decision making. The North and South Alps therefore lend themselves well to a cross-study comparison.

Table S1.

Figure 3 .

 3 Figure 3. Distributions of ecosystem service bundles (ESBs) for the North and South French Alps. Barplots indicate the relative provision of ES within each bundle type. Values are ES z-scores averaged across all municipalities belonging to a specific bundle. Positive z-scores refer to aboveaverage, negative z-scores to below-average values regarding the ES for the regions.

  our case study, potential social-ecological drivers included social and ecological components used in the modelling or quantification of the ES in question (including land cover, elevation, climatic factors), in addition to variables that directly or indirectly drive individual ES and their associations as identified in the literature (biodiversity, NPP) (TableS1). Land cover variables and population density are frequently cited drivers of ES magnitude and distribution(Kienast et al. 2009), including mountainous regions (Grêt-Regamey et al. 2012) and have been widely used as a proxy of ES demand and supply in ES assessments (e.g.[START_REF] Burkhard | Mapping ecosystem service supply, demand and budgets[END_REF]. Protected area coverage relates to an ecosystem's governance and accessibility, has been used as a proxy for spiritual, aesthetic and recreational services[START_REF] Van Jaarsveld | Measuring conditions and trends in ecosystem services at multiple scales: the South African Millennium Ecosystem Assessment (SAfMA) experience[END_REF] and has been shown to be positively correlated with measures of aggregated ecosystem service supply across Europe[START_REF] Maes | Synergies and trade-offs between ecosystem service supply, biodiversity, and habitat conservation status in Europe[END_REF]. Full details are given in Appendix S3.

  of ES indicators within the North and South Alps (p ≤ 0.001): the coverage of grassland, forest, semi-natural, urban land area, protected area coverage, elevation, NPP, plant species richness and population density. The adjusted R 2 values, representing the amount of variance of ES indicators explained by the social-ecological variables were 0.46 for the north and 0.42 for the south. Full methodological details and results are in Appendix S3.2.4 Step 4: Assessing whether ES bundles are associated with different social -ecological systemsIn a call to develop general rules about ES relationships and their implications for management of ES,[START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF] asked whether there exist consistent sets of ES associated with particular socialecological systems. As these systems are not only defined by land cover type,[START_REF] Bennett | Understanding relationships among multiple ecosystem services[END_REF] suggested that the 'anthrome' approach of[START_REF] Ellis | Putting people in the map: Anthropogenic biomes of the world[END_REF] might be useful for identifying a social-ecological system classifications, with distinct systems derived from overlays of social and land use/land cover (LULC) data.[START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF] tested this assertion and quantified the percentage of land area occupied by different anthrome types (derived from overlays of population and LULC data) and bundles of locally derived provisioning ES across South Africa.Hamann et al.

  The crop-dominated bundles in the North and South (ESB1(N) and ESB1(S)) overlapped with SEBs characterised by agricultural land coverage at low elevation and low to intermediate cover of other land uses (Figs.4 and 5; SEB1(N) and SEB1(S)). In the north, the bundle characterised by high provision of forest services (ESB2(N)) broadly overlapped with a bundle characterised by high forest cover (SEB2(N)). The North ES bundle dominated by tourism (ESB3N) did not overlap neatly with

Figure 4 .

 4 Figure 4. Distributions of SEBs for the North and South French Alps. Barplots indicate the relative magnitude of social-ecological variables within each bundle type. Values are variable z-scores averaged across all municipalities belonging to a specific SEB. Positive z-scores refer to aboveaverage, negative z-scores to below-average values regarding the variables for the region. See appendix S2 for variable descriptions).

Figure 5 .

 5 Figure 5. Overlap between ES bundle and SEBs for the north (left) and south (right) of the French Alps, expressed as a percentage of municipalities.

  Figure S3. Principal component analysis biplot for ES in the North French Alps. The first axis (PC1) represents a spatial trade-off between crop and most other services, most strongly with wood and cstock. The second axis (PC2) represents a synergy between high fodd and eros services, and their trade-off with rec. The third axis (PC3) represents a synergy between fodd and tourism, trade-off with wqt. The angles between ES represent the strength of their correlations. The first three components accounted for 61% of the total variation in ES.

Figure S4 .

 S4 Figure S4. Spearman's rank correlation coefficients of pairs of ES across the South French Alps

Table 1 .

 1 Examples of studies that have assessed social-ecological drivers of spatially explicit ES bundles. The studies included here identified and 126 produced maps of bundles of ecosystem services derived from spatially explicit multivariate analyses of ES*. Three studies at the European scale (extent) are not plotted. See Table1.

	127
			Service
	Study	Region	categories (total number	Grain	Method used to obtain bundles	Interpretation of ES bundles
			of variables) #
	Raudsepp-Hearne (2010)	Quebec, Canada	P,C,R(12)	Municipality	k-means clustering	Qualitatively interpreted with regards to coincidence with social-ecological systems as defined by dominant land uses.
	Haines-Young et al. (2011)	Part of Europe	P,C,R(15) (Not just ES)	NUTS-2 regions Unknown	Mean service loadings and marginal impacts of land use and cover change for four services across two time periods were clustered to define groupings of NUTS-2 regions with similar change trajectories.
	Martin-Lopez et al. (2012)	Iberian Spain Peninsula,	P,C,R(14)	Respondents	Hierarchical clustering
	Qiu and Turner (2013)*	Yahara Watershed southern Wisconsin (USA)	P, C, R (10)	30-m grid cells (within 1,336 km 2 watershed)	Factor analysis	Identified three orthogonal axes that represented synergies as well as trade-offs for ES supply. Interpreted interactions by mapping factor scores that represented synergies and trade-offs in ES.
		Southern				Qualitatively interpreted with regards to spatial coincidence with
	Hanspach et al. (2014)	Transylvan ia,	P, C,R,B(9)	Village	Hierarchical clustering	socio-demographic data, derived from commune level statistics, including e.g. total population size, proportions of the main ethnic
		Romania				groups, unemployment, migration levels.

Used redundancy analysis to analyse associations between the relative importance of ecosystem services perceived by people and three types of explanatory variables: stakeholders' characteristics (e.g. education, income), land management strategy (e.g. protection level) and ecosystem type (e.g. presence of mountains). First three axes of the RDA were clustered to obtain bundles.

2. Current approaches to understanding spatially explicit ES associations

	Figure 2. Approach of the spatially explicit analyses of ES associations, organized into four
	conceptual steps.

2.1 Step 1: Assessment, aggregation and harmonisation of ecosystem service indicators Studies that have examined drivers of spatial ES bundles exhibit considerable variation regarding the

	Step 2: Assessment of spatial ES
	associations and delineation of ES
	bundles
	-e.g. PCA, cluster analysis
	Step 3: Identification of social-
	ecological variables important in
	determining or predicting ES bundles
	Step 4: Assessment of congruence
	between ES bundles and social-
	ecological system classifications
	number and types of ES considered, and in how individual ES are quantified (Table 1). Studies have
	typically considered a relatively large number of ES (averaging ~12 ES), encompassing a range of
	provisioning, regulating and cultural ES, and also biodiversity metrics (Table 1). Contrasting large
	numbers of ES within different ES categories can contribute to a better understanding of ES trade-offs
	(Raudsepp-Hearne et al. 2010; Crossman et al. 2013).
	ES maps often vary in the units, range of output values, and spatial resolution. To enable bivariate or
	multivariate analyses, ES datasets have been aggregated to a common resolution. While studies have
	mapped ES at scales ranging from local to global (see Crossman et al. 2013 and Malinga et al. 2015
	for recent reviews), studies mapping ES bundles tend to be conducted for parts of countries at the
	spatial resolution of administrative boundaries, typically the smallest political units such as
	municipalities (Table 1). The use of administrative boundaries has been advocated as relevant for
	multi-ES studies (Raudsepp-Hearne et al. 2010), as municipalities represent the smallest scale of

e.g. qualitative interpretation of maps of ES associations (through maps of PCA site scores or cluster identity) -e.g. quantitative analyses linking drivers to multiple ES or ES bundles such as redundancy analysis

Table 1 )

 1 are classified into groups (bundles) and projected onto maps(Fig 3), allowing the researcher to identify which localities exhibit similar ES associations[START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF][START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF]. ES interactions have also been visualised by mapping the site scores of factor analysis and PCA of multiple ES[START_REF] Qiu | Spatial interactions among ecosystem services in an urbanizing agricultural watershed[END_REF][START_REF] Turner | Bundling ecosystem services in Denmark: Trade-offs and synergies in a cultural landscape[END_REF]. This approach has allowed for the identification of where trade-offs and synergies are the most pronounced in the landscape. Mapping ES associations in these ways has enabled qualitative interpretation of mapped bundles with respect to known distributions of dominant land uses or principal human activities within regions (e.g.

Raudsepp

[START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF] Quieroz et al. 2015;[START_REF] Turner | Bundling ecosystem services in Denmark: Trade-offs and synergies in a cultural landscape[END_REF][START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF]

.

  Whichever quantitative method is used, a critical step is the identification of candidate socialecological variables that are important in explaining or predicting different ES bundles. This initial selection is based either on relationships demonstrated in the primary literature or on expert knowledge, and of course depends on the ES considered in the study.[START_REF] Meacham | Social-ecological drivers of multiple ecosystem services: what variables explain patterns of ecosystem services across the Norrström drainage basin[END_REF] explored four theories of the driving forces behind human impact on ecosystems and tested their relative ability to predict ES bundles. The four models were created by distilling the different driver variables that each theory emphasises. Using random forest analysis, they found that models based on socioeconomic variables performed better than those based on land use.[START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF] used multiple logistic regression to predict the distribution of three ES bundles characterised by low, medium and high levels of direct ES use across South Africa. Drivers were chosen based on variables thought to contribute to the use of natural resources at the household level. They found bundle distribution was determined by social factors, such as household income, gender of the household head, and land tenure, and only partly determined by the supply of natural resources. Qui and Turner (2013) used logistic regression to determine social-ecological determinants of win-win areas, with candidate variables including land use, population density, slope and soil properties. See[START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF] and Table

. Studies are increasingly applying these techniques in ES research to determine how drivers and ES are related to one another, by replacing the latter table with a site × ES table

[START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF][START_REF] Meacham | Social-ecological drivers of multiple ecosystem services: what variables explain patterns of ecosystem services across the Norrström drainage basin[END_REF]

), including, for example, redundancy analysis and canonical correspondence analysis. Other approaches have used regression-based or machine-learning methods with a single response variable, such as ES bundle type (e.g.

[START_REF] Hamann | Mapping social-ecological systems: Identifying 'greenloop' and 'red-loop' dynamics based on characteristic bundles of ecosystem service use[END_REF][START_REF] Meacham | Social-ecological drivers of multiple ecosystem services: what variables explain patterns of ecosystem services across the Norrström drainage basin[END_REF][START_REF] Schulze | Assessing Regional-Scale Impacts of Short Rotation Coppices on Ecosystem Services by Modeling Land-Use Decisions[END_REF]

, or whether a locality represents a win-win or not

(Qui & Turner et al. 2013

).

  Table 1 exhibit considerable variation in the number (mean ~12 ES) and types of ES considered, and in how individual ES are quantified. It is important to distinguish what aspect of a service is being measured by an ES indicator; the potential value provided by an ecosystem, or the service that is actually realised by humans (Jones et al. 2016). Most previous ES bundle analyses, including this study, have mixed indicators ranging from potential supply to actual use values. Two key problems with mixing

  2015) highlighted that positive associations between ES that are actual or potential do not necessarily reflect synergies and can even represent conflicts once the ES are utilised. The selection of which ES are analysed jointly is particularly critical to cross-study comparisons; studies that have analysed associations of different ES, or ES measured or modelled in different ways, are not straightforward to compare. Ultimately, ES bundles delineated by cluster analysis are not generalizable to other regions because a clustering solution is entirely dependent upon the variables used. This issue is already recognised as a limitation for the use of composite indicators of ES (Rodríguez-Loinaz et al. 2012). Raudsepp-Hearne & Peterson (2016) demonstrated that ES bundle spatial patterns were highly dependent on the numbers and types of ES included in the cluster analysis.

. Hamann et al. (2015) focused on bundles of one type of ES, direct use of locally available ES in South Africa (e.g. wood for heating), potentially allowing for a deeper understanding of linkages between ES use and human wellbeing. There is a second difficulty of interpreting bundles of mixed ES indicators: Crouzat et al.

(

  We only considered variables for which continuous spatial data were available in the French Alps, but other unmeasured factors or practices (relating to management history, age of abandonment, or forest age structure) could affect synergies and trade-offs among ES in the regions. This emphasises the need for careful consideration of what actually constitutes a driver of individual ES and ES bundles. Bennett et al.

	(2009) considered many drivers as finer scale management interventions; for example, exogeneous
	drivers (e.g. industrial production) causing environmental change in the social-ecological system, and
	pressures (e.g. use of fertilizers) quantifying the effect of exogenous drivers on a given social-
	ecological system (Mouchet et al. 2014). By using LULC as a determinant, much ES research states
	the obvious about LULC-ES relationships. A danger of circularity exists in such associations, as when
	crop yield is necessarily associated with agricultural lands, and forest-based recreational services can
	only be provided by forests.

Requirement 1: Design studies to test specific hypotheses about specific predictors of key relationships between key ES of interest.

  specific, focused studies should enable a bottom-up understanding of ES bundles in a way that studies that consider all ES simultaneouslylike this case study -cannot.[START_REF] Mitchell | Reframing landscape fragmentation's effects on ecosystem services[END_REF]

	The quantification and mapping of associations between a wide range of ES including provisioning,
	cultural, and regulating services, is thought to enable the identification of a diverse range of trade-offs
	and synergies that might be missed if only individual ES, or a few more commonly quantified ES are
	considered (Lee & Lautenbach 2016). However, as outlined earlier, differences in the distributions
	and types of ES found in different regions mean that determining causal drivers of bundles of all
	available ES is likely impossible.
	Given the diversity and complexity of drivers that affect different ES, a promising approach for
	understanding the degree of generality of different predictors of relationships between ES may be to
	test specific predictions about the importance of specific drivers of relationships of key policy-

relevant ES, based on putative mechanistic relationships. For example, a study might set out to test the relative importance of forest management history and forest age in determining the value of multiple ES across heterogeneous stands (as in

[START_REF] Sutherland | Seeing the forest for its multiple ecosystem services: Indicators for cultural services in heterogeneous forests[END_REF]

. Such 'unpacking' of ES bundles into series of

Table S3 .

 S3 Loadings of ES onto each principal component.

		PC1	PC2	PC3
	cstock	-0.44	-0.06	-0.13
	eros	-0.28	-0.55	-0.21
	fodd	-0.2	-0.44	0.54
	hunt	-0.25	-0.18	-0.05
	rec	-0.28	0.62	-0.01
	tour	-0.33	0.21	0.39
	wood	-0.46	0.18	-0.06
	wqt	-0.23	-0.06	-0.67
	crop	0.42	-0.01	-0.19
	South French Alps

Table S4 .

 S4 Loadings of ES onto each principal component.

		PC1	PC2	PC3
	cstock	-0.50	0.23	-0.14
	eros	-0.26	0.46	-0.41
	fodd	-0.38	-0.35	-0.06
	hunt	-0.07	0.63	0.26
	rec	-0.27	-0.33	-0.15
	tour	-0.28	-0.20	0.31
	wood	-0.45	-0.20	0.06
	wqt	0.13	-0.08	-0.79
	crop	0.40	-0.18	-0.01

Acknowledgements

A University of Southampton IfLS Research Stimulus Fund to RS, MS, KP, JMB and FE supported RS. RS, RL and FE received support from ERC Starting Grant 'SCALEFORES' (Grant number 680176). RS and JMB were funded by CEH project NEC05264. RL, EC and SL received funding from OPERAs (grant number FP7-ENV-2012-two-stage-308393). MS received funding from an ESPA Early Career Fellowship Grant (grant number FELL-2014-104) with support from the Ecosystem Services for Poverty Alleviation (ESPA) programme. The ESPA programme is funded by the Department for International Development (DFID), the Economic and Social Research Council (ESRC) and the Natural Environment Research Council (NERC). PV received support from ERC Grant 'GLOLAND' (grant number 311819). Funding for MGT was provided by the US National Science Foundation (grant numbers DEB-1038759, DEB-1440297, and DEB-1440485). EB received funding from Natural Sciences and Engineering Research Council of Canada (grant number RGPIN 327077-2013) and NSERC EWR Steacie Fellowship. GDP was funded by Social-ecological dynamics of ecosystem services in the Norrström basin (SEEN) project, financed by the Swedish Research Council Formas (grant number 2012-1058) and Swedish Research Council MISTRA (through a core grant to the Stockholm Resilience Centre). We are grateful for constructive feedback from two anonymous reviewers.

Supplementary Information

Appendix S1. Ecosystem service data and spatial characteristics of the case study area

The French Alps study region (Figure S1) covers a total of 52,149 km² (SPCA 1991) and is characterised by high variation in biodiversity, ecosystems and ES provision relative to European averages (Tappeiner et al. 2008[START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF], typical of mountain regions [START_REF] Grêt-Regamey | Muntain ecosystem services: who cares?[END_REF]). The region is dominated by forests and semi-natural areas (67% of the region), with arable lands mainly concentrated in the western broad valleys and piedmonts (27% of the region), while artificial areas cover only 5% of the region. This leads to a clear distinction, typical for mountain regions, of low elevation high-density urban areas surrounded by intensive agriculture in the valleys, and more isolated rural areas (Tappeiner et al. 2008[START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF]. We selected nine ES that have been quantified and mapped in the French Alps previously by [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF]. These are services that were deemed socially, ecologically and economically relevant to the region following consultation with scientists and local collaborators [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF], and include three provisioning (crop, fodd, wood) three cultural (hunt, rec, tour) and three regulating ES (wqt, cstock, eros; see Table 1 for variable codes). All ES are based on either primary data or bespoke modelled surfaces of ES. Full details of these ES are in Crouzat et al. 2015 and Appendix S1. Within the region, elevation, climate and vegetation gradients have had historical consequences on social dynamics and economic activities, resulting in the common separation into the North (Rhône-Alpes) and the South Alps (Provence-Alps-Côte d'Azur; [START_REF] Crouzat | Assessing bundles of ecosystem services from regional to landscape scale: insights from the French Alps[END_REF]. The social-ecological north-south divide is also recognised by an administrative boundary at the NUTS II level, providing a

Appendix S2. Delineation of ES bundles across the study regions

Methods used to delineate ES bundles

Associations among individual ES were quantified using pairwise correlation coefficients based on Spearman's rho (Fig. S2), as is frequently used to inder relationships between ES [START_REF] Mouchet | An interdisciplinary methodological guide for quantifying associations between ecosystem services[END_REF]. We adopted the spatially explicit ES-bundle approach of [START_REF] Raudsepp-Hearne | Ecosystem service bundles for analyzing tradeoffs in diverse landscapes[END_REF] and used cluster analysis to delineate ES bundles across the N and S of the French Alps. A two step clustering approach was adopted (Turner et al. 2015). To minimise skew and make the ES variables dimensionless and comparable in terms of their magnitudes and variability, Box-Cox transformation (Box & Cox, 1964), centring and scaling was applied. First, a PCA was used to quantify the main multivariate relationships between the ES variables to assess whether ES co-occur in spatial bundles.

The number of PCA axes deemed sufficient to characterize the non-random structure in the data in both the N and S ES datasets was selected according to the Kaiser-Guttman criterion, which selects the axes whose eigenvalues are greater than the mean of all eigenvalues (Legendre and Legendre, 1998;Turner et al. 2015). For both regions, the first three components were selected and accounted for 61% and 57% of the total variation in ES in the N and S, respectively. As a precursor to cluster analysis, PCA can serve to separate signal from noise and lead to a more stable clustering solution, with the first axes extracting the essential information while the latter are restricted to noise (Husson et al. 2010). Clustering of PCA axes i.e. uncorrelated components also means that correlated services are not counted more than once or more heavily weighted (Turner et al. 2015). We applied k-means clustering to the relevant PCA axes to delineate ES bundles with 1000 random starts and 10,000 iterations to find a solution with the lowest within-cluster sum of squares according to the relevant PCA axes. K-means clusters municipalities so that the composition of ES values are more alike within than between clusters. Three clusters was deemed appropriate for both the N and S datasets according to a hierarchical cluster analysis using Ward's method and squared euclidean distance (with the FactoMineR package; Lê et al. 2008) and qualitative assessment of ES bundles. Following Renard et al. (2015), we quantified the effective number of ES provided in each bundle using a transformation (H) of the Gini-Simpson's index (S):

for a bundle with N ES with varying proportions p of each service (i).

Results of correlation and PCA analyses for the North and South Alps North French Alps

Agricultural land (agric)

Agricultural land was calculated by taking the percentage of total land area of each municipality occupied by area classed as agricultural by CORINE.

Grazing land (grass)

Grazing land was calculated by taking the percentage of total land area of each municipality occupied by area classed as grasslands or pastures by CORINE.

Forested land (forest)

Forest land was calculated by taking the percentage of total land area of each municipality occupied by area classed as forest by CORINE.

Urban land (urban)

Urban land was calculated by taking the percentage of total land area of each municipality occupied by area classed as urban by CORINE.

Open semi-natural (semi)

Open semi-natural land was calculated by taking the percentage of total land area of each municipality occupied by area classed as semi-natural, other than forest by CORINE.

Elevation (elev)

Elevation values at 30-m resolution were taken from ASTER global digital elevation model raster files (GTOPO30; available at https://asterweb.jpl.nasa.gov/gdem.asp). The median value for each municipality was used.

Mean annual temperature (bio1)

We used the mean annual mean temperature (°C) for the 1950-2000 period using the variable 'bio1' from WorldClim Global Climate Data (Hijmans et al., 2005). The median value for each municipality was used.

Annual precipitation (bio12)

Annual trends of precipitation for the 1950-2000 period in mm. 'bio12' from the WorldClim Global Climate Data (Hijmans et al., 2005). The median value for each municipality was used.

Population density (pop)

We used the log of population density per square kilometre, obtained by dividing the municipality population size by its area. Data were compiled for the most recently available year (2007) (INSEE).

Details of redundancy analysis used to select social-ecological variables important in explaining covariation of ES

We initially considered all of the variables included in Table S2 for explaining covariation in ES. We inspected pairwise correlations between the variables and removed variables with correlation coefficients of >0.80 to reduce multicollinearity (bio1 and bio12). We then computed a global RDA with the remaining potential candidate variables: agric, forest, grass, semi, urban, pop, elev, natura and plant (see above for variables codes). Linear dependencies were further explored by computing variables' variance inflation factors (VIF) for the global model. For both the North and South analyses, the variable agric (proportion of land area covered by agriculture) was not included in the global models to reduce variance inflation factors (all below 5 for the global models without agric).

Forward selection using the packfor package (Dray et al. 2007) was run on the global model to select social-ecological variables important in explaining variation in ES. This procedure selects the model with the combination of variables with the highest R² and p-value [START_REF] Legendre | Numerical ecology, 3rd edition[END_REF].

North French Alps

Appendix S4. Supplementary discussion of case study results

Here we discuss the correspondence between bundles of social-ecological variables and ES bundles, as revealed by RDA and cluster analysis, in order to test the extent to which social-ecological variables could be relevantly used in the French Alps to explain and predict ES bundles. Based on the qualitative interpretation of RDA site scores and cluster distribution, and expert knowledge, we show that results of cluster analysis are contrasted when comparing North to Southern Alps.

In the North Alps, ESB1(N) is characterised by a high level of crop production and a below-average level of supply for all other ES (Fig. 2). Spatially, this bundle is clustered over the broadest lowland valleys and the fringe of the external Alps (Figures 2 andS6). This appears very consistent with its broad overlap with SEB1(N) (Fig. 5), where agricultural and artificial areas are overrepresented (Fig. 4). ESB2(N) presents opposite patterns, as all ES are supplied at above-average levels except crop production which is supplied far below-average (Fig. 2). In particular, this bundle supplies the highest regional levels of forest-related ES (carbon stocks, erosion mitigation and wood stocks). This is coherent with its large overlap with SEB2(N) (Fig. 4), where forested areas are overrepresented compared to the two others (Fig. 4). In between, ESB3(N) supplies an average level of most ES except those specific to bundles ESB1(N) and ESB3(N) (Fig. 2). ESB3(N) is essentially concentrated over areas of intermediate altitude of the external Alps (e.g., the Chartreuse range), which can be captured at municipal scale as mosaic areas containing a mix of forests, grasslands, built-up and semi-natural open areas. This could explain the mixed overlap of ESB3(N) with the three SEB identified for the North Alps (Fig. 4).

In the Southern Alps, ESB2(S) is characteristic of rural mosaic landscapes of the internal Alps,

comprising forested and open areas at generally high altitudes (Figs 2 and S6). It supplies the highest regional levels of numerous ES, in particular fodder production and wood production, recreation and tourism, and carbon stocks (Fig. 2). This is consistent with its large overlap with SEB2(S) (Fig. 4), which over-represents high altitude grasslands and also contains a regional average level of forested areas (Fig. 4). ESB3(S) is a rich and diverse bundle as it supplies an average level of most ES, and the highest regional levels of erosion mitigation and leisure hunting in particular (Fig. 2). This high multifunctionality could be related to its heterogeneous spatial patterns, as inferred by the large combined overlaps of ESB3(S) with both SEB1(S) and SEB2(S) (Fig. 5), i.e. with areas of low to average elevations with contrasting land use features (Fig. 4). ESB1(S) echoes ESB1(N) in the North Alps and is characterised by a high level of crop production and an under-average level of supply for all other ES (Fig. 2). This mono-functional bundle is broadly located in the main intensive agricultura l and urbanized valleys (Rhône and Durance rivers) but poorly overlaps with the Southern SEBs (Fig. 4).

Overall, we highlight a discrepancy in our ability to understand, and further predict, ESBs from SEBs in comparing the results from the North versus the Southern Alps. We hypothesize that this could be linked to discrepancies in the clustering patterns between ES and social-ecological variables in these two regions, as captured at municipal scale. The North Alps can be characterised with easily distinguishable entities in terms of ES and ecological variables at the municipal scale, relying on aggregated land use and biophysical patterns (e.g., internal versus external Alps, large contrasts in terms of elevation and land uses). In contrast, the Southern Alps, and in particular the Southern pre-Alps, have more heterogeneous landscapes. Such landscapes are composed by a fine-grained mosaic of open spaces (pastures) and secondary forests, related to ecological secondary succession after agricultura l abandonment (post World War II) (Hinojosa et al. 2016). We hypothesize that SEBs at municipal scale are too coarse in these heterogeneous landscapes for being predictive of aggregated ES. SEBs in the Southern Alps can be considered as a typology of land cover types, combined with elevation, and these appear not sufficient for predicting ES bundles in municipalities characterized by a high landscape heterogeneityespecially in the southern Pre-Alps (northwestern part of the map). ES were initially modelled at a finer resolution (1-km) and then aggregated at municipal scale to be coherent and comparable with SE variables. But in heterogeneous landscapes, the effect of fine landscape patterns on ES supply might not be negligible [START_REF] Mitchell | Reframing landscape fragmentation's effects on ecosystem services[END_REF], explaining why SEBs at municipal scale cannot be good predictors of ESB. These results point to a discrepancy between the municipal scale required for many SE variables (in particular social ones), and the finer scale required for understanding ecological processes and ES patterns. As a conclusion, while the SEBs identified could be relevant predictors of ESBs at municipal scale in aggregated and contrasted landscape types as in the North Alps, they appear insufficiently comprehensive for heterogeneous areas as in the Southern Alps.