Ecosystem Services Related to Carbon Cycling - Modeling Present and Future Impacts in Boreal Forests.

Forests regulate climate, as carbon, water and nutrient fluxes are modified by physiological processes of vegetation and soil. Forests also provide renewable raw material, food, and recreational possibilities. Rapid climate warming projected for the boreal zone may change the provision of these ecosystem services. We demonstrate model based estimates of present and future ecosystem services related to carbon cycling of boreal forests. The services were derived from biophysical variables calculated by two dynamic models. Future changes in the biophysical variables were driven by climate change scenarios obtained as results of a sample of global climate models downscaled for Finland, assuming three future pathways of radiative forcing. We introduce continuous monitoring on phenology to be used in model parametrization through a webcam network with automated image processing features. In our analysis, climate change impacts on key boreal forest ecosystem services are both beneficial and detrimental. Our results indicate an increase in annual forest growth of about 60% and an increase in annual carbon sink of roughly 40% from the reference period (1981-2010) to the end of the century. The vegetation active period was projected to start about 3 weeks earlier and end ten days later by the end of the century compared to currently. We found a risk for increasing drought, and a decrease in the number of soil frost days. Our results show a considerable uncertainty in future provision of boreal forest ecosystem services.

Introduction

Ecosystem services (ES) are defined as the contributions that ecosystems make to human wellbeing (e.g., Costanza et al., 1997; Daily and Matson, 2008). The ES concept has become widely used and serves to emphasize the dependency of human society’s welfare on natural ecosystems (MEA, 2005; Costanza et al., 2017). A wide variety of ES are found in forests, which cover about 30% of global terrestrial area (FAO, 2018). Forests provide timber, and food, they conserve biodiversity, regulate water resources, and provide recreational opportunities (Shvidenko et al., 2005; Saastamoinen et al., 2014). Forests are essential factors of the global carbon (C) cycle, with an important role in regulating atmospheric concentrations of carbon dioxide (CO2) (Kurz et al., 2013). Boreal forests represent 29% of global forests (Shvidenko et al., 2005), comprising the circumpolar vegetation zone of high northern latitudes that covers one of the world’s largest biogeoclimatic areas (Brandt et al., 2013). The ES provided by boreal forests thus benefit human society both locally and on the global scale (Gauthier et al., 2015). Boreal forests account for about 20% of the global C sinks (Pan et al., 2011). Thereby boreal forests provide a climate regulating ES with bearing on climate change mitigation on the global level, although trade-offs are also recognized (Luyssaert et al., 2018).

Estimates of recent trends in above-ground biomass C indicate an increase for boreal forests (Liu et al., 2015). The stand and landscape level characteristics of the boreal forests of North America, Fennoscandia and Russia vary between the regions because of historical and current regional differences in natural disturbances and management practices (Gauthier et al., 2015; Kuuluvainen and Gauthier, 2018). Statistics from two neighboring boreal countries show increasing trends in roundwood increment in Finland (Vaahtera et al., 2018) and Sweden (Skogsdata, 2018). The future rate of C uptake by forests depends on how ambient temperature, land use and resource management practices evolve (Ahlström et al., 2012; Shao et al., 2013; Lappalainen et al., 2016). Rapid climate warming projected for the boreal zone (IPCC, 2013) has been observed in North America (McKenney et al., 2006; Price et al., 2013), Fennoscandia (Mikkonen et al., 2015; Hedwall and Brunet, 2016), and Russia (Schaphoff et al., 2016), which may have significant consequences for future boreal forest ES. Management practices may interact with climate change with consequences for forest biodiversity (Tremblay et al., 2018; Cadieux et al., 2019). Schaphoff et al. (2016) found that the impacts of climate change on Russia’s boreal forests are often superimposed by other environmental and societal changes, while Kalliokoski et al. (2018) estimated likely increases in future forest productivity in Finland, although with a high uncertainty. Climate change may affect the occurrence and extent of natural disturbances, such as insect outbreaks and fire (Schaphoff et al., 2016), e.g., Navarro et al. (2018) report a northward shift of spruce budworm (Choristoneura fumiferana) during the 20th century. Earlier thermal growing season has been found to be associated with earlier onset of biospheric C uptake, whereas earlier termination of biospheric activity was associated to later termination of thermal growing season (Barichivich et al., 2012). Warmer springs have consequences for bird reproduction (Ludwig et al., 2006; Wegge and Rolstad, 2017) and for moth multivoltinism, in combination with warmer summers (Pöyry et al., 2011). Less severe winter colds may promote the expansion of forest pests (Fält-Nardmann et al., 2018). Warmer winters may also have impacts on the possibilities for winter harvest on drained peatlands (Hökkä et al., 2016) and for winter recreation (Neuvonen et al., 2015).

