1 Introduction

Climate zones characterize regionally similar features of the complex climate system and highlight links between climate and terrestrial systems (Koppen 1936; Thornthwaite 1943; Brovkin et al. 1997; Kottek et al. 2006; Beck et al. 2018). As a result, shifts in climate zones imply substantial climate change and its profound impacts on terrestrial environments, for example, exacerbating aridification, desertification, and ecological degradation (Jylha et al. 2010; Lobell et al. 2011; Huang et al. 2016). Thus, climate zones can efficiently serve as an integrated but simplified climate indicator for understanding climate change and assessing its impacts on our living environments.

The climate zone concept has been widely adopted in assessments of climate change, and its impacts on hydrology, agriculture, and socioeconomy. For example, based on the Koppen climate classification (Kottek et al. 2006), approximately 5.7% of the global land area shifted to warmer and drier types during 1950–2010, and especially, the arid and high-latitude continental climate zones expanded (Chan and Wu 2015). Moreover, across climate zones, the drying and wetting trends and extreme precipitation vary quite heterogeneously (Chen et al. 2017; Ragulina and Reitan 2017). The climate classification scheme is also used to pinpoint biases in model simulations (Gallardo et al. 2013; Phillips and Bonfils 2015). In addition, climate zones, according to various indices, also play an essential role in impact assessments of climate change on hydrology, agriculture/ecology, and socioeconomy (Goodrich and Ellis 2006). For instance, Haines et al. (1988) developed a global river regime classification based on streamflow, assisting with an accurate description of global river features. Hydrologic classification also underpins water resource assessments and planning for the impacts of climate change (Kennard et al. 2010). In agriculture, using temperature and soil moisture, global land was classified into 100 climate zones to evaluate the potential yield of crops around the world (Licker et al. 2010). Moreover, climate zones have an influence on the patterns of global energy transport, trade flows, etc., and they are increasingly used in impact assessments of climate change on the socioeconomy (Landau 1983; Guan et al. 2018).

Because of their indication of climate change and its impacts, geographic shifts in climate zones further gain attention on various scales. In terms of global scale, the large shifts between climate zones were found to lead to expanded arid regions and smaller polar regions (Rubel and Kottek 2010). On regional scales, the variation in the Koppen climate boundary was examined, and no evidence was found to support a northerly migration of the climate boundary in the central United States (Suckling and Mitchell 2000). However, in China, decadal variations in the arid-semiarid boundary show a moving trend toward drier climate zones during the second half of the 20th century (Ma et al. 2005; Kim et al. 2008). A similar expansion of the semiarid zone is also detected in soil moisture variability (Li and Ma 2013). In Europe, climate shifts have strongly affected most regions in Central Europe and Fennoscandia (Gallardo et al. 2013). Furthermore, climate zone shifts leading to dry climate expansion and polar climate contraction are projected to be intensified in the future (Rohli et al. 2015; Chan et al. 2016; Zhang et al. 2017). In particular, drylands will acceleratingly expand under future climate scenarios (Huang et al. 2016).

Previous studies on climate zone shifts were mostly based on atmospheric-only variables (Chan et al. 2016; Huang et al. 2017; Zhang et al. 2017; Defrance et al. 2020); however, from the aspect of soil moisture, global climate zone shifts remain unclear. Soil moisture, on various time scales, directly connects to both atmospheric climate anomalies and almost all aspects of land processes and uniquely modulates the soil-ecosystem-climate relationships (Koster et al. 2004; Seneviratne et al. 2006). For instance, soil moisture anomalies have a substantial effect on summer precipitation known as “hot spots” of land-atmosphere coupling (Koster et al. 2004). By changing albedo, aerodynamic roughness, and evapotranspiration, it also controls surface energy flux partitioning, water balancing, and ecosystem functioning (Shrestha et al. 2018; Schoener and Stone 2019; Liu et al. 2020). Furthermore, for its “memory” of atmospheric anomalies, soil moisture can depict stable dry-wet conditions on a monthly time scale, that is, filter out high-frequency fluctuations (Entin et al. 2000). On such a time scale, vegetation-soil moisture feedbacks are important to the composition and structure of water-limited ecosystems (D’Odorico et al. 2007). As a result, based on soil moisture, the assessment of shifts in climate zones can efficiently capture climate change and its potential impacts on terrestrial ecosystems and environments.

Given these important roles of soil moisture variability in the climate system and ecosystem, increasing efforts have been made in monitoring and modeling to increase our knowledge of soil moisture-climate interactions (Seneviratne et al. 2010). To date, observationally constrained climate model projections can largely be explored for potential future changes in land hydrology and ecosystems (Flato et al. 2013). These efforts gradually make it feasible to investigate climate zone shifts in terms of soil moisture classification (Seneviratne et al. 2013; Li et al. 2020). In the present study, the focus is on potential shifts in climate zones under the CMIP5 RCP8.5 scenario using root-layer soil moisture as an integral indicator. The goal is to provide insights into the potential impacts of soil moisture changes on terrestrial ecological environments under future climate change scenarios. By exploring a perturbed parameter ensemble of climate projections with the latest Hadley Centre model HadGEM3-GC3.05, the ranges of uncertainties are also estimated. Models are first evaluated using observations in terms of land water cycles. Section 2 describes the data and methods in detail, and Sect. 3 describes the main results, followed by a discussion. The entire paper is summarized with conclusions in Sect. 5.

2 Model, data, and methods

2.1 Model and experiments

