## Abstract

This research utilized Bayesian and quantile regression techniques to analyze trends in discharge levels across various seasons for three stations in the Gorganroud basin of northern Iran. The study spanned a period of 50 years (1966–2016). Results indicate a decrease in high discharge rates during springtime for the Arazkouseh and Galikesh stations, with a steep slope of −0.31 m^{3}/s per year for Arazkouseh and −0.19 and −0.17 for Galikesh. Furthermore, Tamar station experienced an increase in very high discharge during summer, with a slope of 0.12 m^{3}/s per year. However, low discharge rates remained relatively unchanged. Arazkouseh station showed a higher rate of decreasing discharge levels and this trend was most prominent during spring. Additionally, the Bayesian quantile regression model proved to be more accurate and reliable than the frequency-oriented quantile regression model. These findings suggest that quantile regression models are a valuable tool for predicting and managing extremely high and low discharge changes, ultimately reducing the risk of flood and drought damage.

## HIGHLIGHTS

For better flood management, high discharge changes are estimated via quantile regression model.

Quantile regression models are suitable for analyzing changes in different ranges of data series.

Bayesian quantile regression is more accurate than frequency-oriented quantile regression models.

High rates of trend variation are visible in extreme values of discharge.

## INTRODUCTION

River discharge is a crucial factor in hydrology and water resources, and it is affected by various climatic factors like temperature and precipitation. Any change in these factors can lead to alterations in the discharge trend (Koutalakis *et al.* 2019; Wei *et al.* 2020). Reducing runoff poses a significant challenge in managing water resources, especially in arid and semi-arid regions such as Iran (Kanani *et al.* 2019). Lower rainfall can decrease runoff and subsequently reduce water resources (Sabbaghi *et al.* 2020). Conversely, heavy rainfall can result in a sudden increase in discharge, leading to economic and social damage in sectors like agriculture, livelihood, and infrastructure. This can seriously disrupt regional development plans (Baba'eeyan *et al.* 2013). Therefore, studying the changes in river discharge and runoff volume due to rainfall can help predict flood time and reduce damages caused by it. Recently, many basins in Iran have undergone significant changes, with a decrease in the amount of runoff and river flow, or in some cases, an increase in the form of floods (Sharafati *et al.* 2020; H. Ahmadi *et al.* 2020a; M. Ahmadi *et al.* 2020; Kousali *et al.* 2022). Trend analysis is a suitable model to evaluate hydro-climatological conditions in basins, as it helps study changes in a variable over time (Binesh *et al.* 2017). Most studies on trend detection use the Mann–Kendall non-parametric test, which was developed by Kendall (1975) and Mann (1945). This method has been used to investigate trends in temperature by Ahmadi *et al.* (2018) and Alhaji *et al.* (2018), precipitation by Basarir *et al*. (2018), Nyikadzino *et al.* (2020) and Salehi *et al.* (2020), flow by Ghorbani *et al*. (2016), Tehrani *et al*. (2019), Galkate *et al.* (2020), Norouzi (2020), Kuriqi *et al.* (2020), and Owolabi *et al.* (2021), and droughts by Kalisa *et al.* (2020) and Ghorbani *et al.* (2020). However, these studies focus on the average and median of the data, which might not capture trends that occur within a specific range of the data. Additionally, the mean alone is not sufficient to fully understand the distribution pattern of the variable being studied at different levels of the independent variable, as noted by Bararkhanpour *et al.* (2020). The pattern of discharge is not consistent and can be affected by changes in weather conditions, resulting in abnormalities and asymmetry. Traditional parametric models may not accurately capture these changes over time. Meanwhile, non-parametric methods can be problematic when dealing with small data samples with short lengths. To strike a balance between efficiency and accuracy, Ker & Coble (2003) proposed a semi-parametric model based on quantile regression. This model can plot regression curves that match different percentile points without the constraints of normal linear regression. Many studies have since applied this model to analyze changes in climatic and hydrological parameters, such as rainfall intensity and trends (Fan & Chen 2016; Dunn *et al*. 2019; Gao & Franzke 2017; Salarijazi 2017; Salman *et al*. 2018). For example, Dhakal & Tharu (2018) observed changes in heavy rainfall distribution and their relationship to tropical storms in the southeastern United States. Abbas *et al.* (2019) found that the quantile regression model can provide patterns for very dry and very wet conditions. Bararkhanpour *et al.* (2020) studied seasonal changes in precipitation in Gorgan, Iran, and found that the lower quantiles in autumn and the middle quantiles in spring and winter have a decreasing trend slope. During summer and winter, there is an upward trend in the upper quantiles. Mohsenipour *et al.* (2020) conducted a study on the distribution of precipitation during the monsoon months in Bangladesh (June to September). The findings revealed that rainfall decreased in many quantiles from June to August and increased in September at most stations. Zhang *et al.* (2020) found that the Arctic Sea ice's quantiles from 0.01 to 0.99 had a statistically significant decrease over time. The study also confirmed the important role that climate patterns play in reducing sea ice, especially in the upper and lower quantiles of sea ice concentrations. Treppiedi *et al.* (2021) researched the trend of changes in precipitation in the Mediterranean region and found an increasing trend in less than an hour of rainfall, particularly in the upper quantiles of precipitation. Despite the quantile regression model's advantages, the quantile curves in this model are often estimated independently. As a result, the quantile curves are evaluated independently, and estimates at extreme quantiles do not consider information from estimates in near-distributed masses (Ramsey 2020). Due to their independent estimation, quantile curves can intersect, which invalidates the distribution of independent variables (Ramsey 2020).

Quantile regression Bayesian models are a solution to the limitations of frequentist models and provide accurate estimates by utilizing information from quantile levels and locations. This method surpasses the shortcomings of the frequentist model and has been proven to be more precise and accurate (Hu *et al.* 2021; Alahmadi *et al.* 2022). Reich *et al.* (2011) proposed Bayesian quantile regression, which is an effective model for measuring variable distribution over time. This model offers a smooth function for the relationship between all possible quantiles, unlike the frequentist model that provides point estimates. Additionally, the Bayesian model borrows information from all quantile levels, leading to more precise estimates. Despite limited research on the subject, some studies have utilized the Bayesian quantile regression model to analyze factors affecting various parameters and their trend over time. Tan *et al.* (2019) modeled the distribution of changes in winter rainfall in Canada and discovered an increase in the upper and lower quantiles in eastern Canada. Tharu & Dhakal (2020) used the Bayesian regression model to analyze trends in monthly and daily maximum precipitation indices across 1,108 locations in the United States over 65 years. The study revealed that the upper quantiles of heavy rainfall distribution changed at a much higher rate than the lower quantiles. Amini *et al.* (2020) investigated the relationship between oceanic–atmospheric indices and standardized precipitation index (SPI) in Iran and found that La Niña intensified drought conditions in several regions of the Caspian Sea, while the North Atlantic Oscillation (NAO) had the weakest effect on drought in different quantiles and regions of Iran. Yong *et al.* (2021) also conducted research on the dependence of heavy precipitation on dew point temperature and surface air temperature at 78 stations on the Tibetan Plateau. The study concluded that there is a stronger correlation between heavy precipitation and dew point temperature. This study aims to examine the trend in flow values, particularly low and high flows, in the Gorganroud River basin. The investigation will be conducted in three stations, namely Arazkouseh, Tamar, and Galikesh, located in Golestan province, Iran, during different seasons. Two statistical methods, frequentist quantile regression and Bayesian quantile regression, will be used to analyze the data. The study will also compare the results obtained from both methods to determine the significance of using the Bayesian approach in quantile regression.

## MATERIALS AND METHODS

The Gorganroud basin, situated in Golestan province in Iran, is a flood-prone area that causes significant damage to its economic and social sectors every year. Despite being located in northern Iran's green and forested belt, the region is susceptible to drought. Therefore, examining the discharge in its rivers at different time scales is crucial. Arazkouseh, Tamar, and Galikesh stations, the oldest hydrometric stations with daily discharge statistics spanning half a century, were selected as a case study in this basin.

The Gorganroud basin covers an area of approximately 11,380 km^{2}, which is equivalent to 48% of Golestan province. Its main branch, running through the Gorgan and Gonbad plains, spans 359.4 km and contains 67% of the surface water in the province, amounting to 828 million cubic metres. The region's average annual rainfall ranges from 200 to 700 mm. Forests dominate the basin's southern and eastern regions, while alluvial plains are used for agriculture and pasture in the north and west.

*et al*. 2021). Figure 1 illustrates the basin's location and stations.

