Baseflow recession characterization and groundwater storage trends in northern Taiwan

Climate change has led to detectable variations in the water cycle (Collins et al. 2013; Gao et al. 2015), which influence the way people utilize water. Therefore, water resources management has become a global issue. Oceans are the primary water resource since they account for 97% of the water distribution on Earth; however, salt water cannot be directly used by humans. Terrestrial freshwater in liquid form includes surface water, groundwater, and soil moisture, and that in solid form includes glaciers, ice sheets, frozen lakes, snow, and permafrost (Collins et al. 2013).

Taiwan is an island in East Asia, located in the subtropical zone along the southeast coast of Asia. Although almost two-thirds of Taiwan is covered by mountains, the island's western side is mostly alluvial plain, where most of its population resides. About 80% of water consumed in Taiwan is supplied by surface runoff in combination with surface storage reservoirs, with the remaining water supplied by groundwater (Yang 2010). Based on the observation data measured by weather stations all over Taiwan, the rainfall distribution decreases gradually from the northeast to the southwest, forming several closed rainfall centers in mountainous areas. Although the long-term mean annual precipitation in Taiwan is up to 2,510 mm, which is approximately 2.5 times more than the global average, water resources available per capita in Taiwan are only one-fifth of the world average owing to uneven precipitation distribution, topography, and high population density (Water Resources Agency 2000). The spatial and temporal heterogeneity of rainfall distribution leads to an apparent distinction in the surface water quantity between wet and dry seasons (Figure 1). Rainfall in the wet season accounts for 70% of annual precipitation. The Intergovernmental Panel on Climate Change assessment report (Collins et al. 2013) indicated that the contrast of mean precipitation between wet and dry seasons will increase over most of the globe as temperatures increase. Increased heating leads to greater evaporation and increases the intensity and duration of drought. At the same time, the water holding capacity of air increases by about 7% per 1 °C warming, resulting in more intense precipitation events (Trenberth 2011).

Despite surface water still acting as the main water resource in Taiwan, there are some problems that need to be dealt with, including reservoir sedimentation, water contamination, and unstable precipitation caused by climate change. An unstable water supply has a significant impact on water resources utilization. Thus, as a relatively stable water source, groundwater storage is a critical component of the terrestrial water budget (Brutsaert 2008). The groundwater monitoring network of Taiwan was established in the 1990s, but most of the observation wells have only recorded water level data during the last decade. Very few reliable records are available of groundwater storage that are long enough to analyze trends or quantity. Even with some direct records from some observation wells, the scale problem has to be considered. The measurements of water level in isolated wells represent the local scale, whereas water resources management requires information at the catchment scale (Zhang et al. 2014). Unlike groundwater level observations, streamflow data have been recorded by the Taiwan Water Resources Agency since the 1950s. These long-period monitoring data are relatively complete and available. Since the hydrometric records represent the integrated behaviors of upstream catchment above the gauge station, understanding the river–groundwater interactions is important for sustainable water resources management (Cai et al. 2016). Recession analysis methods are widely used to investigate the storage–outflow relationship of a catchment area and characterize the features of a specific catchment through the recession curve (Tallaksen 1995; Biswal & Maranni 2010; Mutzner et al. 2013).

In this study, we use flow recession analysis known as the Brutsaert method to parameterize the baseflow recession characterization in northern regions of Taiwan. This method was developed by Brutsaert & Nieber (1977) based on theoretical solutions of the Boussinesq equation and groundwater hydraulic theory. The International Glossary of Hydrology defines low flow as ‘flow of water in a stream during prolonged dry weather’ (Smakhtin 2001). This flow is also referred to as baseflow. During rainless periods, streamflow in natural streams is sustained mostly by the cumulative outflow from upstream riparian aquifers. The Boussinesq equation can be applied to the problem of outflow into a fully penetrating stream channel from an unconfined rectangular aquifer placed on a horizontal impermeable layer (Troch et al. 1993). Therefore, the temporal trends of the annual low flows were then directly calculated by simple linear regression, from which the trends of the annual groundwater storage were estimated (Brutsaert 2008). The long-term groundwater storage changes in our study area can be inferred from streamflow records under several assumptions and data processing procedures.