The main body of this study is based on the perturbed parameter ensemble (PPE) simulations (Murphy et al. 2004) using the coupled Hadley Centre Global Environment Model version 3.05 [a closely related configuration, HadGEM3-GC3.05, is described in Walters et al. (2017a)]. The land surface component of this coupled model is the Joint UK Land Environment Simulator (JULES) (Best et al. 2011; Clark et al. 2011), in which soil water and heat fluxes are represented with a 4-layer scheme, with soil layer thicknesses of 0.1, 0.25, 0.65, and 2.0 m, including the coupling processes between soil water and heat, such as the latent heat exchange caused by the phase change in soil water. The processes of subgrid soil moisture heterogeneity and an interactive water table are represented based on the rainfall-runoff model TOPMODEL (Beven and Kirkby 1979), hence increasing runoff heterogeneity driven by topography. With the river routing model TRIP (Total Runoff Integrating Pathways) (Oki and Sud 1998), the total runoff is routed from grid cells to the sea or inland basins. The water in inland basins can flow back into the soil moisture and the remaining portions flow into seas, which can then contribute to a component of the thermohaline circulation. Based on this coupled model, the PPE experiments were dedicated to achieving a broad range of plausible climate projections and sample uncertainty to assess climate change impacts in the future. To this end, there are 47 parameters selected (Sexton et al. submitted, Table 1), which are associated with the processes of aerosol, land surface, boundary layer, gravity wave drag, convection, cloud, and radiation, according to their ability to control physical processes (Sexton et al. 2019). Then the plausible ranges and prior probability distributions of the parameters were determined based on the analyses of high-resolution process models, such as large-eddy simulations and cloud resolution models, as well as constraints of observations [the distributions of the ranges for the 47 parameters shown in Fig. 1 of Sexton et al. (submitted)]. Subsequently, these 47 parameters were perturbed using a maximin Latin hypercube design (McKay et al. 1979) to generate a large sample of approximately 3000 different combinations. For these combinations, two coarse-resolution experiments were conducted: one is a climate simulation forced with time-varying CO2 concentration as CMIP5 (CRP8.5 from 2005 forward) and daily SST and sea ice during the 1999–2008; the other is a 5-day weather forecast. By calculating the weather forecast errors, mean climate biases against observations and reanalysis data, a subset of 25 parameter combinations were identified in terms of their response to SST and representation of climate outcomes driven by CO2 and aerosol that are as diverse as possible (Murphy et al. 2018; Karmalkar et al. 2019) to run the CMIP5 historical (1900–2005) and RCP8.5 emissions (2006–2100). Details of the choice of parameter combinations are explained in Sexton et al. (submitted), and the experimental configuration for coupled PPE simulations is reported in Yamazaki et al. (submitted).

Fig. 1
figure 1

Comparisons of correlation coefficients between globally averaged annual mean soil moisture from 20 HadGEM3-GC3.05 PPE members, four GLDAS simulations, and two reanalyses over 1979–2008. "Asterisk" on the color bar denotes significance at the p < 0.05 level as n = 30. "Caret" indicates the four PPE members removed from the following analyses

Full size image
Table 1 Information on in situ soil moisture observations and observational networks
Full size table

Since these simulations were designed to be plausible and diverse in their climate response, they have been used as part of a package of information on future climate for the UK Climate Projections [UKCP18; Murphy et al. (2018)]. 5 of the 25 members were not included in UKCP18 because of poor performance at representing the Atlantic meridional overturning circulation and historical global temperature trends. Although this does not necessarily mean that they perform poorly in their representation of climate zones, this study still excluded the five members for a more prudent expectation. Yamazaki et al. (submitted) show that the PPE is diverse in terms of precipitation changes for many worldwide regions, often to a similar degree to CMIP5. HadGEM3-GC3.05 also has a relatively high climate sensitivity [of approximately 5.4 K for the standard variant (Andrews et al. 2019)], so the results here should be considered in the context that these 20 PPE runs are all high warming scenarios. In the present study, we will further evaluate soil moisture and associated land water components in the 20 PPE members in Sect. 3 and identify adequate members for indicating shifts in climate zones. The variables of model data used include soil moisture, evapotranspiration, groundwater table, precipitation, and surface air temperature, with an N216 horizontal resolution (approximately 0.56° × 0.83° of latitude/longitude, 60 km in mid-latitudes), from simulations of CMIP5 historical before 2005 and RCP8.5 emissions for the period of 2005–2100.

2.2 Observation-based data

