The Discharge-Storage Relationship and the Long-Term Storage Changes of Southern Taiwan The Discharge-Storage Relationship and the Long-Term Storage Changes of Southern Taiwan

In this study, the water balance concept is used to understand the relationship between discharge and storage in a basin. Three low flow analysis models developed by Brutsaert, Vogel and Kroll, and Kirchner are used to select recession curves that are parameterized using lower envelope, linear regression, and binning methods to characterize basin hydrological behavior. Furthermore, the lowest annual groundwater storage, which is analyzed according to trend, slope, and a change point test, is used to assess the long-term storage properties of southern Taiwan and is also quantified. The water balance concept is used to assess the impact of discharge on groundwater storage that is affected by the different low flow analysis models. This can lead to a more clear understanding of the relationship between groundwater storage and discharge. After statistical tests of trend, it was determined that Chaozhou Station, which has a significant decreasing trend of the lowest groundwater storage, should implement precautionary measures such as an underground reservoir, an artificial recharge, and a collection gallery in the Donggang River Basin. The method is based on the water balance concept and low flow analysis, where precipitation and evapotranspiration in the natural river system are causing the base flow to be directly controlled by groundwater storage. Therefore, the observed streamflow can be used to assess the discharge-storage relationship and quantify the catchment groundwater storage. The results from this study show that the methods are valid for evaluating the discharge-storage relationship and long-term groundwater storage trends. concept method indicated different low flow analysis models will affect the of and groundwater The Kirchner is the under consideration.


Introduction
Salt water accounts for 97% of the global distribution of water resources, but only 3% of the fresh water, which is approximately 22% of all water on Earth, is stored as groundwater as a main obtainable resource [1]. Nevertheless, over the past half century, global warming has caused a rise in sea level that has resulted in changes in groundwater [2]. Groundwater depletion, which causes land subsidence, salt water intrusion, the necessity to irrigate agriculture, and reductions in food security, is becoming more severe year by year [3][4][5]. Globally, groundwater storage is much greater than that of atmospheric water and soil water. Because of significant differences in the storage of water resources, determining a method by which to assess groundwater storage is a very important issue [6].
Recently, Taiwan has been affected by climate change, which is resulting in distinct wet and dry seasons with stronger rainfall intensity in the wet season and a lack of rain in the dry season [7]. Variations in Taiwan's topography and uneven distribution of rainfall, with more in the mountains than in the plains and more in the northern than in the southern areas, are causing uneven distribution of water resources in time and space. Taiwan's annual average rainfall is up to 2500 mm, which is approximately 2.5 times more than the annual global average rainfall. Nevertheless, Taiwan is still facing a surface water shortage problem in the dry season. In southern Taiwan, which is the most serious area, the ratio of wet season rainfall and dry season rainfall is as high as 9:1 [8]. Therefore, in this study, the southern area is set as the study area to discuss what groundwater resource problems this area will be facing in the future. The water balance concept is used to explore the storage-discharge relationship, and a low flow analysis is used to assess the lowest groundwater storage in southern Taiwan. The assessment results can be provided to water resource agencies to assist with water resource management plans.
Low flow analysis, which is a hydrological method widely used to estimate hydraulic parameters, is used to characterize basin characteristics and assess groundwater storage trends. Brutsaert and Nieber [9] analyzed six basins in the Finger Lakes region of upstate New York using the lower envelope fitting method to characterize recession curves in order to estimate specific yield. Szilagyi et al. [10] pointed out that complexity of basin shape, hydraulic conductivity heterogeneity, and gently sloping impermeable layers do not affect the estimation of hydraulic conductivity and mean aquifer depth. Field-based estimates of these hydraulic parameters using a low flow analysis are a beneficial method for evaluating discharge data that are not easily obtained, mainly in low population density area [11].
Watershed properties such as hydrology, geology, and topography affect streamflow regimes [12]. Thus, through knowing the main factors dominating streamflow recession characteristics, it is possible to understand drainage systems. Zecharias and Brutsaert [13] selected 19 basins in the Allegheny Mountains of the Appalachian Plateaus to investigate streamflow recession curve characteristics, and their results showed that recession time is affected by mean basin slope, drainage density, and the ratio of the hydraulic conductivity and the specific yield. Vogel and Kroll [14] evaluated Massachusetts basin characteristics, and their result showed that recession time is highly correlated with basin area and basin slope. Brutsaert [15] suggested that streamflow recession time is highly correlated with basin area, basin elevation, and stream length. Kirchner [16] used a single-equation rainfall-runoff model to select streamflow data in order to estimate recession time from the headwaters of the Severn and Wye rivers at Plynlimon in mid-Wales. Sánchez-Murillo et al. [17] suggested that geology affects the length of recession time. Metamorphic and sedimentary rocks result in longer recession time. Low flow in flat basalt landscapes recesses rapidly. Knowing the dominating factors related to recession characteristics makes it easier to understand their effects on subsurface properties.
Recently, several studies analyzed annual groundwater storage trends in order to assess water supply availability in the future. Sugita and Brutsaert [18] researched perennial groundwater storage and low flow trends in the Kanto region. If the water demand remained constant, there was no sign of danger to water supply shortages. Hughes et al. [19] found that precipitation and changes in groundwater storage have a close connection. When annual precipitation is above the threshold, groundwater storage will increase and vice versa. Zhang et al. [20] selected 17 catchments to assess perennial groundwater storage trends. The groundwater storage showed different trends over the various shorter periods, and strong consistencies in the trends existed across most catchments. These results suggest that short-term climate variability can significantly influence groundwater storage and that it will be affected by largescale hydro-climate factors. A better insight into how groundwater storage has been changing will be critical in evaluating sustainable water resource management plans. The objectives of this study are (1) to discuss the relationship between discharge and storage and (2) to assess groundwater storage dynamics using low flow analysis to assess what water resource problems southern Taiwan will be facing.

