Butterfly distribution along altitudinal gradients: temporal changes over a short time period

Cristiana Cerrato; Emanuel Rocchia; Massimo Brunetti; Radames Bionda; Bruno Bassano; Antonello Provenzale; Simona Bonelli; Ramona Viterbi

doi:10.3897/natureconservation.34.30728

Abstract

Mountain ecosystems are particularly sensitive to changes in climate and land cover, but at the same time, they can offer important refuges for species on the opposite of the more altered lowlands. To explore the potential role of mountain ecosystems in butterfly conservation and to assess the vulnerability of the alpine species, we analyzed the short-term changes (2006–2008 vs. 2012–2013) of butterflies’ distribution along altitudinal gradients in the NW Italian Alps. We sampled butterfly communities once a month (62 sampling stations, 3 seasonal replicates per year, from June to August) by semi-quantitative sampling techniques. The monitored gradient ranges from the montane to the alpine belt (600–2700 m a.s.l.) within three protected areas: Gran Paradiso National Park (LTER, Sitecode: LTER_EU_IT_109), Orsiera Rocciavrè Natural Park and Veglia Devero Natural Park. We investigated butterflies’ temporal changes in accordance with a hierarchical approach to assess potential relationships between species and community level. As a first step, we characterized each species in terms of habitat requirements, elevational range and temperature preferences and we compared plot occupancy and altitudinal range changes between time periods (2006–2008 vs. 2012–2013). Secondly, we focused on community level, analyzing species richness and community composition temporal changes. The species level analysis highlighted a general increase in mean occupancy level and significant changes at both altitudinal boundaries. Looking at the ecological groups, we observed an increase of generalist and highly mobile species at the expense of the specialist and less mobile ones. For the community level, we noticed a significant increase in species richness, in the community temperature index and a tendency towards homogenization within communities. Besides the short time period considered, butterflies species distribution and communities changed considerably. In light of these results, it is fundamental to continue monitoring activities to understand if we are facing transient changes or first signals of an imminent trend.

Keywords

Butterfly, community composition, mountain ecosystem, LTER, protected area

Introduction

Global warming and land use changes are considered among the main threats to biodiversity (Sala et al. 2000, Lemoine et al. 2007). Other driving factors may also interact with climate and land use changes to impact biodiversity (Brook et al. 2008); however, substantial alterations in the population and distribution of terrestrial species have already been detected worldwide, mainly in response to both of these effects (Mantyka-Pringle et al. 2012).

Mountain ranges are very sensitive to environmental changes and global warming (Huber et al. 2005, Beniston 2006). It is well established that the Alps have experienced a temperature increase over the last century (Beniston 2003, Brunetti et al. 2009) especially at higher altitudes (Acquaotta et al. 2014). In addition to global warming, the alpine chain has suffered (and is still suffering) loss of open habitats as a consequence of forest expansion (Hunziker 1995, Gellrich et al. 2007) due to traditional land use practices being abandoned (Hinojosa et al. 2016). Such a loss of open habitats is concomitant with a noticeable loss of species (Pauli et al. 2007), which are under threat because climatic and land use changes will probably continue into the future (Nogués-Bravo et al. 2007, Chamberlain et al. 2013, Pellissier et al. 2013).

The alpine biodiversity has already responded to these factors. Upward shifts of alpine plants (Walther et al. 2005, Pauli et al. 2012, Vittoz et al. 2013), invertebrates (butterflies, Wilson et al. 2007, Wilson and Gutierrez 2012; carabid beetles, Pizzolotto et al. 2014) and mountain birds (Scridel et al. 2018, Rocchia et al. 2018) have been documented. In addition, community composition has changed at high alpine sites (Keller and Körner 2003), with an accelerating increase of species richness, in many cases due to an expansion of generalist species (Walther et al. 2005, Pauli et al. 2007).

Although some responses are evident, few investigations have focused on alpine biodiversity temporal changes. As mentioned above, alpine environments are under threat, therefore more information is needed to understand how the main adverse factors (global warming and loss of open habitats) have affected, or are affecting, alpine biodiversity over time. Exploring temporal patterns of biodiversity is of great significance because future warming and related environmental changes are expected to cause substantial changes, not only in spatial distribution of species, but also in species turnover over time (Korhonen et al. 2010).

Long-term monitoring programs are fundamental tools to assess and monitor temporal changes of biodiversity (Morecroft et al. 2009; Magurran et al. 2010). However, long-term series are rarely available (Gaston et al. 2008). Therefore, it seems useful in the meantime to focus on short time-scales and investigate short-term biodiversity responses, establishing the basis for disentangling in the future if we are facing directional changes or just short-term fluctuations.

Protected areas are a key part of conservation strategies to mitigate and monitor losses of biological diversity due to changes in climate and land use (Kharouba and Kerr 2010). In those areas, the reduced direct anthropic pressure allows for more easily focusing on and detecting the effect of climate modifications, and the natural evolution of vegetation, on biodiversity. Moreover, thanks to their staff and logistic resources, protected areas have a primary role in carrying out long-term monitoring programs, with the aim of acquiring and sharing knowledge through research networks, such as the LTER network.

In 2006, a multi-taxa monitoring project in the NW Italian Alps was started (Biodiversity Monitoring Project), with the goal of analyzing long-term biodiversity changes in space and time along altitudinal transects in three protected areas (Viterbi et al. 2013).

Following the protocol developed and tested by the LTER site Gran Paradiso National Park (GPNP), seven taxa were monitored (Lepidoptera Rhopalocera, Orthoptera, Aves, Coleoptera Carabidae, Coleoptera Staphylinidae, Araneae, Hymenoptera Formicidae), using standardized, easy-to-apply and cheap sampling techniques. In 2007, the protocol was adopted by two other protected areas (Orsiera Rocciavré Natural Park, ORNP and Veglia Devero Natural Park, VDNP). The protocol has been repeated every 5 years (2 yr monitoring – 4 yr stop), and as a result, data from two sampling periods are now available (1^st season, 2006–2008; 2^nd season, 2012–2013).

In this framework, we focused on butterfly data deriving from the two sampling periods (1^st, 2006–2008; 2^nd, 2012–2013) of the Biodiversity Monitoring Project. Butterflies are an ideal indicator taxon (Dennis 1993; Hellmann 2002; Parmesan and Yohe 2003), strongly influenced by the abiotic environment and climatically sensitive (Hellmann 2002; Parmesan and Yohe 2003). Their ecological sensitivity and ability to disperse (Devictor et al. 2012) make butterflies a good model system to study species responses to climate and habitat changes (Stefanescu et al. 2011, Roth et al. 2014).

The main aim of our research consisted in assessing changes of butterflies’ species distribution and community composition along altitudinal gradients over a short time period. We investigated butterflies’ responses across multiple levels of biological organization following a hierarchical approach, in order to assess heterogeneity or consistency of changes among different levels. As a first step, we analyzed temporal changes of distribution for single species and for species classified into ecological groups (species level). Then, we explored temporal variation of species richness and of community composition (community level). We finally discussed the results in the light of the possible drivers of change within our study areas.

Materials and methods

Study areas, sampling design and data collection

Our study was carried out in three protected areas in the NW Italian Alps: Gran Paradiso National Park (GPNP; 720 km²; 45°31'N, 7°19'E), Orsiera Rocciavrè Natural Park (ORNP; 110 km²; 45°3'N, 7°19'E) and Veglia Devero Natural Park (VDNP; 86.2 km²; 46°12'N, 8°14'E). GPNP is part of the Long-Term Ecological Research network (LTER, Sitecode: LTER_EU_IT_109).

All areas are characterized by mountain (alpine) conditions with vegetation ranging from mixed forest to alpine prairies and glaciers. The three parks have similar mean elevations and vegetation characteristics, but differ slightly in terms of climatic regime (the highest monthly precipitation and lowest annual mean temperature are found in VDNP).

A total of 11 altitudinal transects were set (one for each valley for each of the three parks) ranging from a minimum of 600 to a maximum of 2700 m a.s.l. The transects cover different altitudinal ranges and, overall, they involve three vegetation belts (montane, subalpine, alpine; Suppl. material 1: S1). Each altitudinal transect was composed of 4–7 sampling units (plots) separated by an altitude range of ca. 200 m, to allow independence of sampled data, for a total of 62 plots. Sampling units were circular plots with a radius of 100 m, where monitoring activities were carried out to obtain semi-quantitative data.

We sampled butterflies using linear transects along the diameter of the plot (200 m in length), walked along at uniform speed. According to the protocol by Pollard (1977), we recorded all butterflies observed within an imaginary 5×5 m square. Sampling was limited to sunny conditions, under calm-to-light wind. Each plot was visited once a month from June to August. Individuals were captured and released after specific identification, except for specimens that were difficult to identify, which were retained for further determination. Sampling efforts were characterized by 2 consecutive years of sampling and we collected data during two separate time periods: 1^st sampling period (2006–2007 in GPNP; 2007–2008 in ORNP and VDNP) and 2^nd sampling period (2012–2013 for all the three parks, GPNP, ORNP, VDNP).

