Research ArticleCLIMATOLOGY

Warming homogenizes apparent temperature sensitivity of ecosystem respiration

See allHide authors and affiliations

Science Advances  09 Apr 2021:
Vol. 7, no. 15, eabc7358
DOI: 10.1126/sciadv.abc7358


Warming-induced carbon loss through terrestrial ecosystem respiration (Re) is likely getting stronger in high latitudes and cold regions because of the more rapid warming and higher temperature sensitivity of Re (Q10). However, it is not known whether the spatial relationship between Q10 and temperature also holds temporally under a future warmer climate. Here, we analyzed apparent Q10 values derived from multiyear observations at 74 FLUXNET sites spanning diverse climates and biomes. We found warming-induced decline in Q10 is stronger at colder regions than other locations, which is consistent with a meta-analysis of 54 field warming experiments across the globe. We predict future warming will shrink the global variability of Q10 values to an average of 1.44 across the globe under a high emission trajectory (RCP 8.5) by the end of the century. Therefore, warming-induced carbon loss may be less than previously assumed because of Q10 homogenization in a warming world.


Terrestrial ecosystem respiration (Re), the sum of autotrophic respiration from primary producers and heterotrophic respiration from consumers and detritivores, represents a fundamental biospheric function and plays a major role in the global carbon cycle and growth rate of atmospheric CO2 (1). However, Earth system models currently in use poorly simulate the temporal variability of Re compared with in situ observations (2, 3), with model performance differing substantially among biomes and latitudes (47). Thus, future Re dynamics are poorly constrained and difficult to evaluate (711). Discrepancies between modeled and observed Re are often attributed to inaccuracies in the parameterized temperature sensitivity of ecosystem respiration (Q10) (9, 12). Actually, difference in parameterized Q10 is partially due to the misunderstanding of Q10 from field ecologists and mechanistic modelers (12). Early-generation models interpreted respiration as a temperature-dependent process in which Q10 was often set at a constant value, typically 2 (8, 13), i.e., Re doubles for a 10°C rise in temperature. However, recent studies have suggested that a constant Q10 is inadequate to predict Re because the processes involving plants and soil across time and space are not the same (11, 1416). In contrast, field observations focus more on the apparent Q10, the observable ecosystem response of Re to temperature that arises from the combination of multiple processes, each with their own temperature sensitivities. These studies have shown that Q10 decreases with increasing temperature and is thus higher in ecosystems at higher latitudes and in colder regions (17). Moreover, global warming rates are highly heterogeneous, with higher latitudes warming faster than lower latitudes. The joint consideration of faster warming and higher Q10 at high latitudes and cold regions implies a much stronger carbon efflux–climate warming feedback from these ecosystems (18). Although the spatial relationship between temperature and Q10 has been well established in terrestrial ecosystems, the temporal change in Q10 in response to ongoing global warming has not yet been studied at depth across the globe.

Knowledge of how Q10 will change with future climate warming can be investigated through two different pathways. One is through manipulative field warming experiments (2, 18), which provide direct, quantitative evidence of site-specific ecosystem responses to warming. Unfortunately, these experiments are hard to conduct at a relevant large scale and are therefore rare. The few existing experiments varied in experimental designs and were often applied to specific respiration components or biospheric components, (e.g., soil warming only and soil respiration only) at small temporal and spatial scales, resulting in some important biomes being neglected. These experiments also carried high uncertainty because of the inherently low degree of replication and to inconsistent warming protocols (19). As a result, it is very difficult to quantify the global-scale impact of climate warming on Q10 over a long time frame using aggregated data from diverse field warming experiments.

An alternative approach is through space-for-time substitution, an indirect approach that relates the contemporary spatial pattern of Q10 to the ambient, climatic spatial temperature gradient to infer the likely response of Q10 to climate warming (8). This approach has been widely used as a compromise when observational time series do not cover a large-enough temperature gradient in developing a temperature response function for measured Q10 (8, 20). However, the adequacy of space-for-time substitution was rarely tested (8). The advent and advancement of long-term carbon flux monitoring stations allow us to test this approach on estimation of Q10. In this study, we used standardized long-term eddy covariance (EC) observations from diverse climates and biomes worldwide to (i) identify Q10 changes and their underlying mechanisms along temporal and spatial temperature gradients, (ii) test the adequacy of space-for-time substitution in the analyses of Q10 patterns, and (iii) predict changes in spatial patterns of Q10 under future climate warming scenarios.


