Impacts of climate change on the gross primary production of Italian forests

Abstract

Key message

The impact of climatic change should not be dramatic over Italian forests in terms of GPP, which should increase particularly for evergreen forest types. This positive effect is less marked for deciduous forests. The increasing trend should be reduced by the end of the century for all forest types except mountain conifers because of increasing temperature and decreasing rainfall.

Context

Estimating the spatial and temporal variability of forest gross primary production (GPP) is a major issue of applied ecology, particularly in relation to ongoing and expected climate change.

Aims

The current study proposes a methodological framework for analyzing large-scale forest responses to climate change in terms of GPP.

Methods

The methodology utilizes the GPP estimates of an NDVI-driven model, C-Fix, to assess the performance of a biogeochemical model, BIOME-BGC. The two models were first applied at 1-km pixel scale in Italy over a period of 15 years (1999–2013). The model outputs, aggregated on annual basis for the main Italian forest types, were inter-compared and analyzed in relation to major meteorological drivers (i.e., temperature and water-limiting factors).

Results

C-Fix and BIOME-BGC responded similarly to these major drivers, which supported the application of BIOME-BGC as a prognostic tool to simulate the GPP during three time slices of the RCP4.5 climate scenario.

Conclusion

The results obtained highlight how the importance of spring temperature and water availability is diversified among the forest types in determining changes of forest GPP all over the Italian peninsula in a future climate.

1 Introduction

Long-term climatic observations indicate that the atmospheric concentration of CO₂ and other greenhouse gases has been increasing since the industrial revolution (Meinshausen et al. 2011) and that this increment produces global climate change affecting all earth ecosystems (Campbell et al. 2009; Moriondo et al. 2011). Both terrestrial and aquatic ecosystems, however, can counteract this change by variably behaving as carbon sink (Waring and Running 2007). In particular, forests play an important role within the global carbon cycle, because they store about 45% of the terrestrial carbon (Bonan 2008).

The Mediterranean climate is characterized by mild winters and warm and dry summers, which can give origin to long periods of water shortage. Water availability is therefore the main factor constraining vegetation growth, while temperature limitation is generally more irregular and confined to colder zones (Bolle et al. 2006). The climatic changes expected in the Mediterranean area includes opposite trends for temperature and precipitation, as the former is supposed to increase between 2 and 4 °C, while the latter is supposed to decrease (Giorgi et al. 2004). The rainfall reduction, however, should not be uniform during the year; the number of rainy days per year should decrease, the intensity of rainfall should increase, but both phenomena are not easy to predict (De Dios et al. 2007). The expected responses of Mediterranean forests to climate change are therefore a non-trivial issue, whose analysis is further complicated by the spatial heterogeneity and fragmentation of existing landscapes (Scarascia-Mugnozza et al. 2000).

The approaches currently available to investigate the effects of climate on vegetation can be based on manipulating experiments (e.g., using Free Air CO₂ Enrichment, FACE, experiment or controlled areas) (Norby et al. 2016), on dendrochronological studies (e.g., Misson et al. 2004) and/or on modeling approaches (Potzelsberger et al. 2015), mostly depending on the spatial and temporal scales under investigation. The advantages of using models are that these tools can be applied over large areas both for diagnostic purposes, analyzing the present conditions of vegetation and trying to understand the role of the main climatic drivers on ecosystem functions, and for prognostic aims, therefore trying to evaluate the effect of plausible climate changes. However, the use of models is often complex, since it should be based on the comprehension of all biogeochemical processes and on the selection of those to be simulated. Furthermore, most models require abundant climatic and ecophysiological information plus a proper calibration phase to correctly simulate the behavior of specific ecosystems (Waring and Running 2007).

Numerous investigations have been conducted on the development and application of forest ecosystem models at different spatial and temporal scales (e.g., Jung et al. 2007; Campbell et al. 2009; Misson et al. 2004; Gea-Izquierdo et al. 2015). Some of these studies focused on modeling the main ecosystem processes of Mediterranean forests at regional scale (e.g., Maselli et al. 2009a, 2012). Among the models tested, C-Fix combines meteorological data and remotely sensed observations of the fraction of photosynthetic radiation absorbed by vegetation (fAPAR) for directly estimating the gross primary production (GPP) of terrestrial ecosystems (Veroustraete et al. 2002; Maselli et al. 2009b). The use of fAPAR estimates allows an efficient accounting for the effects of major environmental factors but of course prevents the model applicability for prognostic purposes. A more complex simulation of all forest processes, including respirations and allocations, was obtained by a well-known biogeochemical model, BIOME-BGC (Thornton et al. 2002; Turner et al. 2016). This model has daily meteorological data as main driving variables and is capable of working at high spatial and temporal scales. Consequently, it can be applied for prognostic studies and specifically for simulating the impact of climate change on environmentally heterogeneous and complex ecosystems (Potzelsberger et al. 2015).

