Landscape metrics as functional traits in plants: perspectives from a glacier foreland

Study area

Field work was performed on the Nardis glacier foreland (46°12′14″N, 10° 40′21″E), located in the Adamello-Presanella group (Rhaetian Alps, southern sector of the Italian Central Alps) on the southern slope of Presanella Peak (3,556 m a.s.l.) (Fig. 1A). The glacier has a surface area of approximately 1.67 km² (SAT, 2007) and its tongue extends down to an altitude of 2,720 m a.s.l. (Fig. 1B). The bedrock consists primarily of acidic granitoid material. The geology is characterized by the large Adamello-Presanella-Monte Re di Castello batholith (29.4–41 Ma), consisting of tonalite, an igneous, plutonic intrusive rock. Available climatic data taken from a nearby weather station (46°25′33″N, 10°41′51″E) located at the same altitude indicate a mean summer temperature of 5.7°C and a mean annual precipitation of 897 mm. The study area where the sampling took place was approximately 7 hectares in size and corresponded to the zone in front of the glacier tongue, where the glacier was still present in 1945 (Fig. 1C).

Data collection

We sampled 46 randomly distributed points in the study area. Using the closest individual method (Krebs, 1999), from each of these sampling points we selected the closest 1m × 1m sample plot, which was (i) safe from landslides and flat (<5°) and (ii) without boulders (d > 256 mm). This selection procedure was done to avoid marked differences in site conditions (Vetaas, 1997), even if it were not possible to completely control for a certain amount of variability in topsoil texture. Using historical maps and aerial photographs, frontlines of the glacier tongue, i.e., lines of equal terrain age (isochrones), were established for 1945, 1954, 1970, and 1996. Terrain ages (ta) of each sample plot were then classified as follows: (ta₁) between 15 and 41 years (n = 11); (ta₂) between 41 and 57 years (n = 17); and (ta₃) between 57 and 66 years (n = 18) (Fig. 1C). Each plot was subdivided into 10,000 1 cm × 1 cm-grid cells. In August 2011, vascular plant species distribution intersecting the central axes of each 1 cm × 1 cm-cell were then mapped and digitised using ESRI ArcGIS 9.3, yielding for each species the number of 1 cm² cells occupied per 1 m² (Fig. 2). The same was done to classify the topsoil texture since, amongst the environmental variables, topsoil texture was considered to be a major physical influence on patch development. As a measure of topsoil texture we used the cover of cobbles, which, according to the Wentworth scale (Wentworth, 1922), have a diameter >64 mm, a size (>32 cm²) comparable to the observed mean size of plant patches (47 cm²). Then, we distinguished three classes of cobble cover: (ts₁) <7% (n = 15); (ts₂) between 7% and 14% (n = 15); and (ts₃) >14% (n = 16).

Data analysis

Landscape metrics were calculated according to the procedure of Teixido et al. (2007). For each 1 m × 1 m-plot, we calculated seven metrics (Table 1). Calculations were made using the Patch Analyst 5.0 extension for ArcGIS 9.3 (Rempel, Kaukinen & Carr, 2012), adopting a four-neighbour rule to identify the patches. A synthetic description of each patch metric is reported in Table 1. Two metrics, patch size and shape, can be computed both for every patch in the landscape and as an average for all patch metrics at the plot level (see also McGarigal & Marks, 1994). In the first case, the metrics may be related to each species, i.e., species-specific traits.

We did not sample plant traits in the field because this would have been too time-consuming and too destructive to the vegetation. When available, plant traits were taken from Cerabolini et al. (2010). Data for Salix herbacea L. and, partially, Sedum alpestre Vill. data were taken from Carbognani (2011) and Pierce et al. (2007). Leaf nitrogen and carbon contents of S. alpestre were calculated as a mean of the values of Sedum acre L. and Sedum album L. from Cerabolini et al. (2010). Saxifraga stellaris L., Polygonum viviparum L. and Athyrium distentifolium Tausch ex Opiz, due to their negligible frequencies and cover, were omitted from further analyses. A brief description of each selected trait is provided in Table 1 (see also Cornelissen et al., 2003; Gross, Suding & Lavorel, 2007; Wilson, Thompson & Hodgson, 1999). The plant traits database used here is reported in Table 2.