We extracted long-term time series of half-hourly Re and meteorological data from the global network of FLUXNET sites. All available half-hourly data were gathered, but sites with less than 5 years data were excluded. If more than 25% of 1-year data or if continuous observations for more than 1 month were missing from a site, that particular site-year was excluded from further analysis completely (detailed data filtering is described in Materials and Methods). As a result, a set of 74 sites (604 site-years; table S1) was used to derive Q10 values for each site-year (Q10 year, site) (Eqs. S1 and S2, Materials and Methods). We introduced priori tests for temperature and moisture effects on Q10 year, site. Our results showed that moisture did not significantly influence Q10 year, site when considering the temperature effects on Q10 year, site (fig. S1 and see Materials and Methods). Therefore, in this study, we chose to focus only on the effect of temperature on Q10.

We found that Q10 year, site was negatively correlated with annual mean air temperature at the site (Tyear, site) across almost all sites (slopes of dashed lines in Fig. 1A), implying lower Q10 year, site in warmer years for the majority (95%) of sites. To understand the long-term (years to decades) Q10 year, site response to temperature variation, we identified the temporal temperature sensitivity of Q10 year, sitetemporal), i.e., the change in Q10 year, site per 1°C change in Tyear, site, which we calculated from the regression slope between Q10 year, site and Tyear, site at individual sites (Materials and Methods and Fig. 1A). We found that δtemporal was significantly and positively correlated with annual mean site temperature (Tyear¯,site) (P < 0.0001) (red line in Fig. 1B), exhibiting larger warming-induced declines in Q10 year, site at colder sites than at warmer sites. Specifically, δtemporal at the colder arctic and boreal sites was more negative (average: −0.39; 95% confidence interval (CI): −0.51 to −0.26) than those at temperate (−0.11; CI: −0.14 to −0.08) and (sub-)tropical sites (−0.10; CI: −0.13 to −0.06) (table S2).

Fig. 1 Spatio-temporal patterns of Q10 and its temperature sensitivity against mean annual temperature through time and across space.

In (A), the points show annual mean Q10 values at individual sites (Q10year¯,site). The red dashed lines represent how annual Q10 year, site values vary temporally at individual sites. The blue solid line represents the exponential regressions of spatial Q10year¯,site against temperature. In (B), the points show site-specific δtemporal that are the linear slopes of the temporal regression of Q10 [red dashed lines in (A)] against temperature. The red line is an exponential fitting of δtemporal along temperature. The blue solid line is the first derivative curve of the contemporary spatial regression of Q10 in (A) (blue line) (δspatial). The shadow represents the difference between δspatial and δtemporal. Sites are grouped by ecosystem types according to the International Geosphere Biosphere Programme (IGBP) classification (Materials and Methods). The colored band shows the latitudinal position of the individual sites. Error bars are the standard error of the corresponding mean values.

Our analyses thus revealed strong support for a lower Q10 in a warmer world, with larger declines occurring in colder climate regions where Q10 is presently higher, and smaller reductions where Q10 is presently lower. As a result, the global variability of Q10 is likely to shrink in response to the ongoing global warming, especially given the stronger warming trend in the colder regions. Unfortunately, this analysis did not allow for extrapolation to a future, much warmer climate because the current interannual range of long-term mean annual temperature change at individual sites (typically about 1° to 3°C) did not encompass the projected values expected under future warming. However, the spatial gradient of Q10 did encompass a large enough temperature change range, and for this reason, a space-for-time substitution would provide more insight into how Q10 might change with pronounced future warming. Therefore, we investigated site-specific Q10 year¯,site, calculated as the arithmetic mean of all Q10 year, site values among all study years at individual sites.

Q10 year¯,site varied between 1.1 and 2.5, with an overall arithmetic mean of 1.94 (95% CI: 1.87 to 2.01) (Fig. 1A and tables S1 and S2). No significant differences in Q10 year¯,site were detected among ecosystem types (P > 0.1). However, significant differences in Q10 year¯,site were found among climatic regions (P < 0.0001) (table S2). Arctic and boreal regions generally had greater Q10 year¯,site values than temperate regions, which in turn exhibited higher Q10 year¯,site values than subtropical and tropical regions (table S2). Spatially, Q10 year¯,site significantly and exponentially decreased with increasing Tyear¯,site (R2 = 0.44, P < 0.0001; Fig. 1A), consistent with the negative temperature dependency of the temporal Q10 year, site at individual sites.