In each plot, we also recorded the microclimatic conditions placing data-loggers (Thermochron iButton, DS1922L, Maxim, Sunnyvale, CA, U.S.) which recorded the air temperature every hour throughout the field season. Loggers were located in the center of each plot, at least 1 m above ground and were covered with a white shield.

By using land cover data derived from aerial photos and validated in the field (Agroselviter 2009; Meloni et al. 2009), we classified each plot according to four categories of main habitat cover type: woodland, ecotone (corresponding to transition habitats), grassland (including all grassland types jointly grouped, even if above or below the tree line) and rocks.

Moreover, during each sampling session, we noted, for each plot, land management data (dominant land use) considering three categorical variables: grazing, mowing and no activities.

Statistical analysis

Species level

The first step of our analysis focused on distributional changes of butterflies’ species and groups of species. Following ecological requirements, species could change distribution along altitudinal gradient, either vertically or horizontally with the colonization or the loss of sites at similar altitudes.

To obtain a comprehensive framework on butterflies species distribution changes, we firstly assessed temporal changes in plot occupancy to detect an overall signal of expansion or retraction without considering an altitudinal directionality, then, we investigated altitudinal range changes between sampling periods.

Changes in plot occupancy

We defined occupancy as the number of plots occupied by each species in each sampling session (1^st vs 2^nd) and compared it by using a t-test for paired samples (significance level assessed after 999 randomizations, following Legendre and Legendre 2012). The randomization process allowed to perform the t-test even with non-parametric data (Legendre and Legendre 2012).

To identify which group of species changed the most over time, we analyzed if the occupancy increased/decreased equally between functional groups. We compared the changes in the number of plots per species (Δ plot, 2^nd sampling session minus 1^st sampling session), between the ecological groups of conservation interest by using non-parametric tests (Kruskal-Wallis or Mann-Whitney tests, depending on the number of ecological categories). Following the classification proposed by Balletto and Kudrna (1985) and Balletto et al. (2015), we classified butterflies’ species according to: feeding specialization (from polyphagous to monophagous); altitudinal range (generalist, medium altitude, specialized); high altitude species (species typically found only in the subalpine and alpine vegetation belts); light (“shade loving”, “sun loving”); dispersal capacity; habitat preferences (woodland, ecotone, open areas, screes) (Suppl. material 1: S2). We selected the ecological groups according to their potential responses in the light of possible drivers of change. In particular, feeding specialization and habitat preferences groups could highlight responses in relation to micro habitat changes while altitudinal range and high altitude species groups could have changed the occupancy status because of temperature increase. Moreover, we considered dispersal capacity because we hypothesized that species with high vagility can colonize easily new areas following both drivers of change.

To assess if individual species changed their occupancy of plots over time (if the frequency of presence/absence was statistically different between the two sampling periods), we applied the McNemar’s test for paired nominal data in a 2×2 table (Sokal and Rohlf 1995). We organized data about gains and losses in plots according to the four following combinations: plots with species presence (no change) or species absence in both sampling periods (no change); plots only colonized during the first sampling season (loss) or only during the second one (gains). We tested the 133 species common to both sampling periods, and we considered the species with p-values still significant after Bonferroni correction as “winners” or “losers”, depending on the direction of changes.

Altitudinal range changes

We described the altitudinal range of each species with the following parameters:

– altitudinal optimum (mean and median value);

– higher limit (absolute maximum, 90^th percentile);

– lower limit (absolute minimum, 10^th percentile).

To quantify the amount of change, we compared these parameters between sampling periods using the t-test for paired samples. To compare non-parametric data we assessed significance level after 999 randomizations (Legendre and Legendre 2012). As in the case of occupancy rates, we also compared the changes in altitudinal range among ecological groups.

To test for species responses over time along the altitudinal range, we used the indicator species analysis (IndVal method) proposed by Dufrêne and Legendre (1997). The IndVal method combines the “specificity” of a species (its uniqueness in a specific group of sites) and its “fidelity” (its frequency within the same group). For each species, IndVal can range from 0 (no indication) to 1 (maximum indication). Statistical significances of IndVal are tested by means of a Monte Carlo test, based on 999 randomizations and were performed using the function multipatt of the ‘indicspecies’ package (De Cáceres and Legendre 2009).

For this purpose, a species matrix (150 species) was created using the average abundance data between years for each plot and sampling period. A logarithmic transformation was applied to reduce the influence of extreme abundance values (Legendre and Legendre 2012).

We grouped together sampling sites into four altitudinal bands (1: 550–1200 m; 2: 1250–1650 m; 3: 1700–2150 m; 4: 2200–2700 m) and two time periods (season I, 2006–2008; season II, 2012–2013).

This methodology allowed us to test the association of each species with the altitudinal bands over time, resulting in 255 possible groups, but we restricted our analyses to the 62 groups that could be interpreted ecologically in terms of altitudinal constancy or shift (Suppl. material 1: S3). Depending on the degree of species association (IndVal) to altitudinal bands between sampling periods, we identified the following altitudinal changes. We classified as altitudinal expansion the case in which a species showed, in the second sampling period, new significant associations for higher altitudinal bands, at the same time maintaining the association to the altitudinal bands of the first sampling period. An altitudinal upward shift occurred when the species lost its association to the lower altitudinal bands and gained new associations to higher ones during the second sampling period. At the same time, we identified downward retraction when a species lost its association to the lower bands, without gaining new associations with new higher ones. Then, we pointed out when the species became (spread) IndVal species or lost (disappear) their indicator value for specific altitudinal bands in the second sampling period. When no changes of associations occurred, the species remained stable over time (Suppl. material 1: S3).

Community level

Species’ responses may affect communities’ structure and composition (Wilson and Gutierrez 2012). Therefore, after a species level analysis, it is interesting to assess changes at the community level in order to compare the pattern of changes and to detect potential relationships. In this section, we investigated temporal changes of butterflies’ communities from three different perspectives: species richness, community composition and community temperature index.

Species richness

To evaluate if sampling efficiency remained comparable throughout sampling periods, and if we consequently can carry out comparisons at the assemblage level, we performed a preliminary analysis of species richness estimators and sample completeness (expressed as the proportion of the observed species richness to the averaged estimated species richness). We calculated the averaged estimated species richness, through the mean value of four abundance-based estimators (ACE Abundance-based Coverage Estimator; Chao1; Jacknife1; Jacknife2), and correlated it with the observed species richness. Correlation was high and significant in both sampling periods (Spearman test; 1^st sampling period, ρ = 0.937, p < 0.001; 2^nd sampling period, ρ = 0.918, p < 0.001). Sample completeness was also high in both seasons (mean sample completeness among plots ± standard error; 1^st = 67.59 ± 1.56; 2^nd = 72.17 ± 1.31). Consequently, we considered the sampling efficiency between sampling periods as comparable, allowing us to carry out direct data comparisons.

To analyze how species richness per plot changed over time, we compared species richness between sampling periods using the t-test for paired samples (significance level assessed after 999 randomizations, following Legendre and Legendre 2012).

To understand if changes in species richness were mainly related to characteristics of specific plots, we analyzed it as a function of the following variables: altitude, temperature, geographic location, dominant vegetation cover (habitat type), and dominant land use. We considered the rate-of-change as a dependent variable (hereinafter ROC), defined as the differences in species richness between sampling periods (1^st, 2006–2008; 2^nd, 2012–2013), divided by the species richness of the first sampling period (1^st, 2006–2008). We analyzed the ROC through linear regression, and we compared variables in a multi-model context, according to two criteria: (i) avoiding the simultaneous use in the same model of highly correlated predictors (Spearman’s rs >0.5); (ii) choosing predictors to represent biologically meaningful combinations of predictive variables, and consequently avoiding data dredging. All models were compared with the null model (intercept only) and all continuous variables were standardized, to enable comparisons between variables. The multivariate model selection was performed using Akaike information criteria, in its form corrected for small samples (AICc). As a measure of goodness-of-fit, we calculated the adjusted R². These analyses were performed by the R package ‘MuMIn’ (Barton 2018).

Community composition

We analyzed community compositions by testing for changes both in location (significant changes in community composition per site over time) and dispersion over the years (significant changes in observed differences in community composition between sites, over time). Changes in location were tested by applying non-parametric MANOVA to Bray-Curtis distance matrixes, to test if the multivariate centroids of species composition were similar in the two groups, or not (Anderson 2001; McArdle and Anderson 2001). Non-parametric MANOVA (permutational MANOVA) is an analysis of variance using distance matrixes, which partitions sums of squares of a multivariate dataset. It was performed by the function adonis of the ‘vegan’ package (Oksanen et al. 2012). The significance of the test was assessed by using F-tests based on sequential sums of squares obtained from permutations of the raw data (999 permutations). Since we had to keep in mind the temporal structure and the spatial dependencies of our sampling design (62 sites at two points in time), we applied a restricted randomization, which did not allow for permutations across samples. Changes in dispersion were tested by the betadisper function of the package ‘vegan’, a multivariate analogous of the Levene’s test for comparing group variances (Anderson 2001). To test if one time period is more variable than the other, non-Euclidean distances between objects and group centroids were handled by reducing the original distances to principal coordinates. To test for significance, we applied a similar randomization approach, as previously explained.