### Data

The study utilized a daily time series of discharge data provided by the Regional Water Authority of Golestan Province. The selection of stations was based on their sufficient statistics for over 30 years and not having a dam to collect and store water from upstream. Therefore, Arazkouseh, Tamar, and Galikesh stations were chosen as they had 50-year daily discharge data from 1966 to 2016 (1346–1395). The study calculated the average monthly discharges, formed their seasonal series, and conducted trend tests for analysis.

### Trend detection models

This study utilized two models, namely frequency-oriented quantile regression and Bayesian multiple regression, to analyze the seasonal trend of changes in various monthly discharge values.

### Quantile regression model

*β*

_{0}(

*ρ*) is the intercept and

*β*

_{1}(

*ρ*) is the slope coefficient, both of which change depending on the value of the

*ρ*th quantile being studied;

*ε*represents error with zero expectation and the range of

*ρ*values is from 0 to 1 (Lee

*et al.*2013).

*ρ*th estimation of quantile regression is performed by minimizing the following equation:where

*i*= 1,2,…,

*n*and

*y*(

_{ρ}*x*

_{i}) =

*β*

_{0}(

*ρ*) +

*β*

_{1}(

*ρ*)

*x*. In other words, the absolute value of the difference between a

_{i}*y*observation of the corresponding

_{i}*ρ*th quantile for

*y*(

_{ρ}*x*) is weighted by (1 −

_{i}*ρ*), if the observations are below the quantile line, and weighs

*ρ*if the observations are above the quantile line (Lee

*et al.*2013). The fitted line in this model presents 100%×

*θ*of the points below the line and the rest are above the line. The linear conditional quantile function used in this research is in the form of the following equation:where

*x*is the vector of the independent variables and