Based on this contemporary distribution of Q10 year¯,site across the spatial temperature gradient, we derived the spatial temperature sensitivity of Q10 year¯,sitespatial), i.e., the change in Q10 year¯,site per 1°C change in Tyear¯,site, from the derivative of the exponential relationship between Q10 year¯,site and Tyear¯,site for each site-specific annual mean temperature (Materials and Methods). Compared with δtemporal, δspatial showed a consistent temperature-dependent trend and was more negative at higher latitudes and in colder regions (Fig. 1B and table S2). However, the magnitude of δspatial was significantly lower than that of δtemporal in all climatic regions, especially in colder climates (Fig. 1B and table S2). Therefore, as applied here, the space-for-time substitution substantially underestimated the decline of Q10 year¯,site with climate warming, at least in the current-day period, and especially at high latitudes and in cold regions.

While this systematic difference between Q10 derived through time series and that derived through space-for-time substitution might be real, it could also be an artifact introduced by the use of different biomes or different temperature ranges during their respective fitting processes. We therefore first tested whether the mixing of biomes in the space-for-time substitution influenced the Q10 estimates. Specifically, we estimated Q10 on the basis of a series of specific temperature ranges (STRs) and tested how Q10 varied spatially while controlling for temperature (Materials and Methods). Results showed that the negative spatial correlation between Q10 year¯,site and site temperature (Fig. 1A) is no longer significant when the same temperature ranges are used for all sites (Fig. 2, A to J). However, the mean Q10 values at different specific temperature bins [Q10,str(i-j), 4°C ≤ i, j ≤ 23°C] are significantly correlated with average measurement temperature (P < 0.0001) (Fig. 2), which is consistent with the negative spatial correlation between Q10 year¯,site and site temperature (Fig. 1A). Thus, site-specific temperature was crucial in determining the spatial variation of Q10 year¯,site, while the geographical locations or biome type of the observational sites was not (Fig. 2, fig. S3, and table S2; Materials and Methods).

Fig. 2 The spatial trend of Q10 estimation in STRs (Q10,str(ij)¯) against site measured temperature gradient.

The inset graphs show Q10 estimation in STRs [Q10,str(i-j)] from 4° to 23°C by a 10°C moving window (for more temperature moving window in fig. S3). The central blue dots are the arithmetic mean of the Q10,str(i-j) values (Q10,str(ij)¯) with standard error (bars) for these 10 groups’ STRs (A to J), and the red line is an exponential fit for their relationship with temperature.

We calculated δtemporal using much smaller differences in mean annual temperature (generally <3°C) than were used in estimating δspatial (about 25°C across all sites, from −5° to 20°C; Fig. 1A), which may also have reduced the comparability of these indicators. To test whether the large difference between δtemporal and δspatial originated from temperature ranges, we additionally estimated δspatial values for the smaller temperature ranges equivalent to those used in the temporal estimate, reducing the spatial temperature change window from 25° to 2°C (in 1°C steps; Materials and Methods). Subsequently, we calculated the δspatial values for all sites across all temperature ranges from 25° to 2°C (Fig. 3A, figs. S4 and S5, and table S3). In our dataset, the interannual mean temperature changes were generally <3°C (Fig. 1A); thus, we grouped δtemporal of all sites into three groups of approximately 1°, 2°, and 3°C interannual temperature changes (Materials and Methods). We then tested the agreement between δtemporal and δspatial values obtained from the same temperature range and the same sites. We found no significant differences between δtemporal and δspatial when their associated temperature changes were similar (figs. S4 and S5). These results suggest that the space-for-time substitution required to project Q10 changes under future climate change is valid, but only if the temperature changes through time and across space are at comparable magnitude.

Fig. 3 Potential convergence in temperature sensitivity of ecosystem respiration under climate warming.

In (A), all lines show the spatial patterns in temperature sensitivity of Q10 (δ) against site temperature across different temperature change ranges: blue dashed lines are the spatial δ (δspatial) (fig. S4 and table S3, Materials and Methods), and red solid lines are temporal δ when the interannual temperature change ranges are approximately 1°C (δtemporal,1°C), 2°C (δtemporal, 2°C), and 3°C (δtemporal, 3°C), respectively (fig. S5 and table S3, Materials and Methods). The primary model is shown to directly calculate δ at different site temperature (T, °C) and temperature change ranges (∆T,°C) (fig. S6 and table S3). In (B), the red solid line shows the contemporary Q10year¯,site spatial patterns identical to the blue line in Fig. 1A, and the blue dashed lines are predictions of spatial Q10year¯,site patterns under different Representative Concentration Pathway scenarios (RCPs) until the end of the 21st century (2081–2100) based on the Coupled Model Intercomparison Project (CMIP5) (Materials and Methods). All site-scale predications are shown as dots in fig. S7 (F to I).