The 46 vegetation relevés were subjected to agglomerative cluster analysis using the Bray–Curtis coefficient of dissimilarity and Ward’s clustering method. Next, ANOVA was applied to verify for differences in landscape metrics between (i) the three classes of terrain age and (ii) topsoil texture distinguished, and (iii) the three floristic relevé clusters identified.

To relate plant traits to spatial configuration, taking into account species cover in the plots, we applied RLQ-analysis, a tool to assess how the environment filters certain species traits (Dolédec et al., 1996; Dray, Chessel & Thioulouse, 2003). The RLQ procedure performs a double inertia analysis of an environmental-variables-by-samples (R-table) and a species-by-traits (Q-table) matrix, with a link expressed by a species-cover-by-samples matrix (L-table). RLQ-analysis combines three unconstrained separate ordinations, correspondence analysis of L-table and centred normed principal component analyses of Q- and R- tables, to maximise the covariance between environmental factors and trait data by the use of co-inertia analysis (Bernhardt-Römermann et al., 2008). Here, we studied the joint structure of three data tables, namely (i) a plot-by-landscape metrics data table R-table, (ii) a plot-by-species table containing the abundances of the plant species present in our set of 46 plots table (L-table), and (iii) a species-by-trait data Q-table.

This RLQ analysis (basic-RLQ) was followed by a partial-RLQ, with the aim of checking the effects of the covariates terrain age and cobble cover, i.e., to verify the need to partition the variation related to these factors. This type of analysis is a special case of RLQ, where the covariable represents a partition of samples into groups. If the percentage of co-inertia explained by the most representative axis of partial-RLQ were to be much higher than in the basic-RLQ, this would mean that the influence of the covariate is relevant. The same approach was followed by Wesuls et al. (2012) to partition the response of plant traits to grazing-related environmental parameters from other environmental and temporal variations.

A permutation method was used to compare the hypothesis H₀: X = 0 (trait and landscape are unrelated) against H₁: X ≠ 0 (trait and landscape are related), where X is the fourth corner, a trait-by-landscape table, whose parameters cross the traits (Q-table) to the landscape variables (R-table), via the abundance table L-table (Legendre, Galzin & HarmelinVivien, 1997). The null hypothesis consists of three null joint hypotheses: both R and Q are linked to L (L ↔ Q, L ↔ R), only R is linked to L (L ↮ Q, L ↔ R), only Q is linked to L (L ↔ Q, L ↮ R). The overall null hypothesis is rejected when both null hypotheses are rejected (L ↮ Q and L ↮ R). Dray & Legendre (2008) proposed to set the alpha argument to $\alpha =\sqrt{0.05}$, but recently it has been shown that α should be 0.05 instead (Ter Braak, Cormont & Dray, 2012). Given the limited power of this test with few species (Ter Braak, Cormont & Dray, 2012), as in the present study, we show the results according to both the Dray & Legendre (2008) and the Ter Braak, Cormont & Dray (2012) alpha argument settings. A multivariate permutation test was applied to evaluate the global significance of the traits-spatial configuration relationships, implemented by the function ‘randtest’ of the package ade4 ( Dray & Dufour, 2007). Next, we tested the associations of spatial configuration and trait variables with the axes of the basic-RLQ. The strength of the association of landscape metrics and plant traits was measured with the D2 statistic (Dray et al., 2014). All tests were performed using the combined fourth-corner statistic (Dray et al., 2014) with 49,999 permutations.

All statistical analyses were performed using the open source R software (R Core Team , 2013). We used the library vegan (Oksanen et al., 2011) for the cluster analysis, the library stats (R Core Team , 2013) for ANOVA and the library ade4 ( Dray & Dufour, 2007) for the RLQ analysis.