The application of BIOME-BGC to investigate the medium-long-term evolution of major eco-physiological phenomena over large areas, however, should be preceded by a clear demonstration of the model capacity to correctly respond to major climatic factors in all relevant environmental conditions. The current work postulates that this demonstration can be obtained by inter-comparing the outputs produced over large areas and a known historical period by C-Fix and BIOME-BGC, which work using different logics and input data, thus providing different perspectives for the analysis of forest ecosystem processes. In particular, bearing in mind that the inter-comparison between different models is recognized as an indirect validation approach in case of absence of reference data (Dan et al. 2007), the GPP estimates obtained by C-Fix can be used to assess the functioning of BIOME-BGC and its response to main meteorological drivers over large areas and relatively long time periods.

Based on this premise, the current study is aimed at analyzing the impact of ongoing and future climate change on the GPP of Italian forests. The ecological complexity and spatial heterogeneity of the study area are addressed by the application of BIOME-BGC at high spatial (1 km) and temporal (daily) resolutions. The model capacity to respond to major climate drivers is preliminarily assessed by inter-comparing its GPP estimates with those obtained by C-Fix. Next, BIOME-BGC is applied at the same resolutions in a prognostic mode to simulate the effects of the RCP4.5 climate scenario. The choice of this scenario, which considers a mitigation level intermediate between RCP2.6 (implying a worldwide agreement for aggressive mitigation practices) and RCP8.5 (implying no mitigation practices), responds to our basic research question, i.e., if a reasonably mitigated climate could have significant effects on Italian forest ecosystems.

The study is therefore divided into two main phases: (i) production and analysis of large-area GPP estimates from C-Fix and BIOME-BGC over a reference 15-year period (1999–2013); (ii) application of BIOME-BGC to simulate the impact of the RCP4.5 scenario during three time slices.

2 Materials and methods

2.1 Main features of Italian forests

Italy is characterized by the presence of two main mountain chains, the Alps in the North and the Apennines along the whole peninsula, from the Centre to the South (Fig. 1). The Italian climate therefore varies with the latitudinal and altitudinal gradients and the distance from the sea, ranging from Mediterranean warm to temperate cool. These characteristics influence Italian forests, which, according to the last national forest inventory (www.infc.it, INFC), cover an area of about 8.8 millions of hectares. The most widespread broadleaf forests are dominated by oak species (genus Quercus), some of which are evergreen (both Q. suber and Q. ilex), and by beech (Fagus sylvatica). Concerning conifers, the most abundant species is Norway spruce (Picea abies), together with mountain and Mediterranean pines (genus Pinus). Following autoecological criteria, these forests were grouped into seven forest types: Mediterranean macchia, plain/hilly conifers, evergreen oaks, deciduous oaks, chestnut, beech, and mountain conifers (Chirici et al. 2016) (Fig. 1). The main environmental characteristics of each forest type are listed in Table 1: a gradient of decreasing mean annual temperature and potential evapotranspiration (PET) and increasing rainfall can be observed moving from Mediterranean macchia to mountain conifers.

Fig. 1

Map of the seven forest types (FTs) in Italy: Mediterranean macchia, plain/hilly conifers, holm oak, deciduous oaks, chestnut, beech, and mountain conifers. The map is comprised between 5° 30′–18° 30′ long. E and 36°–47° 30′ lat. N. The upper right box shows a DEM of the country

Table 1 Main environmental characteristics of the seven forest types (averages ± standard deviations). Meteorological data are referred to the 1999–2013 period (Tav, Rain, and PET stand for mean annual air temperature, total rainfall, and total potential evapotranspiration, respectively)

2.2 Models applied

C-Fix is a Monteith-type model, which is capable of directly estimating forest GPP based on the quantity of absorbed photosynthetically active radiation converted into biomass (Veroustaete et al. 2002; Maselli et al. 2009b). Daily GPP is therefore computed as:

$$ {\mathrm{GPP}}_{\mathrm{i}}=\varepsilon \ast {\mathrm{CO}}_2\mathrm{fert}\ast {\mathrm{Tcor}}_{\mathrm{i}}\ast {\mathrm{Cws}}_{\mathrm{i}}\ast {\mathrm{fAPAR}}_{\mathrm{i}}\ast {\mathrm{PAR}}_{\mathrm{i}} $$

(1)

where ε is the maximum radiation use efficiency (1.2 g C MJ⁻¹), CO₂fert is the CO₂ fertilization factor, Tcor is the temperature correction factor, Cws is the water stress factor, and fAPAR is the fraction of absorbed photosynthetically active radiation (PAR), all referred to day i. Tcor is derived from minimum daily temperature, while Cws is obtained from a simplified site water balance, as fully described in Maselli et al. (2009b); fAPAR is derived from satellite observations of the normalized difference vegetation index (NDVI) following the linear transformation of Myneni and Williams (1994). C-Fix has been widely applied and tested in several Italian sites; the model was demonstrated to be accurate in reproducing both annual totals and intra- and inter-annual variations of forest production (Maselli et al. 2009b; Chiesi et al. 2016, 2017).

BIOME-BGC (version 4.2) is capable of simulating the water, carbon, and nitrogen fluxes and storage within terrestrial ecosystems (Thornton et al. 2002). The model works on a daily time step, requiring in input meteorological data (daily minimum and maximum temperature, rainfall, and solar radiation hereafter Tmin, Tmax, Rain, and Rad), site-specific information (i.e., geographical position, altitude, soil texture, and depth), and ecophysiological information about the considered biome types. The simulation of evapotranspiration is based on the Penman-Monteith’s equation, and photosynthesis is based on Farquahar’s approach (see Golinkoff 2010). BIOME-BGC considers the effect of increasing atmospheric CO₂ concentration and nitrogen deposition and has been applied as a prognostic tool over large areas and long temporal periods (e.g., Mu et al. 2008; Hlásny et al. 2011). The current model version was calibrated for the main Italian forest types by Chiesi et al. (2007) and assessed against eddy-covariance flux tower data by Chiesi et al. (2011) and Chirici et al. (2016).

2.3 Study data

A digital elevation model (DEM) of Italy with 1-km spatial resolution was used to characterize the orography of the seven forest types. The spatial distribution of these formations was derived from the CORINE Land Cover 2006 map of Italy (ISPRA 2010). The original dataset divides forests into three broad classes, i.e., broadleaf, conifers, and mixed forests. This third level legend was however expanded to the fourth level enabling to identify the main forest species; these were then regrouped following autoecological criteria into the aforementioned seven forest types (Fig. 1).

Daily meteorological observations (Tmin, Tmax, and Rain) were derived from the Pan-European E-OBS dataset (Haylock et al. 2008). This dataset, which was originally distributed at 0.25° spatial resolution, was downscaled over the Italian territory to 1-km spatial resolution for the years from 1981 to 2013 (Maselli et al. 2012; Fibbi et al. 2016). This operation was performed by applying spatially (or geographically) weighted regressions to a DEM, as fully described in Maselli et al. (2012). The obtained rainfall estimates, which were affected by significant bias, were further corrected as described in Fibbi et al. (2016).

The outputs of the COSMO-CLM regional circulation model (RCM) for RCP4.5 scenario (daily Tmin, Tmax, and Rain) were selected to evaluate the impact of a changing climate on Italian forests. Bucchignani et al. (2016) observed that COSMO-CLM, working at high resolution (8 km as in this study), provides good performances in simulating Italian climate, where biases in predicted temperature and rainfall were generally lower than values usually encountered in RCM simulations working at lower spatial resolution (from 12 to 50 km). In any case, despite these good performances, since biases exist in COSMO-CLM simulations, the outputs of this model cannot be used directly to feed BIOME-BGC model but a downscaling procedure had to be applied. This was carried out by means of a delta-change approach (Maraun and Widmann 2018), according to which the bias is corrected by calculating the monthly average differences between the RCM COSMO-CLM simulations for the baseline period and the simulations under the RCP4.5 future scenario. These differences are then added to the observed meteorological dataset in order to derive daily datasets at a 1-km spatial resolution for the relevant time slices. More specifically, for Tmin and Tmax, the differences were calculated as absolute values between future and present monthly averages; for Rain, the change was calculated as ratio between future and present cumulated monthly rainfall. Rad was then obtained using the MT-Clim algorithm (Thornton et al. 2000), both for the current period and the future climate scenario. The CO₂ increases for the same periods were derived from the RCP4.5 climate scenario (Meinshausen et al. 2011). Since forcing factors related to the RCP4.5 scenario are implemented from 2006 onwards (Bucchignani et al. 2016), the baseline period corresponded to the first 25 years for which downscaled meteorological data were available (1981–2005), and the same length was maintained for the three time slices of the RCP4.5 future scenario (2006–2030, 2036–2060, 2066–2090).