We further analyzed the change in Q10 year¯,site under future warming conditions using space-for-time substitution at comparable temporal and spatial temperature change magnitudes. We extracted site-specific outputs of warming projections obtained from the Coupled Model Inter-comparison Project (CMIP5) under all Representative Concentration Pathway scenarios (RCPs; 2.6, 4.5, 6.0, and 8.5) until the end of the 21st century (2081–2100) (21) as nonequal warming scenarios (table S1; see Materials and Methods also for equal warming scenarios across sites). Using space-for-time substitution and the site-specific nonequal warming scenarios from RCP 2.6 to RCP 8.5, respectively (table S3 and fig. S6, Eq. S5 in Materials and Methods), we found that the curves of the best-fitting exponential regression for Q10 year¯,site and Tyear¯,site became increasingly flatter as climate warming intensified (Fig. 3B). This result is consistent with homogeneous warming scenarios (fig. S7, A to E). Subsequently, Q10 year¯,site gradually converged to a global average of 1.44 (95% CI: 1.38 to 1.50) for the RCP 8.5 warming scenario (Fig. 3B and fig. S7, F to I), a value over 25% lower and less spatially variable than the current global average of 1.94. In addition, the magnitude of δspatial consistently decreased with rising temperature, especially in colder regions (Fig. 3A), suggesting a progressive decline in the spatial difference in δspatial under ongoing warming. Clearly, future warming is likely to homogenize global Q10 year¯,site.


Several mutually nonexclusive mechanisms have been postulated to explain the lower Q10 values at higher temperatures. Manipulative field warming experiments (2, 17, 22), incubation experiments (14), bio-kinetics theory (8), and modeling studies (23) have suggested that Q10 is likely to decline with increasing temperature because of ecosystem thermal acclimation, either through biochemical (e.g., declining increase in molecules that have sufficient energy to overcome enzyme activation energies at increasing temperatures) or through structural changes (e.g., in community composition or reduced microbial biomass) (2224). Other studies have provided evidence for increasing limitations on substrate availability for decomposers as an important constraining factor (25, 26). It is likely that multiple mechanisms work in tandem at varying strengths (8, 2628). In this study, the relationship between Q10 year¯,site and Tyear¯,site was no longer significant when we grouped spatial sites with the same measured temperature ranges (Fig. 2, A to J, and fig. S3), suggesting that spatial variation in Q10 year¯,site was likely caused by intrinsic factors such as the composition, quantity, and metabolic rates of microbes and plants, as well as inherent properties of enzymes under certain temperature levels (23, 2931). Moreover, the positive temperature-dependent trends observed in both δtemporal and δspatial suggest that the temperature response of these intrinsic factors slowed with rising temperature. A meta-analysis of field warming experiments conducted during the past two decades (1994–2017), consistent with our result, also showed similar positive exponential temperature-dependent trends in δtemporal (fig. S8, Materials and Methods). In addition, to test the compounding influences of the seasonal variability, especially the transition periods between growing season and non–growing season, we have reanalyzed our data after limiting to growing season for all sites. We also found that the results from these two methods showed similar trends in Q10 spatial-temporal variations (fig. S9), although the magnitude is different because of different available substrates (25, 26) and temperature change ranges at different time scales.

Our findings, based on global field observational data, provide fresh insights on how the temperature sensitivity of Q10 has changed in recent time (δtemporal) and across space (δspatial), as well as for the future warmer world. These results provide evidence for the potential convergence of global Q10 under ongoing climate warming, a warming response that is currently not reproduced by most Earth system models (six of nine CMIP5 models in Fig. 4 and fig. S10; Materials and Methods). Therefore, warming-induced carbon loss through Re at high-latitude regions may be less than assumptions in Earth system models because these projections underestimate the strong downward-adjustment of Q10 under climate warming, especially at higher latitudes and in colder climates. Thus, this study calls for a careful reassessment of the predicted large carbon loss, especially at higher latitudes and in colder regions.

Fig. 4 Spatial patterns of Q10 and its temperature sensitivity against site mean annual temperature based on the nine model outputs from CMIP5 at the flux tower site.

Solid curves represent those models considering the different temporal variations in Q10 year, site under climate warming, while the dashed curves are not. The shadows are 95% confidence bands. Detailed regression and CMIP5 models’ information are in the table S5, and the points are shown in fig. S10.


EC-based data