Community Temperature Index

The term “Species Temperature Index” (STI), refers to a quantitative description of the realized climatic niche of a species (Tayleur et al. 2016). To obtain such quantification over the Italian territory, we used presence data given by the database CkMap (Ruffo and Stoch 2005). CkMap is an atlas of the Italian fauna, promoted by the Italian Ministry of the Environment, which summarizes data related to the distribution of the Italian fauna on a gridded map (10×10 km). In our case, because we required STI that referred to Alpine populations, we only focused on North Italy, considering 1396 cells. Every cell had a value of “1” (if occupied by the species) or “0” (if species presence was uncertain). Temperature data were obtained by the maps of Metz et al. (2014), using annual mean temperature (BIO1). In this way, we calculated mean temperature values for each species (realized niche optimum).

We used STI to calculate the “Community Temperature Index” (CTI), quantified as the mean STI of all the species present in a given community. We calculated the CTI for each community (plot) and sampling period, and we analyzed the changes in CTI over time.

We compared CTI values between sampling periods by using a t-test for paired samples. As for the case of ROC, we analyzed the temporal changes in CTI (Δ CTI) by a linear regression in a multi-model context regression, as a function of the same environmental variables and models.

Temperature analysis

To quantify changes in weather conditions between sampling periods, identified as a possible driving factor of changes in biodiversity, we compared monthly temperatures (from June to September), as a function of year, altitude and additive or interactive effects. Models were compared with the null one (intercept only) with Linear Mixed Effect Models (using the plot as random effect) and all continuous variables were standardized, to allow comparisons between variables. The multivariate model selection was performed using Akaike information criterion, in its form corrected for small samples (AICc). The analysis has been carried out using the R packages ‘lme4’ (Bates et al. 2015) and ‘MuMIn’ (Barton 2018).

All the analyses were carried out with the R software, version 3.3.3 (R Core Team 2017).

Results

Considering all sampling sites and time periods, we recorded 150 butterfly species, of which 133 were common to both periods, 5 were exclusive to the first sampling session and 12 were exclusive to the second sampling session (Suppl. material 1: S2).

Species level

Changes in plot occupancy

We observed a general increase in mean occupancy levels (n = 150, t = -8.15, p = 0.001; plot/species 1^st = 8.85 ± 0.74, 2^nd = 12.50 ± 0.96, change = 3.65 ± 0.45).

The occupancy did not change equally between ecological groups. Concerning feeding specialization, strictly monophagous species differed from the other feeding groups (KW test, χ² = 9.82, df = 3, p = 0.020), even showing a slight decrease in the number of plots per species (polyphagous = 3.26 ± 0.99, one family = 4.10 ± 0.59, one genus = 3.68 ± 1.05, monophagous = -1.17 ± 0.98). We also recorded significant differences regarding the relationship with altitude. Altitudinal specialists increased less than the generalists (KW test, χ² = 13.13, df = 2, p = 0.001; generalists = 6.12 ± 1.03, medium = 3.04 ± 0.48, specialists = 2.32 ± 1.55) and also high altitude species showed a significantly less marked increase (MW test, W = 1070, p = 0.013, high altitude = 1.76 ± 1.32, others = 4.02 ± 0.46). “Shade loving” species showed the highest increase in mean occupancy levels (MW test, W = 2269, p = 0.041, “shade loving” = 4.53 ± 0.67, “sun loving” = 2.71 ± 0.57).

No significant differences have been recorded considering habitat (KW test, χ² = 2.04, p = 0.564) preferences and dispersal capacity.

We identified 17 species which significantly increased their distribution within the study area, of which 7 were also significant after Bonferroni correction. No species significantly decreased its occupancy level (Table 1).

Table 1.

Download as

CSV

XLSX

Species which significantly changed their area of occupancy, expressed as the number of plots of presence during the first (Plot1) and the second (Plot2) sampling period. The number of plots gained or lost (Δ), the p-value of the McNemar’s test (p-value) and its significance level (* p < 0.05, ** p < 0.001) after Bonferroni correction (p-adj) are also expressed.

Species	Plot 1	Plot 2	Δ	p-value	p-adj
Cyaniris semiargus	8	36	28	<0.001	**
Pieris bryoniae	16	43	27	<0.001	**
Cupido minimus	22	43	21	<0.001	**
Aglais urticae	37	55	18	0.004
Callophrys rubi	9	27	18	0.001	*
Pyrgus malvoides	20	38	18	<0.001	*
Argynnis niobe	20	36	16	<0.001	**
Maculinea arion	10	24	14	0.001
Colias alfacariensis	2	14	12	0.001	*
Lasiommata petropolitana	6	18	12	0.003
Pieris rapae	20	32	12	0.009
Polyommatus coridon	31	43	12	0.001
Eumedonia eumedon	3	14	11	0.001
Pyrgus cacaliae	2	13	11	0.010
Erebia euryale	34	44	10	0.009
Pyrgus carthami	1	11	10	0.003
Argynnis aglaja	33	42	9	0.009

Altitudinal range changes

We observed significant changes both at the minimum and at the maximum boundary of altitudinal range of the species. At the lower altitudinal limit, we observed a significant decrease in the absolute minimum value (t-test, n = 133, t = 3.03, p = 0.004, change = -96.62 m ± 31.85 m). At the higher altitudinal limit, we observed an increase in both the absolute maximum (t-test, n = 133, t = -3.08, p = 0.006, change = 75.19 m ± 24.01 m) and in the 90^th percentile (t-test, n = 133, t = -2.63, p = 0.014, change = 55.15 m ± 20.97 m).

We also observed significant differences in the changes of the altitudinal range between ecological groups. In particular, high altitude species, compared to the others, showed a significant increase in the mean (MW test, W = 1657.5, p = 0.041; high altitude = 59.21 ± 21.24, others = -11.07 ± 18.09), median (MW test, W = 1657, p = 0.040; high altitude = 59.37 ± 19.90, others = -8.72 ± 17.58) and 10^th percentile values (MW test, W = 1737, p = 0.012; high altitude = 68.12 ± 50.24, others = -72.11 ± 28.05).

“Shade loving” species lowered their minimum (MW test, W = 2738.5, p = 0.014; shade = -164.18 ± 45.34, sun = -28.03 ± 43.47) and 10^th percentile values (MW test, W = 2702, p = 0.027; shade = -74.10 ± 36.66, sun = -19.09 ± 34.14) compared to the “sun loving” species.

Strongly vagile species increased their minimum boundary, while weakly and medium vagile species lowered their minimum boundary (KW test, χ² = 8.34, df = 2, p = 0.015; high = 139.29 ± 110.22, medium = -103.49 ± 37.93, low = -178.79 ± 62.12).

We identified 87 species associated with one of the selected group combinations (Fig. 1; Suppl. material 1: S4). In most cases, IndVal species resulted as stable throughout the time period: there were 52 species belonging to 7 different combinations of altitudinal bands. A total of 12 species expanded their altitudinal range, spreading the associations through the higher altitudinal band and remaining stable at the lower ones (altitudinal expansion, 5 combinations of groups). Even if with a low indicator value (Celastrina argiolus, IndVal = 0.510, p-value = 0.018; Fig. 1, Suppl. material 1: S4), only one species showed an altitudinal upward shift, gaining the association to a higher altitudinal band and losing it into the lower one. The general increase in species occupancy also resulted in 19 species becoming indicator species during the second sampling period of specific altitudinal bands (spread). Only three species were indicators of the lowest altitudinal range only during the first time period, showing a decline during the second sampling season (disappear). For no species did we record a downward retraction.

Figure 1.

Changes along altitude. The number of significant IndVal species are shown for each combination of sampling periods (I = 1^st sampling period, II = 2^nd sampling period) and altitudinal bands (1, 550–1200 m; 2, 1250–1650 m; 3, 1700–2150 m; 4, 2200–2700 m). Different colours represent different statuses. Species that lost their indicator value with specific altitudinal bands during the second period (disappear; black, with continuous border); species stable over time (stable; white, with continuous border); species that showed, in the second sampling period, new significant associations for higher altitudinal bands, at the same time maintaining the association to the altitudinal bands of the first sampling period (altitudinal expansion; grey, with continuous border); species that lost their association to the lower altitudinal bands and gained new associations to higher ones during the second sampling period (altitudinal upward shift; grey with dotted border); species that became IndVal of specific altitudinal bands during the second sampling period (spread; white, with dotted border).

Community level

Species richness

Species richness significantly increased from the first to the second sampling period (t-test, n = 62, t = -9.76, p = 0.001, change = 8.82 ± 0.90).

The analysis of the ROC showed a significant effect of both land cover and land use: wooded habitats and managed plots increased the most (Table 2).

Table 2.

Download as

CSV

XLSX