Satellite NDVI images taken from the Spot-Vegetation (VGT) sensor were freely downloaded from the archive of VITO (http://free.vgt.vito.be). These images are distributed in a pre-processed format, which includes the radiometric, geometric, and atmospheric corrections (Maisongrande et al. 2004). Thirty-six 10-day maximum value composite (MVC) images per year, having a 1-km spatial resolution, were downloaded for the whole Italian peninsula, covering the period for which complete annual series were available (1999–2013).

2.4 Data processing

2.4.1 Application of C-Fix and BIOME-BGC in current condition

C-Fix and BIOME-BGC were first applied to simulate the GPP of Italian forests over the 1999–2013 common period. The two models were fed with the meteorological data of the pan-European E-OBS dataset downscaled to 1-km spatial resolution and corrected for rainfall bias (Maselli et al. 2012; Fibbi et al. 2016). NDVI data were pre-processed as fully described in Chirici et al. (2016). The two models were then applied on a per-pixel basis for the common period 1999–2013, yielding daily GPP estimates, which were aggregated on an annual basis for each forest type.

Since the effects of annual meteorology on forest GPP are quite complex, the analysis of this issue was based on the results of recent studies. In particular, Maselli et al. (2014) demonstrated that inter-annual GPP variations of Mediterranean forests are mostly determined by the start of the dry season (SDS), which can be obtained by elaborating spring Cws over the months March–July; the intensity of this control increases along an aridity gradient. The same gradient determines a decreasing intensity of temperature limitation, which is most influential in spring (from February to May), when it affects the beginning of the growing season (Chiesi et al. 2016). The dependence of the GPP simulated by both C-Fix and BIOME-BGC on spring Tmin and SDS was therefore investigated by applying correlation and regression analyses to the relevant annual values averaged for forest types. Next, an inter-comparison of the annual GPP estimates of the two models was carried out to assess the capacity of BIOME-BGC to simulate large-scale inter-annual GPP variations. The results were summarized using the regression coefficients (slope and intercept), the correlation coefficient (r), the root mean square error (RMSE), and the mean bias error (MBE).

2.4.2 Application of BIOME-BGC under the RCP4.5 scenario

The second part of the investigation was devoted to applying BIOME-BGC in a prognostic mode for simulating and analyzing the impact of climate change over Italian forests according to the RCP4.5 scenario. The first step of this phase consisted in the creation of the climate daily dataset for different future periods following the above-mentioned procedure. The model was then initialized over the baseline period (1981–2005) using the spin-up and go mode in order to identify the quasi-equilibrium conditions to be used as a starting point for the following simulations. The simulations for the time slices referred to the years 2006–2030, 2036–2060, and 2066–2090 were performed by using the restart files produced in the previous period. The effects of the climate change on GPP were then assessed on a per-pixel basis all over the Italian forests.

3 Results

3.1 Analysis of C-Fix and BIOME-BGC GPP estimates

During the period 1999–2013, a slight and not significant increase of spring temperatures is observed for all forest types, with the only exception of macchia (Table 2). Contemporaneously, SDS shows an increase, which is more marked and significant for species growing in the driest conditions (i.e., macchia, plain/hilly conifers, and evergreen oaks). These increases are associated with similar variations of NDVI averaged over the growing season (i.e., from April to October). The NDVI trend is particularly evident for the forest types growing in the driest conditions, whose r between NDVI and sequential year is higher than 0.88.

Table 2 Averages and slopes of spring Tmin, SDS, and growing season NDVI over the 1999–2013 years (see text for details) (* = significant slope, P < 0.05; ** = highly significant slope, P < 0.01)

In accordance with the temporal variations observed in both meteorological and satellite data, the annual C-Fix GPP estimates show similar increasing trends (Fig. 2). The increase of annual GPP is again most marked for the three forest types growing in the driest conditions, whose trends are all highly significant (r > 0.85). Figure 3 confirms the effects of spring temperature and SDS as main drivers of C-Fix GPP inter-annual variations; the highest correlations are found with SDS, which is the most important driver for typical Mediterranean species. On the contrary, spring minimum temperature mostly affects species living at the highest altitudes, i.e., chestnut, beech and, at a lesser degree, mountain conifers.

Fig. 2

Evolution of GPP obtained by C-Fix and BIOME-BGC for the seven FTs in Italy during the years 1999–2013. The letters indicate the corresponding forest types: a Mediterranean macchia; b plain/hilly conifers; c evergreen oaks; d deciduous oaks; e, chestnut; f, beech; g, mountain conifers. The statistics of each model are referred to the linear GPP trend on the sequential year, considering 1999 as year 1; the last statistics are referred to the model inter-comparison (** = highly significant correlations, P < 0.01)

Fig. 3

Correlation coefficients found between inter-annual variations of meteorological factors and GPP estimated by C-Fix and BIOME-BGC for the seven forest types during the 1999–2013 period. Panel a refers to spring minimum temperature, while b refers to SDS (** = highly significant correlations, P < 0.01)

These patterns are mostly reproduced well by BIOME-BGC. The GPP averages simulated by this model are in some cases moderately different from those of C-Fix (mainly for macchia and conifers), but the temporal trends from the two models are all similar. Figure 3 indicates that the inter-annual GPP variations estimated by the two models are also similarly dependent on both spring Tmin and SDS. As a consequence, the accordance between C-Fix and BIOME-BGC is generally good in terms of correlation coefficient (Fig. 2); r is higher than 0.85 for Mediterranean macchia, evergreen oaks, and plain/hilly conifers; it is about 0.8 for mountain conifers and higher than 0.65 for the three deciduous species. Mountain conifers show the highest mean difference between the GPP found by C-Fix (1350 g C m⁻² year⁻¹) and BIOME-BGC (1100 g C m⁻² year⁻¹).

3.2 Analysis of future BIOME-BGC GPP simulations

Table 3 shows the expected mean spring Tmin and SDS changes during the three selected RCP4.5 time slices. The spring temperature increases progressively in the three future time slices, up to reach about 2.3 °C for all ecosystems. Contemporaneously, a reduction of spring rainfall is expected for all future periods; SDS occurs progressively earlier for all forest types, up to a maximum of 10–14 days in the last period.

Table 3 Averages of spring Tmin (°C) and SDS (DOY) for the base period and the three future periods of the RCP4.5 scenario. All differences from the base to the future periods are highly significant from the Student’s t test (P < 0.01)

Figure 4 shows the mean forest GPP simulated by BIOME-BGC for the base period and the last slice of the future scenario. Clear GPP increases can be observed for the Mediterranean species (mainly macchia and evergreen oaks) in the southern part of Italy and, above all, for the mountain conifers on the Alps. The effects of the future climate scenario are confirmed by Fig. 5, where the mean annual productions for the base and the three selected periods are shown for each forest type. Forest GPP tends to variably increase from the base to the first and second periods. The highest production increases are observed for evergreen species (i.e., macchia, evergreen oaks, and conifers) during the first period (2006–2030, about + 10%), and, above all, during the second period (about + 15%, with a peak of + 37% for mountain conifers). In the third period (2066–2090), all evergreen species of arid and semi-arid zones show a marginal GPP decrement with respect to the previous scenario, while the mountain conifers go on increasing. The GPP of all deciduous species increases by around 5% in the first and second slices, but goes back to the base values during the last period considered.

Fig. 4

Forest GPP maps simulated by BIOME-BGC in the base period (1981–2005) (a) and in the last slice of the future scenario (2066–2090) (b)

Fig. 5

Average annual GPP simulated by BIOME-BGC for the seven forest types in the base period and in the future scenario, during three time slices. All differences from the base to the future periods are highly significant from the Student’s t test (P < 0.01)

4 Discussion

Italy shows a gradient of decreasing temperature and increasing humidity moving from southern, plain areas characterized by the presence of Mediterranean macchia, towards northern, Alpine zones dominated by mountain conifers. This determines variable and complex patterns of prevailing limiting factors, where forest growth is mainly constrained by water availability in the driest peninsular and insular areas and by temperature in the highest mountain zones.

The Monteith model used to estimate the GPP of Italian forests, C-Fix, has been widely applied and tested in numerous previous investigations (Veroustraete et al. 2002; Yuan et al. 2014). These investigations supported the model capability of reproducing seasonal and annual patterns of forest production in all eco-climatic conditions examined and specifically in Mediterranean areas (Maselli et al. 2009b; Chirici et al. 2016). BIOME-BGC was also successfully applied and tested in various environments all over the globe (Jung et al. 2007). Accordingly, the annual GPP estimates currently obtained by both models are reasonable and consistent with the tower GPP observations reported in Chiesi et al. (2011) and Chirici et al. (2016).

More recent research has demonstrated that C-Fix is also capable of reproducing inter-annual GPP variations, both in water- and temperature-limited areas (Chiesi et al. 2017): high and significant correlations were in fact found between GPP estimates and tree-ring widths collected in a Mediterranean pine forest (Maremma) and in a mountain beech forest (Mugello). This capacity can be ascribed to the combination of remotely sensed fAPAR estimates and ground meteorological data, which both characterize inter-annual production variations. In particular, Maselli et al. (2014) found that the dependence of C-Fix on SDS increases in intensity following an aridity gradient. This is confirmed by the current results, which also indicate the existence of an inverse gradient of influence for Tmin, with the partial exception of mountain conifers.

The responses of BIOME-BGC to both meteorological factors are, as currently demonstrated, similar to those of the NDVI-driven model. Consequently, this biogeochemical model simulates inter-annual forest GPP variations, which are significantly concordant with those of C-Fix. This is less the case for deciduous species, likely due to possible problems of BIOME-BGC phenology subroutine. This subroutine in fact works on the basis of site temperature, radiation, and day length and was tuned for temperate cold areas; consequently, its functioning might be problematic when it is applied in different environments (Trusilova et al. 2009). These findings are also in accordance with those of Jung et al. (2007), who inter-compared the performances of various models of forest processes in different environments: while BIOME-BGC was found to outperform the other models for Mediterranean evergreen ecosystems, the model simulation of GPP was more problematic for deciduous species. However, the use of the phenology subroutine allows a more efficient accounting for inter-annual meteorology variations and was therefore preferred to the alternative possibility of setting a fixed growing season length determined on climatic basis (Golinkoff, 2010).

The similar responses of the two models to the main drivers determine the simulation of similar GPP evolutions during the 1999–2013 study period. In this period, spring Tmin and, above all, SDS, are generally increasing. While this observation is expected for Tmin, due to the well-known tendency to global warming, the SDS trend is, to our knowledge, a novelty. The strength of the trend is particularly high in the driest areas, i.e., macchia, plain/hilly conifers, and evergreen oaks, where spring water availability increases markedly during the study period. The increases of both meteorological drivers are confirmed by the similar, significant increases of growing season NDVI found for all forest types. These NDVI trends are particularly evident for the driest ecosystems, supporting the hypothesis that they are mainly dependent on increasing water availability during the growing season. The combined effect of all these trends is obviously reflected in the GPP estimated by both models, which are similarly increasing for all ecosystems examined. It should be noted, however, that the general significance of these trends might be limited, due to the relatively short time period examined (1999–2013) and to the inherently high inter-annual variability, which particularly affects rainfall patterns in semi-arid zones (Deitch et al. 2017). Consequently, the persistence of these trends over longer time periods should be confirmed by the analysis of additional meteorological and NDVI datasets.

The results of the first research phase have provided a solid basis for applying BIOME-BGC in a prognostic mode and simulating forest GPP on a national scale for different time slices of the RCP4.5 climate scenario. Most regional circulation models converge in simulating a general increase in temperature and a decrease in precipitation in the Mediterranean area (Jacob et al. 2014) and the currently used model fully reproduces this trend for Italy (Bucchignani et al. 2016). More specifically, in Italy, temperature uniformly varies among all forest types, while for rainfall, the situation is more variable: most of the country shows a slight decrease in annual precipitation, which is specially due to a large drop in spring-summer period, while in the northern Alps, rainfalls are expected to be almost stable. This is evident from the spring temperature and SDS values shown in Table 3, which are only partly related to the annual totals. While these patterns do not fully agree with the previously observed SDS trends, it should be kept in mind that all predictions are generally more reliable for temperatures than for precipitation, particularly referred to specific seasons (Giorgi et al. 2004).

Overall, the prediction of annual GPP obtained in this study does not outline a pessimistic picture for Mediterranean forests in the future climate. This is in contrast with the evidences of the European 2003 heat-wave, when drought was the main reducing factor of carbon fluxes over most Mediterranean ecosystems (Reichstein et al. 2013).

As a matter of fact, the expected negative impact due to higher temperatures and reduced rainfall must be evaluated considering their interaction with the expected changes in atmospheric CO₂ concentration and nitrogen deposition relevant to the selected scenario. These higher concentrations should have a fertilizing effect on forest GPP, as long as nitrogen does not become a limiting factor for plant growth, thus overriding the negative effect of warmer climate (Gea-Izquierdo et al. 2017). This is especially true for mountain conifers, which, growing in temperature-limited areas, are the most advantaged from joint increase in temperature and carbon dioxide, yielding the highest forest production variations found in the scenario. The high GPP levels found in the last period for Mediterranean evergreen species could instead be attributed to the high resistance of these species to arid conditions, which use water-saving strategies. On the contrary, the last simulated period seems to have a negative impact on deciduous forests, which, being mostly placed on medium-high elevations, might be poorly adapted to increasing aridity. This last finding, however, could be partly due to the above mentioned model inefficiency in correctly simulating the GPP of deciduous species, particularly in consequence of changes in growing season length.

When examining these results, it should be also kept in mind that the model version currently applied considers only the direct effect of higher CO₂ concentration on photosynthetic efficiency, while the increased water use efficiency in response to the same phenomenon is not accounted for. The expected increase in water use efficiency should further reduce the impact of warmer temperature and drought as simulated in the future scenario (Keenan et al. 2013).

These results generally agree with previous studies obtained with different models for the same region (Gea-Izquierdo et al. 2017), which depict the high adaptability of Mediterranean forests to climatic change, but some elements of uncertainties should be highlighted. While the CO₂ concentration increase in the air is widely recognized to affect vegetation, the existing literature does not clearly converge on the effects of its interaction with climate change on forest adaptation (Baig et al. 2015). In many cases, experiments at high CO₂ concentration were carried out under current climate conditions on temperate forests (e.g., Battipaglia et al. 2013; Norby et al. 2016), while the combined effect of temperature and CO₂ should be assessed on biomes with more extreme conditions as expected in the next future. Consequently, until further evidences become available, the response of Mediterranean forests to the interaction between temperature and CO₂ is still a source of uncertainty for evaluating the adaptability of these ecosystems to a changed climate (Gea-Izquierdo et al. 2017).

The current simulation of the impact of climate change on forest GPP does not take into account the effects of abiotic and biotic factors (such as fires and winds among the former, insects and pathogens among the latter), which are not implemented in BIOME-BGC or show an unknown regional scale spatio-temporal variability. The same is the case for the consequences of possible future variations in forest extension and species distribution, which would require the application of a completely different modeling strategy that is beyond the scope of the present article.

5 Summary and conclusions

The current study proposes a methodological framework for analyzing large-scale forest responses to climate change. The methodology utilizes the GPP estimates of an NDVI-driven model, C-Fix, to assess the performance of a biogeochemical model, BIOME-BGC, which is then applied in a prognostic mode to simulate the impacts of a future scenario. This approach was applied in a typical Mediterranean country, Italy, where forests are particularly fragmented, fragile and vulnerable to climate change.

The first research phase was carried out over a 15-year period over all main Italian forest areas. The results of this phase confirm previous findings on the dependence of annual forest GPP on spring temperature and spring-summer water availability. This dependence varies in intensity from water limited to temperature-limited environments and is correctly accounted for by both models. Consequently, the positive trends of both temperature and water factors observed in the study period determine similar and concordant increases of the GPP simulated by the two models.

These results allowed the application of BIOME-BGC to simulate the impact of climate change as outlined in three different RCP4.5 time slices. This second research phase indicates that in the two intermediate periods GPP should increase particularly for evergreen forest types, both in the driest and most humid areas; the positive effect is instead less marked for deciduous forests, partly due to the uncertainty of the modeling approach. In all cases except mountain conifers, the expected GPP increases are limited by the decreased water availability, which should manifest in the most remote period as a consequence of increasing temperature and decreasing rainfall.

Appendix 1

Table 4 Parameter settings of BIOME-BGC utilized for the seven forest types