The pathway from ecosystems, their biophysical structure, processes and functions to their benefit and value to human society, is illustrated by the ES cascade model, which portrays ES as emerging from the functional and structural properties of the ecosystem (Potschin and Haines-Young, 2011). In the Common International Classification of ES (CICES), services are classified into provisioning, regulating and maintenance, and cultural services (Haines-Young and Potschin, 2018). Several studies have discussed mapping of ES provision potential at different scales (Burkhard et al., 2009; Vihervaara et al., 2010; Maes et al., 2013; Albert et al., 2016a). The European Commission seeks to improve the basis for implementing the EU Biodiversity Strategy for 2020 (European Commission, 2011) by encouraging national ES assessments of the member states of EU through its flagship project MAES (European Commission, 2018). National accounting of natural capital needs information from national ES accounts and ongoing development of natural-capital accounts also supports national ES assessments (Schröter et al., 2016). Saastamoinen et al. (2014) applied the CICES hierarchy on the boreal forest ES in Finland, reporting their results in a conceptual and historical context. Mononen et al. (2016) developed a framework of ES indicators for Finland that complies with both national circumstances and international typologies such as the cascade model (Potschin and Haines-Young, 2011) and the CICES framework (Haines-Young and Potschin, 2018). ES are increasingly used to inform policies, and frameworks including ES analysis for decision support are becoming available (Albert et al., 2016b). Potentially highly promising avenues for underpinning resource allocation decisions are based on detailed place-based analyses of supply and demand of ES (e.g., Kopperoinen et al., 2014; Vauhkonen and Ruotsalainen, 2017). On the other hand, broad, unspecified ES may be more easily adopted by policy actors, according to Van Oudenhoven et al. (2018).

Our aim is to illustrate potential impacts of climate change on boreal forest ES. In this paper we (i) apply two dynamic ecosystem models (JSBACH and PREBAS) driven by a set of climate change scenarios to simulate present day and future values of a set of biophysical variables; (ii) compare simulated present day values with available observations and statistical data; (iii) suggest interpretations of the biophysical variables as ES; (iv) discuss the consequences of changes in these biophysical variables on future ES of boreal forests (Supplementary Figure 1).

Materials and Methods

Ecosystem Models

We applied the land ecosystem model JSBACH (Raddatz et al., 2007; Reick et al., 2013; Goll et al., 2015; Gao et al., 2016) and the stand growth model PREBAS (Valentine and Mäkelä, 2005; Peltoniemi et al., 2015b; Minunno et al., 2016, 2019) for the land area of mainland Finland, using input data in spatial resolution ranging from 16 m (forest resources, PREBAS), 0.1° (land surface characteristics, JSBACH) to a 0.2° × 0.1° longitude-latitude grid (climate data). We report simulated results for a set of biophysical variables for four time periods (1981–2010, 2011–2040, 2041–2070, 2071–2100). Two variables reflect directly ecosystem C fluxes: gross primary productivity (GPP; gCm-2 yr-1) and net ecosystem exchange (NEE; gCm-2 yr-1). Ecosystem phenology is described by three variables: the length of vegetation active period (VAPlength; days), when terrestrial vegetation is assimilating C through photosynthesis, the start (end) of vegetation active period (VAPstart; VAPend; days from January 1st, when the daily level of photosynthesis first exceeds, or returns to below 15% of its summertime value). Stemwood growth (m3 yr-1) was simulated with PREBAS, and the number of summer dry days (soil moisture < 5th percentile of reference period) and number of winter days with soil frost were simulated with JSBACH. The models have been compared earlier and found to produce similar GPP estimates at local and national level (Peltoniemi et al., 2015a). Here we use the two models in parallel to provide further information about the changes in biophysical variables under climate change, and also in supplementing each other for variables predicted by only one of the models. We present the results as boxplots and tabulated percentile values, and cumulative distribution functions aggregated over the whole area, as well as time series and maps for some variables.

JSBACH is a land surface model of an earth system model of Max Planck institute for meteorology (MPI-MET) and describes the biogeophysical processes that regulate the balances of water and CO2. The storage of water and C into the ecosystem as well as their release to the atmosphere is regulated by the climatic variables. Land vegetation is divided into plant functional types (PFT), and for our domain (the Finnish mainland) we based the PFT distribution on Finnish CORINE land cover data (CLC, 2012), which contains information about soil type, thus providing soil characteristics consistent with the vegetation cover data. Seasonal development of leaf area index is regulated by air temperature and soil moisture with PFT specific maximum leaf area index as a limiting value. For the generation of the seasonal cycle, PFTs are divided in summergreen, evergreen, grass and crop phenology types (Böttcher et al., 2016). Photosynthesis is described according to Farquhar et al. (1980), using PFT-specific parameters for the maximum carboxylation (Vmax) and electron transport (Jmax) rates. Global parameter values were used in this study, and as JSBACH was run without explicit nitrogen cycle, the mean values at 25°C (Vmax25) were applied, with Jmax25 = 1.9 ⋅Vmax25 for all PFTs (Wullschleger, 1993; Kattge et al., 2009). The photosynthetic rate is resolved first under non-water-stressed conditions to attain photosynthetic activity (Schulze et al., 1994), and limitations in water availability are accounted for Knorr and Heimann (1995) and Knorr (2000). The radiation absorption within the vegetation canopy is calculated for three layers (Sellers, 1985), accounting for clumping of the leaves in sparse canopies (Knorr, 2000). In addition to the canopy processes the model consists of a 5-layer soil moisture description (Hagemann and Stacke, 2015) and the Yasso soil C module (Goll et al., 2015). We adopted a formulation that delays the beginning of photosynthetic activity of evergreen species in spring (Kolari et al., 2007), as according to Böttcher et al. (2016) the start date of the photosynthetically active season of coniferous evergreens in the model is ahead of the observed. We decreased the threshold of the temperature sum regulating the bud-break from 4 to 2°C in accordance to findings by Böttcher et al. (2016). Furthermore, a condition that reduces stomatal conductance under supersaturation was removed because it falsely prohibits photosynthesis under conditions of very high humidity.