For this study, we defined terrestrial ecosystem respiration (Re) as nighttime net CO2 fluxes. The major dataset used in this study was the database of EC observations from a worldwide network of FLUXNET sites (32). All available half-hourly CO2 flux and meteorological data from the FLUXNET ( were gathered, but sites with less than 5 years’ data were excluded (called Data Set 0 comprising 929 site-years of 113 sites). Gap-filled data were excluded from the selected dataset to avoid spurious effects, and conventional quality control criteria were previously applied to the standard data files (27, 33). If more than 25% of 1-year data or if continuous observations for more than 1 month were missing from a site, that particular site-year was excluded from the analysis completely. If less than five such full years of data were available for a site, we also excluded that site from the analysis. This data filtering was to ensure the continuity and reliability of the data used for further analyses. Air temperature data were extracted from the same databases to match the corresponding time of Re data. Here, we used air temperature (T) rather than soil temperature because both soil and vegetation respiration in Re responded to temperature (34).

Q10 calculation

Results from the statistical analysis (as following the “Statistical analysis” section) showed that temperature dominated seasonal variation of ecosystem respiration (table S1). Several models have been proposed to describe the temperature dependence of respiration (35, 36). Here, we used a Q10-based model (37) (S1), which describes the proportional increase in the rate of respiration per 10°C rise in temperature (38), derived from the Van’t Hoff exponential chemical reaction–temperature equation (39)Re=a×eb×T(S1)where Re and T are half-hour nighttime (local time is later than 20 and earlier than 8) net CO2 fluxes (ecosystem respiration) and corresponding air temperature in a particular year, respectively. a and b are the regression parameters that can be used to estimate Q10 values in a particular site-year from the following formula (37, 40).Q10 year,site=e10×b(S2)

However, recent studies have also revealed that variation in temperature accounts for most of the variation in respiration in a particular site only in the absence of water stress (38, 41). In other words, soil water content or precipitation deficits may also affect the temperature dependence of terrestrial respiration (17, 4244) and can even be the main controlling factor of terrestrial respiration under drier conditions (45, 46). Therefore, we chose only sites for which the Van’t Hoff equation showed a significant fit to the data (P < 0.05), which ensured that the influences from other factors (e.g., moisture or artificial management) on respiration were small in comparison to that of temperature. Last, 74 sites were extracted from Data Set 0 (called Data Set 1) and used to estimate the temperature sensitivity for a particular site-year (Q10 year, site) (604 site-years in total; see table S1). All available nighttime respiration data from within a single year were used to estimate Q10 year, site. The arithmetic means of annual Q10 year, site and Tyear, site were used to estimate the spatial pattern of Q10 (Q10 year¯,site).

Furthermore, to inspect water restrictions, we tested the relationship between Q10 year, site and Tyear, site when binned according to the annual drought index (Dryear,site), the ratio of annual evapotranspiration (ETyear,site), and annual precipitation (PPTyear, site). Also, we tested the relationship between Q10 year, site and Dryear,site when binned by Tyear, site. The relationship between Q10 year, site and Tyear, site are significant in all Dryear,site groups (fig. S1A). However, we did not find any significant relationship between Q10 year, site and Dryear,site (fig. S1B). In other words, the Dryear,site did not change the response of Q10 year, site to temperature despite the fact that a small number of site-years (15.6% of all site-years) have a higher Dryear,site (>1) (red line in fig. S1A). Thus, in this study, we focused on revealing the temperature effect on Q10. Here, we did not use the soil water content because of the large amount of missing data and the inconsistent depth of field observations across the globe.

Temperature sensitivity of Q10 (δ)

We defined a new parameter δ, the rate of Q10 change per 1°C increase in mean annual temperature. This parameter is equivalent to the slope of the simple linear regression between Q10 year, site and Tyear, site from long-term (>5 years per site) EC observational data (S3).Q10 year,site=δtemporal×Tyear,site+εsite(S3)where δtemporal represents the temporal temperature sensitivity of Q10 for an individual site. Q10 year, site is estimated for all available years at individual sites based on Eqs. S1 and S2 , Tyear, site is the annual mean temperature, and εsite is a site-specific regression parameter as the Q10 year, site when the Tyear, site is 0°C.

However, this data-driven approach does not allow the determination of δtemporal under future global warming conditions, because the current-day temperature ranges at the sites typically do not fully encompass the projected future temperature ranges (47). To overcome this issue, space-for-time substitution provides an alternative, widely used approach to hindcast past and forecast future trajectories of ecosystem from the contemporary spatial patterns (8, 48, 49). However, the validity of this space-for-time substitution under future climate scenarios requires verification (8). Therefore, we also calculated the change rate of Q10 per 1°C change in temperature across the spatial gradient. Like δtemporal, we designated the spatial temperature sensitivity of Q10spatial). By definition, δspatial can be calculated as the first derivative of the exponential, across-site relationship between Q10 year¯,site and Tyear¯,site at each site-specific annual mean temperature (Eq. S4 and Fig. 1A).δspatial=dQ10year¯,sitedTyear¯,site(S4)

