Skip to main content

Dramatic increase in water use efficiency with cumulative forest disturbance at the large forested watershed scale



Forest disturbance induced changes in the coupling of forest carbon and water have important implications for ecosystem functioning and sustainable forest management. However, this is rarely investigated at the large watershed scale with cumulative forest disturbance. We used a combination of techniques including modeling, statistical analysis, and machine learning to investigate the effects of cumulative forest disturbance on water use efficiency (WUE, a proxy for carbon and water coupling) in the 19,200 km2 Chilcotin watershed situated in the central interior of British Columbia, Canada. Harvesting, wildfire, and a severe Mountain Pine Beetle (MPB) infestation have gradually cumulated over the 45-year study period, and the watershed reached a cumulative equivalent clear-cut area of 10% in 1999 and then 40% in 2016.


Surprisingly, with the dramatic forest disturbance increase from 2000 to 2016 which was mainly due to MPB, watershed-level carbon stocks and sequestration showed an insignificant reduction. This resilience was mainly due to landscape-level carbon dynamics that saw a balance between a variety of disturbance rates and types, an accumulation of older stand types, and fast growing young regenerated forests. Watershed-level carbon sequestration capacity was sustained, measured by Net Primary Production (NPP). A concurrent significant decrease in annual evapotranspiration (ET), led to a 19% increase in WUE (defined as the ratio of NPP to ET), which is contrary to common findings after disturbance at the forest stand-level. During this period of high disturbance, ET was the dominant driver of the WUE increase.


We conclude that disturbance-driven forest dynamics and the appropriate scale must be considered when investigating carbon and water relationship. In contrast to the stand-level trade-off relationship between carbon and water, forested watersheds may be managed to maintain timber, carbon and water resources across large landscapes.


People throughout the world rely on water supplied from forested watersheds [1], as well as other ecosystem services and values such as timber production or carbon sequestration [2, 3]. Forest carbon and water are coupled through the photosynthetic process and at the ecosystem level through carbon production and associated water consumption. This relationship is frequently altered by both natural and anthropogenic forest disturbance, which has a range of important implications for ecosystem structure, functions and services [4,5,6,7,8,9].

Water use efficiency (WUE), the ratio of carbon uptake to water use, has been commonly used as a proxy to study the relationship between forest carbon and water in forest ecosystems. WUE is used to evaluate the potential impacts of climate change on food production [10, 11], water supply [7], and forest or land use management [12,13,14]. WUE has been studied from the leaf to global scale, so a variety of definitions are found in the literature. Studies looking at leaf-level WUE concentrate on measuring the ratio of net CO2 assimilation to stomatal conductance at sub-daily time scales [15]. Tree-level WUE can be measured through isotope analysis, sap-flow measurements [16], or controlled environment chamber experiments [17] at time scales from daily to seasonal or annual [18]. Stand-level WUE includes the additional effects of other vegetation and evaporation, and water and carbon exchange is often measured using eddy-covariance systems [19, 20]. Watershed or landscape-level studies have involved the use of modeling or analysis based on the water budget [21, 22]. Regional to global studies of WUE often rely on modeling or remote sensing [13, 14, 23], comparing trends across global gradients or investigating the impact of climate change. Among the existing studies on WUE at various spatial scales, large watershed- or landscape-level studies are rare, likely due to the difficulty of conducting field measurements [24].

Disturbance can change stand-level WUE [25]. However, studies are relatively limited with the body of research focusing on how drought or climate change affects WUE [12, 26, 27]. The majority of research into WUE after forest disturbance is at the stand-level using eddy covariance data. These studies found that severe forest disturbance often results in a decrease in WUE, most likely because of greater relative surface evaporation [28], while lower intensity disturbance such as insect attack or prescribed burning can also produce a reduction or no change in WUE [25, 29, 30]. There is an emerging view that the effect of forest disturbance on WUE is variable, depending on tree mortality, site conditions, remaining vegetation recovery, and time since disturbance. In addition to the site-level nature of most WUE measurements, another limitation is that study periods usually cover just the few years following disturbance. In large watersheds or landscapes, forest disturbances of different types cumulate over space and long periods of time, and consequently their effects on forest carbon, water and WUE may be different. To our knowledge there are no measurement-based studies and few modeling-based studies that have examined the effects of long term cumulative forest disturbance on WUE at the large watershed-level.

To advance this topic in the literature, we investigated the forest carbon–water relationship in the Chilcotin watershed situated in the interior of British Columbia, Canada, a large forested watershed that experienced high level of cumulative forest disturbance over the last 45 years. We hypothesise that cumulative forest disturbance is one of the key drivers of carbon, water and WUE. We expect that WUE decreased under significant cumulative forest disturbance in this large forested watershed, in a similar manner found at the forest stand or small watershed scale.


Study area

The Chilcotin watershed has a drainage area of 19,200 km2 and is located 50 km west of Williams Lake on the Fraser plateau in British Columbia, Canada (Fig. 1). The topography is mountainous, ranging from 2,800 m above sea level (m) on the south western-end of the watershed in the coastal mountains down to 500 m in the lower reaches in the east. The large elevation gradient means ecosystems range from non-treed dry Bunchgrass Very Dry Warm Alkali variant (BGxw2) biogeoclimatic (BEC) zone in the valley bottoms up to Coastal Mountain-heather Alpine (CMAunp) BEC zone characterised by alpine tundra [31]. The majority of the watershed is forested in the Interior Douglas-fir (IDF), Sub-boreal Pine–Spruce (SBPS), and Montane Spruce (MS) BEC zones. These are dominated by Pinus contorta Douglas ex Loudon (Lodgepole pine) (78%) with significant components of Pseudotsuga menziesii (Mirb.) Franco (Douglas fir) (9%) and Picea engelmannii Parry ex Engelm (White spruce) (8%) forests, and minor components of other species.

Fig. 1
figure 1

Location of the Chilcotin watershed in the central interior of British Columbia, Canada

Wildfires occur regularly in the study area and ecosystems are adapted to frequent stand-initiating and stand-maintaining fires. The provincial land has supported sustained levels of clear-cut logging since the 1960s. The Mountain Pine Beetle (MPB) is endemic throughout the region and a widespread severe infestation caused the mortality of a large proportion of mature Lodgepole pine trees in British Columbia in the early 2000s.

The climate is continental with an average winter temperature of − 8 °C and 11 °C in the summer. Annual precipitation (P) averages 635 mm which is distributed throughout the year, with a slightly lower proportion in the spring and summer months (19 and 21% respectively). Winter P falls as snow in all but the lowest valley bottoms, leading to a pronounced freshet with the spring melt in May or June.


Climate and streamflow data

Annual timeseries of water, climate, carbon, and forest disturbance data were constructed based on streamflow data availability from 1971 to 2016. The annual climate variables, precipitation (in mm), annual mean, minimum and maximum temperature (Tmean, Tmin, and Tmax) were calculated using the ClimateBC dataset [32] based on a 25 m resolution digital elevation model [32, 33] and were averaged to obtain one estimate for the watershed. Potential evapotranspiration (PET), based on the ClimateBC temperature data was calculated using the Hargreaves method [34] (Fig. 2). Daily discharge was acquired from Water Survey of Canada (station ID: 08MB005) [35], as well as the corresponding drainage area [36]. The average of all daily flows in m3 sec−1 was used to calculate mean annual streamflow (Q), standardized to millimeters per year (mm year−1) based on watershed area. Runoff ratio was calculated as annual Q/P.

Fig. 2
figure 2

Annual climate variables in the Chilcotin watershed from 1971–2016. Where, a minimum (Tmin), mean (Tmean), and maximum (Tmax) daily annual temperature (c), b annual precipitation in millimeters (mm), and c potential evapotranspiration calculated using the Hargreaves equation

Forest disturbance

Forest disturbance information was sourced from two spatial layers maintained by the provincial government. The Vegetation Resources Inventory (VRI) (version 2019) contains the type, year and severity of the most recent logging and non-logging disturbance of each polygon [37]. Additionally, the provincial Consolidated Cutblocks (version 2019) layer was used to capture any recent logging that may have not been updated in the 2019 VRI [38].