*β*is the vector of the parameters related to the

_{ρ}*ρ*th quantile regression. For a set of observations

*i*= 1,2,…,

*n*, (

*x*,

_{i}*y*), the estimation of the

_{i}*β*parameters is in the form of the following equation (Koenker 2005):where the function

_{ρ}*ρ*(.) is defined by the following equation (Koenker 2005):

_{ρ}Koenker (2006) developed Quantreg software in the R programming language for quantile regression problems. The rank inverse model is used to estimate the intercept and slope coefficient of the quantile regression model. The standard errors, confidence intervals, *t*-statistic, and significance for these coefficients are also estimated using the Quantreg software package.

### Bayesian quantile regression

*y*is observed in the spatio-temporal locations (

_{i}*s*,

*t*)

*. Equation (3) can be extended as follows:*

_{i}Ramsey (2020) explains that the covariate effect varies in both quantile level and location. To ensure a valid quantile process, the quantiles must be non-decreasing for all predictor values. Estimating the quantile using this model has the advantage of sharing information or data across quantile levels. However, formulating Bayesian quantile regression is more complex than classical quantile regression because the effects of parameters must be theoretically determined over an infinite number of quantile levels.

*M*is some basis functions, and

*B*(

_{m}*τ*) is a basis function of

*τ*. The coefficients of the quantile process can then be modeled as follows:where

*j*= 1,…,

*k*and

*α*(

_{jm}*s*) are the basis coefficients that determine the shape of the quantile process. Reich

*et al.*(2011) first proposed this semi-parametric approach, which studied trends in summer ozone levels using quantile regression.

This is Bernstein's basis polynomial equation, and Equation (4) is specifically a Bernstein polynomial of degree *M*. Bernstein's polynomials have proven to be commonplace for economic and statistical problems in which shape constraints must be applied (Wang & Ghosh 2012); As Bernstein's polynomial degree tends to infinity, it converges uniformly to any continuous function at intervals [0,1].

*τ*. Using Bernstein's polynomial basis, Reich

*et al.*(2011) showed that this reduces to

*γ*=

_{jm}*α*−

_{jm}*α*

_{jm}_{−1}≥ 0 for

*m*= 2,…,

*M*. Constraints are provided by a hidden infinite variable where and

Determining as a spatial Gaussian process leads to smoothing across spatial locations. Two approaches may be considered for model estimation. The complete model can be estimated with a Monte Carlo Markov Chain Sampling Program (MCMC) in small and moderate-sized samples. Reich *et al.* (2011) developed an approximate model that first fits quantile regressions separately over some quantiles and then uses the Bayesian spatial model to construct the quantile process during these initial estimates. In terms of mean square error, coverage, and power, complete and approximate models perform better than classical quantile regression. MCMC algorithms are described in the Appendix to Reich *et al.* (2011). The Bayesian quantile regression model eventually provides a tool for examining changes in a conditional distribution while accepting a high degree of flexibility in periods of covariate effects. Covariates affect higher sequences and the distribution shape; these effects change across locations and multiple levels. The model is calculated using the BSquare package of R software (Smith *et al.* 2013).

## RESULTS AND DISCUSSION

### Statistical characteristics analysis

Table 1 displays the descriptive statistics for three stations. The Arazkouseh station has a higher average monthly discharge compared with the other two stations. Additionally, the discharge data from this station have a larger deviation from the average data, as indicated by a larger standard deviation. All three stations have positively skewed data with positive kurtosis. The positive skewness is due to the presence of high and very high discharge data in the data series, which is more pronounced in the Arazkouseh station. On the other hand, positive kurtosis implies an excessive concentration of data around the mean, with the highest kurtosis observed in the Tamar station.

Station . | Mean (m^{3}/s)
. | SD . | Skew . | Kurtosis . |
---|---|---|---|---|

Arazkouseh | 5.73 | 6.79 | 21.17 | 6.03 |

Tamar | 1.5 | 1.63 | 3.44 | 18.52 |

Galikesh | 2.5 | 2.44 | 3.36 | 16.01 |

Station . | Mean (m^{3}/s)
. | SD . | Skew . | Kurtosis . |
---|---|---|---|---|

Arazkouseh | 5.73 | 6.79 | 21.17 | 6.03 |

Tamar | 1.5 | 1.63 | 3.44 | 18.52 |

Galikesh | 2.5 | 2.44 | 3.36 | 16.01 |

Station . | Type of quantile regression (QR) . | Trend coefficients in selected quantiles . | |||||||
---|---|---|---|---|---|---|---|---|---|

0.05 . | 0.1 . | 0.25 . | 0.5 . | 0.75 . | 0.9 . | 0.95 . | |||

Arazkouseh | Frequentist QR | Slope | −0.038* | −0.058* | −0.117* | −0.15* | −0.22* | −0.29 | −0.31* |

Intercept | 53.55* | 81.32* | 164.06* | 215.93* | 313.85* | 418.05 | 459.25* | ||

Bayesian QR | Slope | −1.094 | −1.09 | −1.1 | −1.1 | −3.14 | −3.14 | −3.14 | |