Here, we did not apply a linear regression as was done for δtemporal because the spatial temperature change range was much larger than the within-site interannual temperature fluctuations. This, combined with the much larger number of data, revealed a better fit to the data when applying an exponential relation. Therefore, we chose to use the first derivative values to estimate δspatial.

Q10 estimates in STRs

We used STR Q10 estimates to test whether the difference between the spatial and temporal Q10 variations was explained solely by the differences in temperature or possibly attributable to geographical factors. First, we investigated a common temperature range that covered as many site-years in Data Set 1 as possible. After analyzing the annual mean temperature and the maximum and minimum temperatures of the 74 sites, we found a common temperature range (4° to 23°C) across all site-years (about 90% of all site-years) (fig. S2). Second, we created 10 groups of STRs from 4° to 23°C using a 10°C moving window to calculate Q10 values in each specific temperature bin [Q10,str(i-j), 4°C ≤ i, j ≤ 23°C] based on all available years’ date at individual sites. The algorithm for calculating Q10,str(i-j) is identical to the conventional Q10 calculation method (Eqs. S1 and S2 ). Here, we selected 10°C for the moving window in the estimation of Q10,str(i-j) because by definition, Q10 is the proportional change in Re per 10°C rise in temperature. This way, we ensured comparability with most of the other sensitivity analyses of ecosystem respiration. In addition, to incorporate the likely effects of temperature moving windows on Q10,str(i-j) variation, we also selected a series of temperature moving windows with the same mean measurement temperature to estimate Q10,str(i-j) (fig. S3). Third, we investigated the relationships between Q10,str(i-j) and site mean temperature in different temperature windows across all site-years (Fig. 2, A to J). Moreover, Q10,str(i-j) values were respectively aggregated by their arithmetic mean (Q10,str(ij)¯) in different temperature ranges and temperature moving windows (Fig. 2 and fig. S3). Last, we tested the differences of mean Q10,str(ij)¯ at the same mean measurement temperature among these series of temperature moving windows (fig. S3).

Comparison of δtemporal and δspatial across different temperature change gradients

The fact that δspatial was calculated from a large temperature variation (about 25°C in our Data Set 1), while δtemporal only over a much smaller interannual range in mean annual temperature (generally 1° to 3°C) (Fig. 1A), reduces the comparability of both indicators for the temperature dependence of Q10. Therefore, we want to test the agreement between δspatial and δtemporal at both the same temperature and temperature change range. First, we calculated the spatial patterns of Q10 year¯,site for temperature change ranges using moving windows from 25°, 20°, and down to 2°C (by 1°C steps, 20 groups in total). As the spatial temperature gradient across all sites is about 25°C in Data Set 1, we extracted multiple exponential spatial patterns of Q10 year¯,site along each moving temperature gradient window from 24° to 2°C. For example, for a moving window of 20°C, we derived the exponential spatial patterns of Q10 year¯,site against site temperature changes from −5° to 15°C, but also from 0° to 20°C (dots in Fig. 1A). All these patterns were included to ensure robust results across the globe.

Second, using Eq. S4, we calculated δspatial values across different spatial temperature change windows from 25° to 2°C, respectively (black dots in fig. S4). The optimal exponential curves of these δspatial values were designated as δspatial, specific temperature windows (e.g., δspatial, 20°C and δspatial, 2°C) (red lines in fig. S4 and table S3). However, the spatial pattern of Q10 year¯,site across 1°C temperature windows showed a great variation, and it was very difficult to obtain a robust value of δspatial, 1°C. Furthermore, we grouped the δtemporal of all observational sites according to the interannual temperature changes, the difference between the maximum and minimum annual mean temperature of one site, into approximately 1°C (0.5° to 1.5°C; blue dots in fig. S5), 2°C (1.5° to 2.5°C; red dots in fig. S4T), and 3°C (2.5° to 3.5°C; red dots in fig. S4S), respectively. Last, we inspected the agreement between δspatial and δtemporal across similar temperature change window from 2° to 3°C, respectively. If the space-for-time substitution is valid (i.e., it does not show a significant difference between δspatial and δtemporal), then it would imply that we can use the δspatial across different temperature gradients to predict temporal variations of Q10 under climate warming.

Q10 prediction under climate warming