Cumulative equivalent clear-cut area (CECA) was chosen as the watershed-level indicator of forest disturbance [39, 40]. CECA represents the cumulative effect of multiple distinct types of disturbance and recovery over time as vegetation re-grows. CECA has been widely used in British Columbia where winter snowpack and the effects of vegetation alteration on peak flows are a major concern [41,42,43,44]. Stand-level annual equivalent clear-cut area (ECA) is calculated as disturbed area multiplied by a coefficient, which ranges from 0 to 100%. The ECA coefficient is set at 100% after harvesting and wildfire, to reflect changes in hydrological processes such as infiltration and evapotranspiration (ET) and is reduced as the forest recovers [45]. After MPB mortality, an ECA coefficient is applied gradually as tree death and needle drop occur progressively over many years [46, 47]. In this study, hydrological recovery starts around 30 years after MPB mortality (see Additional file 1: Figure S1). In the years after disturbance, each polygon’s ECA score within the watershed was calculated based on the type of disturbance and years since disturbance. All ECA values in the watershed were summed and divided by the gross watershed area to calculate the total percent CECA for the watershed (Fig. 3). The Chilcotin watershed was essentially undisturbed at the start of the study period in 1971 with a CECA of 1%. Low levels of harvesting (an average of ~ 4,500 ha or 0.2% per year) dominated the disturbance profile until 2000 when the CECA reached 10%. From 2000 to 2016, forest mortality caused the CECA to increase to 39.6%, primarily from the MPB epidemic which affected just over 700,000 ha or 36% of the watershed to varying degrees (Fig. 3). Based on this disturbance timeline, we divided the study into two periods—from 1971 to 1999 which was characterised by a low rate of harvesting and a CECA <  = 10%, and from 2000 to 2016 where MPB, wildfire, and harvesting occurred at a higher rate (CECA from > 10–40%).

Fig. 3
figure 3

Chilcotin watershed Cumulative Equivalent Clear-Cut Area (CECA) total annual area disturbed from 1971–2016. Where, CECA is colored by disturbance type (left axis), and grey bars represent total annual area disturbed in hectares (ha) from all types of disturbances (right axis)


Estimation of evapotranspiration

We used the Budyko–based Zhang’s equation (Eq. 1) to calculate annual watershed-level ET [48,49,50]. Zhang’s equation contains the w parameter that reflects the forest’s ability to transpire water at the watershed scale. We varied w to represent the change in ET with forest disturbance. To obtain independent estimates of w, we developed a simple linear regression equation between w and CECA from watersheds with similar topography and climate in the interior of British Columbia. Groups of five to 10 years with low, moderate, and high CECA values were chosen in the Chilcotin, Chilko, Baker, and Moffat watersheds. Annual P, PET, and Q were based on data ClimateBC [32] and Water Survey of Canada [35], and were calculated consistent with the description above. ET was calculated using the water balance equation (Eq. 2), where the change in storage over these five to 10 year grouped periods was assumed to be zero, thereby minimising the error associated with annual variations in storage. We then used the ‘uniroot’ function from the ‘rootSolve’ package [51] in R [52] to solve for w. Finally, using data from all four watersheds, we developed a linear equation to calculate w based on CECA (Additional file 1: Figure S3). This equation was used to calculate annual w for the Chilcotin watershed, which was then used along with annual P and PET in Eq. 1 to calculate annual ET in the Chilcotin watershed. We used this methodology to conceptually represent the change in ET from forest disturbance, however it still dependent on the accuracy of estimating precipitation over large areas. The water balance method showed similar trends, as did using three and five year averages, which gave us more confidence in the ET estimate calculated using Zhang’s equation (Additional file 1: Figure S4, Table S4 and Table S5). In this way, watershed level data has been used to scale point estimates of P partitioning to ET to the watershed [53].

$$\frac{ET}{P}=\frac{1+w\frac{PET}{P} }{1+w\frac{PET}{P}+{\left(\frac{PET}{P}\right)}^{-1}}$$

where, ET = evapotranspiration in mm, P = precipitation in mm, w = plant-available water coefficient (unit-less), and PET = potential evapotranspiration in mm, at the annual timescale.

$$P=Q+ET+\Delta S$$

where, P = annual precipitation, Q = mean annual streamflow, ET = annual evapotranspiration, and ∆S = annual change in storage, which was assumed to be zero.

Carbon modeling and validation

The Carbon Budget Model of the Canadian Forest Service version 3 (CBM-CFS3) [54, 55] was used to simulate annual historical carbon dynamics in the Chilcotin watershed. CBM-CFS3 is a yield curve based model widely used in Canada to simulate terrestrial carbon pools and their fluxes with the atmosphere [56, 57]. Forest information including age, disturbance history, species composition, and the productivity of the forest from the VRI [37], form the main inputs into the carbon modeling.

For stands with no disturbance or non-stand replacing MPB disturbance (mortality < 90%), the forest age in 1971 was calculated as the age in 2019 minus 49 years. For stands with disturbance that resets the stand age (wildfire, MPB with >  = 90% mortality, or harvesting), the VRI doesn’t record the previous stand type, so in order to fill this information gap, we populated these areas with the average of the non-disturbed stands by BEC. Yield curves were calculated using standard provincial growth and yield models. The Variable Density Yield Projection (VDYP) model uses VRI data to project yield curves for stands of natural origin. The Table Interpolation Program for Stand Yields (TIPSY) program was used to calculate yield curves for stands after harvesting. Disturbance types and other CBM-CFS3 parameters are detailed in Additional file 1: Table S2. Historical disturbances were explicitly scheduled in the CBM-CFS3 model [58], and the stand-level outputs were summed up to the watershed scale. We calculated the following watershed-level carbon variables: above ground biomass (AGBIO), dead organic matter (DOM), total ecosystem carbon (TEC), net primary production (NPP), and net biome production (NBP) [54, 55, 58,59,60] as shown in Fig. 4. CBM-CFS3 queries and definitions are shown in Additional file 1: Table S3.

Fig. 4
figure 4

Watershed-level carbon variables in the Chilcotin watershed from 1971–2016. Where carbon stock pools are expressed in grams of carbon per square meter (g C m−2) and a is above ground biomass (AGBIO), b total ecosystem carbon (TEC), c dead organic matter (DOM), and carbon fluxes are in grams of carbon per square meter per year (g C m−2 year−1) d net primary production (NPP), e net biome production (NBP)

CBM-CFS3 is a data driven model, with a key input being gross merchantable yield curves. CBM then applies allometric equations to estimate biomass and turnover parameters for dead organic matter and soil carbon [55]. We used 167 measurements of net merchantable volume from the provincial ground plot database [61] throughout the Chilcotin watershed, to validate our yield curves as a proxy for carbon stocks as in Smiley et al. [62]. We found good agreement between the modelled and measured values, where the modelled volumes average 95% of the plot data (R2 = 0.94, P < 0.001) (Additional file 1: Figure S7).

Calculation of WUE

The average annual watershed-level WUE was calculated as the ratio of NPP to ET in grams of carbon per millimeter of water per year (g C m−2 mm−1 H2O year−1) (Eq. 3, Fig. 5).

Fig. 5
figure 5

Annual water related variables in the Chilcotin watershed from 1971–2016. Where a mean annual streamflow in millimeters (mm), b runoff ratio, c evapotranspiration (ET) in mm, and d watershed scale annual water use efficiency (WUE) in grams of carbon per square meter per millimeter of water (g C m−2 mm−1 H2O), and the blue line is the loess regression smoothed line with 95% confidence limits shown as grey shading

$$WUE= \frac{NPP}{ET}$$

Trend analysis