Species composition

We recorded a total of 19 plant species (including S. stellaris, P. viviparum and A. distentifolium which were omitted from further analyses). Each species frequency and mean patch-level metrics are reported in Table 2. This table shows that the most frequent species were P. alpina, Saxifraga bryoides L., Leucanthemopsis alpina (L.) Heywood, O. digyna and Veronica alpina L., present in more than half of the plots. Mean patch size, total cover and mean shape index of species did not necessarily reflect their frequency. Considering only the species where the number of sampled patches was higher than 30, the most compacted shapes belonged to S. alpestre, the most irregular to L. alpina, with the largest patches to S. bryoides and the smallest to S. alpestre. The dispersion of patch size values around the mean was much larger than for patch shape. The shape of the frequency distribution of patch sizes was usually positively skewed (Fig. S1), while patch shape values showed more frequently symmetrical distributions (Fig. S2).

Three clusters (CLs) resulted from the analysis of the 46 vegetation relevés, characterised by differences in the cover of six plant species (Fig. 3). In the Saxifraga-cluster (CL₁), S. oppositifolia was much more abundant and V. alpina much less abundant than in the other two clusters; while in the Leucanthemopsis-cluster (CL₂), L. alpina was more abundant. Finally, the Oxyria-cluster (CL₃) was characterised by relatively high covers of O. digyna and Cardamine resedifolia L. (Fig. 3). ANOVA showed that the floristically defined relevé clusters differed significantly with regard to certain landscape metrics. The Shannon’s index of diversity and the patch type (=species) richness were significantly higher in the Leucanthemopsis-cluster (CL₂) than in the Saxifraga-cluster (CL₃) and the Oxyria-cluster (CL₁) (Table 3). Also, the mean size of the vegetation patches differed significantly among the three relevé clusters, being greatest in the Saxifraga-cluster (CL₁) and smallest in the Leucanthemopsis-cluster (CL₂).

As expected, none of the landscape metrics differed amongst the three classes of terrain age, that is 15–41 y, 41–57 y and 57–66 y, and of cobble cover, that is <7% , between 7% and 14%, and >15% (results reported in Table S1). Correlations among landscape traits are reported in Table S2.

Relationships between plant traits and landscape metrics

The percentages of total co-inertia explained by the first two axes of the basic-RLQ was 96.1% while those of the two partial-RLQs were 95.3% and 95.7%, respectively adopting terrain age and cobble cover as covariates. The first axis of the basic-RLQ explained 89.8% co-inertia, while the percentage explained by the first axis of the two partial-RLQs was lower or almost identical (86.2% and 90.5%, same order above), meaning that the spatial configuration gradient along the first axis of both the partial-RLQs was not more pronounced than the basic-RLQ (Table 4). Moreover, the ordination diagrams of the basic-RLQ did not show any grouping of plots according either to the terrain age nor to the cobble cover factor (Fig. 4).

The test for the model (H₁: L ↔R) showed that the distribution of species with fixed traits was influenced by the spatial configuration (p = 0.004), while the test for the model (H₁: L ↔Q) showed that species composition of plots with fixed spatial configuration was not influenced by the species traits (p = 0.590). This means that the traits-spatial configuration relationships were not globally significant.

The first basic-RLQ axis was significantly and negatively correlated with mean patch size and patch size coefficient of variation and positively to Shannon’s diversity. The second basic-RLQ axis was negatively correlated to mean shape index. Among the plant traits, canopy height, leaf fresh and dry weight and leaf area showed a positive significant correlation with the first basic-RLQ axis (Table 5).