Best linear regression model for Rate Of Change (ROC). Coefficients (± standard error) of the selected variables are indicated in the cells. Adjusted r squared is indicated as a measure of goodness of fit. alt = altitude; alt² = altitude (second order); park = geographic location; rme = change in mean temperature; rmi = change in minimum temperature; vegetation = dominant cover type (land cover); use = land use; Tmin = mean seasonal minimum temperature during the first season. adj r² = adjusted R². ** p < 0.01; * p < 0.05; ° p = 0.06.

		ROC
alt
alt²
park
rme
rmi
vegetation
	ecotone	-0.394 ± 0.141 **
	meadow	-0.277 ± 0.120 *
	rock	-0.349 ± 0.185 °
use	yes	0.257 ± 0.102 *
Tmin
R² adj		20.64

Community composition

Species composition significantly differed between time periods, but the obtained R² is extremely low, meaning a scarcely relevant change (non-parametric MANOVA, F-value = 5.87, r² = 1.91, p = 0.001). More interestingly, we observed a significant change in dispersion between sampling periods, with a lower dispersion around the median during the second sampling period (distance to centroids 1^st sampling period = 0.584, 2^nd sampling period = 0.528, p = 0.001; indicating a tendency towards homogenization; Fig. 2).

Figure 2.

Homogenisation of community composition. Box-plot of the distances to the centroid of community composition, during each sampling period. A reduction in species heterogeneity at the community level from the first to the second period can be clearly seen.

Community Temperature Index

We observed that the CTI significantly increased from the first to the second sampling period (t = -3.59; p = 0.001), indicating a common trend toward thermophily in all the sampled areas. Moreover, interestingly, we observed that the change in CTI over time was mainly dependent on the geographic position of the sampling plots (Fig. 3; R² = 14.17, p = 0.007; geographic location, p = 0.007). In particular, we recorded the significantly higher increase in CTI in the colder areas, i.e. in the plots located in VDNP, which is characterized by the lowest temperature among our study areas.

Figure 3.

Changes in Community Temperature Index. We observed significant differences in ΔCTI between geographic locations (protected areas). The parks located in colder areas clearly show the highest increase in CTI.

Temperature

The analysis of field-recorded temperatures indicated significant differences between sampling periods, concerning the monthly mean and minimum temperatures observed during July, August and September (Table 3). In particular, we observed a second sampling period significantly warmer than the first one with a significant increase of minimum (+1.22 °C) and mean (+0.83 °C) temperatures. The values of increased temperatures reported refer to the seasonal average (from June to September).

Table 3.

Download as

CSV

XLSX

Results from the best linear model indicating the role of altitude and year (and their interaction) for determining differences between sampling seasons for each month and temperature parameters (Temp).

Temp	Month	Intercept	Altitude	Year	R²
Mean	June	11.43	-2.46		80.2
	July	13.22	-2-72	0.11	91.9
	August	12.65	-2.62	0.2	87.8
	September	8.39	-2.45	0.28	89.2
Max	June	17.08	-1.76		46.8
	July	20.16	-2.03	-0.14	50.5
	August	19.35	-1.92		53.8
	September	15.38	-1.50		42.4
Min	June	7.33	-2.88		83.7
	July	8.3	-3.08	0.25	96
	August	8.29	-2.96	0.26	90.6
	September	4.22	--2.88	0.37	95.9

Discussion

Long-term monitoring, which is periodically and systematically repeated over time (usually decades), is crucial for correctly understanding the variables and mechanisms that determine species’ distribution responses and patterns of community composition as a consequence of climate and land use changes (Magurran et al. 2010, Magurran and Henderson 2010, Legendre and Gauthier 2014). Interestingly, some recent studies in the Alps have shown how significant changes can also occur over short periods, highlighting the importance of verifying distributional and community changes at temporal intervals of less than 10 years (short-term changes).

Roth et al. (2014) described changes in community composition of birds, butterflies and plants in the Swiss Alps along an altitudinal gradient, over just 8 years. Erschbamer et al. (2009) showed how plant species richness in the Dolomites was significantly higher after only 5 years. In this framework, our study on butterflies composition and distribution along altitudinal gradients in the NW Italian Alps represents an important step towards a better comprehension of biodiversity patterns in mountain ecosystems, even if it was restricted to a short time-frame (2006–2007 vs 2012–2013). Indeed, understanding the spatial and temporal dynamics of species-rich communities is critical for understanding how environmental changes can affect biodiversity (McCann 2007).

Driving factors of change

Although our research is based on a short time frame, it is important to consider the role of potential driving factors of change, by analyzing climatic and land cover changes (whenever occurring) as explicative predictors.

Concerning variations in temperature between sampling periods (2006–2008 vs 2012–2013), we observed that the second sampling period was warmer than the first one. These observed differences are coherent with the trend recorded by Beniston (2006) for the Alps and are likely correlated to climate change. Indeed, Beniston (2006) observed an increase in minimum and mean temperatures over time, but no change in the maximum temperature.

Unfortunately, no punctual and detailed data on precipitations were available to introduce them into the analytical framework. However, analyzing data from the official weather stations of the Regional Meteorological Service located in the study areas, we observed an interesting pattern in snow cover, characterized by a reduction in the number of days with snow and by a seasonal shift in snow cover towards the spring, in particular since 2005. Such conditions are indicative of a warmer climate and they correspond to a common pattern for the entire Alpine chain (Scherrer et al. 2016).

Due to the short time frame of analysis, we considered the main habitat type (dominant land cover) as a constant variable because no huge structural changes have potentially occurred. However, because of the natural evolution of vegetation following climatic changes and land abandonment, minimal but effective land cover changes in terms of microhabitat scale (e.g. shrubs encroachment) can be taken into account (Zurlo 2018).

Some sampling sites (plots) were characterized by low-intensity and sustainable grazing activities without any changes between sampling periods.

No other direct land use changes and habitat alterations occurred within our study sites because of the park status and LTER accredited site of our study area.

Species level

Plot occupancy

Despite the short time-frame under analysis, butterflies showed a general increase in mean occupancy levels but species-specific plot occupancy changes differed markedly among species. However, when grouping species according to homogeneous ecological traits, we observed consistent distributional responses, in accordance with several studies that have proved the effect of ecological and life-history traits in shaping species distribution (Forister et al. 2010, Chen et al. 2012, Auer and King 2014). Few studies to date have examined the underlying reasons for distribution changes and so far results have been equivocal (Angert et al. 2011; Buckley and Kingsolver 2012).

Monophagous, altitudinal specialists and high-altitude species appeared to be more limited than other species. These categories comprise species with a high level of specialization, consequently less prone to colonize new environments, even if climatic or environmental constraints are relaxed. In particular, monophagous species are strictly limited by the presence and the quality of their single larval host plant, and are already observed and predicted to be highly vulnerable to climatic/environmental changes (Romo et al. 2014). Our results concerning butterfly specialization are interesting and reflect what has been observed in central Europe concerning habitat specialization, where a decrease of specialized and low vagile species has been observed, along with an increase of generalist species and good dispersers (Habel et al. 2016).

A significantly higher increase in plot occupancy by “shade-loving” species compared to other species can be associated with a tendency towards a higher coverage of shrubs in open areas at low and medium altitude. In the European Alps, the effect of climate change is regionally affected by human activities. Cattle grazing in the alpine pastures has been decreasing throughout the past century, allowing rapid recolonization by trees and shrubs, where the treeline had been artificially lowered (e.g., Vittoz et al. 2008a; EEA 2010; Rocchia et al. 2018).

Interestingly, no species significantly decreased its area of occupancy, but 17 species increased their occupancy area, 7 of which in a particularly marked manner, still significant after Bonferroni correction (Table 1).

Aglais urticae and Argynnis niobe showed pronounced increase in plot occupancy. Our results are in contrast with other researches that reported a decline of these two species in other parts of Europe (Kulfan et al. 1997, Saarinen et al. 2003). Those species are altitudinal generalists with a high dispersal capability like most of the other species which increased their area of occupancy within our study area (e.g. Pieris rapae and Colias alfacariensis). We therefore consider that our observations could be only transient, due to the high dispersal capacity of those species and their tendency to carry out vertical migrations. Particularly interesting was the increase in mean occupancy of C. rubi. Since its low feeding specialization (Balletto et al. 2015), this species could have increased its abundance or spread into new areas, without having limitations due to food availability. However, although we highlighted important plot occupancy changes both at ecological groups and species level, we should also consider that they could be influenced by many factors, such as species detectability, population dynamics and phenology, which could, in turn, affect plot occupancy detection. To disentangle such patterns, long-term monitoring and more detailed study on the population dynamics of target species will be essential.

Altitudinal range changes

As for plot occupancy analysis, we noticed a general process of altitudinal change, considering all butterflies species together.