Intercept | 0.76 | 1.26 | 2.89 | 7.95 | 14.48 | 22.84 | 29.08 | ||

Tamar | Frequentist QR | Slope | −0.009* | −0.012 | −0.016 | −0.013 | −0.049* | −0.045 | −0.1 |

Intercept | 13.12* | 17.3 | 23.3 | 20.77 | 71.18* | 66.3 | 152.36* | ||

Bayesian QR | Slope | −0.33 | −0.33 | −0.33 | −0.496 | −1.12 | −1.13 | −1.17 | |

Intercept | 0.36 | 0.47 | 5.09 | 1.98 | 3.05 | 4.71 | 6.33 | ||

Galikesh | Frequentist QR | Slope | −0.022* | −0.024* | −0.028* | −0.046* | −0.089* | −0.19* | −0.17* |

Intercept | 32.3* | 34* | 40.98* | 66.47* | 115.56* | 275.54* | 245.5* | ||

Bayesian QR | Slope | −0.58 | −0.74 | −0.66 | −0.66 | −0.815 | −0.812 | −0.68 | |

Intercept | 1.15 | 1.48 | 2.1 | 3.4 | 5.47 | 8.4 | 11 |

Station . | Type of quantile regression (QR) . | Trend coefficients in selected quantiles . | |||||||
---|---|---|---|---|---|---|---|---|---|

0.05 . | 0.1 . | 0.25 . | 0.5 . | 0.75 . | 0.9 . | 0.95 . | |||

Arazkouseh | Frequentist QR | Slope | −0.038* | −0.058* | −0.117* | −0.15* | −0.22* | −0.29 | −0.31* |

Intercept | 53.55* | 81.32* | 164.06* | 215.93* | 313.85* | 418.05 | 459.25* | ||

Bayesian QR | Slope | −1.094 | −1.09 | −1.1 | −1.1 | −3.14 | −3.14 | −3.14 | |

Intercept | 0.76 | 1.26 | 2.89 | 7.95 | 14.48 | 22.84 | 29.08 | ||

Tamar | Frequentist QR | Slope | −0.009* | −0.012 | −0.016 | −0.013 | −0.049* | −0.045 | −0.1 |

Intercept | 13.12* | 17.3 | 23.3 | 20.77 | 71.18* | 66.3 | 152.36* | ||

Bayesian QR | Slope | −0.33 | −0.33 | −0.33 | −0.496 | −1.12 | −1.13 | −1.17 | |

Intercept | 0.36 | 0.47 | 5.09 | 1.98 | 3.05 | 4.71 | 6.33 | ||

Galikesh | Frequentist QR | Slope | −0.022* | −0.024* | −0.028* | −0.046* | −0.089* | −0.19* | −0.17* |

Intercept | 32.3* | 34* | 40.98* | 66.47* | 115.56* | 275.54* | 245.5* | ||

Bayesian QR | Slope | −0.58 | −0.74 | −0.66 | −0.66 | −0.815 | −0.812 | −0.68 | |

Intercept | 1.15 | 1.48 | 2.1 | 3.4 | 5.47 | 8.4 | 11 |

**P* value < 0.05 in frequentist quantile regression.

This research is based on seasonal data, and its results for different seasons are as follows.

### Monthly discharge trend detection using regression models during spring

For Tamar and Arazkouseh stations, the slope of changes in the extreme upper and upper quantiles was much higher than in the lower quantiles, indicating a noticeable decrease in very high discharge values. However, the trend slopes of different quantiles in Galikesh station showed values close to each other in almost all studied quantiles.

The covariate effect in different quantiles was measured with a lower degree of uncertainty for the Tamar station and the lower and median quantiles at the Galikesh station. Comparing the results of the frequentist quantile regression and the Bayesian quantile regression models based on the slope and intercept of the quantiles, it can be concluded that in the Bayesian model, the trend slope values in different quantiles are estimated more uniformly compared with the frequentist approach. Additionally, the Bayesian approach shows smoother skewness and kurtosis curves than the frequentist quantile model.

Table 2 shows the regression coefficients (slope and intercept) and statistical significance of the trends. Arazkouseh and Galikesh stations have a significant decreasing trend in most selected quantiles, while the Tamar station has only had a significant trend in a limited number of quantiles. The highest significant trends are found in the upper quantiles of the Arazkouseh and Galikesh stations. For example, the extreme upper quantile of 0.95 decreased with a slope of −0.31 (equivalent to a decrease of 0.31 m^{3}/s per year) at Arazkouseh station and with a slope of −0.17 (equivalent to a decrease of 0.17 m^{3}/s per year) at Galikesh station.

### Studying monthly discharge trends using regression models during the summer

Table 3 provides statistical significance and regression coefficients for summer trends. Most selected quantiles across various stations did not have statistically significant trends, except for extremely lower quantiles (quantile < 0.5). At Arazkouseh, the 0.05, 0.1, and 0.25 quantiles showed a significant increasing trend, indicating that discharge values equivalent to those quantiles in summer have increased by 0.0017, 0.002, and 0.0026 m^{3}/s per year. Conversely, lower quantiles (0.1 and 0.25) at Galikesh and Tamar stations showed a significant decreasing trend with negative slopes. Therefore, it can be stated that the 0.1 and 0.25 quantiles of discharge at Tamar station have decreased by 0.0008 and 0.0032 m^{3}/s per year, respectively, and the 0.25 quantile of discharge at Galikesh station has decreased by 0.0093 m^{3}/s per year. The most significant and noticeable estimated slope was found at the 0.95 quantile in the Tamar station, showing an increasing positive slope of 0.12 m^{3}/s per year.

