Potential role of permafrost thaw on increasing Siberian river discharge

Despite the increasing Siberian river discharge, the sensitivity of streamflow to climate forcing/permafrost thawing is poorly quantified. Based on the Budyko framework and superposition principles, we detected and attributed the changes in streamflow regimes for the three great Siberian rivers (Ob, Yenisei, and Lena) during 1936–2019. Over the past 84 years, streamflow of Ob, Yenisei and Lena has increased by ∼7.7%, 7.4% and 22.0%, respectively. Intensified precipitation induced by a warming climate is a major contributor to increased annual streamflow. However, winter streamflow appears to be particularly sensitive to temperature. Whilst rising temperature can reduce streamflow via evapotranspiration, it can enhance groundwater discharge to rivers due to permafrost thawing. Currently, every 1 °C rise in temperature likely leads to 6.1%–10.5% increase in groundwater discharge, depending on the permafrost condition. For permafrost-developed basins, the contribution to increased streamflow from thawing permafrost will continue to increase in the context of global warming.


Introduction
The Arctic freshwater cycle is changing rapidly in recent decades due to climate warming (Rawlins et al 2010, Morison et al 2012, Fichot et al 2013, Lique et al 2016. Arctic rivers, which contribute approximately 40% of total freshwater into the Arctic Ocean (Serreze et al 2006, Haine et al 2015, are highly integral to the freshwater circulation in the Arctic (Bring et al 2017). Historical observations and climate model simulations demonstrated that the annual discharge of freshwater from the Arctic rivers has increased significantly over the past several decades (Peterson et al 2002, McClelland et al 2006, Overeem and Syvitski 2016. Under climate change conditions, its magnitude was projected to magnify in the 21st century (Haine et al 2015, Brown et al 2019.
The Siberia is recognized as a hot spot with substantial changes in the hydrological cycle across the Arctic (Bring et al 2017, Golubeva et al 2019. Annually, all three great Siberian rivers (i.e. Ob, Yenisei and Lena from west to east, see figure 1) carry roughly 1600 km 3 of freshwater into the Arctic Ocean (Magritsky et al 2018). Typically, these rivers are characterized by a nival hydrological flow regime with peak discharge during spring as a result of high snowmelt runoff (Yang et al 2003(Yang et al , 2007. Over the past few decades, the hydrological regime of Siberian rivers has experienced significant changes (Magritsky et al 2018), and these changes are expected to amplify further in the future (Zhang et al 2018). In addition to increase in annual river discharge, global warming and related changes in permafrost are causing a notable shift toward an earlier peak discharge in spring  (Yang et al 2002, Shiklomanov et al 2007, Tananaev et al 2016, Melnikov et al 2019a. Among the three great Siberian Rivers, the Yenisei has the largest average annual runoff (591 km 3 yr −1 ), followed by the Lena (543 km 3 yr −1 ) and the Ob (407 km 3 yr −1 ) rivers in the period 1936-2019. Evidence of increasing Siberian river discharge to the Arctic Ocean has been presented in pervious publications (Peterson et al 2002, Berezovskaya et al 2005. Earlier works suggested that the increase in precipitation induced by climate warming was the primary driver for the increase in Siberian river discharge (Fukutomi et al 2003, McClelland et al 2004, while permafrost thawing played a minor role in the observed long-term river runoff increase. Conversely, recent studies proposed the importance of meltwater from permafrost degradation (due to climate warming) in increasing the groundwater storage and river discharge across the major Siberian river basins (Lamontagne-Halle et al 2018, Mu et al 2019, Melnikov et al 2019b. Permafrost thaw can cause localized increase in groundwater storage by thickening the active layer (Walvoord andKurylyk 2016, Lamontagne-Halle et al 2018), and enhancement of regional surface-groundwater interactions as a result of extensive permafrost loss (Evans et al 2020). This can lead to more groundwater being discharged into the rivers.
By analyzing hydrological and meteorological data for the period 1936-2019, this study aims to (a) examine streamflow changes of the three largest Siberian rivers using observed daily river discharge of 1936-2019; and (b) reveal different contrasting responses of streamflow from mid-to high-latitude of the Northern Hemisphere to climate warming.

Groundwater-fed baseflow separation
Streamflow commonly originates from surface water (e.g. precipitation, snow/glacial melting, and etc) and groundwater discharge (e.g. lateral flow, permafrost thawing, and etc.). Therefore, we divide the streamflow (Q) into two components, surface-water-fed flow (Q s ) and groundwater-fed flow (Q g ). Different from Q s , Q g as a portion of relatively stable streamflow changes slightly in-season throughout a year (Tan et al 2020). However, Q g may increase slowly over years throughout the northern Eurasia (Evans et al 2020) due to enhanced groundwater discharge and (or) increased shallow groundwater diffusive recharge caused by permafrost thawing (Walvoord andKurylyk 2016, Biskaborn et al 2019). As shown in figure 2(a), when the permafrost warms up, the active layer and talik thicken (O'Donnell et al 2017). As a result, the upper boundary of permafrost starts to decline from State 1 to State 2. With subsequent increase in groundwater storage and groundwater flow (Lamontagne-Halle et al 2018), Q g gradually increases from State 1 to State 2, correspondingly, as shown in figure 2(b).
In cold regions, surface runoff during winter is negligible due to freezing conditions. Thus, winter streamflow is presumably from groundwater (St. Jacques and Sauchyn 2009). Assuming that winter streamflow is the minimum continuous stable flow throughout a year for cold regions (Paznekas and Hayashi 2016), the yearly Q g described in figure 2(b) can be estimated by the observed winter streamflow. Whilst the influence of climate warming on streamflow varies (Tan et al 2020), long-term changes in Q g can be perceived as an indicator of climate change. In this study, the winter period (ice-covered period) in Siberia is from November to April (Wild et al 2019).

Sensitivity of streamflow to climate change
Following the Budyko framework (Budyko 1948), the basin-scale evapotranspiration (E) and streamflow (Q) can be calculated as a function of precipitation (P), potential evapotranspiration (E 0 ) and a parameter that describes basin properties (n) (Roderick and Farquhar 2011). Budyko equation, derived rigorously from a single mathematical assumption concerning the Budyko hypothesis, is widely applied to estimate the climate elasticity of streamflow (Sposito 2017). The generalized form of Budyko equation can be expressed as (Choudhury 1999, Roderick andFarquhar 2011): Assuming a steady state water balance for basins larger than 10 000 km 2 , an analytical expression using first order approximation can be derived from equation (1) to quantify the sensitivity coefficients of streamflow to climate (P, E 0 ) and basin properties (n) (Roderick and Farquhar 2011): where ε P , ε E0 , and ε n are the sensitivity coefficients of P, E 0 , and n to Q, respectively.
For cold regions, climate warming induced permafrost thawing can increase winter streamflow (St. Jacques and Sauchyn 2009), which eventually contribute to the streamflow. We assume that the winter streamflow is completely contributed by groundwater and it remains constant throughout a year. With this assumption, we calculate the annual groundwater-fed baseflow (Q g ) in figure 2 from the measured winter streamflow. As mentioned earlier, changes in Q g is expected to be a result of variation in temperature, which is confirmed for the three great Siberian rivers (section 3.3). Therefore, the relationship between changes in Q g and T can be described as: where dQ g = Q g − Q g is the groundwater-fed baseflow departure, mm; dT = T − T is the temperature departure, • C; α T is the change rate of dQ g /Q g as a function of dT, • C −1 ; Q g and T are long-term mean annual groundwater-fed baseflow and temperature, respectively. Next, the changes in streamflow induced by temperature can be written as: For a large basin, Q g /Q in equation (4) represents the percentage of groundwater-fed streamflow of the total streamflow, and it is highly dependent on the surfacewater-fed streamflow (Q s ) (see figure 2(b)). Q s is a partial product of P and E 0 . Assuming other hydrological and hydrogeological conditions remain constant, Q s increases when P increases and E 0 decreases. Therefore, we propose an empirical dependence of Q g /Q on the (P − E 0 ) as follows (see example of figure S5 in supplementary files (available online at stacks.iop.org/ERL/16/034046/mmedia)): where β and γ are the empirical parameters.   (2019)) with 0.36 • C/decade, 0.35 • C/decade and 0.36 • C/decade (all p < 0.001), respectively, which are all greater than the trends in the annual mean temperature.
The multiyear mean (1936-2019) annual precipitation (P) of the Ob, Yenisei and Lena river basins were 423 mm, 418 mm, and 348 mm, respectively (figure 3(b), table S1). During the past 84 years, precipitation had positive trends in these three basins, i.e. 9.18 mm/decade, 5.30 mm/decade and 7.67 mm/decade (all p < 0.001), respectively. Increased precipitation in this region could be related to enhanced regional evapotranspiration (McClelland et al 2004) and atmospheric moisture transport (Bintanja 2018), and also the Arctic Oscillation (Frey and Smith 2003) under climate warming.
The mean annual potential evaporation (E 0 ) in the Ob River basin (620 mm) is much larger than that in the Yenisei (472 mm) and Lena (433 mm) river basins (figure 3(c)). The E 0 slightly trended upward across the Ob, Yenisei and Lena river basins with 5.23 mm/decade, 3.38 mm/decade and 2.51 mm/decade (all p < 0.001), respectively. In the period 1936-2019, the P-Q (i.e. difference between P and Q) of the Ob, Yenisei and Lena river basins were about 41%, 38% and 29% of E 0 , respectively. Similar to the multi-year trend of E 0 from west to east, the P-Q in these three basins increased with 7.7 mm/decade, 3.2 mm/decade, 2.3 mm/decade, respectively.
Overall, among these three basins, the Ob River Basin is relatively warm, has ample precipitation and the highest potential evaporation, while the Lena River Basin is characterized by the coldest climate, less precipitation and the lowest potential evaporation (figure 3(e)).

Impacts of climate changes on streamflow
Following the methods of Risbey and Entekhabi (1996) and Fu et al (2007), we calculate the annual percentage departures for streamflow (∆Q = Q−Q Q × 100%), surface-water-fed streamflow (∆Q s = Qs−Qs Qs × 100%), groundwater-fed baseflow (∆Q g = Qg−Qg Qg × 100%), precipitation (∆P = P−P P × 100%), and temperature (dT = T − T) for a specific basin and plot the results on a precipitation-temperature plane (figure 4). To better compare streamflowprecipitation-temperature relationship among different basins, we divide the precipitation and temperature departures by their corresponding standard deviations.
As shown in figure 4, it is clear that changes in Q are highly sensitive to changes in P but they are less sensitive to changes in T for the Ob and Lena river basins, which is consistent with a previous study (Xu et al 2020). By contrast, changes in Q are sensitive to both changes in P and changes in T for the Yenisei river basin. As a major component of streamflow, Q s is even more sensitive to P for all three river basins, which is likely due to direct influence of P on Q s (Ficklin et al 2016). To differ from Q and Q s , changes in Q g are positive response to changes in T, probably due to the enhanced groundwater discharge induced by permafrost thawing under a warming climate (Liu et al 2003). Additionally, Q g also appears sensitive to P for Ob and Lena river basins. This is because Q g reacts slowly to other delayed sources (such as P), although Q g usually comes from groundwater storage (Hall 1968, Eckhardt 2008. Specifically, the response of changes in Q and Q g to changes in climate variables in the Yenisei river basin is quite different from that in the other two river basins. This is probably because of reservoir regulations (Yang et al 2004, Stuefer et al 2011 and a lack of continuous discharge data for the Yenisei River during 1963-1979. Similar to Roderick and Farquhar (2011), we assume that changes in river basin storage are less relative to the magnitude of fluxes (P, E, Q) over the prebreakpoint  and post-breakpoint (1988-2019) periods of temperature, i.e. dQ = dP − dE.
Typically, values of n in equation (1) range between 0.6 and 3.6, and lower values of n denote a lower estimate of E for a given P and E 0 (Roderick and Farquhar 2011). As shown in table S3, for the Ob basin, the parameter n changes slightly from 0.99 to 1.08 from one state (pre-breakpoint) to another state (post-breakpoint). For the other two basins, the parameter n does not change between two states (n = 0.75 for the Yenisei basin, and n = 0.61 for the Lena basin), indicating the basin properties generally remain unchanged over the past 84 years (figure S3).
The decreasing values of n from the Ob, Yenisei to Lena river basins suggest that relative value of E for given P and E 0 is getting smaller from the west to east Siberian basins, which consistent with the air temperature gradient across the Siberian. Previous studies (Li et al 2013, Shi et al 2019) estimated the values of n as between 1.1 and 1.2 for all three basins, which are higher than our calculated values (especially for the Yenisei and Lena river basins). These differences could be due to the different reanalysis data of P and E 0 used in analysis.
On this basis, theoretical result of Budyko framework (table S3) quantitatively predicts a positive response of streamflow to precipitation, and a 10% increase in P would result in increases in Q by approximately 16.2%, 13.5% and 12.6% for the Ob, Yenisei and Lena rivers, respectively. However, a 10% increase in E 0 would reduce Q by approximately 6.2%, 3.5% and 2.6% for the Ob, Yenisei and Lena rivers, respectively. These results indicate that the relative impact of P on the Q is much larger than E 0 for all three basins. Besides, the response of Q to climatic factors P and E 0 decreases from West to East Siberian, indicating that streamflow in colder regions with developed permafrost is less sensitive to P and E 0 .
As mentioned earlier, changes in Q g reflects the influence of temperature changes on streamflow. Analysis of the 'naturalized' daily discharge data by removing the dam/reservoir impact shows that Q g of the Ob (1950Ob ( -2007, Yenisei  and Lena (1959Lena ( -2007 river basins changed notably with increasing trends of 1.57 mm, 4.17 mm and  (Shiklomanov 2010) for the Ob (1950Ob ( -2007, Yenisei  and Lena (1959Lena ( -2007 rivers. dT and std(dT) are the temperature departure from the average annual temperature and its standard deviation; △P and std(△P) are the relative changes in annual precipitation to the mean annual precipitation and its standard deviation; △Q, △Qs and △Qg are the relative changes in annual Q, Qs, and Qg to their mean annual values, respectively. 0.90 mm per decade (p < 0.05), respectively. Analysis of relationship between annual mean temperature departure (T − T) and changes in Q g (dQ g /Q g ) indicates that an increase of 1 • C annual mean temperature would increase Q g of the Ob, Yenisei and Lena river basins by approximately 6.1%, 10.5% and 7.8% (all p < 0.05), respectively. Equivalently, 1 • C warming would increase the annual streamflow of the Ob, Yenisei and Lena river basins by approximately 2.2%, 2.7%, and 1.0%, respectively. This suggests that a warmer climate would cause the baseflow to rise, which could lead to an increase in total streamflow in cold regions. Nonetheless, the influence of rising temperature on the streamflow is rather complex (Tan et al 2020).

Predicting streamflow for future climates
By neglecting the changes in basin properties (n), we apply equations (2), (4) and (5) and consider the effect of T, P and E 0 on streamflow to predict the streamflow of the three great Siberian Rivers (listed in table S4). Figure 5 shows the relationship between observed and simulated dQ/Q during the period 1936-2019. In general, the simulated dQ/Q is in good agreement with the observed, especially for the Lena River. The correlation between simulated and observed dQ/Q is quite low for the Yenisei River, probably related to the intense reservoir operations which altered natural changes and variations in river regime (Yang et al 2004, Stuefer et al 2011. Additionally, the lack of continuous daily streamflow data during the period 1963-1979 for the Yenisei River weakens the reliability and accuracy in streamflow simulation. Relative to the period 1980-2000, the annual mean air temperature in Siberia is likely to increase by ∼3 • C-5 • C by the end of the 21st century (Groisman et al 2013). An increase in temperature would be accompanied by increasing precipitation, intense and frequent extreme precipitation events (Burt et al 2016). To understand the further changes in Q that corresponds to climate changes, we adopt the modeling results for the highest emission scenario  average temperatures of the Ob, Yenisei, and Lena river basins would increase by 4.54 • C, 4.64 • C, and 4.71 • C, respectively; and their annual precipitations would increase by 17%, 21%, and 24%, respectively (Martynova and Krupchatnikov 2018). We use the equations described in table S4 to estimate streamflow in 2081-2100. We correlate changes in E 0 with change T (by using the correlation between annual E 0 and annual average temperature during the period 1936-2019), assuming the changes in E 0 only depend on the T changes ( figure S4).
As shown in table S5, by comparing 2000-2019 with 1981-2000, the mean annual temperatures of the Ob, Yenisei and Lena river basins have increased by 0.53 • C, 0.62 • C, and 0.71 • C, respectively. This confirms that climate warming rate increased from west to east Siberia. This is similar to the predicting scenario trend in 2081-2100 and also shows that warming effect in colder regions is much greater than the warmer regions. Observations demonstrate that over the past 20 years, streamflows of the Ob, Yenisei, and Lena rivers have increased by 5%, 2% and 9%, respectively. It is expected that the intensified climate warming to the end of 21st century would significantly increase the annual streamflow of the Lena (∼28%) and Yenisei (∼20%) rivers, and some increase in the streamflow of the Ob River (∼9%) (table S5). Our predictions are in agreement with the previous estimates, which suggested that an increase of annual mean streamflow would be up to 15% for the Ob River, 20% for the Yenisei River, and 25% for the Lena River by 2081-2100 (Groisman et al 2013).
This suggests that in the colder region (such as Lena River), where most of the basin area is covered by continuous permafrost, streamflow is more sensitive and responsive to a warmer climate (Liu et al 2003, Zhao et al 2020. Considering an empirical relationship between Q g /Q and P-E 0 (figure S5), we further estimated the changes in Q g under a warming climate for Lena River, a drainage basin with larger extent of continuous permafrost (table S1). From 1981-2000 to 2000-2019, a 0.71 • C warming resulted in an increase in Q g of approximately 2 mm, accounting for 10% of the total increase in Q (∼20 mm). When the temperature increases by 4.71 • C by the end of 21st century, a percentage of the increase in Q g (∼12 mm) to the total increase in Q (∼62 mm) is expected to hit 19% (table  S5). This suggests that the contribution of groundwater discharge to annual streamflow due to permafrost thawing in permafrost-dominated basins would continue to increase under global warming.  (Schaner et al 2012, Chesnokova et al 2020 in cold regions induced by climate warming also caused the streamflow to increase, which is evidenced by the increase in winter streamflow (Streletskiy et al 2015, Evans et al 2020. Historical data analysis demonstrated that the rate of warming is faster and the increase in Q is greater in basins with large permafrost distribution (e.g. Lena River). Our projection also indicated that the contribution of groundwater discharge to the streamflow due to permafrost thawing would continue to increase in the context of global warming, especially for the permafrost-dominated river basins. This implies that river basins with significant permafrost coverage are more sensitive to a warming climate (Liu et al 2003). Although temperature plays an important role in interannual streamflow variability (Milly et al 2018), the influence of rising temperature on the streamflow in permafrost region is complicated (Woo et al 2008, Tan et al 2020. For a relatively warm river basin with low permafrost coverage (e.g. Ob), rising temperature can increase the basin evapotranspiration, which can in turn reduce the streamflow (Shi et al 2020).

Discussion
It is worth noting that for ice-rich permafrost in the Siberian river basins, the excess ground ice provides plenty of water for subsurface flow. Meanwhile, ground ice may influence the rate and timing of permafrost thaw (Lee et al 2014). More importantly, ice-rich permafrost degradation drives widespread landscape collapse, forming numerous thermokarst (Farquharson et al 2019, Nitzbon et al 2020) with increasing surface water storage (Fedorov et al 2014). In a sense, the whole hydrological systems are dynamic and changeable (Liljedahl et al 2016). Additionally, the interaction between surface water and groundwater is largely restricted to taliks (Bense et al 2012), and open talik breakthrough most likely leads to enhanced exchange between surface water bodies and sub-permafrost aquifers (Walvoord and Kurylyk 2016).
From hydrographic perspective, this study presented the various responses of streamflow to warming climate under spatially heterogeneous climate, land-cover and hydrogeological conditions. Following a space-for-time substitution approach, streamflow in a colder basin (e.g. Lena) would evolve into that of a warmer basin (e.g. Ob). Permafrost thaw due to rising temperatures would increase streamflow (Zhao et al 2019). However, with further permafrost degradation, such positive effect on streamflow would probably be offset by the negative effect of increased basin evapotranspiration over time. This would result in a situation where runoff hit a threshold level and then decline. This situation clearly appeared in the Ob River basin, which is characterized by the highest precipitation but the lowest streamflow among the three river basins. Regarding the Yenisei and Lena rivers, the climatic conditions for the threshold of maximum streamflow to occur remain unclear, thus further research would be needed.
Our knowledge about the response of streamflow changes in large basins to climate warming is partly limited by the data shortage and data quality (Woo et al 2008), and partly constrained by the complex feedbacks that precede the permafrost-vegetation-hydrology-climate system. Given the complications to estimate global trends in evapotranspiration and uncertainties in precipitation projections, future changes of surface water budget over the Northern Eurasia remain uncertain (Ohmura andWild 2002, Liu et al 2014). Therefore, uncertainties in streamflow predictions still exist in a warm climate over Siberia, which are limited by future climate projections and our understanding of feedback between land surface hydrology and climate systems.

Conclusions
In this study, we attempted to gain insight into the response of streamflow to climate changes in the cold regions from mid-to high-latitude of the Northern Hemisphere by analyzing the long-term meteorological and hydrological data  in the three largest Siberian river basins. To evaluate the impact of warming temperatures on streamflow, a groundwater-fed baseflow (Q g ) associated with the observed winter streamflow was proposed. We assumed the changes in Q g to be the results of enhanced groundwater discharge to river due to permafrost thawing, although increases in precipitation (Bintanja 2018) could also lead to Q g to increase.
This study made an attempt to explain how climate warming leads to permafrost thawing, which further affects streamflow, ignoring the changes in topography, vadose zone structure, vegetation patterns, and surface-groundwater interactions resulting from permafrost thawing. In our view, closely linked permafrost and climate change is the key driver that alters hydrological processes in cold regions. Hence, it is necessary to strengthen critical zone observations and mathematical modeling on nonlinear freezethaw processes (Woo et al 2008, Lamontagne-Hallé et al 2020 at sub-basin scales under different climate and permafrost conditions.

Data availability statement
The data that support the findings of this study are openly available at the following DOI: https:// doi.org/10. 5281/zenodo.4430254. Science Foundation (No. 2020M670434). The authors gratefully acknowledge the anonymous associate editor and reviewers for their valuable comments and suggestions, which have led to substantial improvements over an earlier version of the manuscript.