In addition, simulations were performed with PREBAS, a C-balance based stand growth and gas exchange model (Valentine and Mäkelä, 2005; Peltoniemi et al., 2015b; Minunno et al., 2019). In PREBAS, photosynthesis (GPP) and evapotranspiration are calculated using a light-use-efficiency approach linked to soil moisture and driven by daily environmental inputs of radiation, temperature, vapor pressure deficit, precipitation and ambient CO2 concentration (Peltoniemi et al., 2015b; Minunno et al., 2016; Kalliokoski et al., 2018). GPP is allocated to mean-tree growth and respiration at an annual time step, and allocation of growth to different tree components is based on conservation of structural constraints, e.g., the pipe model (Valentine and Mäkelä, 2005). Tree mortality due to crowding is included in the model. The growth module updates canopy leaf area index which is used as input to the gas exchange module in the following year. To calculate NEE, PREBAS has been linked with the soil C model YASSO07 through annual litter inputs (Liski and Westman, 1995; Tuomi et al., 2009). In addition to weather data, PREBAS requires information on the initial state of the simulated forest, and forest management actions, including thinning, clearcut and regeneration, to obtain a realistic dynamic pattern of forest development. The relatively simple structure of the model allows us to apply it at a large regional scale and make simulations of C and water fluxes that are both climate and management sensitive. Here, we defined forest management on the basis of current management recommendations in Finland (Sved and Koistinen, 2015; Minunno et al., 2019). The simulations were carried out by linking PREBAS to information on soil fertility, stand basal area, mean height and mean diameter, derived from 16 m resolution multisource forest inventory data maps (Tomppo et al., 2014; Mäkisara et al., 2016; LUKE, 2018a,b). The forest resource maps were divided into 16 km grid cells. Within each grid cell, simulations were carried out for 50 forest classes with the highest proportional areas, accounting also for the remaining classes by distributing their area to these dominant classes. The forest classes represented different combinations of dominant tree species, age, basal area and mean height. Categories of soil fertility were herb-rich heath, mesic heath, sub-xeric heath, xeric heath, and poorer types (Cajander, 1949). To account for mineral and peat soil types, we used two sets of soil parameters, one describing typical mineral soil, and one for a soil layer with high soil water retention capacity (depth = 430 and 1000 mm, respectively), as PREBAS does not simulate organic soils as such. The simulated annual values for each class were aggregated to grid cell level based on 16 m resolution forest maps. Climate data from the nearest climate grid point were used for each 16 km forest simulation grid cell.

We compared the simulated annual NEE for the period 1981–2010 with reported values of the C sink of Finnish forests (Official Statistics of Finland, 2019). The simulated stemwood growth values for the period 2009–2012 were compared to observed forest growth for the same period. The simulated VAPstart days were compared to estimates based on webcam images (Arslan et al., 2017; Peltoniemi et al., 2018a,b), as well as to those estimated from satellite observations (SYKE Vegetation Phenology 2001–2017).

Climate Change Scenarios

The climate change scenario data through years from 1981 to 2100 were obtained as the output from five global climate models (GCMs; CanESM2, CNRM-CM5, GFDL-CM3, HadGEM2-ES and MIROC5) from the fifth phase of the Coupled Model Intercomparison Project (CMIP5; Meehl et al., 2009; Taylor et al., 2012). We used output for climate models driven by three representative concentration pathways (RCPs), that lead to radiative forcing levels of 2.6, 4.5 and 8.5 W m-2 by the end of the century (Moss et al., 2010; Van Vuuren et al., 2011). These pathways include one mitigation scenario (RCP2.6), one stabilization scenario (RCP4.5) and one high emission scenario (RCP8.5) (Van Vuuren et al., 2011). The resulting climate variables were down-scaled to a 0.2° × 0.1° longitude-latitude grid by a quantile-quantile type bias correction algorithm for daily mean temperature, relative humidity, shortwave radiation and wind speed (Räisänen and Räty, 2013) and parametric quantile mapping for daily precipitation (Räty et al., 2014). Gridded harmonized FMI meteorological data by Aalto et al. (2013) were used. Long-wave radiation data were interpolated by a bilinear interpolation method. Global mean CO2 concentrations from the RCPs 2.6, 4.5 and 8.5 were linearly interpolated to monotonously increase through the calendar years. The changes of temperature and precipitation from a baseline period 1981–2010 to periods 2011–2040, 2041–2070, and 2071–2099 in Finland varied considerably between data from five down-scaled CMIP5 climate models (Figure 1). PREBAS was run with all three RCPs, while JSBACH input was taken from RCP4.5 and RCP8.5. The RCP8.5 output of HadGEM2-ES was not included in JSBACH runs.