Station . | Type of quantile regression (QR) . | Trend coefficients in selected quantiles . | |||||||
---|---|---|---|---|---|---|---|---|---|

0.05 . | 0.1 . | 0.25 . | 0.5 . | 0.75 . | 0.9 . | 0.95 . | |||

Arazkouseh | Frequentist QR | Slope | 0.0017* | 0.002* | 0.0026* | −0.0008 | −0.012 | −0.023 | −0.032 |

Intercept | −2.26* | −2.85* | −3.59* | 1.5 | 18.13 | 34.7 | 48.16 | ||

Bayesian QR | Slope | 0.72 | 0.59 | 0.0147 | 0.002 | −0.22 | −0.5 | −0.65 | |

Intercept | −2.69 | −1.64 | −0.5 | 0.1 | 0.6 | 2.044 | 3.14 | ||

Tamar | Frequentist QR | Slope | −0.0002 | −0.0008* | −0.0032* | −0.0024 | −0.0008 | 0.028 | 0.12* |

Intercept | 0.3 | 1.17* | 4.43* | 3.63 | −0.28 | −37.36 | −168.07* | ||

Bayesian QR | Slope | 0.016 | 0.016 | 0.045 | 0.0006 | 0.0006 | 0.399 | 1.64 | |

Intercept | −4 | 2.56 | −0.814 | 0.0476 | 0.35 | 1.6 | 3.15 | ||

Galikesh | Frequentist QR | Slope | −0.0044 | −0.0045 | −0.0093* | −0.008 | −0.0023 | 0.0039 | 0.029 |

Intercept | 6.67 | 6.8 | 13.7* | 12.85* | 5.08 | −2.72 | −36.9 | ||

Bayesian QR | Slope | −0.072 | −0.11 | −0.122 | −0.122 | −0.122 | 0.16 | 0.93 | |

Intercept | 0.46 | 0.63 | 0.93 | 1.35 | 1.98 | 3.2 | 4.62 |

Station . | Type of quantile regression (QR) . | Trend coefficients in selected quantiles . | |||||||
---|---|---|---|---|---|---|---|---|---|

0.05 . | 0.1 . | 0.25 . | 0.5 . | 0.75 . | 0.9 . | 0.95 . | |||

Arazkouseh | Frequentist QR | Slope | 0.0017* | 0.002* | 0.0026* | −0.0008 | −0.012 | −0.023 | −0.032 |

Intercept | −2.26* | −2.85* | −3.59* | 1.5 | 18.13 | 34.7 | 48.16 | ||

Bayesian QR | Slope | 0.72 | 0.59 | 0.0147 | 0.002 | −0.22 | −0.5 | −0.65 | |

Intercept | −2.69 | −1.64 | −0.5 | 0.1 | 0.6 | 2.044 | 3.14 | ||

Tamar | Frequentist QR | Slope | −0.0002 | −0.0008* | −0.0032* | −0.0024 | −0.0008 | 0.028 | 0.12* |

Intercept | 0.3 | 1.17* | 4.43* | 3.63 | −0.28 | −37.36 | −168.07* | ||

Bayesian QR | Slope | 0.016 | 0.016 | 0.045 | 0.0006 | 0.0006 | 0.399 | 1.64 | |

Intercept | −4 | 2.56 | −0.814 | 0.0476 | 0.35 | 1.6 | 3.15 | ||

Galikesh | Frequentist QR | Slope | −0.0044 | −0.0045 | −0.0093* | −0.008 | −0.0023 | 0.0039 | 0.029 |

Intercept | 6.67 | 6.8 | 13.7* | 12.85* | 5.08 | −2.72 | −36.9 | ||

Bayesian QR | Slope | −0.072 | −0.11 | −0.122 | −0.122 | −0.122 | 0.16 | 0.93 | |

Intercept | 0.46 | 0.63 | 0.93 | 1.35 | 1.98 | 3.2 | 4.62 |

**P* value < 0.05 in frequentist quantile regression.

### Studying monthly discharge trends using regression models during autumn

Table 4 displays the statistical significance of the trend during autumn. The findings indicate that most of the quantiles for all three stations showed no significant trend. Specifically, the low discharge values in the Araskouseh (quantile 0.05), Tamar (quantile 0.1), and Galikash (quantile 0.1) stations have significantly decreased by 0.0045, 0.0055, and 0.0063 m^{3}/s per year, respectively. On the other hand, the high discharge values in Arazkouseh (quantile 0.9) and Galikesh (quantile 0.75) stations have significantly decreased by 0.069 and 0.008 m^{3}/s per year, respectively, while it has increased by 0.038 m^{3}/s per year in Tamar station (quantile 0.95).

Station . | Type of quantile regression (QR) . | Trend coefficients in selected quantiles . | |||||||
---|---|---|---|---|---|---|---|---|---|