The non-parametric Mann–Kendall test [63] was implemented on pre-whitened data, following the process recommended in Yue et al. [64]. Although often used in hydrologic analysis [65,66,67], some recent studies have debated its applicability to non-stationary hydrologic data [68, 69]. However, we consider its use appropriate as it is one of multiple statistical tests used to help interpret the data. Trend analysis was applied separately to the periods of low and high cumulative disturbance to enable the detection of opposite monotonic trends in each period. All statistical tests are at a significance level of 5%.

Correlation analysis

Two rank-based non-parametric tests (Kendall correlation and Spearman correlation) [70] were used to test the direction and strength of the relationship between forest disturbance and carbon or hydrological variables. As in the trend analysis, they were applied separately to each period.

Gradient boosting machine

The tree-based machine learning Gradient Boosting Machine (GBM) method was used to investigate the relative importance of possible drivers in water and WUE variables. The ‘caret’ package in R [52, 71] was used to fit a GBM model for each variable. Data was pre-processed to scale predictor variables between 0 and 1. Highly correlated variables were removed based on the ‘cor’ function in the ‘caret’ package with a cut off of 90%, resulting in Tmean being dropped from the list of potential GBM model predictor variables. The optimal number of variables for each model was selected using the ‘rfeControl’ function in the ‘caret’ package that minimized the Root Mean Square Error (RMSE). 15 fold cross validation repeated three times was used to tune the GBM models, select the one with the lowest RMSE, and estimate model performance in terms of RMSE and R2. While not a strictly independent estimate of model accuracy, repeated cross validation is commonly used for machine learning models, particularly when an independent dataset is not available [72, 73]. Consistent with other studies, we reported on the relative importance of each predictor variable, the RMSE, and R2 of each final model. We ran four GBM models in total, one each for Q, ET, runoff ratio, and WUE with independent climate and forest disturbance variables (CECA, Tmax, Tmin, P, and PET) being the predictor variables in all cases.

Drivers of WUE variation

We calculated the magnitude of change in WUE and its components (NPP and ET) using Sen’s slope and Eq. 4 following El Masri [74]. The percentages of change in WUE, NPP, and ET were calculated for the low (1971–1999) and high (2000–2016) disturbance periods using the smoothed line from the loess regression (Figs. 4 and 5).

$$Change (\%)=100 \times \frac{y \times s}{m}$$

where, y = number of years in period, s = slope from the Sen’s slope test, m = mean value of NPP, ET, or WUE over the period.


Trends in climate, carbon, and water

Forest disturbance, as measured by CECA significantly increased in the Chilcotin watershed during both the period of low disturbance and high disturbance, increasing from 1% in 1971 to 39.6% in 2016 (Fig. 3, Table 1). Contrary to our expectations based on stand-level dynamics, with increasing disturbance carbon stocks (AGBIO, TEC, and DOM) all significantly (P < 0.001) increased over the period from 1971 to 1999. Also over this first period, NPP increased significantly while NBP decreased or became more negative (P < 0.001). In the period of high disturbance from 2000–2016, AGBIO decreased significantly (P < 0.001) and DOM increased significantly (P < 0.001). Interestingly, TEC, NPP and NBP did not have any statistically significant trends in the high disturbance period.

Table 1 Time-trend analysis of annual variables in the Chilcotin watershed by period

Temperature increased significantly in the Chilcotin watershed from 1971 to 1999, but not in the period of high disturbance from 2000 to 2016. Tmean increased from an average of 0.86 °C to 1.82 °C in the first and last decade of the study period. P and PET did not have any statistically significant trend in either period. Of the hydrological variables, only ET had a significant negative trend (P = 0.027) from 2000 to 2016. WUE significantly increased from 2000 to 2016 (P = 0.043) (Table 1).

Forest disturbance and carbon

Unexpectedly, the correlation tests for the period of high disturbance from 2000–2016 showed that total carbon stocks (TEC) and carbon sequestration (NPP) were not significantly related to CECA (Table 2). The lack of the significant correlations between CECA and carbon stock and sequestration variables, along with their insignificant trends clearly suggests that the severe cumulative disturbance from 2000 to 2016 did not cause significant reduction in carbon stock and sequestration. However, AGBIO was significantly negatively correlated with CECA (Table 2), while DOM had a significant and positive correlation with CECA from 2000–2016.

Table 2 Kendall and Spearman correlation results between forest disturbance and selected carbon and water variables

Forest disturbance and water

Correlation analysis from the period of high disturbance found that all water variables had a significant correlations with forest disturbance (CECA). ET had a significant (P = 0.010) negative correlation with CECA from 2000–2016 with correlation coefficients of – 0.46 and – 0.62 for the Kendall and Spearman tests, respectively (Table 2). Both Q and runoff ratio showed significant (P < 0.05) positive correlations with CECA in the high disturbance period.

P was the most important variable in the runoff ratio and ET GBM models (R2 >  = 0.80), with relative importances of 41.6%, and 70.9% respectively (Table 3). The runoff ratio GBM model did not depend on CECA to explain variation. In contrast, P contributed 70.9% towards the prediction of ET in the GBM model, with CECA contributing the rest of the 29.1%. For the Q GBM model, the predictor variable with the highest importance was Tmin with 30.5%. P was a close second at 25.7%, while CECA accounted for 16.8%. In the Q GBM model, all climate variables contributed, reflecting the complex climatic and biological interactions and their influences on Q.

Table 3 Gradient Boosting Machine model results for selected variables

Forest disturbance and water use efficiency

WUE averaged 0.53 g C m−2 mm−1 H2O over the study period and was significantly (P = 0.017, 0.018) positively correlated to CECA from 2000—2016, with correlation coefficients of 0.43 and 0.57 for the Kendall and Spearman tests, respectively. The WUE GBM model (R2 = 0.93, RMSE = 0.04, n = 45) placed importance on P (55.6%) and CECA (44.4%), reflecting the predictor variables that were relatively important to NPP and ET. In the period of low disturbance (1971—1999), WUE increased by 7% (Sen's slope P = 0.0001). This was largely attributable to increases in NPP which increased by 11%, compared to a 2% increase in ET. In the period of high disturbance (2000–2016) WUE increased by 19% (Sen's slope P < 0.001) driven primarily by a sharp reduction in ET. ET was reduced at a higher rate (− 19%) than NPP (− 2%) causing WUE to increase during this period.


The effects of cumulative forest disturbance on carbon stocks and sequestration

Our study showed that the severe cumulative disturbance from 2000 to 2016 did not cause significant reduction in TEC stocks, but it significantly changed the AGBIO and DOM carbon pools. The significant reduction in AGBIO in the high disturbance period is consistent with both stand-level and landscape-level studies of carbon stock after disturbance [75,76,77]. Watershed-level DOM increased significantly after the high levels of MPB mortality in the early 2000s, as MPB killed living biomass in trees, converting it to DOM [55]. This suggests that there is a relationship between forest disturbance and some carbon sub-pools, and the conversion from AGBIO to DOM caused by MPB infestation was the key reason for insignificant reduction in the total amount of carbon stored in the watershed (TEC) despite severe cumulative disturbance. The CBM-CFS3 model does not explicitly account for the transfer of dissolved organic matter from the soil into rivers, which may be an important physical process to investigate further in relation to increasing disturbance [78, 79].

Studies have shown that forest disturbance can affect and drive changes in ecosystem services such as carbon sequestration and water provisioning [5]. Timber production and forest carbon are often presented as a trade-off relationship, with harvest driving lower carbon stocks or vice versa [2, 77]. At the stand-level, this is most certainly the case. However, this study showed that this relationship may not hold at the large watershed-level. Both the initial state of the watershed’s forest, and the slow rate of accumulation of forest disturbance until the early 2000s, allowed substantial TEC carbon stocks to persist despite disturbance. In areas with stand replacing natural disturbance regimes, harvesting may produce desired outcomes for carbon storage [80], especially if the carbon stored in harvested wood products are considered [81].