In summary, the first basic-RLQ axis represented a gradient of increasing diversity as a response to the presence of larger plant species patches (Fig. 5A). Moreover, as can be seen in Fig. 5B and Fig. 5C, the first basic-RLQ axis represented a gradient of increasing cover of taller species with larger and softer leaves, on the right hand side, like Luzula alpino-pilosa (Chaix) Breistr., O. digyna, Geum reptans L. and Adenostyles leucophylla (Willd.) Rchb., to tiny species with smaller, heavier and harder leaves on the left hand side, like saxifrages and Gnaphalium supinum L.. The second basic-RLQ axis, even if less pronounced than the first axis, represented a gradient of increasing conservation of acquired resources, that is from species with higher specific leaf area, leaf nitrogen content and tendency for lateral spread, like Cerastium uniflorum Clairv. and O. digyna, to species with larger dry matter content, like saxifrages (Fig. 5C). Interestingly, this axis also represents a gradient of increasing patch compactness.

There were eleven significant associations between plant traits and landscape metrics. Mean patch size was negatively correlated with canopy height, leaf dry and fresh weight, and leaf area, the latter two of which were also negatively correlated with the coefficient of variation of patch size. Among the spatial metrics, the Shannon’s index of diversity was most frequently and positively associated with plant traits, namely with canopy height, leaf dry and fresh weight, and leaf area. Finally, patch richness was positively associated with canopy height (Fig. 6). We should stress that the majority of the tests were significant only when using a significance level of $\sqrt{0.05}$. The complete list of correlations between landscape and plant traits is reported in Table S3.

Species composition and landscape metrics

The terrain age and plant composition observed here correspond to early-successional stages of glacier foreland colonization; it does not include any mid- or late-successional species according to Caccianiga et al. (2006), apart from S. bryoides. However, this species, according to same authors, is present in the same proportion in early- and mid-successional sites. The composition also includes several ubiquitous species, like P. alpina, L. alpino-pilosa and Hieracium alpinum L. These species are typical wind-dispersed ruderals. The spectrum of plant traits encountered is limited and only some species, namely A. leucophylla, P. alpina, Ranunculus glacialis L., G. reptans, L. alpino-pilosa, and O. digyna exhibited relatively strong departures from the average characteristics of the main group of species.

We have shown that relevé clusters differ in landscape metrics. For instance, S. bryoides a species with small leaves but relatively high nitrogen content, showed the highest patch sizes amongst the frequent species and characterised a cluster with a correspondingly higher mean patch size. V. alpina and L. alpina, both medium-sized species, were present with many patches and, accordingly, characterised a cluster with a higher Shannon’s diversity. Finally, O. digyna, a relatively low compact species, with medium lateral spread and nitrogen content, was the most frequent in the cluster with the highest mean patch shape.

Spatial configuration as a correlate of plant traits

The RLQ-analysis showed that plant species with efficient conservation of nutrients, which have a higher dry matter content in leaves (Pierce et al., 2013; Wilson, Thompson & Hodgson, 1999), also have lower values of the shape index, or, in other words, their patches are more compact. This pattern, maintained due to processes of intra- and interspecific competition for space and nutrients, ensures efficient acquisition-conservation trade-offs in plants characterized by slow growth, as in P. alpina, which, among the common species, is second in patch compactness only to S. alpestre. Second, we found that Shannon diversity increased with increasing cover of upright-growing plant species, characterised by larger and heavier leaves. A possible explanation for this lies in the way plants with these traits compete directly and/or indirectly, and how they modify one another’s biotic and abiotic environment, thereby generating a more equitable distribution of patch sizes, combined with a higher number of species.