0.05 . | 0.1 . | 0.25 . | 0.5 . | 0.75 . | 0.9 . | 0.95 . | |||

Arazkouseh | Frequentist QR | Slope | −0.0045* | 0.0002 | −0.0055 | −0.014 | −0.03 | −0.069* | −0.056 |

Intercept | 6.98* | 0.74 | 9.27 | 21.49 | 45.1 | 100* | 85 | ||

Bayesian QR | Slope | −0.17 | −0.22 | −0.12 | −0.5 | −0.94 | −0.94 | −0.94 | |

Intercept | 0.9 | 1.17 | 1.79 | 2.94 | 4.44 | 6.24 | 7.92 | ||

Tamar | Frequentist QR | Slope | −0.004 | −0.0055* | −0.003 | −0.005 | −0.0001 | 0.02 | 0.038* |

Intercept | 6.16 | 8.24* | 4.74 | 7.97 | 1.53 | −26.04 | −50.87* | ||

Bayesian QR | Slope | −0.18 | −0.16 | −0.096 | −0.163 | −0.14 | −0.088 | −0.056 | |

Intercept | 0.5 | 0.59 | 0.73 | 1 | 1.41 | 2.095 | 2.8 | ||

Galikesh | Frequentist QR | Slope | −0.004 | −0.0063* | −0.004 | −0.0027 | −0.008* | −0.0017 | 0.01 |

Intercept | 6.24 | 9.46* | 6.57 | 5.08 | 12.9* | 4.46 | −11.8 | ||

Bayesian QR | Slope | −0.18 | −0.139 | −0.124 | −0.173 | −0.17 | −0.125 | −0.125 | |

Intercept | 0.62 | 0.76 | 0.99 | 1.35 | 1.75 | 2.1 | 2.36 |

Station . | Type of quantile regression (QR) . | Trend coefficients in selected quantiles . | |||||||
---|---|---|---|---|---|---|---|---|---|

0.05 . | 0.1 . | 0.25 . | 0.5 . | 0.75 . | 0.9 . | 0.95 . | |||

Arazkouseh | Frequentist QR | Slope | −0.0045* | 0.0002 | −0.0055 | −0.014 | −0.03 | −0.069* | −0.056 |

Intercept | 6.98* | 0.74 | 9.27 | 21.49 | 45.1 | 100* | 85 | ||

Bayesian QR | Slope | −0.17 | −0.22 | −0.12 | −0.5 | −0.94 | −0.94 | −0.94 | |

Intercept | 0.9 | 1.17 | 1.79 | 2.94 | 4.44 | 6.24 | 7.92 | ||

Tamar | Frequentist QR | Slope | −0.004 | −0.0055* | −0.003 | −0.005 | −0.0001 | 0.02 | 0.038* |

Intercept | 6.16 | 8.24* | 4.74 | 7.97 | 1.53 | −26.04 | −50.87* | ||

Bayesian QR | Slope | −0.18 | −0.16 | −0.096 | −0.163 | −0.14 | −0.088 | −0.056 | |

Intercept | 0.5 | 0.59 | 0.73 | 1 | 1.41 | 2.095 | 2.8 | ||

Galikesh | Frequentist QR | Slope | −0.004 | −0.0063* | −0.004 | −0.0027 | −0.008* | −0.0017 | 0.01 |

Intercept | 6.24 | 9.46* | 6.57 | 5.08 | 12.9* | 4.46 | −11.8 | ||

Bayesian QR | Slope | −0.18 | −0.139 | −0.124 | −0.173 | −0.17 | −0.125 | −0.125 | |

Intercept | 0.62 | 0.76 | 0.99 | 1.35 | 1.75 | 2.1 | 2.36 |

**P* value < 0.05 in frequentist quantile regression.

### Studying monthly discharge trends using regression models during winter

Based on the data presented in Table 5, it appears that there is no significant trend in any of the selected quantiles at the Galikesh station. However, at the Tamar station, there was a notable decrease in the lower values of discharge, specifically in quantiles 0.05 and 0.1, with trend slopes of 0.0048 and 0.0064 m^{3}/s per year, respectively. Similarly, at the Arazkouseh station, the upper 0.75 quantile exhibited a significant decreasing trend with a slope of −0.09 m^{3}/s per year. It is worth noting that most of the slopes observed were relatively minor and insignificant, particularly during the autumn season, suggesting limited changes in discharge during this period.

Station . | Type of quantile regression (QR) . | Trend coefficients in selected quantiles . | |||||||
---|---|---|---|---|---|---|---|---|---|

0.05 . | 0.1 . | 0.25 . | 0.5 . | 0.75 . | 0.9 . | 0.95 . | |||

Arazkouseh | Frequentist QR | Slope | 0.001 | −0.02 | −0.024 | −0.042 | −0.09* | −0.07 | −0.14 |

Intercept | 0.57 | 30.49 | 37.13 | 65.8 | 135* | 112.8 | 218 | ||

Bayesian QR | Slope | −0.53 | −0.53 | −0.53 | −0.18 | −0.8 | −0.8 | −0.8 | |