Looking at the ecological groups, we observed that, even along altitude, generalist and highly mobile species were more prone than specialist and less mobile ones to change their distributional ranges. Particularly, high altitude species showed a significantly higher increase in their mean, median and 10^th percentile altitudinal parameters, corresponding to a reduction in their lower altitudinal boundaries and in the surface available. Those species are already limited in their distribution. Their presence is, in many cases, limited by minimum temperatures (Pellissier et al. 2013) and, consequently, they cannot lower their altitudinal range nor, in many cases, raise it, due to the drastic decrease of vegetation cover at higher altitudes (a high occurrence of rock cover and a strong reduction of the availability of herbs and grasses). Moreover, the reduction in permafrost, which is a relatively new and rapid phenomenon, makes high altitude rock and screes unstable, preventing colonization by plant species (Cannone et al. 2007). However, we should also consider a limit of our sampling design, as we did not monitor plots over 2700 m a.s.l., consequently reducing our possibility to observe an expansion towards higher altitudes and a colonization of new plots by high altitude species.

Concerning altitudinal patterns, we used the IndVal test to identify species associated with specific altitudinal zones. Our analysis highlighted 52 species associated with one or more altitudinal zones at both time periods: these species are widely spread within their altitudinal zone, and stable over the investigated period, consequently representing species which – in the future – could be monitored to investigate important altitudinal changes. For 12 species, we observed new significant associations for higher altitudinal bands (altitudinal expansion). We recorded a change in association from a lower to a higher altitudinal band (upward shift) only for one species, Celastrina argiolus, which was associated with the first altitudinal zone (550–1000 m) during the first sampling period, while during the second sampling period, it resulted characteristic for the second zone (1250–1650 m) only (see statistical analyses – altitudinal range changes for the description of the relationships between IndVal association and altitudinal changes). It is a species that favors woodland and shade, and probably an increase in minimum temperature with an increase in the shrub land coverage had an important effect in the shift of this species.

In contrast, other studies in European countries observed on the opposite a higher amount of species shifting their altitudinal distribution, disappearing from low elevation sites and colonizing higher altitudes (e.g., Hill et al. 2002; Konvička et al. 2003). The differences from our work, where most of the species changing their altitudinal association did not disappear from the lowest altitude sites, could be related both to the shorter temporal frame we analyzed and to the lower altitudinal gradient we monitored. Indeed, our low altitude sites are located in the montane belt, where the effects of climate and environmental changes could be less severe than in the lowland plains.

Finally for 19 species, the changes in their distribution determined new associations for specific altitudinal bands during the second sampling period (spread). These species represent important taxa to be followed over time, to understand if this pattern represents the first signal of an upward expansion.

Community level

Species richness

We observed significant changes in species richness per plot within the analyzed period. Butterfly communities are known to quickly change their arrangement due to environmental changes (Thomas 2005), and previous studies have indicated that butterflies could respond even faster than other taxa (Devictor et al. 2012).

We observed a clear and significant increase in species richness within our temporal frame, mainly related to land cover and land management. The highest rate of change was clearly evident in the wooded areas (woodland clearings), while ecotones (transitional areas dominated by shrubs, and mainly located inside the subalpine belts) showed the lowest rate of change. Other authors also observed similar results, although they were mainly related to individual species abundance. For example, Sgardeli et al. (2016) observed a higher increase in species abundance simultaneously with an increase of temperature within forest areas with respect to other areas. Indeed, on days with high temperatures and solar radiation (which are exacerbated in open areas), wooded areas can exert a buffering effect, protecting individuals from extreme temperatures and reducing temperature leaps (Oliver and Morecroft 2014). Even thermophilous and “sun-loving” butterfly species can be inhibited in their flight activity by high temperatures and solar radiation, with potentially dangerous consequences that have not yet been fully understood (Cerrato et al. 2016).

In grazed areas, we observed an increase in species richness that was twice as much as in the unmanaged area. This can be probably explained by the low-intensity, sustainable grazing of these sites. Although intense grazing activities can be detrimental for butterflies (Scalercio et al. 2014), it has been previously observed that a balanced grazing pressure can increase the presence of plant species belonging to Poaceae and Fabaceae (Fischer and Wipf 2002), which represent the most exploited plant families that host the larvae of many butterfly species. Moreover, grazing maintains woodland clearings and open herbaceous areas below the tree line, that would be rapidly colonized by shrubs and trees, without such management activities (Nagy and Grabherr 2009). Indeed, balanced grazing activities have the fundamental role of maintaining micro-habitat heterogeneity, which is, in turn, of high importance in reducing the potential detrimental effects of climate change to butterflies in mountain ecosystems (Klečková et al. 2014; Konvička et al. 2016).

Community composition

Although we observed major differences in species distribution between our sampling periods, we did not observe significant differences in community composition. These results were probably due to a consequence of the short time-frame of analysis.

Studies that have demonstrated variation in butterfly community composition clearly take into account longer time-frames (e.g., Habel et al. 2016) and, to our knowledge, no investigations have reported significant changes in community composition in such a short time-frame. Interestingly, however, we noticed a tendency towards biotic homogenization in butterfly community composition. Biotic homogenization refers to the increase in biological similarity between communities, a replacement process leading to a decrease in distinctiveness in community composition over time, as a result of the replacement of some specialist species with other generalists, which become more uniformly distributed across previously different assemblages (Olden and Rooney 2006). Indeed, species respond individually to changing environmental conditions, depending mainly on their physiological characteristics and habitat requirements (e.g., Wilson and Gutierrez 2012). This determines new species assemblages, which can only be appreciated by examining the entire communities throughout time (e.g., Wilson et al. 2007).

For example, a similar change in community composition over time, accompanied by an increase in community similarity, has been observed in the data analysis from the UK Butterfly Monitoring Scheme over a period of 20 years (Gonzàlez-Megìas et al. 2008). This tendency towards biotic homogenization has been observed in recent decades in different taxa, following land cover and climatic changes (e.g., Bühler and Roth 2011; Eskildsen et al. 2015). These phenomena often determine an increase of generalist and of highly vagile species, to the detriment of other taxa (e.g. Menéndez et al. 2006; Bonelli et al. 2011).

Community Temperature Index

The observed increase in CTI tended to derive from the increase in mean occupancy levels recorded by a large number of generalist, more termophilous species in the second sampling period compared to the first one. This phenomenon was particularly marked in the coldest areas, representing a potential threat for high altitude and microthermic species. Indeed, the generalist species successfully colonizing a higher altitude from a lower one can outcompete the more specialized and vulnerable species of the alpine belt, challenging the survival of their local populations in the long term. Our results are consistent with an increase in CTI already observed in other geographic areas (e.g., in Greece by Zografou et al. 2014 over 13 years; in the Swiss Alps by Roth et al. 2014 over 8 years; at European scale by Devictor et al. 2012 over two decades). In any case, our trend was observed over a shorter time-frame, and, if confirmed during the next monitoring sessions, could represent a warning signal for alpine butterfly fauna.

Conclusions

In recent decades, changes in species distribution and composition of communities have already been analyzed in some parts of the Alps. Nevertheless, data usually concern plants (e.g. Grabherr et al. 1994; Pauli et al. 2001; Vittoz et al. 2008a, 2008b), or birds (Popy et al. 2010; Scridel et al. 2017), with very little information regarding other taxonomic groups, and mainly concerning changes in the distribution of single species (e.g., Battisti et al. 2005).

Responses to climate and habitat changes vary widely between species of the same communities (Walther 2010). Consequently, it is of crucial importance to understand the heterogeneity of species responses and the implications of the changes in community composition that follow (Wilson and Gutierrez 2012). Moreover, many previous studies that aimed at understanding changes over time in community composition have relied on the comparisons of contemporary data with historical data (atlases, collection specimens), which were often collected in a non-standardized way and/or referred to a much coarser spatial grain (Wilson and Gutierrez 2012).

Sampling units placed in well-specified areas represent a more appropriate tool for investigative purposes (e.g., Archaux 2004; Viterbi et al. 2013; Brandmayr and Pizzolotto 2016), and our data from the Italian Alps, taken from exactly the same sites at 4-year intervals, represent a first attempt to fill these gaps, focusing on changes in mountain ecosystems and mainly considering community level responses. Moreover, thanks to the LTER network, such standardized and temporally replicated data are available for broader spatial scale analysis and comparison with data coming from other sites of the same network, greatly improving our capacity to understand biodiversity patterns.

Interestingly, and somewhat alarmingly, our results suggest that, even if the time-frame under analysis is relatively short, we already observed considerable changes, in particular considering that the research was carried out in protected areas and LTER sites where habitat alteration by direct human effects is strongly reduced.

To summarize, in just 5 years, we observed:

– a general increase in mean occupancy level and in species richness of butterflies;

– no significant changes in the mean altitudinal optimum, but significant changes at both altitudinal limits.

Moreover, the observed changes differed across species, determining:

– an increase in shared species (a tendency to homogenization) within communities, even if the overall community composition did not change;

– a significant increase in Community Temperature Index (CTI).

Considering these results, it is now even more important to continue our monitoring in order to understand in the near future if our observed patterns represent only transient changes or are the first signals of an imminent trend.

Acknowledgements