Water balance concept
In a catchment, a hydrological system can be represented with the mass conservation equation: where S [LT À1 ] is the unit area water stored in the catchment; P [LT À1 ] is the rate of precipitation; E [LT À1 ] is the evapotranspiration rate; and Q [LT À1 ] is the unit area discharge.
Assuming that discharge is only related to storage, f(S) can be expressed as the storagedischarge relationship Assuming that Q is an increasing single-valued function of S since storage-discharge function is invertible The discharge change rate through time is yielded by combining Eqs. (1) and (2) Q is assumed to be an increasing single-valued function of S. Additionally, S cannot be directly measured, and thus dQ dS can be defined as a function of Q g(Q) refers to the streamflow sensitivity function because it expresses the sensitivity of discharge to changes in storage [16]. Because S cannot be directly measured in the catchment, combining Eqs. (4) and (5) to yield Eq. (6) solves the measurement problem of S, where discharge only has to be expressed as changes in storage In a catchment, when the precipitation and evapotranspiration are much smaller than discharge, the relationship between discharge and groundwater storage can be most accurately expressed Finally, the relationship between discharge and storage is derived from Eq.

Low flow analysis
In the water balance conceptual method, Eq. (8) represents the relationship between discharge and storage. Brutsaert [21], Vogel and Kroll [14], and Kirchner [16] used low flow analysis models to select recession curves and the same method has been used in this study to characterize southern Taiwan's hydrological behaviors in order to quantify the lowest groundwater storage and to assess southern Taiwan whether facing a groundwater shortage problem in the future.