The observed correlation of landscape metrics with species composition together with the correlations with specific life-form traits seem to indicate some level of life-form or species-based spatial self-organisation. Self-organization does not imply specific causalities between vegetation patterns and the environment, but is induced by internal variation independent of external drivers (Bolliger, Sprott & Mladenoff, 2003). Initial establishment of any particular species in this microhabitat depends on successful seed establishment. The first selection of species is mainly trait-driven because only species with a wind dispersal strategy are selected from the geographical species pool (Keddy, 1992). However, the later establishment is dependent on random abiotic factors such as wind-aided dispersal and by small-scale variation of the soil surface characteristics, such as texture, aiding germination. Stochastic factors thus potentially lead to a high degree of heterogeneity in seedling distribution due to the variability of seed rain, the soil seed bank, germination, mortality rates of the seedlings as well as other factors (Erschbamer, Kneringer & Schlag, 2001; Marcante, Winkler & Erschbamer, 2009). Large-scale successional stages usually contain a wide array of different microhabitats where species composition is mainly driven by habitat conditions and local disturbances (e.g., floods, rock falls, and avalanches) (Burga et al., 2010).

Here we tried to maintain, at a minimum level, these potential environmental filters by selecting a restricted range of terrain ages, slope, and topsoil texture. This has probably favoured us in finding the observed relationships between landscape patch metrics and functional traits.

Landscape metrics as functional traits

Landscape-scale variation, usually surveyed at larger resolutions than those used here, is often seen as a filter of plant traits rather than a plant trait itself. This can be partially explained by the need to study how environmental change and disturbance, especially anthropogenic, may affect the landscape mosaic and in turn filter specific traits (Gámez-Virués et al., 2015). Another possible explanation is the difficulty to select the appropriate grain and extent to correctly detect each species’ patch boundaries. According to O’Neill et al. (in Turner, Gardner & O’Neill, 2001) the grain size of a map should be two to five times smaller than the spatial features being analysed, and map extent should be two to five times larger than the largest patches. The mean, upper quartile and maximum patch size was 47 cm², 35 cm², and 3,104 cm², respectively. This means that the extent (1 m²) we chose was appropriate and that the grain (1 cm²) could be apparently increased. However, some species presented very low mean patch sizes, e.g., S. alpestre (only 7 cm²). This means that to perform a multi-species patch-level calculation of landscape metrics on a glacier foreland with terrain ages <70 y, the grain should be maintained at least lower than 3.5–1.4 cm², in practice, units of 3.0–1.0 cm². Finally, we have to stress again that the identification of each patch could be difficult in late successional stages, where their edges could not be easily distinguished.

Databases of plant functional traits such as LEDA (Knevel et al., 2003), TRY (Kattge et al., 2011) or BIOPOP (Poschlod et al., 2003) are very useful since they provide possible determinants of the response of primary producers to environmental factors, and how they can affect other trophic levels, ecosystem processes, services and diversity (Kattge et al., 2011). Existing databases of functional traits do not include data about patch-level landscape metrics of plant species. We suggest that these databases should incorporate such data.

One possible criticism of this proposal might be that some landscape traits might already be incorporated in existing plant species attributes. However, this should be treated with caution. For example, while patch size might be associated with the Hodgson et al.’s (1999) lateral spread measure, included as a class parameter into the attribution of Grime’s (1977) CSR strategy, here we observed that a reduction in patch size is more clearly associated with canopy height and leaf traits. Also the capability to perform clonal growth, resulting in large and compact patches, is an attribute largely included into existing databases and is well known for many alpine grasses and sedges (Groenendael et al., 1996). On the other hand, patch size and shape could be linked to ontogenetic features of the individual (i.e., individual age) for species with indefinite growth such as many woody species, like many cushion species. Individual age, in turn, could be related to stochastic station parameters such as terrain age, which may be the case in a glacier foreland, and should not be considered as a plant trait. Compactness of the patch itself could be the result of many plant attributes and thus requires care if treated as a plant trait. Moreover, non-clonal single individuals are rather compact and we have shown that the most compact patches are also the smallest. This could be due to the situation where the species occurs in single individuals. In general, we have shown that the distribution of species-specific landscape trait values, in particular patch size, might be rather skewed, meaning the mean might not be the best summarising statistical indicator and that care should be taken in the selection of the individuals where measuring the traits (Cornelissen et al., 2003). All these aspects call for further research to collect these data and to link them with physiological attributes directly measured in the field.