The Brutsaert method has been successfully applied to many basin-scale studies with different climatic and physiographic conditions. The initial derivation of this method began from the problem of outflow from an unconfined aquifer placed on a horizontal impermeable layer into a stream channel. Brutsaert & Nieber (1977) obtained the baseflow coefficient and related it to the geomorphological characteristics, including drainage area and stream length, in the Finger Lakes region in upstate New York. Basins with arid and mountainous conditions, which seem not to meet the theoretical equations, still successfully predicted groundwater hydraulic parameters that are consistent with field measurements (Mendoza et al. 2003). In the cold regions, permafrost thawing causes the thickness of active groundwater layer change. The rate of permafrost thawing can be quantified by measuring low flows during the open water season in the upstream river basin. Thus, Lyon et al. (2009) and Brutsaert & Hiyama (2012) utilized this characteristic to determine the permafrost thawing trends. The groundwater storage changes inferred from the baseflow were compared with trends of related hydrologic variables, such as average groundwater level changes measured in an observation well network (Brutsaert 2008; Zhang et al. 2014), changes in sea surface temperature in the tropical Pacific Ocean (Gao et al. 2015), annual average water inflow rates (Sugita & Brutsaert 2009), and precipitation over the same area (Sugita & Brutsaert 2009; Hughes et al. 2012). Several studies have discussed the variability of groundwater storage and catchment characteristic parameters caused by climate change and land use change (Gao et al. 2015; Bogaart et al. 2016).

The climate of the present study area is much more humid than those in the research mentioned above. Precipitation with significant seasonal differences greatly affects the stability of a water supply. In 2015, Taiwan had the most severe drought damage in nearly 67 years (Yang 2015). The northern region is the major megalopolis in Taiwan and accounts for 45% of the population. Therefore, investigating the situations and trends of basin water resources is important. The novelty in this study is the application of the Brutsaert method to explore groundwater resources and catchment characteristics in a sub-tropical climate zone with high population density. The specific objectives of this study are: (1) to evaluate the main topography factors controlling the baseflow coefficients and characteristic drainage time scale and (2) to estimate long-term groundwater storage trends in northern Taiwan from streamflow.

Taiwan is an island on the Pacific Ocean with the geographic coordinates of 25 °03′N latitude and 121 °30′W longitude. The terrain in Taiwan is divided into two parts, namely, flat plains and rugged forest-covered mountains. The Central Mountain Range extends from the northeast to the southern tip of the island, forming a ridge of mountains that acts as the principal watershed.

The northern region of Taiwan has four major river basins, namely Lanyang River, Danshui River, Fengshan River, and Touqian River (Figure 2). The Water Resources Agency (WRA) has been monitoring river discharge in these areas since the 1950s. In order to analyze the long-term groundwater storage trends and robust characteristic drainage time scale, eight gauge stations were selected for which daily streamflow records were available. The four river basins have at least 24 years of monitoring daily streamflow and rainfall records. Table 1 shows the characteristics of the study basins. The area of the subbasins ranges from 115.93 to 820.69 km². The length of river channels ranges from 18.42 to 74.95 km. These regions have varied terrain characteristics; the terrains range from less than 100 m in elevation (plains), to less than 500 m in elevation (hills), to more than 1,000 m in elevation (mountains) (Yeh et al. 2015). The variability of topography results in significant main channel slope distinctions between upstream and downstream subbasins. The climate of this area is warm and humid; monthly mean temperatures are typically 15–18 °C in winter and 27–29 °C in summer. The precipitation data were recorded by automatic weather stations near the gauge stations mentioned in Table 1. Rainfall in the northern regions of Taiwan during the summer and winter seasons primarily comes from frontal rain, typhoons, and northeast monsoons. The average annual rainfall in Lan-Yang Bridge is up to 2,932 mm. Due to the orographic rain caused by the northeast monsoons, long-period mean rainfall values in the northeast region are higher than those in the northwest.

As suggested by Brutsaert (2008), data points that meet the following criteria were excluded from the analysis:

1. Data points with positive or zero values of dy/dt.

2. Two data points before the changes in the flow rate that are zero and positive; 3 data points after the changes in the flow rate that are zero or positive.

3. Four data points after major events (flow duration curve is used to determine the major event (Kingsford et al. 2002)).

4. Suddenly anomalous data points with −dy/dt that is much larger than the previous one.

5. Runoff on days of rainfall, and a number of days after the rainfall event.