Interestingly, our analysis also showed that the cumulative disturbance from 2000 to 2016 did not cause significant reduction in carbon sequestration capacity of the watershed. This result is contradictory to the simple expectation that more disturbance equals lower NPP, a rule that generally holds at the site level, but may not be linearly extrapolated to large areas with significant variation or for all types of disturbance. Eddy covariance studies after MPB mortality have shown a resilience in NPP attributable to residual treed and non-treed vegetation [82, 83]. Kurz et al. [75] predicted an 11% drop in NPP and that the MPB epidemic would cause British Columbia’s central interior forests to be a carbon source at least until the year 2020 (NBP < 0). In the Chilcotin watershed however, the impact might not be severe enough to cause the watershed to be a substantial carbon source (average NBP from 2000–2016 was 4.43 g C m2 year−1), except for 2009 when wildfires in the central part of the watershed emitted carbon (Fig. 4). Nevertheless, the average NPP in the disturbance period of 247.5 g C m2 year−1 is within the range of what other modeling studies in the interior of British Columbia have found (NPP ranging from 118 to 660 g C m2 year−1) [75, 84, 85].

There are two main reasons for non-reduction of carbon sequestration in the watershed as described above. Firstly, based on our knowledge of stand-level carbon dynamics after disturbance, we would expect NPP to decrease with increasing levels of disturbance [76, 77, 83, 86,87,88,89,90]. However, this stand-level relationship cannot simply be extrapolated to the watershed-scale over long periods of time because of forest recovery from regenerated forests or understory vegetation of the MPB stands. The type of disturbance and stand productivity determine the rate of recovery, with NPP recovering to pre-disturbance levels by an average of 20 years after harvesting (Additional file 1: Figure S2). Secondly, the low and consistent levels of harvesting that occurred in the watershed converted a portion of the watershed from older forest to young faster growing stands. In British Columbia, all harvested stands must be replanted after harvesting, often with control of competing vegetation, seedling spacing, and seed selection to enhance growth rates [91]. The cohort of young, fast growing forest has the effect of increasing the overall growth rate (or NPP) of the watershed. This finding is similar to other retrospective studies in British Columbia that found carbon uptake in recovering young forests was able to partially offset disturbance emissions [57] and that younger forests had higher CO2 uptake [80].

Cumulative forest disturbance and water processes

Our results showed that the partitioning of P into ET and Q were significantly affected by forest disturbance. ET significantly decreased over the high disturbance period from 447 mm in 2000 to 427 mm in 2016, despite a 12% increase in PET over the same period. Correlation analysis showed a significant (P < 0.010) negative relationship between ET and CECA. The negative correlation was expected for ET because forest disturbance causes the mortality of forests and consequently reduced transpiration. This result is consistent with other studies that found ET was reduced after disturbance such as harvesting, fire, and insect disturbances [29, 30, 92,93,94]. The ET GBM model showed that after P, CECA had the next highest relative importance of 29.1%.

Table 2 also showed cumulative forest disturbance (CECA) has a positive correlation with Q and runoff ratio in the high disturbance period (2000 to 2016). This is expected as the higher levels of forest disturbance reduced ET and consequently increased the partitioning of P into Q. The significant increase of Q and runoff ratio is consistent with other studies [95,96,97,98].

It is important to mention that in this study, we used a novel approach to varying the watershed parameter, w, in Zhang’s equation (Eq. 1) to account for changes in ET due to forest disturbance. This approach improves our estimation of watershed-based ET which is critical for quantifying WUE. It is widely recognized that accurate estimation of watershed-based ET is challenging mainly due to the lack of widely applicable methodology [99, 100]. Many studies have investigated the long term relationship between climate and water partitioning at the watershed-level utilizing the Budyko framework [48, 101,102,103,104]. Often watershed properties are represented as a constant, such as w representing plant-available water capacity in Zhang’s Eq. (50) or the watershed characteristic m in Fuh’s Eq. (105). Increasingly, it is recognized that this static representation is not sufficient in mixed cover watersheds [106] or as forests change through time [107]. In this study, we developed a relationship between w and CECA based on long-term data in some comparable watersheds, and then applied this relationship to assign w values according to various CECA levels. The estimates of ET from this approach were checked by those calculated from the water balance method. The high consistence in ET estimates between those two methods demonstrated that our ET estimation was reliable. Additional analysis across three and five year windows was carried out to test if ET was sensitive to the assumption that the annual change in storage is zero. Results did not show significant change in either the trend or magnitude of ET, indicating that the analysis and conclusions are robust to this assumption (Additional file 1: Table S4 and Table S5). Our long term average estimated ET (471 mm) compares closely to those from similar watersheds such as Moffatt and Baker at 468 and 439 mm respectively (Additional file 1: Figure S6) [44].

Landscape-level WUE increase under cumulative disturbance

Over the study period, we found that the average annual WUE was 0.53 g C m2 mm−1 H2O in the Chilcotin watershed. Direct comparison with NPP values from other studies is difficult, as many are either based on other measures of carbon uptake such as Gross Primary Production (GPP), or Net Ecosystem Production (NEP), rather than NPP, measurements are from the growing season only, or taken in different climates or land covers than this study. The values of WUE that we found (~ 0.5 g C m2 mm−1 H2O) are lower than those calculated for forests across the continental United States of America (USA)(0.93 g C kg−1 H2O) [13], which are likely more productive biomes and have very different climates. Our WUE is lower than those calculated just for the growing seasons, such as ~ 0.95 g C m2 kg−1 H2O in European temperate evergreen conifer forest [108], but is similar to a study in North Eastern China that found annual WUE calculated as NPP/ET was in the range of 0.46–1.10 5 g C m2 mm−1 H2O [109]. We are not aware of any other studies on annual WUE using NPP/ET in the same region as our study area.

During the high disturbance period (2000–2016), annual watershed-level WUE significantly increased in the Chilcotin watershed and correlation analysis found that WUE was significantly and positively correlated to CECA (P = 0.018). Cumulative forest disturbance played an important role in both carbon and water variables and consequently WUE changes. Additional evidence for the relationship of WUE and forest disturbance was also shown by the WUE GBM model, which found that of the independent variables we tested it with, P (55.6%) and CECA (44.4%) had the highest relative importance for predicting WUE. Further analysis using Sen’s slope, revealed that the dramatic increase in WUE (19%, P < 0.001) during the high disturbance period was mainly due to the significantly larger decrease in ET (19%) than the more resilient NPP (-2%). Clearly, this unexpected result was not only related to significant decrease in ET, but also the resilience of NPP at the large watershed scale. Increasingly, WUE is being used as a cross biome indicator of ecosystem resilience to drought, at the global scale [110] and at the basin scale across India [111] and being used to support irrigation or forest management strategies in arid regions [112]. Our finding that WUE was resilient to disturbance in the Chilcotin watershed adds to this body of work.

Studies have used WUE across large climatic gradients to assess trade offs between water use and terrestrial carbon [113, 114]. Gao et al. [113] found differences between semi-humid and arid zones in China and used this to infer optimal strategies about the location of afforestation. This suggests that our result may be area or regionally specific and underlines the importance of future studies in other geographic, climatic, or ecological types.

The majority of WUE studies are at shorter time intervals and either at the site-level or at continental or global scales [12, 13, 18, 113, 115, 116]. Stand-level studies showed WUE decreasing after some types of disturbance [25, 84, 117], although areas with substantial understory vegetation, or a low severity of disturbance often showed no change in WUE [25, 30]. Those that showed a decrease in WUE concluded that the reduction in carbon uptake was larger than the reduction in ET [117]. Contrary to this, we found that watershed-level NPP was largely sustained under high disturbance, while ET was decreased, which led to a significant increase in WUE. This result is the opposite of expected based on stand-level WUE studies, and was because NPP was sustained at the watershed-level as discussed above. This unique finding adds to the body of literature on WUE to improve our knowledge of how WUE changes with forest disturbance at the large spatial scale.