Intercept | 2 | 2.7 | 4 | 6.84 | 11.48 | 16.54 | 20.76 | ||

Tamar | Frequentist QR | Slope | −0.0048* | −0.0064* | −0.005 | −0.0009 | 0.004 | 0.019 | 0.02 |

Intercept | 7.36* | 9.64* | 8.29* | 2.59 | −3.8 | −24.21 | −30.9 | ||

Bayesian QR | Slope | −0.15 | −0.143 | −0.149* | −0.14 | −0.14 | −0.14 | −0.14 | |

Intercept | 0.69 | 0.76 | 0.92 | 1.3 | 1.9 | 2.67 | 3.36 | ||

Galikesh | Frequentist QR | Slope | −0.0012 | 0.0036 | 0.0002 | 0.0005 | −0.03 | −0.042 | −0.008 |

Intercept | 2.6 | −3.82 | 1.076 | 1.17 | 46.95 | 62.95 | 18.54 | ||

Bayesian QR | Slope | 0.024 | −0.0055 | −0.015 | −0.1 | 3.55 | −0.675 | −0.79 | |

Intercept | 0.89 | 1.05 | 1.42 | 2.02 | −0.78 | 5.3 | 8.35 |

Station . | Type of quantile regression (QR) . | Trend coefficients in selected quantiles . | |||||||
---|---|---|---|---|---|---|---|---|---|

0.05 . | 0.1 . | 0.25 . | 0.5 . | 0.75 . | 0.9 . | 0.95 . | |||

Arazkouseh | Frequentist QR | Slope | 0.001 | −0.02 | −0.024 | −0.042 | −0.09* | −0.07 | −0.14 |

Intercept | 0.57 | 30.49 | 37.13 | 65.8 | 135* | 112.8 | 218 | ||

Bayesian QR | Slope | −0.53 | −0.53 | −0.53 | −0.18 | −0.8 | −0.8 | −0.8 | |

Intercept | 2 | 2.7 | 4 | 6.84 | 11.48 | 16.54 | 20.76 | ||

Tamar | Frequentist QR | Slope | −0.0048* | −0.0064* | −0.005 | −0.0009 | 0.004 | 0.019 | 0.02 |

Intercept | 7.36* | 9.64* | 8.29* | 2.59 | −3.8 | −24.21 | −30.9 | ||

Bayesian QR | Slope | −0.15 | −0.143 | −0.149* | −0.14 | −0.14 | −0.14 | −0.14 | |

Intercept | 0.69 | 0.76 | 0.92 | 1.3 | 1.9 | 2.67 | 3.36 | ||

Galikesh | Frequentist QR | Slope | −0.0012 | 0.0036 | 0.0002 | 0.0005 | −0.03 | −0.042 | −0.008 |

Intercept | 2.6 | −3.82 | 1.076 | 1.17 | 46.95 | 62.95 | 18.54 | ||

Bayesian QR | Slope | 0.024 | −0.0055 | −0.015 | −0.1 | 3.55 | −0.675 | −0.79 | |

Intercept | 0.89 | 1.05 | 1.42 | 2.02 | −0.78 | 5.3 | 8.35 |

**P* value < 0.05 in frequentist quantile regression.

## DISCUSSION

The research findings suggest that quantile regression models offer a more comprehensive understanding of hydrological and meteorological variables than non-parametric tests. In particular, the Bayesian quantile regression model has been found to detect trends that were previously unnoticed by other trend detection methods (Tharu & Dhakal 2020).

A comparison of the results obtained from quantile regression models and Bayesian quantile regression reveals that all studied stations exhibited a decreasing trend across all quantiles during the spring season. This trend was more pronounced in the upper and extreme upper quantiles. Conversely, during the summer season, lower quantiles showed almost negligible slopes for all stations, while upper quantiles, particularly the extreme upper ones, exhibited a significant negative slope at the Arazkouseh station and a positive slope at Tamar and Galikesh stations. In the autumn season, quantile regression showed a negative trend in the upper quantiles for the Arazkouseh station and a significant positive trend for the Tamar station. However, the Galikesh station only showed a noteworthy positive trend in extreme upper quantiles. Bayesian quantile regression showed similar trends in all quantiles of Tamar and Galikesh stations. Finally, during the winter season, the upper quantiles in Arazkouseh and Galikesh stations showed a negative trend, while the upper quantiles of the Tamar station showed a significant positive trend. However, Bayesian quantile regression indicated that the trends in all three stations were approximately equal and negative.

Recent studies have highlighted that the construction of the Golestan Dam has resulted in a decrease in average monthly discharge and an increase in maximum flow indicators in Gonbad and Qazaqali stations in the Gorganroud basin (Daiechini *et al.* 2020). Further, the maximum flow indicators in the Aq Qhala station have increased after the construction of the Voshmgir Dam. H. Ahmadi *et al.* (2020) and (M. Ahmadi *et al.* 2020) suggest that rainfall-runoff in the Gorganroud region, under the SSP8.5 scenario, displays an increasing trend for the future, with maximum discharge occurring in April and May. They further found that an increase in extreme rainfall could lead to more frequent flood events in the region. Additionally, Norouzi Nazar *et al.* (2023) have projected that base flow in the Gorganroud watershed would increase by 166.7% and 77.2% during the periods of 2021–2040 and 2040–2060, respectively, due to climate and land use changes. They also note that runoff is expected to increase by 54.2% and 41% during the same periods, respectively.