Evapotranspiration during recession flow tends to accelerate the decline of the streamflow rate, and leads to higher values of dy/dt. Also, the rate of decline of groundwater storage from a basin phreatic aquifer is lower than that displayed by the other outflow components, which result from overland flow, interflow, or channel drainage. This means that low flows during periods of minimal evapotranspiration can be regarded as the smallest dy/dt for a given flow rate y, corresponding to the slowest recession rate, which is the lower envelope of the data cloud in the recession plot. Because the data points are scattered, the determination of the exact position of the lower envelope is subjective. Hence, the position of the lower envelope was determined such that 5% of the flow data points are below it to make some allowance for the unavoidable error suggested by Troch et al. (1993). In our study, with fixed slope b, we gradually calibrated to allow 5% of data points to lie below the fitted envelope line, then obtained a.

The lower envelope can be analyzed using fixed slope b (equal to 1), representing a linearized outflow from a rectangular aquifer into a fully penetrated channel. Thus, the reciprocal of constant a is equal to the characteristic drainage time scale value.

Raw streamflow records were filtered and treated using the procedure mentioned in the previous section. Figure 3 shows the logarithmic recession plot of data points −dy/dt against the corresponding y data points for eight gauge stations in four major basins of northern Taiwan. With the lower envelope fitting lines with slope 1 (i.e., b = 1), the baseflow coefficient a for each subbasin can be determined. As shown in Equation (6), the reciprocal of constant a can be considered as a characteristic drainage time scale value in the long-time outflow rate. The values of coefficient K for the eight selected gauge stations varied from 34 days to 84 days (Table 2). The mean value is 54 days and standard deviation is 16 days.

Brutsaert (2008) inferred that the e-folding drainage time scale is relatively insensitive for large basins and that the K value is approximating 45 ± 15 days, which is a typical value for a working assumption. Compared with previous studies (Brutsaert & Sugita 2008; Zhang et al. 2014; Sánchez-Murillo et al. 2015), the results of the present study show a higher average magnitude of K with a similar uncertainty. Power law coefficient a depends on the initial storage in the basin. If a basin is more wet, coefficient a will be smaller (Biswal & Kumar 2014; Sánchez-Murillo et al. 2015). The humidity of the basin is related to baseflow coefficient, which could be the reason that northern Taiwan shows a higher average magnitude of K. The long-term mean annual precipitation in Taiwan is approximately 2.5 times more than the global average. Larger mean annual precipitation makes a river basin more wet, and the inverse relationship of a (i.e., K) tends to be higher than in other regions of the world.

Baseflow recession characteristics show apparent regional differences in our study area. The basin on the northeastern side (Lanyang River) displayed a longer drainage time scale than those of basins on the northwestern side (Danshui River, Fengshan River, and Touqian River). Since the mean annual rainfall of Lanyang River is about 690 mm/year higher than that in the other three basins in the northwestern region (Table 1), we believed that initial status of the basin is more wet, and caused a higher drainage time scale. The results also show that even in a given major basin, the physiographic differences between upstream subbasins and downstream subbasins can produce distinct a and K values. The characteristic drainage time scale mean value of San-Hsia, Hsin-Pu, and Jein-Kuo Bridge gauge stations (located in the downstream region of the study area) is 36 ± 2 days; for gauge stations located in the upstream region, namely Xiu-Luan, Shang-Ping, and Nei-Wan, the mean value is 62 ± 1 days. We utilized the results of K to correlate with physiographic factors of each subbasin (excluding the two stations at Lanyang River basin which show significantly different magnitude of K value with other stations) by Pearson's correlation. We observed that K values are strongly correlated with the main channel slope (r = 0.90) (Figure 4). The resulting r value might be largely uncertain due to fewer points in the regression. However, this result could still be interpreted for the catchment characteristics of northern Taiwan. Based on the calculation of the K average value, the uncertainty of parameters is limited in a narrow range of less than 2 days. These data showed strong consistency. This means that the storage coefficient has a high consistency under similar geomorphic conditions. The smaller drainage time scales were found in flatter plain subbasins. Gauge stations situated in mountainous areas had larger K values. Unlike past studies (Brutsaert & Sugita 2008), K had a poor correlation with drainage area (r= −0.44) and length of channel (r= −0.45) for the northern Taiwan basins.