FIGURE 1

Annual and seasonal projected change in precipitation sum (%) versus change in mean temperature (°C) in Finland. Changes calculated from the reference period 1981–2010 to periods 2011–2040 (bottom), 2041–2070 (middle) and 2071–2100 (top) from regionally down-scaled data from simulations with three levels of climate forcing (RCP) driving five down-scaled CMIP5 global climate models (GCM). Seasons ANN annual (left), DFJ winter, MAM spring, JJA summer, and SON autumn (right).

Ecosystem Services

On the basis of the analysis of Saastamoinen et al. (2014) and Mononen et al. (2016), we identified a set of key boreal forest ES. In the provisioning services section of the CICES classification, we selected productivity and supply of harvestable wood as key benefits to human society (Table 1). Forestry land covers 86% of the total land area of Finland, and forest industry products accounted for 20% of the total value of Finnish goods exports in 2017 (Vaahtera et al., 2018). Gross primary production, or the rate at which terrestrial vegetation assimilates CO2 in photosynthesis, (GPP gC m-2 yr-1) was simulated by JSBACH and PREBAS, and forest stemwood growth (m3 yr-1) by PREBAS.

Table 1

Ecosystem services with corresponding biophysical proxies.

Section	Class	Function	Benefit (B), or Loss (L)	Modeled results of biophysical variables
Provisioning services	Wood	Growing stock increment	Productivity (B), Supply of harvestable wood (B)	GPP (gCm-2 yr-1), Stemwood growth (m3 ha-1 yr-1)
Regulation and maintenance services	Climate regulation	Carbon sequestration	CO2 forcing (L)	NEE (gCm-2 yr-1)
	Transpiration	Evapotranspiration	Summer drought (L)	Number of dry summer days
	Soil structure	Bearing capacity	Opportunities for winter harvest (B)	Number of soil frost days
	Maintaining habitats	Reproduction	Viable populations (B/L)	VAPstart, VAPend, VAPlength (days)
Cultural services	Nature tourism	Natural events and phenology	Opportunities for recreational experiences (B/L)	VAPstart, VAPend, VAPlength (days)

For the regulation and maintenance section of the CICES classification, the key service studied was related to regulating climate through C sequestered by growing vegetation. The corresponding benefit to society was avoided CO2 forcing (Table 1). The proxy for CO2 forcing, NEE (gC m-2 yr-1), is the difference between the total respiration rate of the ecosystem, and GPP. When the flux is from the atmosphere to the ecosystem, NEE is negative, and when the flux is from the ecosystem to the atmosphere, NEE is positive. In 2017, the Finnish total greenhouse gas emissions were 55.5 million metric tons CO2 eq, including energy, industrial, agriculture, and waste sectors. In contrast, the net removal from the atmosphere by forest land was 27.0 million metric tons CO2 eq. Altogether, the net sink of the land use, land use change, and forestry sector (LULUCF) was 20.4 million metric tons CO2 eq (Official Statistics of Finland, 2019). Climate change impacts on future values of NEE were simulated by JSBACH and PREBAS.

Risk of productivity losses related to summer drought, and decreasing opportunities for winter harvest on drained peatlands were considered as examples of regulating services. JSBACH simulations of the number of dry summer days, and soil frost days, were used for the climate impact projections. We used the length and the timing of the vegetation active period as proxies for both a regulating service (the maintenance of habitats), and a cultural service (recreational opportunities). The impact of climate change on the length (VAPlength, number of days) and timing of the vegetation active period (VAPstart, VAPend, day of year) was simulated with JSBACH and PREBAS.

Current potential ES provision was estimated as the median values of the simulated biophysical variables for the reference period (1981–2010). The impact of climate change on ES provision was interpreted as future changes in simulated biophysical variables compared to reference period median values. Positive changes in GPP, stemwood growth and number of soil frost days are in our analysis associated with benefits for human society. Positive values of NEE represent a flux from the ecosystem to the atmosphere (CO2 emissions), while negative values of NEE (CO2 removals) represent avoided CO2 forcing. Thus, in our analysis, positive changes in the simulated values of NEE are considered losses, and negative changes benefits to human society. Increasing frequencies of dry summer days are similarly considered losses. An increase in the length of the vegetation active period, however, may lead to both benefits and losses.

Results

Provisioning Services: GPP and Supply of Harvestable Wood