The Brutsaert (2008) model
Low flow is a period without precipitation. In these periods, the hydrologic water flow system in the river is only related to the groundwater drainage from the aquifer in the catchment. As exponential decay type function is more commonly used to describe low flows where _ Q is defined as volume of flow rate in the river [L 3 T À1 ]; _ Q 0 is the value of flow rate at arbitrarily chosen time such as t = 0 [L 3 T À1 ]; and K, which is defined as the characteristic of time scale in the catchment, is also referred as the storage coefficient [T].
A power law function is used to describe the relation between the change in the flow rate and the discharge characterized hydrologic behavior in this research (e.g., Brutsaert and Nieber [9] in which a and b are constants [À]. These constants can be characterized to display catchment hydrologic behavior [22].
Precipitation, evapotranspiration, and artificial input events are excluded in the low flow conditions when groundwater storage is only related to the water flow in the river. Furthermore, mass conservation law of hydrologic system can display as the following integral function: where S and Q are defined as Eq. (1). Substitution of Eq. (9) in Eq. (11) derives following relationship between groundwater storage and drainage from the catchment According to Eq. (12), the temporal trend of catchment groundwater storage can be obtained from the trend of discharge in the river as the following function: Because daily flow varies, in this study, the annual lowest 7-day daily mean flow is used to replace daily flow to reduce uncertainty where Q L7 is the annual lowest 7-day daily mean flows [LT À1 ].
The long-term temporal trend of catchment groundwater storage is significant because it reveals whether the catchment will face a water shortage problem in the future.
Obtaining the trend in the catchment groundwater storage requires the characteristics of time scale K. In this study, a lower envelope was adopted, which is defined as roughly 5% data points below the logarithmic plot of change of flow rate data against the flow rate. In order to delete the effectiveness of precipitation, the criteria for the selection of streamflow data are as follows [21]: a. Changes of flow rate are zero and positive, and sudden anomalous data should be eliminated.
b. Three data points after the changes in the flow rate are zero and positive, and four data points after major events should be eliminated.
c. Two data points before the changes in the flow rate are zero and positive should be eliminated.
d. When there is suddenly a higher value of change in the flow rate in a dry period, it should be eliminated.

The Vogel and Kroll (1992) model
Vogel and Kroll [14] assumed discharge to be linearly related to storage. Eq. (15) is a wellknown base flow equation where K b is base flow recession constant [À]. Eq. (16) is derived by Eq. (9) by comparing the Brutsaert (2008) model Vogel and Kroll's model used a linear regression to define the characteristics of time scale K to represent average basin hydrological behavior. The Vogel and Kroll model has criteria for selection of the streamflow data as follows: a. A recession range is between a 3-day moving average beginning to decrease and a 3-day moving average beginning to increase.
b. Delete 30% of the front total recession length.

The Kirchner (2009) model
Kirchner [16] suggested that when evapotranspiration, precipitation, aquifer recharge, and leakage are much smaller than discharge, they can be neglected. At this time, discharge is directly related to changes in storage as shown in Eq. (17) below: The power law equation proposed by Brutsaert and Nieber [9] is used to represent hydrological low flow analysis. Kirchner's model for streamflow criteria is simpler than the other available alternatives. It only requires the streamflow recession segment.

Mann-Kendall test
Mann-Kendall (MK) [23,24] test is nonparametric method developed from Kendall's tau (τ). It can be used to test the relationship between two sets of data. The advantages of this method are that the extreme values and missing data problems will not seriously affect the certification value. The MK test assesses the trend in a series via comparing the value of the series before and after to determine whether the series exhibits a specific degree of trend. The null hypothesis given that if there is not significant trend in the series, test statistic S is defined as: where {X 1 , X 2 , X 3 ,…, X n } is streamflow data, which are arranged in accordance with time {T 1 , T 2 , T 3 ,…, T n }, n is the number of data. When n is close to infinity, the probability of the S distribution curve will present as a normal distribution with a mean of 0. In addition, when n is more than 10, the variance of S can be substituted into the following approximate solution: In this study, long-term streamflow data are likely to be repeated in the data series; thus, Kendall modified the approximated solution Eq. (19) to Eq. (20) where u is the duplicate value number of the data series.
Finally, the normalized statistical test S values become the Z value, as follows: When Z is a positive value, this indicates that the series is exhibiting an increasing trend; in contrast, when the value is negative, it indicates that the series has decreased. At this time, the obtained Z value should be tested by a significance test to assess whether the series is significant. Assuming, a significance level of α, if |Z|≧Z α , the null hypothesis is rejected, which represents that the series has a significant trend; whereas, the series has no significant trend. In the study, the significance level is set as α = 0.05. When |Z|≧1.96, the series has a significant trend. When it is below this level, there is no significant trend.

Theil-Sen slope
The Theil-Sen slope [25] is used to estimate the magnitude of the trend slope. The Theil-Sen slope estimation method is different from the slope values calculated using a linear regression, because it selects the median value, and therefore, the properties are less affected by extreme values. Thus, it is often used with the MK test. Slope β is defined as follows: The  [26] can be used to search for significant change points in a data series. The definition of a change point is when a data series {X 1 , X 2 ,…, X n } has a change point X t , order {X 1 , X 2 ,…, X t } is F 1 (X) and {X t + 1 , X t + 2 ,…, X n } is F 2 (X), then F 1 (X) 6 ¼ F 2 (X). The definition of U t, n is as shown in Eq. (24). If there is not a change point in the data series, U t, n j jon the function of time, t will continue to rise, and there will be no turning point. On the contrary, if there is a change point, U t, n j jon the function of time t, there will be a decreasing turning point. In the same data series, the turning point may occur several times on behalf of this data series, and there may be more than one change point To confirm that change points exist, Eq. (25) is used to calculate the extreme value of U t, n j jthat is turning point as K n . Eq. (26) is used to calculate the probability of a change point. In this study, P = 0.95 is set as the confidence level, where P > 0.95 judges that the time is a significant changing point However, in some data series, a change point may not exist by itself; thus, Eq. (27) is used to calculate each year's P(t) value. The P(t) value is identified when it is greater than the confidence level

Study area
The study area is set in southern Taiwan, including Chiayi, Tainan, Kaohsiung, and Pingtung area. The main basin includes the Bazhang River, the Jishuei River, the Cengwen River, the Yanshui River, the Erren River, the Gaoping River, the Donggang River, and the Linbian River (see Figure 1). The annual average temperature is between 22 and 26 C, and the annual average precipitation is about 2000 mm, with distinct wet and dry seasons. The precipitation is concentrated in the May to October high flow periods, and on the contrary, November to April is the low flow period, which causes the rivers to dry up in southern Taiwan [8]. Two standards are used for the selection of the streamflow station. First, a streamflow station must provide long-term continuous measurement data; thus, selected stations have at least 40 years of daily flow data. Second, streamflow stations are unaffected by artificial structures such as dams. After screening, Taiwan's southern basin was determined to have six stations complying with the above two criteria. These stations are listed in Table 1 and are located in Figure 1.  In this study, the water balance method is used to discuss the relationship between discharge and storage in six catchments in southern Taiwan using three low flow models from Brutsaert [21], Vogel and Kroll [14], and Kirchner [16] to select streamflow data. However, the groundwater storage cannot be measured directly at the catchment; thus, the streamflow sensitivity function (Eq. (6)) is used to assess the sensitivity of discharge to changes in storage.
First, the storage-discharge relationship in southern Taiwan according to the Brutsaert model is discussed. The results show that the most sensitivity is at the Chungde Bridge Station, and the least is at Chaozhou Station. The order of sensitivity is as follows: Chungde Bridge at the Erren River, Changpan Bridge at the Bazhang River, Laonong at the Gaoping River, Xinbei at the Linbian River, and Chaozhou at the Donggang River. However, the Xinshi catchment cannot exhibit the relationship between storage-discharge because the process of selecting streamflow data using the Brutsaert model is more complex than other alternatives. Thus, the streamflow data for Xinshi Station, which cannot reveal the storage-discharge relationship, are less than that for the other stations (Figure 2(a)). According to Eq. (8), when the discharge sensitivity to changes in groundwater storage is high as indicated by the slope of the storage-discharge relationship, the recession time is shorter. The time of the order of recession is Chaozhou at the Donggang River, Xinbei at the Linbian River, Laonong at the Gaoping River, Changpan Bridge at the Bazhang River, and Chungde Bridge at the Erren River as shown in Figure 3(a). However, Xinshi Station cannot reveal the storage-discharge relationship, which causes the relationship between discharge and time to be inaccurate as well.
The results for the Vogel and Kroll model show that the Chungde Bridge station cannot display the storage-discharge relationship because the Q t ≥ 0.7Q t-1 selection criteria resulted in discontinuous selected data. Figure 2 shows that for the same changes in discharge, the Kirchner model has higher storage variations and shorter recession periods.

Low flow analysis
In this study, low flow analyses of the Brutsaert [21], Vogel and Kroll [14], and Kirchner [16] models are fitted by lower envelope, linear regression, and binning to parameterize streamflow recession curve and discuss the relationship between parameters (e.g., hydrological characteristic constants and recession time) and topography factors. Stream order, average elevation, length of main stream, average slope, slope of main streamflow, and basin area are selected and shown in Table 1 to analyze the correlation between hydrological characteristic constants and recession time (K).
In order to ensure that low flow is a period without precipitation, thus, the Brutsaert [21], Vogel and Kroll [14], and Kirchner [16] models have their own selection procedure criteria for streamflow data. Figure 4 shows that different models will affect the selection of streamflow  data points; thus, the streamflow data points of the Kirchner model are approximately 40 times than the Brutsaert model. Those results in assessing the basin hydrology characteristics that are more subjective. Table 2 shows that the characteristic constant a order is binning > linear regression > lower envelope. Generally, the characteristic constant a using the binning and linear regression has higher values due to the fact that the regression line is moved upward when compared to the lower envelope fit line. Table 3 shows that the recession time order is lower envelope > linear regression > binning. This indicates the inverse relationship between the characteristic constant a and recession time. The six basins with the distribution of recession time show that the binning fitting method exhibits the most quickly recessions, with ranges from 3 to 22 days. The lower envelope fitting method had the longest recessions, with ranges from 32 to 127 days. Because the average value will be influenced by the outlier median is more suitable for representing recession time. It can be noted in Figure 5 that lower envelope fitting with the   Table 2. Characteristic constant a for three extraction procedures (BRU, VOG, and KIR) and three different fitted models (lower envelope, linear, and binning).
Brutsaert model obtained a median of 55.5 days; the Vogel and Kroll model obtained a median of 60 days, and the Kirchner model obtained a median of 45.5. Figure 6 shows a characteristic constant a range from 0.008 up to 0.355. Figure 7 shows that the characteristic constant b varies from 0.678 to 1.432 without outliers.
In the results of the correlation between recession curve and topography factors, the Vogel and Kroll model, which is fitted using a linear regression, is the most suitable result for characterizing southern Taiwan basin hydrological behavior. The recession time is highly correlated with average elevation, length of main stream, and basin area as shown in Table 4. According to Tague and Grant [27], a higher characteristic hydrology constant a represents rapid recession. In this study, only the characteristic hydrology constant a from the Vogel and Kroll model as fitted using a linear regression is found to be highly correlated with the slope of the main stream as shown in Table 5. There is a high correlation the between characteristic hydrology constant b and average slope according to the Vogel and Kroll model, which is fitted using a linear regression, as shown in Table 6. According to above results indicating that the Vogel and Kroll model as fitted using a linear regression is the most optimal result to represent  Changpan Bridge  61  57  39  16  10  5  15  13  3   Xinshi  39  63  46  16  10  5  15  9  3   Chungde Bridge  35  75  32  12  9  4  12  10  3   Laonong  127  102  45  30  90  13  21  15  14   Chaozhou  85  57  83  23  29  14  22  17  12   Xinbei  50  53  60  17  19  6  16  13  4   Table 3. Recession time for three extraction procedures (BRU, VOG, and KIR) and three different fitted models (lower envelope, linear, and binning) (day). southern Taiwan basin characteristics, the Vogel and Kroll model is used to quantify the lowest groundwater storage in southern Taiwan.

The trend test of lowest groundwater storage
The optimal results of the Vogel and Kroll model as fitted using a linear regression is shown in Section 4.2. It is used to quantify the lowest groundwater storage in southern Taiwan. In the trend test results shown in Table 7, Changpan Bridge, Chungde Bridge, and Xinbei Station exhibit the most significant increasing trend. Chaozhou Station exhibits a significant decreasing trend. Xinshi Station has an increasing trend. Laonong Station has a decreasing trend. Figure 8 shows the perennial basin slope test indicating that the increasing trend range for  changes in the lowest groundwater storage ranges from 33.5% to 2.7 times, and the decreasing trend ranges from À16.3% to À96.9%.
In this study, the change point is primarily set from 1983 to 2004 as shown in Figure 9. The lowest groundwater storage is smaller than before the change point at Chaozhou Station and greater than before the change point at Changpan Station, Chungde Station, and Xinbei Station. The lowest groundwater storage before the change point is 0.65 mm, and after the change point, it is 4.07 mm. The difference between before and after is up to 5.3 times than that of the Xinbei Station. 1980s, the average global temperature increased rapidly, resulting in abnormal weather phenomena. Additionally, Lu et al. [29] pointed out that after the 1980s, the temperature increased significantly by about twice than the century before. These studies show that there was a change in Taiwan's climate in the 1980s, thus suggesting that the impact of climate change on Taiwan's hydrological conditions caused a significant change point to occur in the 1980s. Fan et al. [30] indicated that the frequency of extreme rainfall events related to typhoons increased from 2000 to 2009, occurring once on an average 3-4 years between 1970 and 1999 and an average once every year after 2000. According to the research of Tsuang et al. [31], an annual average of 3.3 typhoons occurred in the twentieth century. However, in the research of Tu et al. [32], typhoon frequency increased up to an annual average of 5.7, causing significant increases in rainfall and affecting changes in groundwater storage.
In this study, the Brutsaert [21], Vogel and Kroll [14], and Kirchner [16] models are fitted using lower envelope, linear regression, and binning methods to parameterize the recession curve.  The results show that the Vogel and Kroll model is the most optimal model to describe hydrological behavior in southern Taiwan. Thus, it is selected to quantify the lowest groundwater storage in the six basins. According to the trend, slope, and change point tests, the lowest groundwater storage exhibits a decreasing trend at the Laonong and Chaozhou stations. Most notably, Chaozhou Station exhibits a significant decreasing trend. In the future, shortage problem of Chaozhou station has a greater probability of facing low groundwater storage than the other basins under consideration. Moreover, in the study, the storage-discharge relationships in southern Taiwan are assessed using the water balance concept. All of the above hydrological studies in the long-term in southern Taiwan can be provided to water resources authority as a reference for water resource management planning.

Conclusions
This study examined the discharge-storage relationship and the changes in the lowest groundwater storage in six selected southern Taiwan basins over the last 40 years using observed daily streamflow. The method is based on the water balance concept and low flow analysis, where precipitation and evapotranspiration in the natural river system are ignored, causing the base flow to be directly controlled by groundwater storage. Therefore, the observed streamflow can be used to assess the discharge-storage relationship and quantify the catchment groundwater storage. The results from this study show that the methods are valid for evaluating the discharge-storage relationship and long-term groundwater storage trends.
The water balance concept method indicated that different low flow analysis models will affect the sensitivity of discharge and groundwater storage. The Kirchner model is the most sensitive model. Thus, the recession time in this model was shown to be significantly shorter than that of the others under consideration. To summarize, the Vogel and Kroll model as fitted using a linear regression is most the optimal model by which to represent the hydrological characteristics in southern Taiwan. Overall, the lowest groundwater storage exhibited a significant decreasing trend at Chaozhou Station. The assessment results can be provided to water resource agencies to assist with water resource management plans.