On the basis of contemporary spatial Q10 year¯,site pattern (Fig. 1A) and site-specific δtemporal values, we can directly predict a global spatial Q10 year¯,site pattern under ongoing 1°C equilibrium warming scenario (Eq. S5). However, only quantifying global Q10 year¯,site pattern under a 1°C warming scenario is not explicit enough to illuminate the further Q10 year¯,site patterns under the projected climate warming (50, 51). Thus, we further examined what the global spatial Q10 year¯,site pattern would be under future warming scenarios in 1°C steps based on the δspatial, i °C derived from space-for-time substitution (Eq. S5)Q10year¯,site,+i°C=Q10year¯,site,+(i1)°C+{δtemporal(i=1)δspatial,i°C(i>1)(i=1,2,,n,nZ+)(S5)where i is the projected warming scenario, and Q10year¯,site,+i°C is the projected spatial Q10 year¯,site under the i °C warming scenario. Here, δtemporal is the temporal temperature sensitivity of Q10 year, site, and δspatial, i °C is the spatial temperature sensitivity of Q10 year¯,site with an i °C temperature change window. For example, if one site, with T °C site temperature, is projected to have a 3°C warming scenario (i = 3) and current Q10 is Q10 year¯,site, then we will use the next three steps to estimate the Q10 year¯,site variation: (i) We need to calculate the δtemporal of this site (T °C temperature inputs, δtemporal, site) according to δtemporal variation against site temperature (Fig. 1B and table S3). The sum of δtemporal, site and Q10 year¯,site is the Q10 estimation of this site under 1°C warming scenario (Q10year¯,site,+1°C). (ii) Similar to the first step, we need to check the δspatial, 2°C variation against site temperature to estimate δspatial at this site (note here T + 1°C temperature inputs, δspatial, 2°C, site) (table S3). The sum of Q10year¯,site,+1°C and δspatial, 2°C, site is the Q10 estimation of this site under 2°C warming scenario (Q10year¯,site,+2°C). (iii) In the same way, we can calculate δspatial, 3°C, site at T + 2°C temperature inputs and Q10year¯,site,+3°C for this site under 3°C warming scenario. In addition, we also provided a general model for parameters of δ patterns (a and b in δ = a × eb × T) at different temperature change ranges (Fig. 3A, table S3, and fig. S6). Together with Eq. S5, we can directly estimate Q10 variation response to warming for a certain site.

In this study, we used two kinds of warming scenarios to describe future Q10 year¯,site patterns. One was the spatial equilibrium warming scenario by 1°C warming step, which provided an intuitive result in how Q10 year¯,site is likely to change under future warming conditions from 1° to 5°C (fig. S7, A to E). The other was the site-specific warming scenario, a nonequal warming trend across sites, evaluated by the CMIP5 under all RCP scenarios (2.6, 4.5, 6.0, and 8.5) until the end of the 21st century (2081–2100) based on the Intergovernmental Panel on Climate Change (IPCC) Fifth Assessment Report (50) (fig. S7, F to I, and table S2). The method of Q10 year¯,site projections is identical to the spatial-equilibrium warming scenario, as Eq. S5 shows.

Statistical analysis

To explore the main driving factors of seasonal Re, we developed correlation analysis (COR) and partial correlation analysis (PCOR) between ecosystem respiration (Re) and climate factors at a particular site-year, including temperature (T), vapor pressure deficit (VPD; considering the data missing, VPD could be replaced by soil water content or relative humidity), and PPT (table S1). Both COR (table S1) and PCOR (similar to COR, thus are not included here) showed a greater correlation of temperature at most sites (68 of 72 sites, the remaining 2 sites do not have enough VPD data for COR and PCOR), which indicated that temperature dominated the seasonal variation of ecosystem respiration, not VPD or PPT. To examine the role of ecosystem types or climatic regions in determining the patterns of Q10 year¯,site and δ, we grouped all sites into different components using two schemes. One classified sites by ecosystem types, including croplands (CRO), forests [FOR; including evergreen needle leaf forests (ENF), evergreen broadleaf forests (EBF), mixed deciduous and evergreen forests (MF), and deciduous broadleaf forests (DBF)], shrubs [SHR; including closed shrublands (CSH), open shrub (OSH), and woody savannas (WSA)], grasslands (GRA), and wetlands (WET) (table S1). The other classification grouped sites into three climatic regions: arctic and boreal, temperate, and (sub-)tropical (table S2). These data passed the normality test (Shapiro-Wilk test) and homogeneity of variance test (Bartlett test) (P > 0.05). Next, we used a one-way analysis of variance and Tukey’s post hoc to test for multiple comparisons across different climate biomes or among different ecosystem types to investigate the spatial variability of Q10 year¯,site(δ) at α = 0.05 (table S2). From the results given in table S2, we discovered that only the climatic regions, not the ecosystem types, showed significant Q10 year¯,site(and δ) spatial differences. Therefore, we can ignore the small confounding effects of ecosystem types on Q10 year¯,site(and δ) spatial variation with temperature fluctuations in our next analysis (52). In addition, we tested the differences between the δtemporal and δspatial using paired t test at α = 0.05 among the groups as mentioned above. All statistical and modeling procedures were performed using the R statistical computing packages (Version 3.5.0) (


We searched peer-reviewed journal articles including “experimental warming, temperature sensitivity (Q10), ecosystem respiration, or soil respiration” using Web of Science in the past two decades (1994–2017). Then, we only selected the studies that (i) manipulated warming in the field, (ii) had clear control and experimental groups, (iii) had calculated Q10 for two groups or for which it could be calculated from complete observations of respiration and temperature, and (iv) included at least one whole growing season. Last, we extracted 54 field warming studies across the globe in our meta-analysis (table S4). We synthesized the Q10 values and temperature measurements for the experimental group (warming: Q10, W and TW) and the control group (nonwarming: Q10, C and TC) at individual sites. Using the similar method with EC-based δtemporal calculation (Eq. S3), we also calculated the δtemporal based on these data from field-warming experiments (WE-based δtemporal) (Eq. S6).δtemporal=Q10,WQ10,CTWTC(S6)

Then, we compared the WE-based δtemporal with EC-based δtemporal along a similar temperature gradient.

CMIP5 models comparison

Using the 74 flux tower positions (table S1) as the center pixel, we extracted nine model outputs, including autotrophic respiration (Ra), heterotrophic respiration (Rh), and temperature, from CMIP5 to estimate their Q10 year, site and trend in time (δtemporal) over two decades (1990–2012). The detailed model information is shown in the table S5, and the estimation methods of Q10 year, site and δtemporal are identical to Eqs. S1 to S3 , in which the Re was regarded as the sum of Ra and Rh. Overall, this analysis showed that most Earth system models (six of nine, in this study) ignored variation of δtemporal under climate warming despite considering the spatial variation of Q10 year¯,site, except for CCSM4, BCC-CSM, and INM-CM4, whose temporal trends of Q10 year, site across sites were similar to the EC-based result (Figs. 1B and 4 and fig. S10). Thus, although Q10 input in these models is typically insensitive to temperature and generally well approximated (Q10 = 2) according to specific documentations, some Earth system models also take into account the temperature-dependent Q10 using some specific functions and parameters. For example, the Carbon Exchange between Vegetation, Soil, and Atmosphere (CEVSA) model in the BCC-CSM introduced a temperature-dependent piecewise function to adjust the Q10 inputs (5355). Our results suggested that a series of temperature-dependent decline rate in Q10temporal) across the globe respond to different climate warming scenarios (Fig. 3 and table S3), which can be applicable to adjust the Q10 inputs for many Earth system models in future carbon efflux estimation. However, future warming may enhance water stress in arid regions, which is likely to affect the climate-carbon feedback with more complex situations and make it more difficult to extrapolate to future climates (56). Therefore, future continuous water scaling researches, especially in the arid and semiarid regions, are still needed for an adequate mechanism explanation in carbon-climate feedback.


Supplementary material for this article is available at

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We thank the FLUXNET dataset, the ESGF Search Portal, and the KNMI Climate Explorer for providing model outputs of warming scenarios. Funding: This work was financially supported by the Second Tibetan Plateau Scientific Expedition and Research Program (STEP) (no. 2019QZKK1002); the Strategic Priority Research Program of the Chinese Academy of Sciences (nos. XDA20050102 and XDA19070303), the National Natural Science Foundation of China (no. 41807331), and the National Key Research and Development Program (nos. 2016YFC0502001 and 2017YFA0604801). I.A.J. is supported by the European Research Council Synergy (no. ERC-2013-SyG-610028 IMBALANCE-P). A.R.D. acknowledges support from the U.S. Department of Energy Ameriflux Network Management Project. I.M. acknowledges funding from the Academy of Finland Flagship (grant no. 337549) and ICOS-Finland by the University of Helsinki. Author contributions: X. Zhang, B.N., I.A.J., and S.P. conceived the method. B.N. and X. Zhang. synthesized data, performed the analysis, and prepared figures. B.N., X. Zhang, S.P., and I.A.J. drafted the manuscript. All authors provided data for the analysis, discussed the results and implications, and commented on the manuscript at all stages. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Specifically, the EC-based measurements that support the findings of this study are available from the FLUXNET dataset ( The site-specific warming scenarios that support the findings of this study are also publicly available at the KNMI Climate Explorer ( The CMIP5 model outputs are available at the CEDA ESGF Search Portal ( Additional data related to this paper may be requested from the authors (niub{at} and zhangxz{at}

Stay Connected to Science Advances

Navigate This Article