GPP was projected to increase with climate change, more so with higher radiative forcing (RCP8.5), and more toward the end of the century (Figure 2). For the low emission scenario (RCP2.6) hardly any difference between the time periods can be seen in the simulated GPP values. Differences between the medium and the high emission scenario simulations (RCP4.5 and RCP8.5) are evident only from the mid-century onward (Figure 2 and Supplementary Figure 2). The estimated median increase in GPP was 34% for the mid-century time period (2041–2070), and 46% for the end of the century (2071–2100), over all radiative forcing levels (Supplementary Table 1). Current observed mean annual growth in 15 forestry regions of Finland compared well (R2 = 0.95) with PREBAS simulated results for the same period (2.0–7.6 m3 ha-1 yr-1; Supplementary Figure 3). In the climate change projections for mean annual forest growth, all scenarios give similar results for the first time period (2011–2040, Figure 2 and Supplementary Figure 4), whereas the differences become clear later, especially toward the end of the century. The highest increase in mean annual growth was projected for the north of Finland (Figure 3). The simulated country-wide median forest growth increased from 5.6 m3 ha-1 yr-1 in 1981–2010 by 2.7 to 8.3 m3 ha-1 yr-1 in 2041–2070, and by 3.4 to 8.7 m3 ha-1 yr-1 in 2071–2100 (Supplementary Table 1), corresponding to an increase of 48% by mid-century, and 60% by the end of the century.

FIGURE 2

Summary of simulated values for (A) GPP (gC m-2 yr-1); (B) NEE (gC m-2 yr-1); (C) forest growth (m3 ha-1 yr-1); (D) VAPlength (days); (E) VAPstart (days from January 1st); (F) VAPend (days from January 1st) for three time periods (2011–2040; 2041–2070; 2071–2100). The results are obtained from simulations with PREBAS and JSBACH using future climate data from five climate models with driving scenarios RCP2.6, RPC4.5 and RCP8.5 (PREBAS), and RCP4.5 and RCP8.5 (JSBACH). Forest growth based on PREBAS results only.

FIGURE 3

Mean modeled change in annual growth (m3 ha-1 yr-1) from 1981–2010 to 2041–2070, calculated over all climate change scenarios. Forestry regions and their borders indicated on the map.

Regulating Services

Avoided CO2 Forcing

Interpreted as an ES, a more negative value of NEE means that the ecosystems of the region have a potential to mitigate the radiative forcing by atmospheric CO2, which is the most important contributor of human induced climate change. Across Finland, the simulated median annual NEE was –50 gC m-2 yr-1 for the period 1981–2010 (Supplementary Table 1). Reported current national total annual C sink of forest biomass and soil varied in the range –35 to –18 Mt CO2 eq yr-1 in the period 1990 to 2017 (Official Statistics of Finland, 2019). Averaged over total forestry land (262,000 km2), the reported range corresponds to annual removals of –36 to –19 gC m-2 yr-1. Our simulated reference period annual C sink thus clearly exceeds the range of reported current values (Official Statistics of Finland, 2019).

Simulations with the low emission scenario (RCP2.6), which include only PREBAS results, represent clearly larger C sink values (more negative NEE) than simulations with the other scenarios, for all time periods. The high emission scenario (RCP8.5) shows the largest range in simulated NEE results (Figure 2). In projections across Finland from both ecosystem models and all climate scenarios, the median removal of CO2 increased with 48 and 39% by mid-century and the end of the century, respectively. The corresponding simulated median values were NEE –74.1 and –66 gC m-2 yr-1 for the period 2041–2070 and 2071–2100 (Supplementary Table 1).

Increasing Risk for Drought and Decreasing Opportunities for Winter Harvest

Warming temperatures, changing precipitation patterns (Figure 1) and increasing forest growth (Figure 2, 3) altered the water availability in the JSBACH simulations. The risk for drought was slightly increasing for all scenarios. Much of the change is due to earlier snow melt and thus possibility for earlier drought events. Decreasing precipitation was projected for some climate models in winter for all time periods, and also in summer for the end of the century (Figure 1). The mean number of dry summer days was approximated to around 4 days for the reference period. By the end of the century (2071–2100), the mean simulated number of dry summer days increased to about 15 days in the south and to 10 days in the north (Supplementary Table 2). In southern Finland, simulations for the mid-century gave even higher number of dry days (23). In the JSBACH simulations, the number of soil frost days was decreasing in all simulations (Supplementary Table 2). The changes were largest in southern Finland under the high forcing scenario (RCP8.5), from 134 days in 1981–2010 to 25 days in 2071–2100. This means the opportunities for winter harvest on drained peatlands are decreasing with all scenarios, more in the south than in the north.

Regulating and Cultural Services: Maintaining Habitats, Nature Tourism

Estimates of the start dates of the vegetation active period on the basis of webcam images were similar to those simulated by the JSBACH model (see section “Webcam Network” in Supplementary Material and Supplementary Figure 6). The vegetation active period ended slightly later according to the JSBACH results than based on the webcam images. VAPstart in Finland occurred in median on 28 April in all simulations for the present period (1980–2010) (Supplementary Table 1). In comparison, for the period 2001–2017, VAPstart in Finland was determined from satellite observations (SYKE, 2018). The satellite-derived median VAPstart occurred on 13 April in evergreen coniferous forest. JSBACH simulations included all PFTs whereas PREBAS simulations included forested areas only. Thus, the simulated present day VAPstart seems to be slightly ahead of satellite observations with a median VAPstart on 5 May for all forest types (Figure 4). In comparison to reference conditions, the vegetation active period started (VAPstart) about 2 weeks earlier in the mid-century, and 3 weeks earlier toward the end of the century. This means spring would start earlier in April than at present. The vegetation active period ended (VAPend) by the end of the century in early October, compared to late September in the reference conditions (Figure 2 and Supplementary Table 1). In the climate change simulations, the length of the vegetation active period (VAPlength) increased from 162 days with about 20 days in the mid-century and with 30 days by the end of the century (Figure 2 and Supplementary Table 1).