In situ observations of soil moisture were taken from the international soil moisture network (ISMN, https://ismn.geo.tuwien.ac.at), which is a global soil moisture database maintained by international cooperation. The ISMN collects available observations from around the world through coordination by the Global Energy and Water Exchanges Project (Dorigo et al. 2011). To augment spatial coverage, we use observations with records longer than 36 months (frozen period excluded based on soil temperature), which can represent the local climate and meet the statistical requirement of the large sample theory (n > 30). The selected observations cover 20 networks with 644 sites from low to high latitudes. With respect to these sites, we filtered out partial records of less than “Good” quality: soil moisture (1) values exceeding the range between 0 and 0.6 m3 m− 3 or larger than saturation points; (2) records of soil temperature below zero; (3) peaks without precipitation events; (4) spikes and negative/positive jumps in the soil moisture spectrum; (5) low constant values; and (6) saturated plateaus in the soil moisture spectrum, according to Dorigo et al. (2013). Finally, all the records of the original time intervals (Table 1), which have “Good” quality ratios of over 85% over selected periods, are used to estimate monthly means. The measurements were conducted at various depths ranging from 5 to 50 or to 100 cm below the soil surface with sensors for high-frequency measurement and the oven-drying method for those of lower frequencies in China, Mongolia, and the former Soviet Union regions. These in situ observations have played an important role in the evaluations of satellite and model products (e.g., Robock et al. 2000; Dorigo et al. 2011; De Lannoy et al. 2014; Wang et al. 2016). Further details are presented in Table 1, and the information about instruments and data quality control processes can be found in the reports of networks and the references therein (available from https://ismn.geo.tuwien.ac.at).

Two reanalyses and four offline land model soil moisture products are used to evaluate global spatiotemporal patterns of soil moisture produced by all ensemble members. These are the National Centers for Environmental Prediction Climate Forecast System Reanalysis (CFSR) (Saha et al. 2010), the fifth-generation atmospheric reanalysis by the European Centre for Medium-Range Weather Forecasts (ERA5) (Hersbach et al. 2020), and four members of the Global Land Data Assimilation System (GLDAS) (Rodell et al. 2004). CFSR has a high horizontal resolution (0.31° × 0.31° over 1979–2010, 0.2° × 0.2° over 2011–2017), and ERA5 has a 0.28° × 0.28° latitude-longitude resolution. GLDAS simulations were driven by atmospheric forcing from the combination of models, satellite retrievals, radar measurements, and observations. We use monthly outputs from 1979–2009 at a 1° × 1° latitude-longitude resolution.

To indirectly examine soil moisture simulations, land surface evapotranspiration and terrestrial water storage simulations are also validated with observation-based data. The evapotranspiration (ET) is estimated using latent heat flux (LE) data given that ET = LE/2.501, where the units are mm and 106 Joule m− 2 s− 1, respectively. The observation-based LE dataset is derived from a global network of eddy-covariance towers (FLUXNET) (Jung et al. 2009), covering the period of 1982–2011 with a 0.5° × 0.5° latitude-longitude resolution. Another satellite-based ET data set, the Moderate Resolution Imaging Spectroradiometer (MODIS), is also used for cross-examination. The MODIS data set has a 0.5° × 0.5° resolution, covering the 2000–2014 period (Mu et al. 2011).

The terrestrial water storage anomalies (integrated from the anomalies of the groundwater table, soil moisture, snow water, canopy water, and river water in the HadGEM3-GC3.05 model) are assessed based on satellite retrievals from the Gravity Recovery and Climate Experiment (GRACE). Herein, the mean anomalies of GRACE products are used based on the regularized mass concentration solutions from two data centers (NASA’s Jet Propulsion Laboratory, the University of Texas at Austin Center for Space Research), according to Watkins et al. (2015) and Save et al. (2016). The data products have a 0.5°× 0.5° resolution, covering 2002–2017. Details about GRACE data products can be found in GRACE Tellus (https://grace.jpl.nasa.gov) and the publications therein.

2.3 Methods

2.3.1 Comparisons of PPE mean soil moisture with in-situ observations

To validate the soil moisture of PPE simulations by HadGEM3-GC3.05, we conducted site-grid cell comparisons between monthly observations and model outputs for the grid cells that contain the observational sites. The consistency measured by the Pearson correlation coefficient includes the annual cycle because it is still a challenge for land models to regenerate soil moisture in phase with observed annual cycles. Furthermore, the relative variations in RCP8.5 to historical soil moisture largely determine the results of the present study. Owing to the sparsity of in situ observations and despite discontinuous monthly data, the sites that still have over 36 samples after the frozen periods were excluded were included in the comparisons by scoring out the model outputs at the corresponding time. In addition, their errors and variance are also compared at the site-grid scale. To address the mismatches of soil layer depths in HadGEM3-GC3.05 with in situ observations and the other reanalysis and land models, the HadGEM3-GC3.05 soil moisture was linearly interpolated to the corresponding observation or model depths, analogous to the reanalysis or observations in Li et al. (2005) and Albergel et al. (2013). Because of the various and short durations of the individual in situ observations, we integrated them into a time series by averaging soil moisture across sites along the timeline from 1979 to 2017, serving as a subset of the statistical population of climate during the whole period, to assess model fidelity on multiple time scales, especially the long-term trend. The model time series is correspondingly generated at the same time using outputs with observations. The resulting time series may not represent a real climatology, but it is fair and more robust to assess consistency between observations and simulations considering the sparsity of in situ observations. Then, we decomposed the resulting time series into long-term, seasonal, and subseasonal and noise signals using a method of seasonal trend decomposition based on Loess (Cleveland et al. 1990) in an additive model, which is described as follows:

$${Y}_{t}={T}_{t}+{S}_{t}+{R}_{t},$$
(1)

where Yt denotes monthly time series of observations or simulations; and Tt, St, and Rt denote decomposed long-term, seasonal, and the remaining components during the period of t, respectively. The long-term component includes linear and nonlinear signals. The nonlinear part reflects repeated but nonperiodic fluctuations. The seasonal component occurs over a fixed period (e.g., annually). The residual portion represents subseasonal signals and noise.

2.3.2 Classification of climate zones

The climate zones, in the present study, are classified using root-layer (0–1 m) soil moisture content, which includes water in both liquid and solid phases. As a result, the total soil water budgets and water resource potential can be represented in the context of a warming climate, as well as water-heat interactions. This representation enables the detection of climate impacts from both precipitation variability and rising temperatures in future scenarios, especially in the mid-high latitudes. Soil moisture-based climate zones contain changes not only in atmospheric climate change but also in interactions with land processes. Based on the spatial patterns of well-established climate classifications, such as the AI index, Koppen classification, and precipitation (Kottek et al. 2006; Huang et al. 2016), as well as the features of soil moisture, we classified the global land areas into four basic dry-wet climate zones: arid, semiarid, subhumid, and humid zones by the mean volumetric soil moisture of 18, 23, and 33% in the root layer in HadGEM3-GC3.05 during 1976–2005 (30 years) as a climatic reference, relative to which their changes are considered to be the potential shifts in climate zones under a future global warming scenario. The thresholds are chosen such that the classification gives nearly consistent results to classification based on the established indices during 1976–2005 (for more details, see Sect. 3.2). The former three types comprise the areas of global drylands. The drylands are also detected using an atmospheric-only aridity index to indicate the land process effects of climate change. The aridity index is a ratio of annual total precipitation to potential evapotranspiration (P/PET). Values less than 0.65 denote drylands, which further fall into four subtypes of hyperarid (< 0.05), arid (0.05–0.2), semiarid (0.2–0.5), and dry subhumid (0.5–0.65) areas (Middleton and Thomas 1992).

2.3.3 Relationships of areal changes in climate zones to the climatic driver and other statistical analysis methods

The relationships of areal changes in climate zones to changes in climate factors (precipitation and surface air temperature) are investigated using ordinary least squares linear regression. Then, coefficients are tested using a t-test, with the t statistic estimated as follows:

$$t=r\sqrt{\frac{n-2}{1-{r}^{2}}},$$
(2)
$$r=\frac{{\sum }_{i=1}^{n}({x}_{i}-\overline{x})({y}_{i}-\overline{y})}{\sqrt{{\sum }_{i=1}^{n}{({x}_{i}-\overline{x})}^{2}}\sqrt{{\sum }_{i=1}^{n}{({y}_{i}-\overline{y})}^{2}}},$$
(3)

where r denotes the correlation coefficient, n is the sample size, and xi and yi are the soil moisture and climatic driver, with corresponding means of \(\overline{x}\) and \(\overline{y}\), respectively. Then, the probability value p is estimated using a t-distribution with n–2 degrees of freedom. When p is less than the significance level (e.g., 0.05), the regression coefficient is statistically significant. The significance of Pearson correlation coefficients (formula 3) and the difference in means is tested using the same method. The absolute errors are computed as the differences between observations and simulations, which are divided by the means of observations as the relative errors. The root mean square error (RMSE) is computed as follows:

$$RMSE=\sqrt{\frac{\sum _{i=1}^{n}{\left({x}_{HadGEM,i}-{x}_{obs, i}\right)}^{2}}{n}},$$
(4)

where xHadGEM,i and xobs,i are modeled and observational soil moisture time series, respectively, with the sample size n. Standard deviations (SD) are estimated using the following formula.

$$SD=\sqrt{\frac{\sum _{i=1}^{n}{\left({x}_{i}-\overline{x}\right)}^{2}}{n-1}},$$
(5)

3 Results

3.1 Evaluation of PPE soil moisture simulations

Given that the PPE experiment aims to provide possible plausible pathways of future climate change, from the perspective of soil moisture based on spatiotemporal consistency to offline simulations and reanalyses over the 1979–2008 period (Fig. 2), we selected 16 of the 20 PPE members to investigate future soil moisture climate zones. Four members (E04, E14, E16, and E19 in Fig. 2) were excluded due to their negative correlations with offline simulation and reanalysis of soil moisture. The selected 16 members are spatiotemporally comparable with offline land model simulations and two high-resolution reanalyses at global and annual scales (p < 0.05).

Fig. 2
figure 2

Correlation coefficients between monthly soil moisture of 16-member means and in situ observations. Sample numbers vary with individual sites; 30 is the minimum sample size, and the responding significant coefficient is 0.349 at the p < 0.05 level

Full size image

Then, we further evaluated the fidelity of the selected 16 PPE members in terms of land water cycles using observations on various spatiotemporal scales. At the site scale, the monthly soil moisture of 16 PPE mean was compared with in situ observations at 644 sites around the world (Fig. 1). The simulated soil moisture (soil frozen times excluded) shows reasonable consistency with observations, with significant correlation coefficients for 67% site-grid cell pairs (431 out of 644, p < 0.05) and 53% (20%) of coefficients larger than 0.5 (0.7). Due to the discontinuity of observations, the time series (for the nonfrozen period sites) partly retain the annual signals. Thus, the significant correlation coefficients indicate the consistency of annual dynamics with the observed, and then the adequacy of soil moisture transmission-related representations in the model. The lower correlations appear regionally in Mongolia, China, and the middle USA. The cause is largely associated with the relatively low sampling frequency of the observations over these networks (once per 10 days, Table 1) and low observational quality and spatial representativeness in the mountainous areas, e.g., the USCAN network, USA (Dorigo et al. 2011). In terms of the errors and amplitudes (Online Resource, Figs. S1–5), overall, the PPE mean soil moisture is systematically wetter than the observations, with smaller variation amplitudes. Specifically, the absolute errors range from − 28.6 to 23.6%, with 97.4% of sites being less than 10 percent; the RMSEs of less than 10% account for up to 65.4% (421) of the total 644 sites; the sites of less than 50% relative errors are up to 73.6% (474 sites). The differences in SD vary from − 4 to 3.8% and less than zero at 77.5% of sites. Although site-grid cell comparisons provide valuable information on model performance, the spatiotemporal mismatches between in situ observations and grid-average simulations might make for a less favorable evaluation (Entekhabi et al. 2010; Zeng et al. 2015; Bi et al. 2016). For filtering out parts of site-specific signals in observations by regional mean, the regionally integrated data hence enhance the regional characteristics, which are more appropriate for evaluations of modeled soil moisture.

To this end, the network-averaged soil moisture is compared separately for the 20 networks. The monthly time series over various periods (Online Resource, Fig. S6) show that model soil moisture can regenerate the significant variations in 19 of 20 networks, with correlations ranging from 0.3 to 0.87 (p < 0.05), larger than 0.5 in 14 networks. The annual cycles for the 20 networks (Fig. 3) also show good agreement between modeled and observational soil moisture despite mismatches in peak months in parts of networks (e.g., Fig. 3g, h, o). Note that these annual cycles might not characterize the real climatological cycles due to some soil moisture being removed for frozen soil at various depths and sites. On the global scale, all available observations were then integrated into a time series along the timeline from 1979 to 2017, also representing a subset for the statistical population of global soil moisture, to verify the model’s ability to generate long-term soil moisture variability and trend. The globally integrated monthly simulations agree well with the observations, with a significant correlation coefficient of 0.81 (0.68 after annual cycles were removed, p < 0.01, Fig. 4a). Although simulations are 16% larger than the observations on average, the systematic bias in the means has only a small effect on the assessment of relative changes in soil moisture between two periods. According to the seasonal trend decomposition, the main features of long-term changes also agree well with those observed (r = 0.85, p < 0.01). However, the linear term of simulations is less than that of the observations (3.0 and 5.7 volumetric percentage points, pp, for simulations and observations over the 39 years, Fig. 4b). The seasonal terms highlight the remarkable ability of the model to regenerate annual cycles in observations (r = 0.98, p < 0.01), despite larger annual ranges (10.6 and 7.6 pp for simulations and observations, Fig. 4c). Regarding the higher frequency variability (residual term), the climate model captures the substantial fluctuations observed (r = 0.38, p < 0.01, Fig. 4d). Generally, on the regional scale, the climate model generated reasonable soil moisture variability on various time scales.

Fig. 3
figure 3

Annual cycles of network-averaged PPEs mean soil moisture and observations. The serial numbers (a–t) in panels correspond to those for networks in Table 1

Full size image
Fig. 4
figure 4

Seasonal trend decomposition of integrated monthly soil moisture for observations and the corresponding PPE mean. Here, r denotes the correlation coefficient, and "double asterisk" denotes the significance level p < 0.01. Linear trends are 3.0 and 5.7 volumetric percentage points for simulations and observations, respectively, over the 39 years

Full size image

To compensate for the scarcity of available soil moisture observations, we further compared another two components of the land water cycle continuum: evapotranspiration, ET, and terrestrial water storage anomalies. In terms of ET, the comparisons with two observation-based products (Fig. 5) showed that simulations are generally larger than observations, particularly over the regions with high ET totals (> 20 mm per month, Fig. 5a, b). Moreover, the biases vary with months (or ET monthly totals), and smaller simulations occur in months with higher ET totals (spatial averages from 0° to 60° N) compared with FLUXNET (Fig. 5c, d). Against the MODIS product, smaller biases are observed in months with low ET totals (Fig. 5c, e). Nonetheless, such systematic biases have less influence on the adequacy of the spatial pattern of the model ET. Furthermore, the comparisons between annual cycles show that the model largely captures the annual cycles in the two observation-based products (Fig. 5d, e). Given the complexity of observations, the comparisons of ET indicate reasonable adequacy of water-heat related representations in the model.

Fig. 5
figure 5

Comparisons of monthly ET between HadGEM3-GC3.05 simulations and observation-based products. Spatial patterns for the ET differences between HadGEM3-GC3.05 and FLUXNET (1985–2005), and MODIS (2000–2014) in a, b. Scatters denote monthly ET averaged over 0° to 60° N during 1985–2005 for FLUXNET and HadGEM3-GC3.05 and 2000–2014 for MODIS and HadGEM3-GC3.05 in c and corresponding annual cycles d, e. Black dots denote no significant differences between two ET data at the p < 0.05 level by the paired sample t test in a, b. Dashed lines in c, e denote the standard deviations for 16 HadGEM3-GC3.05 PPE members

Full size image

From the aspect of terrestrial water storage anomalies, the model also performs reasonably well (Fig. 6). The spatial pattern of simulated anomalies largely characterizes that of GRACE, despite the generally smaller amplitudes of variations (Fig. 6a, b). In addition, the time series of globally averaged simulations captured the temporal patterns in the annual dynamics of observations, with a significant correlation coefficient of 0.71 (p < 0.01). The linear trend in normalized simulations also regenerates the decreases in observed anomalies over 2002–2017, despite a 57% slower speed (Fig. 6c). The discrepancies largely result from the absence of representations of anthropogenic effects on terrestrial water storage anomalies in the HadGEM3-GC3.05 model. For instance, groundwater withdrawals and water diversion in China have significant impacts on the natural patterns of terrestrial water storage (Long et al. 2020). However, the exact mechanism for the difference requires further study. In summary, the terrestrial water storage anomalies integrated from multiple components imply that the model can represent the main processes of land water cycles.

Fig. 6
figure 6

Spatial patterns of averaged terrestrial water storage (TWS) anomalies (relative to the 2004–2009 time-mean baseline) during the period of 2002–2017 for GRACE satellite data (a, in the form of liquid water equivalent in millimeters, lwe mm) and HadGEM3-GC3.05 PPE mean (b) and monthly evolutions of global averaged time series with the regions covered by constant ice removed (c)

Full size image

These evaluations of soil moisture, evapotranspiration, and terrestrial water storage, along with atmospheric studies using HadGEM3-GC3.05 (e.g., Walters et al. 2017a, b), suggest that the selected PPE members can reasonably encompass the observed variability in land water cycles. This result gives us confidence that the model can underpin climate projections under the RCP8.5 scenario.

3.2 Climate zone classification based on soil moisture

We classified the global land into four dry-wet climate zones that are arid, semiarid, subhumid, and humid climate zones (Fig. 7) using 30-year average soil moisture in the root layer (0–1 m) with the reference period of 1976–2005. Herein, soil moisture is the average soil moisture of 16 PPE members (Sect. 2.3.2). Among the four climate zones, the former three-type zones comprise drylands. The areal percentages of arid, semiarid, subhumid, and humid climate zones are 20.9, 11.3, 36.6, and 31.2%, respectively, and 68.8% of drylands, relative to the global land areal total (excluding the Greenland and Antarctic areas). The global pattern of climate zones is generally consistent with those of classifications using atmospheric-only indicators (Belda et al. 2014; Huang et al. 2016; Zhang et al. 2017). However, the boundaries between climate zones can vary substantially with different classification methods (Belda et al. 2014). In the present study, soil moisture is constrained not only by atmospheric processes but also by land processes, such as soil properties, vegetation, topography, and water-heat feedbacks. For example, at mid-high latitudes, soil temperature is a dominant factor controlling soil moisture dynamics (Niu and Yang 2006). Consequently, this classification appropriately represents the combined responses of land processes to global warming.

Fig. 7
figure 7

Geographical distributions of soil moisture-based wet-dry climate zones using 30-year average root layer soil temperature (16 PPE mean from HadGEM3-GC3.05 historical simulations) with the reference period of 1976–2005

Full size image

3.3 Response of soil moisture climate zones to climate change in the RCP8.5 scenario

3.3.1 Spatial shifts in climate zones

Shifts in climate zones indicate changes in regional climate regimes and hence profound climate change. In the RCP8.5 scenario, extensive areal expansions or contractions of climate zones (Fig. 8) suggest that shifts between climate zones will occur throughout the 2030s, 2060s, and 2090s in comparison with the baseline period (1976–2005). In the arid zone, along its peripheries, new regions shift from semiarid zones. Meanwhile, in southern America, mid-eastern Africa, Central Asia to northern China, and northern Australia, parts of arid areas shift into wetter zones, with increased probability after the 2030s (Fig. 8a, e, i). The expansions of semiarid regions mainly emerge in middle Siberia, northern America, the Tibetan Plateau and eastern adjoining areas and parts of the semiarid zone peripheries. Over the contracted areas, the semiarid climate zone mainly shifts into the drier arid zone, while in parts of areas in eastern Africa, Central Asia, and northern China, the zone shifts into the subhumid zone. The newly emerged semiarid areas are remarkably enlarged by the 2070s in comparison with those in the 2030s (Fig. 8b, f, j). The subhumid zone shows a marked northward expansion. This zone can reach approximately 67° N, then 73° N, and further north by the 2030s, 2070s, and 2090s, respectively (Figs. 8, 9). The subhumid regions also spread in most humid lower latitudes, except southern China and most of the Maritime Continent. In particular, the peripheries of the Amazon and Congo basins gradually change into the subhumid zone from the humid climate zone (Fig. 8c, g, k, d, h, l). By the 2090s, the shifts in climate zones are characterized by northward expansions of the drier climate zones. As a result, most humid regions over mid-high latitudes, as well as parts of the Amazon and Congo basins, shift to the subhumid zone (Fig. 9). Overall, the dominant shifts between the four climate zones highlight the tendency of drying soil moisture in the root layer along the boundaries of climate zones throughout the 21st century in the RCP8.5 scenario.

Fig. 8
figure 8

Variations in the four climate zones in the 2030s, 2060s, and 2090s relative to the baseline period of 1976–2005. The red end shows the probability of expansion and the blue end shows the probability of contraction of different climate zones in the future (2030s, 2060s, and 2090s) relative to 1976–2005. The probability is estimated using 16 PPE members

Full size image
Fig. 9
figure 9

Variations in the four climate zones in the 2090s relative to the baseline period 1976–2005. The increased regions indicate the expansion of each climate zone, with colors denoting the shifts from other climate zones based on the baseline period; the decreased regions indicate the contraction of each climate zone, with colors denoting the shifts to others in the 2090s

Full size image

In terms of the temporal evolutions of the areas covered by each climate zone (Fig. 10), in the arid climate zone, it will decrease by approximately 0.63 million km2 by 2.3% in 2099 in comparison with the mean area over 1976–2005. The semiarid and subhumid zones will expand by 7.7% and 9.5%, which comprise 1.15 and 4.59 million km2, respectively. Accordingly, the humid zone contracts 12.6%, which is approximately 5.12 million km2. Furthermore, according to individual projections among the 16 PPE members, the strongest areal expansions are up to 23.6 and 20.0% by 2099 for the semiarid, and subhumid zones, at the expense of a 10.4% and 29.5% contraction in the arid and humid zones, respectively (Fig. 10). For the mean annual cycles during the RCP8.5 scenario period (figures not shown), the strongest expansions toward drier zones occur from July to September, while the most remarkable contractions appear from April to May. In addition, a probability density function analysis (figures not shown) indicates that in the humid zone interannual areal changes tend to intensify during the second half of the 21st century compared with the first half.

Fig. 10
figure 10

Variations in the areas of four wet–dry climate zones according to the soil moisture index during historical (1976–2005) and RCP8.5 scenario (2006–2099) periods relative to the means of the period of 1976–2005. The colored lines denote area anomalies based on the soil moisture of 16 PPEs, and the solid red line is the mean of 16 members. The percentages are minimum, mean, and maximum area changes during the 2090s relative to the mean over 1976–2005. S denotes the linear regression trend by ordinary least square regressions over the scenario period with a p value of the level of significance

Full size image

3.3.2 Areal changes in global drylands

Since drylands are susceptible to variations in regional water budgets, areal changes in global drylands can be expected if there are any shifts in climate zones. The projected areal expansion is remarkable in the RCP8.5 scenario throughout the 21st century in light of the soil moisture index (Fig. 11). Compared with the baseline period (1976–2005), the expanded areas mainly emerge over the mid-high latitudes, such as northern Canada and mid-eastern Siberia, by the end of the 2090s according to the PPE mean. Regional expansions also appear in the European continent, the Tibetan Plateau, the eastern adjoining areas, and the northeastern USA. In addition, the Amazon and Congo basins, along with parts of the maritime continent, namely, the major tropical rainforest regions, undergo the impact of drying expansions from their peripheries. In contrast, dryland contractions are observed in the scattered areas of the western USA, mid-eastern Asia, and middle Africa (Fig. 11a). Regarding temporal evolution, drylands expand to approximately 9 million km2 by the end of the 2090s, 10.5% in the ensemble mean compared with the baseline period (Fig. 11b). The range of expansion across the PPE is from 0.6–19.0%. This diverse range occurs despite the narrow range in global temperature changes explored by the PPE, suggesting that this is not the main driver of the PPE spread. However, the increasing speed tends to slow after the 2050s, accompanied by a decrease in leveling off in the humid zone areas. Subsequently, most humid zone regions will be subjected to soil drying under the RCP8.5 scenario, indicating a substantial decrease in soil water storage.

Fig. 11
figure 11

Spatial pattern of areal changes (over 50% probability based on 16 PPEs) in drylands by soil moisture index between 2090s and 1976–2005 means under RCP8.5 (a) and temporal evolutions of dryland areas from historical (1976–2005) to RCP8.5 scenario (2006–2099) periods (b) relative to the mean of the period of 1976–2005

Full size image

3.3.3 Drivers of soil moisture climate zone shifts

To detect the drivers of areal changes in the four soil moisture climate zones under the RCP8.5 scenario, we conducted trend and regression analyses using PPE projections. The results (Fig. 12) show that the drier/humid climate zones significantly expand/contract during the 21st century, although precipitation is expected to significantly increase over most regions, except northeastern South America, southern Europe, and parts of western Africa (Fig. 12a). In the context of global warming (Fig. 12b), the trends of root-layer soil moisture vary with regions in contrast to those of precipitation due to the regulation of land-atmosphere interactions. For instance, over the mid-high latitudes, the regional mean water balances (north of 60° N, Northern Hemisphere, Online Resource, Fig. S7) show a cumulative decrement in soil moisture accompanied by increases in precipitation, evapotranspiration, and runoff. Since HadGem3-GC3.05’s land model JULES incorporated the groundwater scheme, from the perspective of water budgets, it can be speculated that soil moisture-groundwater exchanges will play an important role in land water balance as the temperature rises under the RCP8.5 scenario during the 2090s. Moreover, decreasing trends in soil moisture can be seen in Mongolia and eastern adjoining regions, the Tibetan Plateau, the Mediterranean regions, and parts of western Africa, eastern North America, and northern South America. Conversely, soil moisture increases significantly in Central Asia, northern China, parts of South Asia and northern Australia, eastern Africa, western North America, and southeastern South America (Fig. 12c). Taking Central Asia as an example (latitudes from 40° N–50° N, longitudes from 55° E–80° E, Online Resource, Fig. S8), the soil moisture increase is mainly contributed by increases in precipitation. Subtracting runoff and evapotranspiration losses and precipitation cumulatively increase soil moisture from the start year, with less increments in 10 years (during the 2090s). However, the mechanisms behind water cycle changes under RCP8.5 require specific observations and model experiments in further studies.

Fig. 12
figure 12

Patterns of the linear trends for annual total precipitation, surface air temperature, and soil moisture under RCP8.5 during 2006–2099, with black dots denoting significance at the p < 0.05 level

Full size image

4 Discussion

Two particular aspects must be addressed before we draw our conclusions: (1) the differentiation in responses to global warming among different components of the climate system and (2) the difference between dryland expansion rates defined by the atmospheric aridity index and by soil moisture.

Figure 13 compares the model spreads among the 16 PPE members between precipitation (Pr), actual evapotranspiration (AET), runoff (RO), surface air temperature (Tas) and soil moisture anomalies in RCP8.5 relative to the baseline period of 1976–2005. The long-term evolution of water and heat processes (Fig. 13) shows that global land-averaged surface air temperature is expected to rise 59.7% (47.4–74.5% across 16 PPEs) averaged in the 2090s against that in the baseline period (1976–2005) and the increase is up to 8.7 K (7.3–10.3 K). Precipitation is expected to increase by 8.7% (4.8–13.2%), and runoff is expected to increase by 26.1% (20.0–34.4%), up to 95.5 mm (13.6–188.8 mm) and 60.9 mm (42.4–82.7 mm), respectively. In contrast, actual evapotranspiration and soil moisture decrease 7.7% and 5.1%, respectively, up to 36.9 mm (− 129.3–9 mm) and 1.5 volumetric percentage points, pp, (− 0.5 to − 3.2 pp), respectively. The increased precipitation at low-middle latitudes is mainly balanced by the increased runoff and slightly increased soil moisture. However, the soil across mid-high latitudes loses water and hence leads to the drying tendencies of soil and drying shifts in climate zones, despite the increases in precipitation (Fig. 12a–c, Online Resource, Fig. S9).

Fig. 13
figure 13

Global land average annual mean precipitation (Pr), actual evapotranspiration (AET), runoff (RO), surface air temperature (Tas) and soil moisture anomalies (in the form of percentage points, pp) relative to their means over the baseline period of 1976–2005, with changed minimum, mean, and maximum percentages (right axes) of the means during the 2090s relative to the baseline period

Full size image

In terms of soil moisture, the correlation analysis (Fig. 14) over the last 30 years (2070–2099) between soil moisture and all parameters perturbed in PPE runs indicates that variations in soil moisture are most sensitive to those in two parameters: the dependence of photosynthesis on temperature (tupp_io) and mixing detrainment rate (amdet_fac), with correlation coefficients of − 0.82 and − 0.67 (p < 0.01), respectively. The parameter tupp_io denotes the optimal temperature for photosynthesis determining the turnover point for temperature, above which further increases in temperature will drive a decline in photosynthesis, and the optimal temperature would be expected to have a large impact on plant functioning. The higher tupp_io is, the broader the temperature range that plants hold active photosynthesis, and the more soil moisture that would be lost by plant transpiration.

Fig. 14
figure 14

The relationship of changes in soil moisture with two key parameters in 16 PPE simulations during 2070–2099. Here, tupp_io denotes the optimal temperature for photosynthesis, whose low values lead to greater temperature dependence of photosynthesis, and amdet_fac is mixing detrainment, whose larger values increase large-scale humidity and temperature profiles. The correlation coefficients, − 0.67 and − 0.82, are both significant at the p < 0.01 level with a sample size of 30 (last 30 years, 2070–2099) and a degree of freedom of 28 (sample size − 2)

Full size image

The other parameter is mixing detrainment, amdet_fac, incorporated with “forced detrainment” (related to buoyancy loss in the core updraft), in the Gregory-Rowntree convection schemes in the HadGEM3-GC3.05 model (Derbyshire et al. 2011). Mixing detrainment regulates the outflow of cloudy air into the environment; as a result, a larger value leads to stronger air mixing between clouds and the ambient atmosphere and then convection inhibition followed by decreases in convective precipitation. On the other hand, the larger values of the parameter tend to enhance large-scale humidity profiles favorable to the formation of large-scale precipitation. However, for the large ratio of convective precipitation to the total precipitation amount over land (Online Resource, Fig. S10), the decreased convective precipitation leads to a decrease in the total precipitation and soil moisture. This mechanism largely accounts for the negative correlation between variations in mixing detrainment and soil moisture (corr = − 0.67, Fig. 14). Nevertheless, it must be noted that the relationships of mixing detrainment with convective and large-scale precipitation are not straightforward, much less those with soil moisture, due to the resulting nonlinear effects on cloud radiation, cloud radiation-precipitation processes, and land-atmosphere feedbacks. In summary, from the water input and output aspects, the two parameters, among the all parameters, mostly affect the soil water budget in the HadGEM3-GC3.05 model.

Regarding the difference in dryland expansion rates from the perspective of the atmospheric aridity index and soil moisture, both dryland definitions, using the two indicators, show similar areal expansions in Europe and lower latitudes. However, more extensive expansion in the mid-high latitudes is revealed by soil moisture in comparison with that by the aridity index (Figs. 11, 15). One of the main reasons behind this result is that warming involves more frozen soil water in the land-atmosphere water cycles, which is further enhanced by boosted vegetation. There has been an increasing body of evidence showing observed Arctic greening trends (Potter et al. 2013; Vickers et al. 2016; Edwards and Treitz 2017). Additionally, with global warming, freeze/thaw-induced transitions in soil hydrological conditions open up pathways for soil water to the groundwater flow system, resulting in substantial soil drying, even though water inputs (precipitation minus evaporation) tend to rise in the mid-high latitudes (Lawrence et al. 2015). Compelling increases in river base flow owing to increased groundwater have been observed across Canadian and Alaskan Arctic rivers (Walvoord and Striegl 2007; Wu et al. 2008, 2010; Bense et al. 2012). However, the atmospheric-only indictor (i.e., P/PET, aridity index) cannot efficiently depict such geophysical-chemical processes. Additionally, PET does not directly correlate with temperature variations (Mckenney and Rosenberg 1993; Kafle and Bruins 2009). As a result, soil moisture may exacerbate the drying effects of atmospheric processes in mid-high latitudes. Unfortunately, atmospheric drying trends are also projected to rise in these regions according to multimodel atmospheric aridity indices (Huang et al. 2016; Zhang et al. 2017). Taken together, in the mid-high latitudes, climate zones shift to drier zones under the RCP8.5 scenario. The unsustainability of the soil water supply to enhanced evapotranspiration and groundwater discharge loss will determine the fate of Arctic greening and the transformation of ecosystems.

Fig. 15
figure 15

Spatial pattern of areal changes (over 50% probability based on 16 PPEs) in drylands by aridity index between the 2090s and 1976–2005 means under RCP8.5 (a) and temporal evolutions of dryland areas from historical (1976–2005) to RCP8.5 scenario (2006–2099) periods (b) relative to the means of the 1976–2005 period

Full size image

However, the shifts from humid to sub-humid zones in the Amazon and Congo basins mainly result from the decreases in precipitation (Fig. 12), which are controlled by atmospheric processes. Although the projections of precipitation in these regions vary across models, a strong consensus with 13 CMIP5 models highlights a drying and lengthening of the dry season (from June to September) (Joetzjer et al. 2013), and similar trends are also seen according to 36 GCMs (Boisier et al. 2015). Based on these findings, climate drying can be propagated into root-layer soil in these regions, resulting in substantial soil moisture decreases, hence increasing the water-stress risk to rainforest ecosystems, which is the most important carbon sink. On the other hand, dryland contractions at middle latitudes according to the atmospheric aridity index are not observed in light of the soil moisture index (Figs. 8, 9, 11, 15), for example, in the western USA, northern Mexico, China, and India, South America, eastern Africa, southern India, and northeastern Australia. The discrepancies suggest that the limited enhancement of precipitation cannot convert the dry regimes of the soil moisture climate and thus has a slight influence on the dry terrestrial environments in these regions under RCP8.5. Furthermore, the shifts in the peripheries of the four subtype zones may substantially alter local climate regimes, which highlight the accompanying vulnerability of the local ecosystems.

5 Summary and conclusions

Climate change is expected to reshape the pattern of climate zones as the temperature rises and precipitation changes. Many indicators are used to define climate zones and examine their shifts, but soil moisture is a metric directly linked to ecosystem functions and a manifestation of the climatic balance among the atmosphere, biosphere, hydrosphere, and lithosphere. Using soil moisture classification and PPE simulations, this study has investigated potential changes in terrestrial climate zones, particularly those associated with snow/glacier melting, frozen soil thawing, and groundwater changes. The PPE allows us to estimate a range of potential future changes due to modeling uncertainty.

An evaluation against available observations suggests that the perturbed physics ensemble produces realistic simulations of soil moisture and its spatiotemporal variability. According to root-layer soil moisture, the global land surface is classified into four dry-wet climate zones: arid, semiarid, subhumid, and humid zones, with a reference period of 1976–2005. The first three zones make up the global drylands (68.8% of land areal total excluding the Greenland and Antarctic areas), including 20.9% of the arid zone, 11.3% of the semiarid zone, and 36.6% of the subhumid zone. The humid zone accounts for 31.2% of areal total areas globally.

Under the RCP8.5 scenario, the global drylands are projected to experience a significant expansion (10.5%) at the expense of the humid zone by the end of the 21st century relative to the 1976–2005 climatology. This result is qualitatively consistent with previous studies using different metrics. The areal expansion consists of 7.7% (–8.3 to 23.6%) and 9.5% (3.1–20.0%) semiarid and subhumid zone growth at the expense of 2.3% (–10.4 to 7.4%) and 12.6% (–29.5 to 2.0%) arid and humid zone contraction. The 10.5% ensemble mean dryland expansion rate globally is slower than the 19.2% estimate using the precipitation-potential evapotranspiration ratio, but regional expansion over the mid-high latitudes can be much greater. This difference is mainly because of frozen soil thawing and accelerated evapotranspiration with Arctic greening and accelerated polar warming, which can be detected in soil moisture but not from atmosphere-only indices.

The projected potential shifts in climate zones highlight areas of ecosystem vulnerability and climate risks. Clearly, high risks mainly exist in the peripheries of subtype zones, especially in the mid-high latitudes and the Amazon and Congo basins. The mechanisms behind the shifts toward drier regimes in the mid-high latitudes involve land processes that cannot be fully described by atmospheric-only indictors. The soil moisture-based assessment presented in this study facilitates a better understanding of climate change from the land process perspective and its direct impacts on terrestrial ecosystems.

This study highlights the importance of land processes in climate change impact assessments using a classification based on soil moisture, which would not be possible if only the atmospheric aridity index is used. The current analysis is based on one particular upper-end projection, the RCP8.5 scenario, which may not be absolutely realistic. Because of limitations from current land surface models and the lack of reliable long-term soil moisture observations, there must be unavoidable uncertainties in quantitative projections. With the involvement of land water balance and vegetation distribution, in principle, it is desirable to use soil moisture classification to define climate zones, in comparison with the atmospheric aridity index. The results are qualitatively valuable for ecosystem and environment impact assessments in the context of future climate change.