We are grateful to the Park Directors for logistic support and data availability and to all Park Wardens, students and researchers for their essential help during the field work. The project leading to this research was coordinated by Gran Paradiso National Park and partially funded through the “Monitoraggio della biodiversità in ambiente alpino” grant (Progetti di Sistema ex cap.1551), provided by the Italian Ministry of the Environment (Ministero dell’Ambiente e della Tutela del Territorio e del Mare). This research also received funding from the European Union’s Horizon 2020 Research and Innovation Programme under grant agreement no. 641762 “Improving Future Ecosystem Benefits through Earth Observations” (Ecopotential) and the Project of Interest “NextData” of the Italian Ministry of Education, University and Research. We thank Sönke Hardersen, Pietro Brandmayr and Martin Konvička for their careful revisions on earlier versions of this paper.

References

Acquaotta F, Fratianni S, Garzena D (2014) Temperature changes in the North-Western Italian Alps from 1961 to 2010. Theoretical and Applied Climatology 122(3–4): 619–634. https://doi.org/10.1007/s00704-014-1316-7

Agroselviter (2009) Attività 1.1. Descrizione della vegetazione alpina in relazione alla conservazione della biodiversità. In: Progetto Interreg IIIA ALCOTRA – GestAlp – Modelli di gestione per la valorizzazione della biodiversità e del pastoralismo dei territori alpini e transfrontalieri.

Anderson MJ (2001) A new method for non-parametric multivariate analysis of variance. Austral Ecology 26: 32–46. https://doi.org/10.1111/j.1442-9993.2001.01070.pp.x

Angert AL, Crozier LG, Rissler LJ, Gilman SE, Tewksbury JJ, Chunco AJ (2011) Do species’ traits predict recent shifts at expanding range edges? Ecology Letters 14(7): 677–689. https://doi.org/10.1111/j.1461-0248.2011.01620.x

Archaux F (2004) Breeding upwards when climate is becoming warmer: No bird response in the French Alps. The Ibis 146(1): 138–144. https://doi.org/10.1111/j.1474-919X.2004.00246.x

Auer SK, King DI (2014) Ecological and life-history traits explain recent boundary shifts in elevation and latitude of western North American songbirds. Global Ecology and Biogeography 23(8): 867–875. https://doi.org/10.1111/geb.12174

Balletto E, Kudrna O (1985) Some aspects of the conservation of butterflies in Italy, with recommendations for a future strategy (Lepidoptera, Hesperiidae and Papilionoidea). Bollettino della Società Entomologica Italiana 117: 39–59.

Balletto E, Bonelli S, Barbero F, Casacci LP, Sbordoni V, Dapporto L, Scalercio S, Zilli A, Battistoni A, Teofili C, Rondinini C (2015) Lista Rossa IUCN delle Farfalle Italiane – Ropaloceri. Comitato Italiano IUCN e Ministero dell’Ambiente e della Tutela del Territorio e del Mare, Roma.

Barton K (2018) MuMIn: Multi-Model Inference. R package version 1.40.4. https://CRAN.R-project.org/package=MuMIn

Bates D, Maechler M, Bolker B, Walker S (2015) Fitting Linear Mixed-Effects Models Using lme4. Journal of Statistical Software 67: 1–48. https://doi: 10.18637/jss.v067.i01

Battisti A, Stastny M, Netherer S, Robinet C, Schopf A, Roques A, Larsson S (2005) Expansion of geographic range in the pine processionary moth caused by increased winter temperatures. Ecological Applications 15(6): 2084–2096. https://doi.org/10.1890/04-1903

Beniston M (2003) Climatic change in mountain regions: A review of possible impacts. Climatic Change 59(1/2): 5–31. https://doi.org/10.1023/A:1024458411589

Beniston M (2006) Mountain weather and climate: A general overview and a focus on climatic change in the Alps. Hydrobiologia 562(1): 3–16. https://doi.org/10.1007/s10750-005-1802-0

Bonelli S, Cerrato C, Loglisci N, Balletto E (2011) Population extinctions in the Italian diurnal lepidoptera: An analysis of possible causes. Journal of Insect Conservation 15(6): 879–890. https://doi.org/10.1007/s10841-011-9387-6

Brandmayr P, Pizzolotto R (2016) Climate change and its impact on epigean and hypogean carabid beetles. Periodicum Biologorum 118(3): 147–162. https://doi.org/10.18054/pb.2016.118.3.4062

Brook BW, Sodhi NS, Bradshaw CJA (2008) Synergies among extinction drivers under global change. Trends in Ecology & Evolution 23(8): 453–460. https://doi.org/10.1016/j.tree.2008.03.011

Brunetti M, Lentini G, Maugeri M, Nanni T, Auer I, Böhm R, Schöner W (2009) Climate variability and change in the greater alpine region over the last two centuries based on multi-variable analysis. International Journal of Climatology 29: 2197–2225. https://doi.org/10.1002/joc.1857

Buckley LB, Kingsolver JG (2012) Functional and Phylogenetic Approaches to Forecasting Species’ Responses to Climate Change. Annual Review of Ecology Evolution and Systematics 43(1): 205–226. https://doi.org/10.1146/annurev-ecolsys-110411-160516

Bühler C, Roth T (2011) Spread of common species results in local-scale floristic homogenization in grassland of Switzerland. Diversity & Distributions 17(6): 1089–1098. https://doi.org/10.1111/j.1472-4642.2011.00799.x

Cannone N, Sgorbati S, Guglielmin M (2007) Unexpected impacts of climate change on alpine vegetation. Frontiers in Ecology and the Environment 5(7): 360–364. https://doi.org/10.1890/1540-9295(2007)5[360:UIOCCO]2.0.CO;2

Cerrato C, Lai V, Balletto E, Bonellil S (2016) Direct and indirect effects of weather variability in a specialist butterfly. Ecological Entomology 41(3): 263–275. https://doi.org/10.1111/een.12296

Chamberlain DE, Negro M, Caprio E, Rolando A (2013) Assessing the sensitivity of alpine birds to potential future changes in habitat and climate to inform management strategies. Biological Conservation 167: 127–135. https://doi.org/10.1016/j.biocon.2013.07.036

Chen IC, Hill JK, Ohlemüller R, Roy D, Thomas CD (2012) Rapid Range Shifts of Species. Science 1024: 17–20. https://doi.org/10.1126/science.1206432

De Cáceres M, Legendre P (2009) Associations between species and groups of sites: Indices and statistical inference. Ecology 90(12): 3566–3574. https://doi.org/10.1890/08-1823.1

Dennis RLH (1993) Butterflies and climate change.Manchester University Press.

Devictor V, van Swaay C, Brereton T, Brotons L, Chamberlain D, Heliölä J, Herrando S, Julliard R, Kuussaari M, Lindström Å, Reif J, Roy DB, Schweiger O, Settele J, Stefanescu C, Van Strien A, Van Turnhout C, Vermouzek Z, WallisDeVries M, Wynhoff I, Jiguet F (2012) Differences in the climatic debts of birds and butterflies at a continental scale. Nature Climate Change 2(2): 121–124. https://doi.org/10.1038/nclimate1347

Dufrêne M, Legendre P (1997) Species assemblages and indicator species: The need for a flexible asymmetrical approach. Ecological Monographs 67: 345–366. https://doi.org/10.1890/0012-9615(1997)067[0345:SAAIST]2.0.CO;2

EEA (2010) Europe’s ecological backbone: recognising the true value of our mountains. Office for Official Publ. of the European Communities.

Erschbamer B, Kiebacher T, Mallaun M, Unterluggauer P (2009) Short-term signals of climate change along an altitudinal gradient in the South Alps. Plant Ecology 202(1): 79–89. https://doi.org/10.1007/s11258-008-9556-1

Eskildsen A, Carvalheiro LG, Kissling WD, Biesmeijer JC, Schweiger O, Høye TT (2015) Ecological specialization matters: Long-term trends in butterfly species richness and assemblage composition depend on multiple functional traits. Diversity & Distributions 21(7): 792–802. https://doi.org/10.1111/ddi.12340

Fischer M, Wipf S (2002) Effect of low-intensity grazing on the species-rich vegetation of traditionally mown subalpine meadows. Biological Conservation 104(1): 1–11. https://doi.org/10.1016/S0006-3207(01)00149-5

Forister ML, McCall AC, Sanders NJ, Fordycec JA, Thorned JH, O’Brien J, Waetjen DP, Shapiro AM (2010) Compounded effects of climate change and habitat alteration shift patterns of butterfly diversity. Proceedings of the National Academy of Sciences of the United States of America 107(5): 2088–2092. https://doi.org/10.1073/pnas.0909686107

Gaston KJ, Jackson SF, Cantú-Salazar L, Cruz-Piñón G (2008) The Ecological Performance of Protected Areas. Annual Review of Ecology Evolution and Systematics 39(1): 93–113. https://doi.org/10.1146/annurev.ecolsys.39.110707.173529

Gellrich M, Baur P, Koch B, Zimmermann NE (2007) Agricultural land abandonment and natural forest re-growth in the Swiss mountains: A spatially explicit economic analysis. Agriculture, Ecosystems & Environment 118(1–4): 93–108. https://doi.org/10.1016/j.agee.2006.05.001