FIGURE 4

Start of vegetation active period (A) from satellite observations 2001–2017, calculated as weighted average for all forest types based on their proportions in the satellite grid – days from January 1st expressed as dates; (B) as mean modeled results from JSBACH and PREBAS 1981–2010 – days from January 1st expressed as dates; and (C) mean modeled number of days difference in 2040–2070 compared to 1981–2010.

Discussion

We aggregated the results from the two ecosystem models in order to account for the uncertainty inherent in the use of different assumptions, e.g., JSBACH simulations do not account for forest management. In this paper we focus on presenting the results on the whole-country level. Data on grid-based results of the biophysical variables for each ecosystem model and each climate scenario separately are available online (Metadata, 2018a,b). These online metadata descriptions include links to the Climateguide.fi service, where maps of the grid-based results can be found by selecting the category “Terrestrial ecosystems,” and the variable of interest “GPP,” “NEE,” etc. The maps are displaying results from each ecosystem model (JSBACH, PREBAS), climate model, and climate forcing level RCP separately. Seasonal patterns of GPP simulated by PREBAS and JSBACH aggregated to the national level were almost identical, and temporal trends of annual totals were also parallel. This is remarkable in light of earlier results, as ecosystem models tend to vary quite much in their GPP estimates globally and in the northern latitudes (Anav et al., 2013; Shao et al., 2013). However, it has to be noted that the current version of the JSBACH photosynthesis model has been modified for boreal conditions using partly the same empirical evidence as PREBAS (e.g., Mäkelä et al., 2004; Kolari et al., 2007).

The climate change simulations of forest growth did not account for growth restrictions due to limited access to nutrients, mainly nitrogen or phosphorus. Forest growth in higher CO2 concentrations may become nutrient limited in nutrient deficient sites (e.g., Wieder et al., 2015). As PREBAS projections of future growth were made assuming current forest practices continue in the future (Sved and Koistinen, 2015; Minunno et al., 2019), here no provision was made for potential adaptive responses by the forest owners, e.g., aiming for more variation in stand structure, stand age or silvicultural management (Blennow, 2012; Montoro Girona et al., 2017), introducing new tree species or varieties (Laakkonen et al., 2018; Sousa-Silva et al., 2018), or striving to maximize positive synergies between ecosystem services (Kuuluvainen and Gauthier, 2018). Current total volume of growing tree stock in Finland is 2.5 billion m3 with a total annual increment of 107 million m3, which corresponds to annual average 6.8 and 3.2 m3 ha-1 in southern and northern Finland, respectively. Total volume of roundwood harvested from forests was a record high 72.4 million m3 in 2017, about 17% higher than the average of the preceding ten-year period (Vaahtera et al., 2018). Current forest policies in Finland aim to increase the annual cut to over 80 million m3 in the next decades. Our simulations indicate that these plans could be feasible in a changing climate, however, increased cutting levels would also imply a lower growing stock than that in our simulations which assumed continuing the current harvest intensity. This would have further implications on C sequestration that need to be explored separately. We did not study the impact of climate change and cutting regimes on biodiversity indicators, such as the amount of dead wood, but others have found significant effects of bioenergy extraction on dead wood levels (e.g., Hof et al., 2018) and conflicting objectives of wood production and biodiversity conservation (Mönkkönen et al., 2014; Angelstam et al., 2019). In their comprehensive review on the retention approach for forestry, Lindenmayer et al. (2012) argued that the practice of permanently retaining significant elements of the original forest is crucial for maintaining multiple forest values, such as biodiversity and carbon stocks. Recently, partial harvest methods were found adequate for regeneration for sustainable forest management in Canadian boreal forests (Montoro Girona et al., 2017, 2018).

Regarding NEE, JSBACH and PREBAS models produced fairly different patterns especially under climate change – PREBAS predicted more negative NEE values than JSBACH (MONIMET, 2017; Peltoniemi et al., unpublished results). We assume that the differences are caused by different approaches to describing tree growth and management in these models, JSBACH assuming steady-state conditions and PREBAS accounting for forest growth and management explicitly. The harvests lead to a drain of stemwood from the forest which reduces the amount of C entering the soil as coarse woody debris, and therefore reduces the rate of heterotroph respiration from the soil. In our presentation, the NEE projections aggregate results from these models with different assumptions concerning vegetation and forest management, and the variance of the results reflect uncertainties in process descriptions and parameterizations. Forests remain a C sink in our simulation results.