In this analysis, the groundwater storage trends were determined to have a linear relationship with discharge in the river. The running averages of the annual series of minimum 7-day flows are regarded as dry weather flows, denoted as y_L7. Compared with the lowest 1-day flow value, y_L7 is less affected by measurement errors (Smakhtin 2001). The running average of minimum 7-day flows was applied because it indicates the lowest groundwater storage level in each year, which can be carried over to the next year (Gao et al. 2015). This means that the quantity value can be used to find out the groundwater storage long-term evolution in a basin-scale study. In the application of Equation (8), we used the results of the characteristic drainage time scale calculated from previous sections. With the product of K and y_L7, the groundwater storage for each subbasin can be obtained. Then, trends of annual groundwater storage can be estimated via the Theil–Sen slope method. The results are listed in Table 3 and Figure 5. In the present study, two different forms of temporal length records were analyzed, namely the period of the entire record of each individual station and the period of the 21st century (i.e., most recent 15 years) of all selected stations. In the period of the entire record time, the results are relatively invariable. The trends of groundwater storage for most subbasins showed merely slight changes, with an average value of 0.103 ± 0.030 mm/year, with exceptions being Chia-Yuang Bridge (0.950 mm/year) and Nei-Wan (0.646 mm/year). However, the variable amplitude of groundwater storage significantly increased in the most recent 15 years. Compared with results from all recording time analysis, the amount of changing values are 1.4 (Nei-Wan) to 20 (Lan-Yang Bridge) times higher at all subbasins. The Chia-Yuang Bridge subbasin at the Lanyang River basin in the northeast of the study area exhibited the strongest downward trend (4.392 mm/year). Three subbasins, namely, Nei-Wan, Shang-Ping, and San-Hsia, had opposing trends between different analysis periods. This phenomenon is commonly seen in areas around the world (Brutsaert 2008; Sugita & Brutsaert 2009). Different land use or regional climate fluctuation influences the groundwater recharge and storage state.

Over the period of 2000–2014, the groundwater storage trends in northern Taiwan had an obvious spatial distribution (Figure 6). River basins with significant negative trends in groundwater storage are located in the northeast of the study area (Chia-Yuang Bridge and Lan-Yang Bridge subbasins). In contrast, the subbasins in the northwest, namely Xiu-Luan (0.245 mm/year), San-Hsia (1.454 mm/year), Hsin-Pu (1.519 mm/year), Shang-Ping (1.193 mm/year), Nei-Wan (0.929 mm/year), and Jein-Kuo Bridge (0.326 mm/year), all showed positive trends in groundwater storage. In the present study, we compared these results with trends of related hydrologic variables, i.e., annual precipitation and annual lowest groundwater table level. Groundwater storage quantity also displays spatial differentiation in this area. The groundwater storage we mentioned here is the water capacity over field capacity in the phreatic aquifer. The values of groundwater storage for the Lanyang River basin are larger than those for the other three major basins in the northwest of Taiwan. Based on a groundwater geological map in a WRA report (Water Resources Agency 2014), the Lanyang River basin is dominated by alluvial deposits, which contain abundant groundwater. The Danshui River basin is dominated by limestone and some alluvium. The Fengshan River basin and the Touqian River basin are dominated by conglomerate and clasolite. The results show that the lithological characteristics of river basins is the main factor that determines the groundwater storage quantity.

The Brutsaert method, based on the concept that baseflow in a natural drainage system is directly controlled by groundwater storage, can provide estimates of basin-scale recession characteristics and evolution of groundwater storage. The advantage of this method is that it uses streamflow data, which are more easy to obtain for longer periods of record than groundwater level data in Taiwan. In this study, we calculated the catchment characteristics and groundwater storage trends in northern Taiwan using streamflow measurements in a number of major river basins. With a high-density population, the study area is under the influence of human activity. Eight selected subbasins with gauge stations that recorded at least 24 years of complete daily streamflow data were used to represent stable basin features. The analysis showed that the characteristic drainage time scale value K varied between 34 and 84 days, with a mean value of 54 days and a standard deviation of 16 days. Compared with previous studies, the results show a higher average magnitude of K, but the uncertainty is consistent with estimates of roughly 15 days by Brutsaert. K values in upstream subbasins and downstream subbasins reflect an obvious spatial pattern. The characteristic drainage time scale mean value for subbasins located in the downstream regions of the study area is 36 ± 2 days, and that of those located in the upstream region is 62 ± 1 days. This finding means that the storage coefficient has a high consistency under similar geomorphic conditions. Our study shows that the K value is strongly correlated with the main channel slope. Over the period of record, groundwater storage of five subbasins showed downward trends ranging from −0.053 to −0.950 mm/year; another three subbasins exhibited upward trends ranging from 0.111 to 0.141 mm/year. Over the period of 2000–2014, the groundwater storage trends in northern Taiwan had an obvious spatial distribution. River basins with significant negative trends in groundwater storage are located in the northeast of the study area. In contrast, subbasins in the northwest all showed positive trends in groundwater storage. Since groundwater is a stable water source in Taiwan, the trend analysis results can be used as reference for water resources management.