Gonzàlez-Megìas A, Menéndez R, Roy D, Brereton T, Thomas CD (2008) Changes in the composition of British butterfly assemblages over two decades. Global Change Biology 14(7): 1464–1474. https://doi.org/10.1111/j.1365-2486.2008.01592.x

Grabherr G, Gottfried M, Pauli H (1994) Climate effects on mountain plants. Nature 369(6480): 448–448. https://doi.org/10.1038/369448a0

Habel JC, Segerer A, Ulrich W, Torchyk O, Weisser WW, Schmitt T (2016) Butterfly community shifts over two centuries. Conservation Biology 30(4): 754–762. https://doi.org/10.1111/cobi.12656

Hellmann JJ (2002) Butterflies as model systems for understanding and predicting climate change. In: Schneider SH, Root TL (Eds) Wildlife Responses to Climate Change. Island Press, Washington, DC 93–126.

Hill JK, Thomas CD, Fox R, Telfer MG, Willis SG, Asher J, Huntley B (2002) Responses of butterflies to twentieth century climate warming: Implications for future ranges. Proceedings. Biological Sciences 269(1505): 2163–2171. https://doi.org/10.1098/rspb.2002.2134

Hinojosa L, Napoléone C, Moulery M, Lambin EF (2016) The “mountain effect” in the abandonment of grasslands: Insights from the French Southern Alps. Agriculture, Ecosystems & Environment 221: 115–124. https://doi.org/10.1016/j.agee.2016.01.032

Huber U, Bugmann H, Reasoner MA (2005) Global Change and Mountain Regions. Page (Huber UM, Bugmann HKM, Reasoner MA, editors) Advances in global change research. Springer Netherlands, Dordrecht. https://doi.org/10.1007/1-4020-3508-X

Hunziker M (1995) The spontaneous reafforestation in abandoned agricultural lands – perception and aesthetic assessment by locals and tourists. Landscape and Urban Planning 31(1–3): 399–410. https://doi.org/10.1016/0169-2046(95)93251-J

Keller F, Körner C (2003) The Role of Photoperiodism in Alpine Plant Development. Arctic, Antarctic, and Alpine Research 35(3): 361–368. https://doi.org/10.1657/1523-0430(2003)035[0361:TROPIA]2.0.CO;2

Kharouba HM, Kerr JT (2010) Just passing through: Global change and the conservation of biodiversity in protected areas. Biological Conservation 143(5): 1094–1101. https://doi.org/10.1016/j.biocon.2010.02.002

Klečková I, Konvička M, Klecka J (2014) Thermoregulation and microhabitat use in mountain butterflies of the genus Erebia: Importance of fine-scale habitat heterogeneity. Journal of Thermal Biology 41: 50–58. https://doi.org/10.1016/j.jtherbio.2014.02.002

Konvička M, Maradova M, Beneš J, Fric Z, Kepka P (2003) Uphill shifts in distribution of butterflies in the Czech Republic: Effects of chaning climate detected on a regionale scale. Global Ecology and Biogeography 12(5): 403–410. https://doi.org/10.1046/j.1466-822X.2003.00053.x

Konvička M, Beneš J, Čížek O, Kuras T, Klečková I (2016) Has the currently warming climate affected populations of the mountain ringlet butterfly, Erebia epiphron (Lepidoptera: Nymphalidae), in low-elevation mountains? European Journal of Entomology 113: 295–301. https://doi.org/10.14411/eje.2016.036

Korhonen JJ, Soininen J, Hillebrand H (2010) A quantitative analysis of temporal turnover in aquatic species assemblages across ecosystems. Ecology 91(2): 508–517. https://doi.org/10.1890/09-0392.1

Kulfan J, Kulfan M, Zach P, Topp W (1997) Ist der Kleine Fuchs, Aglais urticae (Nymphalidae), in Zukunft gefährdet? Nota Lepidopterologica 20: 330–334.

Legendre P, Gauthier O (2014) Statistical methods for temporal and space-time analysis of community composition data. Proceedings. Biological Sciences 281(1778): 20132728. https://doi.org/10.1098/rspb.2013.2728

Legendre P, Legendre L (2012) Numerical ecology. 3^rd edition. Elsevier Science BV (Amsterdam).

Lemoine N, Bauer H-G, Peintinger M, Böhning-Gaese K (2007) Effects of Climate and Land-Use Change on Species Abundance in a Central European Bird Community. Conservation Biology 21(2): 495–503. https://doi.org/10.1111/j.1523-1739.2006.00633.x

Magurran AE, Henderson PA (2010) Temporal turnover and the maintenance of diversity in ecological assemblages. Philosophical Transactions of the Royal Society of London. Series B, Biological Sciences 365(1558): 3611–3620. https://doi.org/10.1098/rstb.2010.0285

Magurran AE, Baillie SR, Buckland ST, Dick JM, Elston DA, Scott EM, Smith RI, Somerfield PJ, Watt AD (2010) Long-term datasets in biodiversity research and monitoring: Assessing change in ecological communities through time. Trends in Ecology & Evolution 25(10): 574–582. https://doi.org/10.1016/j.tree.2010.06.016

Mantyka-Pringle CS, Martin TG, Rhodes JR (2012) Interactions between climate and habitat loss effects on biodiversity: A systematic review and meta-analysis. Global Change Biology 18(4): 1239–1252. https://doi.org/10.1111/j.1365-2486.2011.02593.x

McArdle BH, Anderson MJ (2001) Fitting multivariate models to community data: A comment on distance-based redundancy analysis. Ecology 82(1): 290–297. https://doi.org/10.1890/0012-9658(2001)082[0290:FMMTCD]2.0.CO;2

McCann K (2007) Protecting biostructure. Nature 446(7131): 29–29. https://doi.org/10.1038/446029a

Meloni F, Aronica L, Odasso M (2009) Studio della vegetazione nell’ambito del Progetto di Monitoraggio della Biodiversità nel Parco Orsiera Rocciavrè e nelle Riserve di Chianocco e Foresto.

Menéndez R, Megías AG, Hill JK, Braschler B, Willis SG, Collingham Y, Fox R, Roy DB, Thomas CD (2006) Species richness changes lag behind climate change. Proceedings. Biological Sciences 273(1593): 1465–1470. https://doi.org/10.1098/rspb.2006.3484

Metz M, Rocchini D, Neteler M (2014) Surface temperatures at the continental scale: Tracking changes with remote sensing at unprecedented detail. Remote Sensing 6(5): 3822–3840. https://doi.org/10.3390/rs6053822

Morecroft MD, Bealey CE, Beaumont DA, Benham S, Brooks DR, Burt TP, Critchley CNR, Dick J, Littlewood NA, Monteith DT, Scott WA, Smith RI, Walmsley C, Watson H (2009) The UK Environmental Change Network: Emerging trends in the composition of plant and animal communities and the physical environment. Biological Conservation 142(12): 2814–2832. https://doi.org/10.1016/j.biocon.2009.07.004

Nagy L, Grabherr G (2009) The Biology of Alpine Habitats. Oxford University Press, New York.

Nogués-Bravo D, Araújo MB, Errea MP, Martínez-Rica JP (2007) Exposure of global mountain systems to climate warming during the 21^st Century. Global Environmental Change 17(3–4): 420–428. https://doi.org/10.1016/j.gloenvcha.2006.11.007

Oksanen J, Kindt R, Legendre P, O’Hara B, Simpson GL, Stevens MHH, Wagner H (2012) Vegan: Community Ecology Package.

Olden JD, Rooney TP (2006) On defining and quantifying biotic homogenization. Global Ecology and Biogeography 15(2): 113–120. https://doi.org/10.1111/j.1466-822X.2006.00214.x

Oliver TH, Morecroft MD (2014) Interactions between climate change and land use change on biodiversity: Attribution problems, risks, and opportunities. Wiley Interdisciplinary Reviews: Climate Change 5(3): 317–335. https://doi.org/10.1002/wcc.271

Parmesan C, Yohe G (2003) A globally coherent fingerprint of climate change impacts across natural systems. Nature 421(6918): 37–42. https://doi.org/10.1038/nature01286

Pauli H, Gottfried M, Grabherr G (2001) High summits of the Alps in a changing climate: the oldest observation series on high mountain plant diversity in Europe. In: Walther GR, Burga CA, Edwards PJ (Eds) ‘Fingerprints’ of climate change: adapted behaviour and shifting species ranges. Kluwer Academic/Plenum Publishers, New York. https://doi.org/10.1007/978-1-4419-8692-4_9

Pauli H, Gottfried M, Reiter K, Klettner C, Grabherr G (2007) Signals of range expansions and contractions of vascular plants in the high Alps: Observations (1994–2004) at the GLORIA master site Schrankogel, Tyrol, Austria. Global Change Biology 13(1): 147–156. https://doi.org/10.1111/j.1365-2486.2006.01282.x