The analysis of flow data has revealed that changes in discharge distribution have affected several quantiles above river flow, particularly during monsoon months, in stations such as Tamar and Galikesh, leading to more flooding. Therefore, quantile regression is preferred as it can reveal the temporal dependencies of the variable in question for the mean value and its quantiles (Abbas *et al.* 2019). Saadi *et al.* (2017) used the quantile regression model to examine the distributional changes of annual monsoon rainfall and river flow for 31 stations in Sarawak and found that upward trends were mainly in the lower quantiles and in January and December in the upper quantiles. Similarly, Yves *et al.* (2019) found that severe floods had an increasing trend only in a limited number of basins, while most areas exhibited a decreasing trend in their study of flood trends and their drivers in 171 basins in the south of France. According to Pasch & Engelke (2022), the flood risk assessment requires an examination of changes in high quantiles. To predict floods and very high discharges more effectively, authorities can use the quantile regression method. In some cases, there has been a decreasing trend in high discharge values for certain stations, such as Arazkouseh station in all seasons, which may lead to problems with available water supply. Additionally, changes in flow patterns can pose significant threats to animal and plant species and have adverse effects on the environment. Nourali (2023) suggested that accurate flood forecasting and flood control planning must be considered in response to sudden floods in the Tamar basin area located in the Gorganrood basin in Iran. The results of Donyaii & Sarraf (2021) showed that discharge and temperature values in the Gorganroud basin will experience decreasing and increasing trends, respectively, at the end of the 21st century. Meanwhile, precipitation will not change much. Yves *et al*. (2019) indicated that the decrease in flood events could be due to the decrease in soil moisture caused by an increase in temperature, evaporation, and transpiration, and a decrease in rainfall in the region. Daiechini *et al.* (2020) stated that the construction of the dam and changes in climatic variables mostly influenced the change in hydrological indicators of the flow. Khoie *et al*. (2023) demonstrated that more than 60% of the changes in streamflow in the Gorganroud basin were due to land use changes.

When comparing how frequentist quantile and Bayesian quantile curves are fitted, we find that the estimates of the slope of the quantile curve in the frequentist and Bayesian models are different in some cases. These differences mostly appear when the slope of the quantile curve in the desired quantile is estimated with frequentist quantile regression, causing the crossing of adjacent quantile lines or fitting different slopes with adjacent polynomials. Such estimates have lower accuracy and validity. Tan & Shao (2017) showed that the Bayesian quantile regression model produces smoother confidence intervals of trend coefficients at extremely high and low quantile levels. It also estimates a parameter of uncertainty more accurately than the frequentist quantile regression approach (Kozumi & Kobayashi 2011). Therefore, the Bayesian regression model modifies the invalid measurements of frequentist quantile regression and provides valid estimates for all possible quantiles.

## CONCLUSION

This study investigated changes in discharge values over 50 years for the Arazkouseh, Tamar, and Galikesh stations of the Gorganroud basin on a seasonal scale. Two new regression models, frequentist quantile regression and Bayesian quantile regression, were used to determine trends. The Bayesian quantile regression model provided more accurate results as it considered information in all discharge quantiles. These models are important for flood and drought management as they show the status of changes in different discharges.

The study found that high discharge values had different decreasing and increasing trends, with most cases showing negative and decreasing trends, especially in the higher extreme values of discharge. Tamar and Galikesh stations had an increasing trend in summer, which could result in summer floods. Although the discharge is decreasing in most quantiles, the flood months (upper quantiles) are increasing. This finding is consistent with recent reports of floods in the region. The decreasing trend in the Gorganroud basin, which will become more severe in the future due to global climate changes, along with the increase in the water demand of the region for various purposes due to population growth and expansion of industries, is a serious warning for politicians, planners, and local managers to prevent the possible occurrence of water crisis in the region in the future with appropriate planning. Future research should quantify the climatic changes in land use in the decreasing and increasing trends of discharge in the Gorganroud basin.

## ACKNOWLEDGEMENTS

This research project has been supported by the Gorgan University of Agricultural Sciences and Natural Resources (grant no. 98-412-8). The authors are grateful to Gorgan University of Agricultural Sciences and Natural Resources for providing financial support for conducting this research.

## DATA AVAILABILITY STATEMENT

Data cannot be made publicly available; readers should contact the corresponding author for details.

## CONFLICT OF INTEREST

The authors declare there is no conflict.

## REFERENCES

*Proceedings of 3rd*

*(ISBS 2017), Volume 2*(Fırat, S., Kinuthia, J. & Abu-Tair, A., eds),

*Journal of Water and Soil Resources Conservation*

**5**(4), 18–34.

*Quantile Regression in R: A Vignette*. Available from: http://www.econ.uiuc.edu/∼roger/research/rq/vig.pdf.

**46**(1),

**13**(3),

**37**(11), 4211–4231.

*Asia-Pacific Journal of Atmospheric Sciences*

**53**(4), 489–500.