Future studies on forest disturbance should also consider the effects of climate change on WUE. Studies investigating long term trends in WUE under the effect of climate change have mixed results. Those looking at leaf or tree intrinsic WUE have often found an increase that is attributed to CO2 fertilization or other climate factors [26, 27]. Other studies found a decreasing WUE trend in forest biomes, or over the few last decades, compared to the decades before [14, 74, 118]. There are important regional differences to the general trends, and studies using process based models sometimes showed divergent responses [27, 74]. Although climatic variables were included in our calculations of ET and thus WUE, we did not explicitly separate out the relative role of climate change and forest disturbance in WUE variation. This will be our research priority, particularly with regard to the potential impacts of future climate change. Incorporation of full ecohydrologic modeling and independent sources of watershed-scale NPP and ET validation are also areas of high priority.


Cumulative forest disturbance played an important role in both carbon and water dynamics, and consequently WUE in the Chilcotin watershed over the last 45 years. To our surprise, watershed-level forest carbon stocks and sequestration capacity were resilient after the significant cumulative disturbance, which was likely due to the initial state of the forests and varied types of disturbance at levels low enough to balance an aging carbon store with young fast growing regenerating forest. The resilient carbon variables, along with a substantial decrease in annual ET following the disturbance, drove WUE to increase by 19%. These findings have important implications for managing forest carbon and water in the large watershed or landscape scale.

Availability of data and materials

The timeseries data for the Chilcotin watershed that supports the conclusions of this article is included within the articles additional files as Additional file 1: Table S1.



Above Ground Biomass




Cumulative Equivalent Clear-cut Area


Carbon Budget Model of the Canadian Forest Service version 3


Dead Organic Matter


Equivalent Clear-cut Area



g C:

Grams of carbon


Gradient Boosting Machine


Gross Primary Production






Interior Douglas-fir


Square meter




Mountain Pine Beetle


Montane Spruce


Net Biome Production


Net Ecosystem Production


Net Primary Production




Potential Evapotranspiration


Mean Annual Streamflow


The R Programming Language


Root Mean Square Error


Annual Storage Change


Sub-Boreal Pine–Spruce


Total Ecosystem Carbon


Table Interpolation Program for Stand Yields


Annual Maximum Daily Temperature


Annual Mean Daily Temperature


Annual Minimum Daily Temperature


United States of America


Variable Density Yield Projection


Vegetation Resources Inventory


Plant-Available Water Coefficient