It has been shown that future development of NEE varies much according to ecosystem model and applied climate forcing. According to Ahlström et al. (2012), simulations with a dynamic global vegetation model (LPJ-GUESS), forced by 18 CMIP5 climate projections, showed net release by 2100 for the boreal region. However, the projections were not in full agreement with each other. The magnitude of the release varied, and even turned to uptake in some simulations. Further, according to an ensemble of CMIP5 models both uptake and release of C by land ecosystems is accelerated in the twenty-first century (Shao et al., 2013). The boreal latitudes could become a major C sink by 2100 in many models. It should also be recognized that there is a trade-off situation between the ES supply of “harvestable wood” and “avoided radiative forcing.” According to calculations by the Finnish Climate Change Panel, meeting the Finnish targets for implementing the Paris climate change agreement would, in addition to stringent GHG emission reductions, also require a large increase in the C sinks (Finnish Climate Change Panel, 2018). Similar measures would be required also at the global level (Rockström et al., 2017; Rogelj et al., 2018). Intensification of biomass removals from forests may invoke harmful impacts on forest productivity, biodiversity, soil quality, and climate change mitigation potential (Aherne et al., 2012; Forsius et al., 2016; Mäkelä et al., 2016; Soimakallio et al., 2016; Vanhala et al., 2016; Eyvindson et al., 2018).

In boreal and northern conditions, dry and warm summers have been found to increase drought-related symptoms in trees, such as defoliation or discoloration (Muukkonen et al., 2015), limit tree radial growth (Aakala and Kuuluvainen, 2010) and contribute to severe disturbances (Aakala et al., 2011). Satellite-derived estimates of vegetation growth trends indicate a clear negative response to recent dry summer conditions (Piao et al., 2011). Shorter periods of frozen soil increase wind damages to forest, as the trees are less firmly anchored in non-frozen soil. Wind-thrown trees may also serve as starting points for bark beetle outbreaks (Pukkala, 2018). In boreal forests, timber is mainly harvested with heavy machinery, which can operate only if the soil has a sufficient bearing capacity (Suvinen, 2006). Higher and more variable air temperatures in winter, coupled with longer and more frequent rain events, may cause pronounced edge effects by skidding trails (Montoro Girona et al., 2016), or lead to periods during which the sites on drained peatlands are not accessible (Kokkila, 2013).

For regulating ES, earlier start of vegetation active period means earlier opportunities for birds, insects and other species, but is also linked to the risk for coincident occurrences of cold spells that may be detrimental. For cultural ES, climate warming is expected to lead to earlier opportunities for nature tourism linked to the coming of spring, such as bird watching and observing shoot growth, bud bursting and flowering. However, potential detrimental long-term impacts of climate change on species diversity and conservation are also well recognized (e.g., Thomas et al., 2013; Virkkala et al., 2013; Heikkinen et al., 2015; Virkkala, 2016; Cadieux et al., 2019). For regulating ES, later end of vegetation activity may mean decreased leaching of N in autumn. For provisioning and cultural ES, later end of vegetation may be associated with improved opportunities for the growth and picking of mushrooms and berries, although the success of these species is mostly regulated by local weather conditions.

In terms of regulating ES, this means that there is a positive impact on the reproduction and survival of birds, insects and other species that are dependent on forest vegetation activity. There might be negative impacts due to mismatch of animal production and food supply, or better match between insects and host plants and resulting damages that could alter ecosystem composition (Foster et al., 2013; Pureswaran et al., 2015). The spruce bark beetle (Ips typographus), a serious pest in boreal forests may benefit from warmer temperatures (Annila, 1969) and earlier occurrence of spring (Lange et al., 2006; Zang et al., 2015), increasing the risks for declining forest health (Blomqvist et al., 2018). With the lengthening of the vegetation active period, there is also the option for spatial migration of plant and animals to more northern areas (e.g., Jeganathan et al., 2014).

Summertime nature tourism in the Nordic countries, e.g., for fishing, biking, hiking, kayaking, bird and animal watching, was perceived as interesting activities by respondents in an international survey (MEK, 2010). According to the Finnish national outdoor recreation inventory, almost 70% of annual nature visits were made during the 5 months from May to September, the most popular outdoor recreation activities being walking, swimming, spending time in nature and at a recreational home, picking wild berries, cycling, boating, fishing, cross-country skiing and mushroom picking (Neuvonen and Sievänen, 2011). According to Lankia et al. (2015), the overall value of recreational nature visits was considerable in comparison to other land uses. It is likely that an increase in VAPlength would be beneficial to human society. It should, however, be noted that as the increase in VAPlength may be accompanied by increased probability of spring and autumn rains (Ruosteenoja et al., 2016), the recreational value may not be realized. Here, we do not report simulated future changes in number of snow cover days in response to climate change. It is clear, however, that the occurrence of earlier start of spring and later beginning of autumn season are coupled to shorter winter periods, leading to losses as regards the opportunities for winter tourism. Cross-country skiing is an important recreational activity in Finland (Neuvonen and Sievänen, 2011) and Neuvonen et al. (2015) presented an interactive tool to study the implications of climate change for cross-country skiing. A considerable proportion of recreational nature visits are made to nature conservation areas, and in the case these visits constitute an increased pressure on the nature conservation areas, the goals of conservation and recreation may be conflicting. It is more likely, however, that the goals of forestry and conservation are colliding. The recreational experience is not independent of the maintenance of habitats, and we are aware of the risk of double attribution in estimating quantitative values for these services, e.g., in monetary terms. It was, however, beyond the scope of our work to assess quantitative values for the services we have studied. In the case of the qualitative analysis carried out here, we think we are justified to use the vegetation active period characteristics to illustrate both the potential for habitat maintenance and recreation opportunities.

For the near decades, even in mid-century, the three emission levels (RCP) gave broadly similar projected temperature increase and change in precipitation on the annual and country-wide level (Figure 1). Different climate models, however, gave varying responses especially for seasonal changes. Extreme warming (>5°C) was projected by some climate models by the end of the century for the high emission scenario. Ruosteenoja et al. (2016) call the late century response to RCP8.5 “an alarm signal, demonstrating the furious climatic changes that would be expected if the mitigation of greenhouse gas emissions were totally neglected.” Although spatially explicit data have been used to drive the ecosystem models, the results we report here are mainly estimates aggregated over the whole country. We realize that important detailed information is lost in such an aggregation, especially concerning tradeoffs between ES, and the spatial distribution and balance of ES supply and demand (Potschin and Haines-Young, 2013; Vauhkonen and Ruotsalainen, 2017). The value of boreal forest future ES has also not been addressed in this paper. We have qualitatively described future ES on the basis of simulated biophysical quantities such as forest growth and CO2 exchange with the atmosphere. Analyzing current and future valuation of ES by forest owners and other stakeholders would be needed in order to introduce boreal forest ES into a natural capital framework (Schröter et al., 2016; Lai et al., 2018), or to assess forest carbon policies (Pohjola et al., 2018). Such considerations are, however, beyond the scope of our work.

Conclusion

In our analysis, climate change impacts on key boreal forest ES are both beneficial and detrimental. Increased GPP is beneficial, leading to increased C sequestration, as well as increased forest growth and improved supply of harvestable wood. Warmer temperatures lead, however, also to higher ecosystem respiration, which in some instances cause increasing fluxes of CO2 from vegetation to atmosphere (NEE > 0). Although GPP is steadily increasing in our simulations, there are times in early and late 21st century when GPP is smaller than ecosystem respiration and annual NEE becomes positive. Here, the assumptions concerning future conditions included only different climate change scenarios, not considering any adaptation measures, e.g., regarding changing forestry practices. Thus, future forest growth simulations with PREBAS reflect a continuing business as usual response of the forestry sector. This means the current forestry practices are assumed to continue, and possible shifts toward higher harvest rates to accommodate increased demand for bioenergy have not been considered here. Possible increased risk of forest damage was also not accounted for, although increasing probabilities of more frequent insect outbreaks and wind throw damages have been recognized. There is still a large uncertainty in predicting future forest growth, since the models disregard natural disturbances, and are not constrained by potential nutrient limitations (N, P and base cations). These limitations call for further model developments as well as biomass and flux data to inform process based modeling. Our current simulation results likely represent the high end of the future growth estimates. We also noted the need to consider a balanced approach between forest harvesting, C sequestration potential and nature based ecosystem services in forested ecosystems.

While some of the studied ES are expected to improve with climate change (biomass growth, nature tourism), the projected decrease in number of soil frost days is expected to decrease the opportunities for winter harvest in forestry especially in the south. Future levels of net ecosystem exchange of CO2 are uncertain; some simulations indicate decreasing radiative forcing levels, while other simulations lead to increasing levels. Forests remain C sinks, however, in all simulations. The risk for summer drought is slightly increasing in the whole country. Climate warming is expected to lead to earlier opportunities for spring-time nature tourism, and may improve opportunities for autumn mushroom and berry picking. The opportunities for winter sports are decreasing while hiking opportunities are increasing with earlier occurrences of spring.

Our results indicate the range of uncertainty in future provision of boreal forest ES. The estimates described in this paper may provide input to comprehensive systems for integrated natural capital accounting, as well as sharing such information via hubs such as National Clearing House Mechanism for the Convention on Biological Diversity (Finnish Ecosystem Service Indicators, 2015). Furthermore, the examples in this paper may highlight the potential for advancing ES monitoring, such as improving information on ecosystem variables for instance by high-resolution Earth Observation or novel in situ observation networks such as webcams.

Data Availability

The datasets generated for this study are available on request to the corresponding author.

Author Contributions

MH planned and drafted the manuscript. MH, TA, AnA, AlA, TM, AM, FM, MP, KB, and KR planned the framework. TA and TM provided and interpreted JSBACH model results. AM, MP, and FM provided and interpreted PREBAS model results. IL provided the web interface of the model results. KB contributed revised material on phenology and performed GIS analysis. AnA, IB, KR, PV, and MF contributed and revised material on ecosystem services.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.