Pauli H, Gottfried M, Dullinger S, Abdaladze O, Akhalkatsi M, Alonso JLB, Coldea G, Dick J, Erschbamer B, Calzado RF, Ghosn D, Holten JI, Kanka R, Kazakis G, Kollar J, Larsson P, Moiseev P, Moiseev D, Molau U, Mesa JM, Nagy L, Pelino G, Puscas M, Rossi G, Stanisci A, Syverhuset AO, Theurillat JP, Tomaselli M, Unterluggauer P, Villar L, Vittoz P, Grabherr G (2012) Recent Plant Diversity Changes on Europe’s Mountain Summits. Science 336(6079): 353–355. https://doi.org/10.1126/science.1219033

Pellissier L, Alvarez N, Espíndola A, Pottier J, Dubuis A, Pradervand J-N, Guisan A (2013) Phylogenetic alpha and beta diversities of butterfly communities correlate with climate in the western Swiss Alps. Ecography 36(5): 541–550. https://doi.org/10.1111/j.1600-0587.2012.07716.x

Pizzolotto R, Gobbi M, Brandmayr P (2014) Changes in ground beetle assemblages above and below the treeline of the Dolomites after almost thirty years (1980/2009). Ecology and Evolution 4(8): 1284–1294. https://doi.org/10.1002/ece3.927

Pollard E (1977) A method for assessing changes in the abundance of butterflies. Biological Conservation 12(2): 115–134. https://doi.org/10.1016/0006-3207(77)90065-9

Popy S, Bordignon L, Prodon R (2010) A weak upward elevational shift in the distributions of breeding birds in the Italian Alps. Journal of Biogeography 37(1): 57–67. https://doi.org/10.1111/j.1365-2699.2009.02197.x

R Core Team (2017) R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/

Rocchia E, Luppi M, Dondina O, Orioli V, Bani L (2018) Can the effect of species ecological traits on birds’ altitudinal changes differ between geographic areas? Acta Oecologica 92: 26–34. https://doi.org/10.1016/j.actao.2018.08.001

Romo H, García-Barros E, Márquez AL, Moreno JC, Real R (2014) Effects of climate change on the distribution of ecologically interacting species: butterflies and their main food plants in Spain. Ecography 37:1063/1072. https://doi.org/10.1111/ecog.00706

Roth T, Plattner M, Amrhein V (2014) Plants, birds and butterflies: Short-term responses of species communities to climate warming vary by taxon and with altitude. PLoS One 9(1): e82490. https://doi.org/10.1371/journal.pone.0082490

Ruffo S, Stoch F (2005) Checklist e distribuzione della fauna italiana. Memorie del Museo Civico di Storia Naturale di Verona, 2.serie, Sezione Scienze Della Vita 16.

Saarinen K, Lahti T, Marttila O (2003) Population trends of Finnish butterflies (Lepidoptera: Hesperioidea, Papilionoidea) in 1991–2000. Biodiversity and Conservation 12(10): 2147–2159. https://doi.org/10.1023/A:1024189828387

Sala OE, Chapin FS, Armesto JJ, Berlow E, Bloomfield J, Dirzo R, Huber-Sanwald E, Huenneke LF, Jackson RB, Kinzig A, Leemans R, Lodge DM, Mooney HA, Oesterheld M, Poff NL, Sykes MT, Walker BH, Walker M, Wall DH (2000) Global biodiversity scenarios for the year 2100. Science 287(5459): 1770–1774. https://doi.org/10.1126/science.287.5459.1770

Scalercio S, Bonacci T, Mazzei A, Pizzolotto R, Brandmayr P (2014) Better up, worst down: Bidirectional consequences of three decades of changes on a relict population of Erebia cassioides. Journal of Insect Conservation 18(4): 643–650. https://doi.org/10.1007/s10841-014-9669-x

Scherrer SC, Begert M, Croci-Maspoli M, Appenzeller C (2016) Long series of Swiss seasonal precipitation: Regionalization, trends and influence of large-scale flow. International Journal of Climatology 36(11): 3673–3689. https://doi.org/10.1002/joc.4584

Scridel D, Bogliani G, Pedrini P, Iemma A, von Hardenberg A, Brambilla M (2017) Thermal niche predicts recent changes in range size for bird species. Climate Research 73(3): 207–216. https://doi.org/10.3354/cr01477

Scridel D, Brambilla M, Martin K, Lehikoinen A, Iemma A, Matteo A, Jähnig S, Caprio E, Bogliani G, Pedrini P, Rolando A, Arlettaz R, Chamberlain D (2018) A review and meta-analysis of the effects of climate change on Holarctic mountain and upland bird populations. The Ibis 160(3): 489–515. http://doi.wiley.com/10.1111/ibi.12585. https://doi.org/10.1111/ibi.12585

Sgardeli V, Zografou K, Halley JM (2016) Climate change versus ecological drift: Assessing 13 years of turnover in a butterfly community. Basic and Applied Ecology 17(4): 283–290. https://doi.org/10.1016/j.baae.2015.12.008

Sokal RR, Rohlf FJ (1995) Biometry: the principles and practice of statistics in biological research, 3^rd edn. WH Freeman and company, New York.

Stefanescu C, Carnicer J, Peñuelas J (2011) Determinants of species richness in generalist and specialist Mediterranean butterflies: The negative synergistic forces of climate and habitat change. Ecography 34(3): 353–363. https://doi.org/10.1111/j.1600-0587.2010.06264.x

Tayleur CM, Devictor V, Gaüzère P, Jonzén N, Smith HG, Lindström Å (2016) Regional variation in climate change winners and losers highlights the rapid loss of cold-dwelling species. Diversity & Distributions 22(4): 468–480. http://doi.wiley.com/10.1111/ddi.12412. https://doi.org/10.1111/ddi.12412

Thomas JA (2005) Monitoring change in the abundance and distribution of insects using butterflies and other indicator groups. Philosophical Transactions of the Royal Society of London. Series B, Biological Sciences 360(1454): 339–357. https://doi.org/10.1098/rstb.2004.1585

Viterbi R, Cerrato C, Bassano B, Bionda R, Hardenberg A, Provenzale A, Bogliani G (2013) Patterns of biodiversity in the northwestern Italian Alps: A multi-taxa approach. Community Ecology 14(1): 18–30. https://doi.org/10.1556/ComEc.14.2013.1.3

Vittoz P, Rulence B, Largey T, Freléchoux F (2008a) Effects of Climate and Land-Use Change on the Establishment and Growth of Cembran Pine (Pinus cembra L.) over the Altitudinal Treeline Ecotone in the Central Swiss Alps. Arctic, Antarctic, and Alpine Research 40(1): 225–232. https://doi.org/10.1657/1523-0430(06-010)[VITTOZ]2.0.CO;2

Vittoz P, Bodin J, Ungricht S, Burga CA, Walther GR (2008b) One century of vegetation change on Isla Persa, a nunatak in the Bernina massif in the Swiss Alps. Journal of Vegetation Science 19(5): 671–680. https://doi.org/10.3170/2008-8-18434

Vittoz P, Cherix D, Gonseth Y, Lubini V, Maggini R, Zbinden N, Zumbach S (2013) Climate change impacts on biodiversity in Switzerland: A review. Journal for Nature Conservation 21(3): 154–162. https://doi.org/10.1016/j.jnc.2012.12.002

Walther GR (2010) Community and ecosystem responses to recent climate change. Philosophical Transactions of the Royal Society of London. Series B, Biological Sciences 365(1549): 2019–2024. https://doi.org/10.1098/rstb.2010.0021

Walther GR, Beißner S, Burga CA (2005) Trends in the upward shift of alpine plants. Journal of Vegetation Science 16(5): 541–548. https://doi.org/10.1111/j.1654-1103.2005.tb02394.x

Wilson R, Gutierrez D (2012) Effects of climate change on the elevational limits of species ranges. In: Beever EA, Belant JL (Eds) Ecological consequences of climate change. CRC Press, Taylor and Francis Group (NW, USA)

Wilson RJ, Gutiérrez D, Gutiérrez J, Monserrat VJ (2007) An elevational shift in butterfly species richness and composition accompanying recent climate change. Global Change Biology 13(9): 1873–1887. https://doi.org/10.1111/j.1365-2486.2007.01418.x

Zografou K, Kati V, Grill A, Wilson RJ, Tzirkalli E, Pamperis LN, Halley JM (2014) Signals of Climate Change in Butterfly Communities in a Mediterranean Protected Area. PLoS One 9(1): e87245. https://doi.org/10.1371/journal.pone.0087245

Zurlo M (2018) Application of multi-scale remote sensing imagery in monitoring of a mountain protected area. Doctoral Thesis. https://doi.org/10.13130/zurlo-michele_phd2018-02-06

Supplementary material

Supplementary material 1

Supplementary data

Cristiana Cerrato, Emanuel Rocchia, Massimo Brunetti, Radames Bionda, Bruno Bassano, Antonello Provenzale, Simona Bonelli, Ramona Viterbi

This dataset is made available under the Open Database License (http://opendatacommons.org/licenses/odbl/1.0/). The Open Database License (ODbL) is a license agreement intended to allow users to freely share, modify, and use this Dataset while maintaining this same freedom for others, provided that the original source and author(s) are credited.

Download file (76.17 kb)