Water Use Efficiency


  1. Costanza R, d’Arge R, de Groot R, Farber S, Grasso M, Hannon B, et al. The value of the world’s ecosystem services and natural capital. Nature. 1997;387(6630):253–60.

    Article  CAS  Google Scholar 

  2. Cademus R, Escobedo F, McLaughlin D, Abd-Elrahman A. Analyzing trade-offs, synergies, and drivers among timber production, carbon sequestration, and water yield in Pinus elliotii forests in Southeastern USA. Forests. 2014;5(6):1409–31.

    Article  Google Scholar 

  3. Canadell JG, Raupach MR. Managing Forests for Climate Change Mitigation. Science. 2008;320(5882):1456.

    Article  CAS  Google Scholar 

  4. Jackson RB, Carpenter SR, Dahm CN, McKnight DM, Naiman RJ, Postel SL, et al. Water in a Changing World. Ecol Appl. 2001;11(4):1027–45.

    Article  Google Scholar 

  5. Roces-Díaz JV, Vayreda J, De Cáceres M, García-Valdés R, Banqué-Casanovas M, Morán-Ordóñez A, et al. Temporal changes in Mediterranean forest ecosystem services are driven by stand development, rather than by climate-related disturbances. For Ecol Manage. 2021;480:118623.

    Article  Google Scholar 

  6. Farley KA, Jobbágy EG, Jackson RB. Effects of afforestation on water yield: a global synthesis with implications for policy. Glob Change Biol. 2005;11(10):1565–76.

    Article  Google Scholar 

  7. Ferguson PR, Veizer J. Coupling of water and carbon fluxes via the terrestrial biosphere and its significance to the Earth's climate system. Journal of Geophysical Research: Atmospheres. 2007;112(D24).

  8. Gentine P, Green JK, Guérin M, Humphrey V, Seneviratne SI, Zhang Y, et al. Coupling between the terrestrial carbon and water cycles—a review. Environ Res Lett. 2019;14(8):083003.

    Article  CAS  Google Scholar 

  9. Chen B, Coops NC. Understanding of coupled terrestrial carbon, nitrogen and water dynamics- an overview. Sensors. 2009;9(11):8624–57.

    Article  CAS  Google Scholar 

  10. Hatfield JL, Dold C. Water-Use Efficiency: Advances and Challenges in a Changing Climate. Front Plant Sci. 2019;10:103.

    Article  Google Scholar 

  11. Wallace JS. Increasing agricultural water use efficiency to meet future food production. Agr Ecosyst Environ. 2000;82(1):105–19.

    Article  Google Scholar 

  12. Sun Y, Piao S, Huang M, Ciais P, Zeng Z, Cheng L, et al. Global patterns and climate drivers of water-use efficiency in terrestrial ecosystems deduced from satellite-based datasets and carbon cycle models. Glob Ecol Biogeogr. 2016;25(3):311–23.

    Article  Google Scholar 

  13. Tian H, Chen G, Liu M, Zhang C, Sun G, Lu C, et al. Model estimates of net primary productivity, evapotranspiration, and water use efficiency in the terrestrial ecosystems of the southern United States during 1895–2007. Managing landscapes at multiple scales for sustainability of ecosystem functions. 2010;259(7):1311–27.

    Google Scholar 

  14. Tian H, Lu C, Chen G, Xu X, Liu M, Ren W, et al. Climate and land use controls over terrestrial water use efficiency in monsoon Asia. Ecohydrology. 2011;4(2):322–40.

    Article  Google Scholar 

  15. Zhou S, Yu B, Huang Y, Wang G. The effect of vapor pressure deficit on water use efficiency at the subdaily time scale. Geophys Res Lett. 2014;41(14):5005–13.

    Article  Google Scholar 

  16. Fernandes TJG, Del Campo AD, Herrera R, Molina AJ. Simultaneous assessment, through sap flow and stable isotopes, of water use efficiency (WUE) in thinned pines shows improvement in growth, tree-climate sensitivity and WUE, but not in WUEi. For Ecol Manage. 2016;361:298–308.

    Article  Google Scholar 

  17. Kaminski KP, Kørup K, Nielsen KL, Liu F, Topbjerg HB, Kirk HG, et al. Gas-exchange, water use efficiency and yield responses of elite potato (Solanum tuberosum L.) cultivars to changes in atmospheric carbon dioxide concentration, temperature and relative humidity. Agricultural Forest Meteorol. 2014;187:36–45.

    Article  Google Scholar 

  18. Vickers D, Thomas C, Pettijohn C, Martin J, Law B. Five years of carbon fluxes and inherent water-use efficiency at two semi-arid pine forests with different disturbance histories. Tellus B Chemical Physical Meteorol. 2012;64(1):17159.

    Article  Google Scholar 

  19. Baldocchi D. Breathing of the terrestrial biosphere: lessons learned from a global network of carbon dioxide flux measurement systems. Aust J Bot. 2008;56(1):1–26.

    Article  CAS  Google Scholar 

  20. Beer C, Ciais P, Reichstein M, Baldocchi D, Law BE, Papale D, et al. Temporal and among-site variability of inherent water use efficiency at the ecosystem level. Global Biogeochemical Cycles. 2009;23(2)

  21. Hwang T, Band L, Hales TC. Ecosystem processes at the watershed scale: Extending optimality theory from plot to catchment. Water Resources Research. 2009;45(11)

  22. Troch PA, Martinez GF, Pauwels VRN, Durcik M, Sivapalan M, Harman C, et al. Climate and vegetation water use efficiency at catchment scales. Hydrol Process. 2009;23(16):2409–14.

    Article  Google Scholar 

  23. Kang S, Running SW, Kimball JS, Fagre DB, Michaelis A, Peterson DL, et al. Effects of spatial and temporal climatic variability on terrestrial carbon and water fluxes in the Pacific Northwest, USA. Environmental Modelling Software. 2014;51:228–39.

    Article  Google Scholar 

  24. Jasechko S, Sharp ZD, Gibson JJ, Birks SJ, Yi Y, Fawcett PJ. Terrestrial water fluxes dominated by transpiration. Nature. 2013;496(7445):347–50.

    Article  CAS  Google Scholar 

  25. Clark KL, Skowronski NS, Gallagher MR, Renninger H, Schafer KVR. Contrasting effects of invasive insects and fire on ecosystem water use efficiency. Biogeosciences. 2014;11(23):6509–23.

    Article  Google Scholar 

  26. Keenan TF, Hollinger DY, Bohrer G, Dragoni D, Munger JW, Schmid HP, et al. Increase in forest water-use efficiency as atmospheric carbon dioxide concentrations rise. Nature. 2013;499(7458):324–7.

    Article  CAS  Google Scholar 

  27. Huang M, Piao S, Sun Y, Ciais P, Cheng L, Mao J, et al. Change in terrestrial ecosystem water-use efficiency over the last three decades. Glob Change Biol. 2015;21(6):2366–78.

    Article  Google Scholar 

  28. Mkhabela MS, Amiro BD, Barr AG, Black TA, Hawthorne I, Kidston J, et al. Comparison of carbon dynamics and water use efficiency following fire and harvesting in Canadian boreal forests. Agric For Meteorol. 2009;149(5):783–94.

    Article  Google Scholar 

  29. Dore S, Montes-Helu M, Hart SC, Hungate BA, Koch GW, Moon JB, et al. Recovery of ponderosa pine ecosystem carbon and water fluxes from thinning and stand-replacing fire. Glob Change Biol. 2012;18(10):3171–85.

    Article  Google Scholar 

  30. Reed DE, Ewers BE, Pendall E. Impact of mountain pine beetle induced mortality on forest carbon and water fluxes. Environ Res Lett. 2014;9(10):105004.

    Article  CAS  Google Scholar 

  31. Meidinger D, Pojar J. Ecosystems of British Columbia. Victoria, British Columbia, Canada: BC Ministry of Forests; 1991.

    Google Scholar 

  32. Wang T, Hamann A, Spittlehouse D, Carroll C. Locally Downscaled and Spatially Customizable Climate Data for Historical and Future Periods for North America. PLoS ONE. 2016;11(6):e0156720.

    Article  CAS  Google Scholar 

  33. Natural Resources Canada. Digital Elevation Model of Canada - Canada3D, 2001. 2019–10–03 ed2019.

  34. Hargreaves GH, Samani ZA. Reference crop evapotranspiration from temperature. Appl Eng Agric. 1985;1:96–9.

    Article  Google Scholar 

  35. Environment Canada. Extracted from Environment and Climate Change Canada's HYDAT.mdb, released on 2020–01–29. In: Canada E, editor. 2020.

  36. Environment and Climate Change Canada. National hydrometric network basin polygons. 2020–01–15 ed2020.

  37. Ministry of Forests L, Natural Resource Operations and Rural Development, . Forest Vegetation Composite Rank 1 Layer. 2020–04–14 ed2019.

  38. Ministry of Forests L, Natural Resource Operations and Rural Development,. Harvested Areas of BC (Consolidated Cutblocks). 2020–03–31 ed2020.

  39. Winkler R, Boon S. Revised Snow Recovery Estimates for Pine-dominated Forests in Interior British Columbia. Report. Kamloops, British Columbia: Ministry of Forests, Lands and Natural Resource Operations; 2015. Report No.: Extension note 116.

  40. Winkler R, Spittlehouse D, Boon S. Streamflow Response to Clearcut Logging on British Columbia's Okanagan Plateau. Ecohydrology. 2017.

  41. Lin Y, Wei X. The impact of large-scale forest harvesting on hydrology in the Willow watershed of Central British Columbia. J Hydrol. 2008;359(1–2):141–9.

    Article  Google Scholar 

  42. Wei X, Zhang M. Quantifying streamflow change caused by forest disturbance at a large spatial scale: A single watershed study. Water Resources Research. 2010;46(12).

  43. Zhang M, Wei X, Li Q. A quantitative assessment on the response of flow regimes to cumulative forest disturbances in large snow-dominated watersheds in the interior of British Columbia. Canada Ecohydrol. 2016;9(5):843–59.

    Article  Google Scholar 

  44. Zhang M, Wei X, Li Q. Do the hydrological responses to forest disturbances in large watersheds vary along climatic gradients in the interior of British Columbia, Canada? Ecohydrology. 2017;10(2):e1840.

    Article  Google Scholar 

  45. British Columbia Ministry of Forests. Forest Practices Code of British Columbia. Coastal Watershed Assessment Procedure Guidebook (CWAP) and Interior Watershed Assessment Procedure Guidebook (IWAP). Victoria, British Columbia, Canada: British Columbia Ministry of Forests; 1999.

  46. Forest Practices B. Lodgepole Pine Stand Structure 25 Years after Mountain Pine Beetle Attack. Report. Forest Practices Board; 2007. Report No.: FPB/SR/32.

  47. Lewis D, Huggard D. A Model to Quantify Effects of Mountain Pine Beetle on Equivalent Clearcut Area. Streamline Watershed Management Bulletin. 2010;13(2).

  48. Budyko MI. Climate and Life: Academic Press; 1974.

  49. Zhang L, Dawes WR, Walker GR. Response of mean annual evapotranspiration to vegetation changes at catchment scale. Water Resour Res. 2001;37(3):701–8.

    Article  Google Scholar 

  50. Zhang L, Hickel K, Dawes WR, Chiew FHS, Western AW, Briggs PR. A rational function approach for estimating mean annual evapotranspiration. Water Resources Research. 2004;40(2).

  51. Soetaert K. rootSolve: Nonlinear root finding, equilibrium and steady-state analysis of ordinary differential equations. 2009.

  52. R Core Team. R: A Language and Environment for Statistical Computing. Vienna, Austria: R Foundation for Statistical Computing; 2016.

  53. Brooks PD, Troch PA, Durcik M, Gallo E, Schlegel M. Quantifying regional scale ecosystem response to changes in precipitation: Not all rain is created equal. Water Resources Research. 2011;47(10).

  54. Kull SJ, Rampley GJ, Morken S, Metsaranta J, Neilson ET, Kurz WA. Operational-scale Carbon Budget Model of the Canadian Forest Sector (CBM-CFS3) version 1.2: user’s guide. Edmonton, AB: Natural Resources Canada, Canadian Forest Service, Northern Forestry Center 2016.

  55. Kurz WA, Dymond CC, White TM, Stinson G, Shaw CH, Rampley GJ, et al. CBM-CFS3: A model of carbon-dynamics in forestry and land-use change implementing IPCC standards. Ecol Model. 2009;220(4):480–504.

    Article  Google Scholar 

  56. Bernier PY, Guindon L, Kurz WA, Stinson G. Reconstructing and modelling 71 years of forest growth in a Canadian boreal landscape: a test of the CBM-CFS3 carbon accounting model. Can J For Res. 2010;40(1):109–18.

    Article  Google Scholar 

  57. Trofymow JA, Stinson G, Kurz WA. Derivation of a spatially explicit 86-year retrospective carbon budget for a landscape undergoing conversion from old-growth to managed forests on Vancouver Island. BC Forest Ecology Management. 2008;256(10):1677–91.

    Article  Google Scholar 

  58. Boisvenue C, Smiley BP, White JC, Kurz WA, Wulder MA. Improving carbon monitoring and reporting in forests using spatially-explicit information. Carbon Balance Manage. 2016;11(1):23.

    Article  CAS  Google Scholar 

  59. Shaw CH, Hilger AB, Metsaranta J, Kurz WA, Russo G, Eichel F, et al. Evaluation of simulated estimates of forest ecosystem carbon stocks using ground plot data from Canada’s National Forest Inventory. Ecol Model. 2014;272:323–47.

    Article  Google Scholar 

  60. Smyth CE, Trofymow JA, Kurz WA. Decreasing uncertainty in CBM-CFS3 estimates of forest soil C sources and sinks through use of long-term data from the Canadian Intersite Decomposition Experiment. Report. Victoria, BC.: CIDET Working Group Natural Resources Canada, Canadian Forest Service, Pacific Forestry Centre; 2009. Report No.: Information Report BC-X-422.

  61. BC Gov. Forest Inventory Sample Plots. In: Forest Analysis and Inventory Branch, editor. Victoria, British Columbia2020.

  62. Smiley BP, Trofymow JA, Niemann KO. Spatially-explicit reconstruction of 100 years of forest land use and disturbance on a coastal British Columbia Douglas-fir-dominated landscape: Implications for future watershed-scale carbon stock recovery. Appl Geogr. 2016;74:109–22.

    Article  Google Scholar 

  63. Mann HB. Nonparametric Tests Against Trend. Econometrica. 1945;13(3):245–59.

    Article  Google Scholar 

  64. Yue S, Pilon P, Phinney B, Cavadias G. The influence of autocorrelation on the ability to detect trend in hydrological series. Hydrol Process. 2002;16(9):1807–29.

    Article  Google Scholar 

  65. Déry SJ, Wood EF. Decreasing river discharge in northern Canada. Geophysical Research Letters. 2005;32(10).

  66. Li Y, He D, Li X, Zhang Y, Yang L. Contributions of Climate Variability and Human Activities to Runoff Changes in the Upper Catchment of the Red River Basin, China. Water. 2016;8(9).

  67. Zhang M, Wei X. Contrasted hydrological responses to forest harvesting in two large neighbouring watersheds in snow hydrology dominant environment: implications for forest management and future forest hydrology studies. Hydrol Process. 2014;28(26):6183–95.

    Article  Google Scholar 

  68. Razavi S, Vogel R. Prewhitening of hydroclimatic time series? Implications for inferred change and variability across time scales. J Hydrol. 2018;557:109–15.

    Article  Google Scholar 

  69. Serinaldi F, Kilsby CG, Lombardo F. Untenable nonstationarity: An assessment of the fitness for purpose of trend tests in hydrology. Adv Water Resour. 2018;111:132–55.

    Article  Google Scholar 

  70. Kendall MG. Rank Correlation Methods. New York: Oxford University Press; 1975. p. 202.

    Google Scholar 

  71. Kuhn M. Building predictive models in R using the caret package. Journal of Statistical Software. 2008.

  72. Assal TJ, Sibold J, Reich R. Modeling a Historical Mountain Pine Beetle Outbreak Using Landsat MSS and Multiple Lines of Evidence. Remote Sens Environ. 2014;155:275–88.

    Article  Google Scholar 

  73. Rice JS, Emanuel RE. Ecohydrology of Interannual Changes in Watershed Storage. Water Resour Res. 2019;55(10):8238–51.

    Article  Google Scholar 

  74. El Masri B, Schwalm C, Huntzinger DN, Mao J, Shi X, Peng C, et al. Carbon and Water Use Efficiencies: A Comparative Analysis of Ten Terrestrial Ecosystem Models under Changing Climate. Scientific Reports. 2019;9(1):14680.

    Article  CAS  Google Scholar 

  75. Kurz WA, Dymond CC, Stinson G, Rampley GJ, Neilson ET, Carroll AL, et al. Mountain pine beetle and forest carbon feedback to climate change. Nature. 2008;452(7190):987–90.

    Article  CAS  Google Scholar 

  76. Pfeifer EM, Hicke JA, Meddens AJH. Observations and modeling of aboveground tree carbon stocks and fluxes following a bark beetle outbreak in the western United States. Glob Change Biol. 2011;17(1):339–50.

    Article  Google Scholar 

  77. Seely B, Welham C, Kimmins H. Carbon sequestration in a boreal forest ecosystem: results from the ecosystem simulation model, FORECAST. The Role of Boreal Forests and Forestry in the Global Carbon Budget. 2002;169(1–2):123–35.

    Google Scholar 

  78. Drake TW, Van Oost K, Barthel M, Bauters M, Hoyt AM, Podgorski DC, et al. Mobilization of aged and biolabile soil carbon by tropical deforestation. Nat Geosci. 2019;12(7):541–6.

    Article  CAS  Google Scholar 

  79. Zhao G, Gao Y, Wang L, Hao Z, Wen X, Song X. Isotopically-tracked hydrological changes in carbon cycling and its sources in a Chinese subtropical forested watershed. J Hydrol. 2019;575:1041–51.

    Article  CAS  Google Scholar 

  80. Sharma T, Kurz WA, Stinson G, Pellatt MG, Li Q. A 100-year conservation experiment: Impacts on forest carbon stocks and fluxes. For Ecol Manage. 2013;310:242–55.

    Article  Google Scholar 

  81. Kurz WA, Shaw CH, Boisvenue C, Stinson G, Metsaranta J, Leckie D, et al. Carbon in Canada’s boreal forest—A synthesis. Environmental Rev. 2013;21(4):260–92.

    Article  CAS  Google Scholar 

  82. Bowler R, Fredeen AL, Brown M, Andrew BT. Residual vegetation importance to net CO2 uptake in pine-dominated stands following mountain pine beetle attack in British Columbia. Canada Forest Ecology Management. 2012;269:82–91.

    Article  Google Scholar 

  83. Brown MG, Black TA, Nesic Z, Fredeen AL, Foord VN, Spittlehouse DL, et al. The carbon balance of two lodgepole pine stands recovering from mountain pine beetle attack in British Columbia. Land-Atmosphere Interactions: Advances in Measurement, Analysis, and Modeling – A Tribute to T Andrew Black. 2012;153:82–93.

  84. Meyer G, Black TA, Jassal RS, Nesic Z, Grant NJ, Spittlehouse DL, et al. Measurements and simulations using the 3-PG model of the water balance and water use efficiency of a lodgepole pine stand following mountain pine beetle attack. For Ecol Manage. 2017;393:89–104.

    Article  Google Scholar 

  85. Hof AR, Dymond CC, Mladenoff DJ. Climate change mitigation through adaptation: the effectiveness of forest diversification by novel tree planting regimes. Ecosphere. 2017;8(11):e01981.

    Article  Google Scholar 

  86. Amiro BD, Barr AG, Barr JG, Black TA, Bracho R, Brown M, et al. Ecosystem carbon dioxide fluxes after disturbance in forests of North America. Journal of Geophysical Research: Biogeosciences (2005-2012). 2010;115.

  87. Brown M, Black TA, Nesic Z, Foord VN, Spittlehouse DL, Fredeen AL, et al. Impact of mountain pine beetle on the net ecosystem production of lodgepole pine stands in British Columbia. Agric For Meteorol. 2010;150(2):254–64.

    Article  Google Scholar 

  88. Mathys A, Black TA, Nesic Z, Nishio G, Brown M, Spittlehouse DL, et al. Carbon balance of a partially harvested mixed conifer forest following mountain pine beetle attack and its comparison to a clear-cut. Biogeosciences. 2013;10(8):5451–63.

    Article  CAS  Google Scholar 

  89. Humphreys ER, Black TA, Morgenstern K, Cai T, Drewitt GB, Nesic Z, et al. Carbon dioxide fluxes in coastal Douglas-fir stands at different stages of development after clearcut harvesting. The Fluxnet-Canada Research Network: Influence of Climate and Disturbance on Carbon Cycling in Forests and Peatlands. 2006;140(1–4):6–22.

    Google Scholar 

  90. Meyer G, Black TA, Jassal RS, Nesic Z, Coops NC, Christen A, et al. Simulation of net ecosystem productivity of a lodgepole pine forest after mountain pine beetle attack using a modified version of 3-PG. For Ecol Manage. 2018;412:41–52.

    Article  Google Scholar 

  91. Woods J, Mahoney C. Climate-Based Seed Transfer: Guiding British Columbia’s reforestation investments in a changing climate. In: BC FGCo, editor. 2016.

  92. Biederman JA, Harpold AA, Gochis DJ, Ewers BE, Reed DE, Papuga SA, et al. Increased evaporation following widespread tree mortality limits streamflow response. Water Resour Res. 2014;50(7):5395–409.

    Article  Google Scholar 

  93. Kim J, Hwang T, Schaaf CL, Orwig DA, Boose E, Munger JW. Increased water yield due to the hemlock woolly adelgid infestation in New England. Geophys Res Lett. 2017;44(5):2327–35.

    Article  Google Scholar 

  94. Millar DJ, Ewers BE, Mackay DS, Peckham S, Reed DE, Sekoni A. Improving ecosystem-scale modeling of evapotranspiration using ecological mechanisms that account for compensatory responses following disturbance. Water Resour Res. 2017;53(9):7853–68.

    Article  Google Scholar 

  95. Bosch JM, Hewlett JD. A review of catchment experiments to determine the effect of vegetation changes on water yield and evapotranspiration. J Hydrol. 1982;55:3–23.

    Article  Google Scholar 

  96. Brown AE, Zhang L, McMahon TA, Western AW, Vertessy RA. A review of paired catchment studies for determining changes in water yield resulting from alterations in vegetation. J Hydrol. 2005;310(1):28–61.

    Article  Google Scholar 

  97. Livneh B, Deems JS, Buma B, Barsugli JJ, Schneider D, Molotch NP, et al. Catchment response to bark beetle outbreak and dust-on-snow in the Colorado Rocky Mountains. J Hydrol. 2015;523:196–210.

    Article  Google Scholar 

  98. Zhang M, Wei X, Sun P, Liu S. The effect of forest harvesting and climatic variability on runoff in a large watershed: The case study in the Upper Minjiang River of Yangtze River basin. J Hydrol. 2012;464–465:1–11.

    Article  Google Scholar 

  99. Ford CR, Hubbard RM, Kloeppel BD, Vose JM. A comparison of sap flux-based evapotranspiration estimates with catchment-scale water balance. Agric For Meteorol. 2007;145(3–4):176–85.

    Article  Google Scholar 

  100. Lu J, Sun G, McNulty SG, Amatya D. A comparison of six potential evapotranspiration methods for regional use in the Southeastern United States. J American Water Resources Association. 2005;41(3):621–33.

    Article  Google Scholar 

  101. Berghuijs WR, Larsen JR, van Emmerik THM, Woods RA. A Global Assessment of Runoff Sensitivity to Changes in Precipitation, Potential Evaporation, and Other Factors. Water Resour Res. 2017;53(10):8475–86.

    Article  Google Scholar 

  102. Lee C-H, Yeh H-F. Impact of Climate Change and Human Activities on Streamflow Variations Based on the Budyko Framework. Water. 2019;11.

  103. Shen Q, Cong Z, Lei H. Evaluating the impact of climate and underlying surface change on runoff within the Budyko framework: A study across 224 catchments in China. J Hydrol. 2017;554:251–62.

    Article  Google Scholar 

  104. Zhao J, Huang S, Huang Q, Wang H, Leng G. Detecting the Dominant Cause of Streamflow Decline in the Loess Plateau of China Based onthe Latest Budyko Equation. Water. 2018;10:1277.

    Article  Google Scholar 

  105. Wei X, Li Q, Zhang M, Giles-Hansen K, Liu W, Fan H, et al. Vegetation cover—another dominant factor in determining global water resources in forested regions. Glob Change Biol. 2018;24(2):786–95.

    Article  Google Scholar 

  106. van Dijk AIJM, Peña-Arancibia JL, Bruijnzeel LA. Land cover and water yield: inference problems when comparing catchments with mixed land cover. Hydrol Earth Syst Sci. 2012;16(9):3461–73.

    Article  Google Scholar 

  107. Zhang S, Yang H, Yang D, Jayawardena AW. Quantifying the effect of vegetation change on the regional water balance within the Budyko framework. Geophys Res Lett. 2016;43(3):1140–8.

    Article  Google Scholar 

  108. Kuglitsch FG, Reichstein M, Beer C, Carrara A, Ceulemans R, Granier A, et al. Characterisation of ecosystem water-use efficiency of european forests from eddy covariance measurements. Biogeosciences Discussions. 2008;5:4481–519.

    Google Scholar 

  109. Li X, Farooqi TJA, Jiang C, Liu S, Sun OJ. Spatiotemporal variations in productivity and water use efficiency across a temperate forest landscape of Northeast China. Forest Ecosystems. 2019;6(1):22.

    Article  Google Scholar 

  110. Ponce-Campos GE, Moran MS, Huete A, Zhang Y, Bresloff C, Huxman TE, et al. Ecosystem resilience despite large-scale altered hydroclimatic conditions. Nature. 2013;494(7437):349–52.

    Article  CAS  Google Scholar 

  111. Sharma A, Goyal MK. Assessment of ecosystem resilience to hydroclimatic disturbances in India. Glob Change Biol. 2018;24(2):e432–41.

    Article  Google Scholar 

  112. Zhou S, Yu B, Zhang Y, Huang Y, Wang G. Water use efficiency and evapotranspiration partitioning for three typical ecosystems in the Heihe River Basin, northwestern China. Agric For Meteorol. 2018;253–254:261–73.

    Article  Google Scholar 

  113. Gao Y, Zhu X, Yu G, He N, Wang Q, Tian J. Water use efficiency threshold for terrestrial ecosystem carbon sequestration in China under afforestation. Agric For Meteorol. 2014;195–196:32–7.

    Article  Google Scholar 

  114. Xiao J, Sun G, Chen J, Chen H, Chen S, Dong G, et al. Carbon fluxes, evapotranspiration, and water use efficiency of terrestrial ecosystems in China. Agric For Meteorol. 2013;182–183:76–90.

    Article  Google Scholar 

  115. Baldocchi DD, Verma SB, Anderson DE. Canopy photosynthesis and water-use efficiency in a deciduous forest. J Appl Ecol. 1987;24(1):251–60.

    Article  Google Scholar 

  116. Mikkelson KM, Bearup LA, Maxwell RM, Stednick JD, McCray JE, Sharp JO. Bark beetle infestation impacts on nutrient cycling, water quality and interdependent hydrological effects. Biogeochemistry. 2013;115(1):1–21.

    Article  CAS  Google Scholar 

  117. Jassal RS, Black TA, Spittlehouse DL, Brümmer C, Nesic Z. Evapotranspiration and water use efficiency in different-aged Pacific Northwest Douglas-fir stands. Agric For Meteorol. 2009;149(6–7):1168–78.

    Article  Google Scholar 

  118. Jiang Y, Still CJ, Rastogi B, Page GFM, Wharton S, Meinzer FC, et al. Trends and controls on water-use efficiency of an old-growth coniferous forest in the Pacific Northwest. Environ Res Lett. 2019;14(7):074029.

    Article  CAS  Google Scholar 

Download references


The authors would like to acknowledge Tongli Wang, Faculty of Forestry, University of British Columbia, Vancouver who helped with ClimateBC data.


Funding for this project was made available through the Natural Sciences and Engineering Research Council of Canada (NSERC) Collaborative Research and Development Grant (CRDG) for the project titled “The effects of reforestation on forest carbon and water coupling at multiple spatial scales” (CRDPJ 485176–15).

Author information

Authors and Affiliations



KGH carried out carbon and evapotranspiration modeling, statistical analysis and was the primary author of the manuscript. YH calculated forest disturbance variables and calculated input climate, forest and hydrologic data. YH interpreted results and helped to review and write the manuscript. XW reviewed and interpreted results and contributed to the overall design and message of the manuscript. XW was a major contributor in writing the manuscript. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Xiaohua Wei.

Ethics declarations

Competing interests

The authors declare that they have no competing interests

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1.

Figures and Tables.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Giles-Hansen, K., Wei, X. & Hou, Y. Dramatic increase in water use efficiency with cumulative forest disturbance at the large forested watershed scale. Carbon Balance Manage 16, 6 (2021).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: