Methodological Approaches to Calculating Maximum River Discharge in the Upper Irtysh Basin

INTRODUCTION

This study investigates flood hazards in the East Kazakhstan Region (EKR), the easternmost administrative division of Kazakhstan covering 97,800 km² (Fig. 1). The region shares borders with Russia’s Altai Territory and Altai Republic to the north, China to the east and southeast, and Kazakhstan’s Abay Region to the west. Ust-Kamenogorsk serves as the regional administrative centre.

The East Kazakhstan region’s river network, as documented in the Republic of Kazakhstan’s water and energy cadastre, comprises over 800 rivers exceeding 10 km in length, including 48 rivers longer than 50 km and 20 rivers surpassing 100 km. All waterways in the region constitute tributaries of varying orders within the Irtysh River system. Based on hydrological regime characteristics (Fig. 1), these tributaries can be classified into three distinct groups:

— Group 1: right-bank tributaries of the South-Western Altai (e.g., Bukhtarma, Uba, Ulba, Kurshim, Qalzhyr, Naryn): Characterized by perennial flow and high discharge capacity;

— Group 2: left-bank tributaries of the Kalbinsky Range (e.g., Ulken-Boken, Qaiyndy, Ablaketka, Ulanka, Dresvyanka, Kyzylsu): Exhibiting reduced but generally sustained flow;

— Group 3: southern Zaisan Lake basin rivers (e.g., Qandysu, Uydene, Kenderlyk): Typically, ephemeral systems that frequently terminate in alluvial sand deposits or experience complete desiccation.

The territory of EKR is characterized by several types of floods of different origins, including spring freshets, rain-induced floods, ice-jam and debris-jam floods, and wind-driven surges. The main channel of the Irtysh River is controlled by Bukhtarma Reservoir and rarely inundates coastal areas. Flooding during freshets and flash floods is typical for all tributaries of the Irtysh, particularly in sections with more uniform riverbeds, whereas ice and debris jams occur in mountainous areas where the channel narrows. Inundation of coastal areas due to wind-driven surges is primarily observed in Lake Zaisan and the Bukhtarma Reservoir.

An analysis of materials from the Ministry of Emergency Situations of the Republic of Kazakhstan (MES RK) has identified a number of key areas regularly at risk of flooding. The highest concentration of such areas is observed along rivers of the first group, particularly at the confluences of tributaries of the Irtysh River, including the Ulba River near Ridder, as well as the Krasnoyarka, Glubochanka, and other rivers. Several flood-prone zones are also located along the Bukhtarma River and its tributaries. Certain settlements in the lower reaches of the Kalzhir, Kurchum, and Naryn rivers are also susceptible to flooding. The most flood-prone areas of the first group of rivers include Chapaevo village on the Krestovka River (Altai District), Ust-Talovka settlement on the Talovka River, Ubinka village on the Oba River (Shemonaikha District), Karatogai village on the Kalgutty River (Kurchum District).

Among the most significant floods was the spring freshet that occurred in March 2018 near Ust-Kamenogorsk on the Ulba River. According to media reports, more than ten villages were affected, with over 700 residents losing their homes due to the inundation of more than 480 houses (Vostochny Potop 2018). The direct damage to the region’s infrastructure was estimated at 3.2 billion tenge (Year 2018 Is Not Over Yet 2018). The disaster was triggered by an abrupt and unusually rapid temperature rise, combined with heavy precipitation and frozen ground, which prevented water absorption and intensified surface runoff.

The coastal areas of second-group rivers, due to their low water flows, are significantly less prone to flooding. However, isolated flooding events have occasionally been recorded along the Ulken-Naryn, Kayyndy, Lailinka, Ulanka, and Tainty rivers. Among the most problematic areas are Samarskoye village on the Lailinka River, Mirolyubovka village on the Kayyndy River (Samar District), Ulanskoye village and Zhanuzak village on the Ulanka River, Asubulak village on the Ungyrdy River, Besterek village on the Kolbala Stream (Ulan District).

The tributaries of the third group are generally characterized by rapid water level rises during snowmelt floods, which often leads to inundation, particularly in their lowland sections. Certain areas along these rivers experience recurrent flooding, including The Kandysu River and rivers near Zaisan City, Kensai, Zharsu, Bakasu, and Saryzhira villages on the Uidene River, Tugyl settlement on the Kabyrgatal River (Zaisan District).

In Soviet and subsequently Russian engineering practice, the development of technological and informational tools for flood risk management has been based on approaches established in the mid-20th century (Rozhdestvensky, Chebotarev 1974). These approaches primarily focus on determining design hydrological characteristics (water discharge or levels) with low probability (i.e., rare recurrence intervals) using available observational data series of river flow.

The objective of this study was to conduct a statistical analysis of flood hazards based on actual maximum flow data from Irtysh River tributaries. The analysis employed various methods of stochastic hydrology recommended by both current regulatory documents in force in Russia and Kazakhstan and authoritative literary sources.

The article evaluates the adequacy and effectiveness of various engineering calculation methods for determining maximum water discharges with specified exceedance probabilities, as applied to rivers in the study region. Each computational approach is based on an appropriate probabilistic model incorporating: (1) a theoretical probability distribution for the studied variable, and (2) a parameter estimation method for limited observational datasets. The research employs a comprehensive range of stochastic hydrology methods for developing flow frequency curves, including:

— The standard calculation framework recommended by current regulatory documents in Russia (() SP 529.1325800.2023 (2023)) and Kazakhstan (MSP 3.04-101-2005 (2006));

— An alternative computational scheme developed by Yu.B. Vinogradov (1988);

— Several probability distributions validated through international flood frequency analysis practice (Gubareva 2010, 2011);

— The state-of-the-art L-moments method for distribution parameter estimation, predominantly used in international studies (Hosking & Wallis 1997; Gubareva & Gartsman 2010).

The majority of EKR lies within the mountain system of the Southwestern Altai, except for its southern portion which partially encompasses the Saur-Tarbagatai mountain range. Progressing westward and southwestward, the mountains gradually transition into the more subdued topography of the Kazakh Uplands (Kazakhskiy Melkosopochnik). A prominent river valley, formed by the Irtysh River and its numerous tributaries, cuts through the mountainous terrain from the southwest to the northeast. This valley includes the intermontane basin of Lake Zaisan. Overall, the predominantly mountainous landscape of the region exhibits a wide range of elevations, from 200 to 4,500 meters above sea level, with a general slope trending northwestward and westward (Belyanin et al. 2013; Egorina, Zinchenko, Zinchenko 2000; Egorina et al. 2015).

The climate of the region, situated in the central part of the Eurasian continent, is classified as harsh continental and further complicated by mountainous terrain, following the principles of altitudinal zonation. These characteristics significantly influence the distribution of most meteorological parameters. Winter in the region is cold and prolonged, with mean January temperatures ranging from -12°C to -17°C in lowland areas to -23°C to -27°C in high-altitude zones. Absolute minima in some years can drop to -51°C to -54°C. Summer is hot, with July averages between 15°C and 24°C, while absolute maxima reach 35°C to 45°C. The number of days with temperatures above 0°C varies from fewer than 200 in mountainous areas to 230 in the southern lowlands of EKR. Precipitation is highly unevenly distributed, ranging from 400–650 mm in mountainous regions to less than 200 mm in the Zaisan Depression. Mountainous zones typically experience sufficient or excessive moisture, whereas lowland areas face moisture deficits. Mean annual wind speeds across the oblast generally range from 2–5 m/s, though in some areas, they can exceed 15 m/s (Egorina & Popova 1989; Egorina et al. 2015).

The mountainous landscapes of EKR exhibit distinct altitudinal zonation with four characteristic elevation belts. The lower belt, extending to 500-600 meters above sea level, encompasses plains and foothills. In the northwestern foothills, chernozem soils support feather grass-forb steppe communities, while the left bank of the Irtysh River valley features feather grass-fescue vegetation on dark kastanozems (dark chestnut soils). The Zaisan Depression displays unique arid-environment vegetation including wormwood-fescue communities on light kastanozems (light chestnut soils) and wormwood-anabasis associations on brown soils, with widespread occurrence of solonchaks, solonetz soils, and dune sands. At middle elevations (up to 1900-2000 meters), the forest belt dominates with mixed woody vegetation growing on brown forest soils. Higher elevations (up to 2800-3000 meters) are occupied by the subalpine-alpine belt characterized by meadow communities developing on mountain meadow soils. The uppermost elevations form the nival belt, where mountain peaks contain permanent snowfields, glaciers, and exposed bedrock surfaces. This vertical zonation reflects the transition from steppe ecosystems through forested middle elevations to alpine and ultimately glacial environments (Belyanin et al. 2013; Egorina et al., 2015).

The primary watercourse in EKR is the Irtysh River segment flowing from the border with China to the administrative boundary between EKR and Abai Region. This approximately 800 km long stretch incorporates several major hydrological features: the Kara Ertis (Black Irtysh) river section, the through-flow Zaisan Lake, and three major reservoirs (Bukhtarma, Ust-Kamenogorsk, and Shulba). Within Kazakhstan, the Irtysh River basin covers approximately 545,000 km². The river’s flow regime in EKR exhibits mixed feeding sources, predominantly snowmelt and glacial meltwater. Discharge patterns are influenced by tributary inflows and water withdrawals from China, where Kara Ertis waters are extensively used for agricultural irrigation. Three major hydroelectric dams along the main channel create a fully regulated flow regime that virtually eliminates flood risks in most reaches. Exceptions occur in braided river sections and during wind-driven surge events in lacustrine portions of the channel. The Irtysh maintains year-round flow without stable ice cover due to its broad channel and sustained velocities (1.2-1.4 m/s) even in downstream reaches. The spring freshet period (March-June) is significantly attenuated by dam operations. Long-term monitoring (1961-2022) at the Ust-Kamenogorsk gauging station records mean annual discharge of 559 m³/s at corresponding water level of 287.25 m (Belyanin et al. 2013; Egorina et al. 2015; Pavlenko et al. 2024).

MATERIALS AND METHODS

The study utilized maximum annual discharge records from gauging stations spanning their entire observation periods through 2022. The primary data sources included annual publications and reference materials from the State Water Cadastre maintained by RSE “Kazhydromet” (State Water Cadastre of RSE Kazhydromet 2025). The analysis focused on gauging stations with the most extensive observation histories, though nearly all available series contained significant gaps in their records. Selected stations provided continuous maximum annual discharge data series of at least 30 years, with the shortest complete series covering 33 years and the longest extending to 82 years. Approximately 85% of the analyzed data came from rivers of the first group. The complete list of monitoring stations and their respective observation periods is presented in Table 1.

The research methodology included, first, testing the observation time-series for stationarity and, second, constructing a probabilistic model—a frequency curve of maximum annual water discharges—to obtain design values with low exceedance probabilities required for engineering design and flood prevention planning. The current methodological framework for probabilistic modelling in engineering hydrology is based on the assumption of stationarity in observation series and can only be correctly applied to such series. At present, clearly evident climate changes, as well as significant anthropogenic transformations of river catchments, have led to non-stationarity in hydrological series, including maximum flow data. This necessitates the development of a methodological approach to identify various forms of non-stationarity (threshold changes in mean and variance parameters, linear trends, cyclical variations, etc.) using statistical criteria.

In this study, the maximum discharge series were tested for the presence of a linear trend — that is, a monotonic change in mean values over time — using two statistical criteria described below. For series where no significant linear trend was detected, a probabilistic model was developed by selecting the optimal combination of an analytical distribution curve and a parameter estimation method that provided the closest approximation to the empirical distribution curve constructed from the sample of measured values.

For time series demonstrating significant trends, there exists no generally accepted methodological framework for constructing probabilistic models. Therefore, estimates of design discharges with low exceedance probabilities can only be obtained through special studies extending beyond the analysis of sample data. In this research, we limited ourselves to approximate estimates of design discharges by eliminating the identified trends — that is, by adjusting the sample to stationary conditions corresponding to the maximum mean value.

For series with downward trends, adjustment was made to the first year of observations, while for series with upward trends — to the last year. Comparing the calculation results for the adjusted series with the results of conventional calculations performed on the same series without accounting for the presence of trends allowed for a quantitative assessment of their influence on the design values.

For time series with downward trends, adjusting to the beginning of the series allows for calculations with a safety margin, which is generally beneficial when estimating maximum discharges. However, additional analysis is required to prevent excessive safety margins that could lead to unjustified costs for flood protection measures. In the case of upward trends, adjusting the sample to the end of the observation period still does not provide design values of the required reliability, since the trend will continue into the future to unknown extents. This situation demands particular caution in engineering decision-making. Thus, the results presented in this article for non-stationary series represent preliminary estimates that require special case-by-case verification studies in each specific instance.

Analysis of Linear Trends. To assess the linear trend for each original data series, an equation of the following form is constructed:

(1)

where i is the year number in the multi-year series, counted from the beginning of the series; is the moving mathematical expectation of discharge for the i-th year; q is the linear trend coefficient, calculated based on the correlation coefficient between the values of the characteristic and their chronological sequence numbers; p is the regression intercept (constant term).

Significance of the Linear Trend is assessed using a specialized modification of Student’s t-test (Handbook of Hydrology…, 1993):

(2)

where: T_c is the test statistic, T_{1-a/2, n-2} is the quantile of Student’s t-distribution with (n-2) degrees of freedom of probability (1-α/2), α - is the significance level of the estimate, r is the correlation coefficient characterizing the trend, n is the length of the observation series (sample size). To assess the significance of the trend, the root mean square error of the correlation coefficient may also be used (Rozhdestvensky, Chebotarev 1974), as defined by the formula:

(3)

The trend is considered statistically significant if the correlation coefficient is at least twice as large as the root mean square error. The adjustment of individual series values to achieve stationary conditions is performed using the formula:

(4)

where: and are the original and adjusted (stationary) values of the characteristic in the i-th year of the multi-year series, n* is the duration of the observation period from the first to the last year. It should be emphasized that the duration of the observation period from the first to the last year exceeds the length of the observation series if there are missing years in the record. In this case, the year index i corresponds strictly to the sequential numbering of years in the continuous period from the first observation year to the last.

Probability Distribution Functions and Parameter Estimation Methods. The guidelines of the current Russian regulatory document SP 529.1325800.2023 (2023), which align with the applicable Kazakh standards MSP 3.04-101-2005 (2006), prescribe the use of the following distributions in hydrological calculations: primarily curves derived from the gamma distribution, including the Pearson Type III (binomial) distribution and the three-parameter Kritsky-Menkel gamma distribution. The method of moments is recommended as the primary approach for estimating the parameters of analytical curves based on sample data. For the Pearson Type III distribution, an additional graphical-analytical method is suggested, while the Kritsky-Menkel distribution calls for the approximate maximum likelihood method. Notably, the regulations do not prohibit the use of alternative calculation techniques, provided they are properly justified. Hereafter, we will refer to the normative-recommended computational framework as the baseline approach.

An alternative framework by Yu.B. Vinogradov (Vinogradov 1988) employs a family of functionally normal curves and nonparametric methods for estimating their parameters from samples. In the present study, the three-parameter lognormal distribution and the C3 distribution (described via transformation) were utilized:

(5)

here z is the normally distributed random variable, x – the initial variable, a – the transformation parameter.

The study employed direct numerical fitting of analytical curves to sample points (calibration) using various convergence measures. For sample points, unbiased estimates of the empirical distribution function coordinates were adopted. The parameters of the analytical distribution function were computed based on minimizing the convergence measure.

The total relative divergence between the empirical and analytical curves in terms of probability was used as the convergence measure:

(6)

where p* = m/n is the empirical probability of the order statistic members; p** is the analytical probability; p_a* and p_b* are the confidence interval bounds for probability p given m and n at a specified significance level; m and n represent the rank of the i-th value and the total sample size, respectively (Vinogradov 1988, p. 251). The minimization of the ω metric ensures probability-based convergence between the analytical and empirical curves and characterizes the reliability of the adopted solution.

An additional convergence measure based on absolute magnitude was employed, defined as the root mean square deviation between the ordinates of the empirical and analytical curves. This metric characterizes the precision of the adopted solution.

(7)

where: k_p* and k_p** are equally probable quantiles of the empirical and analytical distributions, expressed as modular coefficients. Thus, Vinogradov’s alternative approach enables the construction of four distribution curve variants:

— Lognormal with ω-based approximation

— Lognormal with s-based approximation

— C3 with ω-based approximation

— C3 with s-based approximation

The metrics ω (reliability) and s (precision) subsequently serve as criteria for comparative evaluation of all computed probability curves. A key adequacy criterion for the analytical probability curve is its containment within the confidence interval bounds (p_a, p_b) of the empirical curve.

Additionally, among the widely used and internationally recommended probability distributions for peak flow calculations, the following were selected:

— Generalized Extreme Value (GEV) distribution

— Generalized Logistic (GLO) distribution

— Generalized Pareto (GPO) distribution

— Two-parameter Gumbel (GM2) distribution

The selection of these distributions is supported by established research (Pisarenko et al., 2002; Naydenov, 2004; Gubareva, 2010), which demonstrates that heavy-tailed distributions are most appropriate for probabilistic modeling of extreme hydrological characteristics, particularly in the upper tail section. The application of these distributions also leads to alternative computational approaches (Gartsman, Gubareva, Kichigina, 2020).

The L-moments method, proposed relatively recently as an alternative approach for characterizing probability distribution shapes (Hosking & Wallis, 1997), offers several key advantages over conventional methods and has gained widespread international adoption. Theoretically, L-moments represent a modification of the probability-weighted moments (PWMs) originally introduced by Greenwood et al. (1979).

The procedure begins by calculating unbiased sample estimates of probability-weighted moments (PWMs) from an ascendingly ordered sample x_1:n≤x_2:n≤...≤x_n:n of size n that can be described as follows:

(8)

The generalized expression for b_r is

(9)

The sample L-moments are derived as

(10)

Or in general

(11)

where coefficients are derived as

(12)

The sample L-moments of r-th order are derived as

(13)

where t – sample L-coefficient of variation, t₃ – sample L-coefficient of skewness, t₄ – sample L-kurtosis. In the study by (Gubareva & Gartsman 2010), algorithms are provided for the mutual computation of L-moments and parameters for several three-parameter distribution laws. These distributions are characterized by location (shift), shape, and scale parameters.

Table 2 presents the analytical probability distributions and parameter estimation methods employed in this study to develop probability curves for peak water discharge values. Each computational approach combines a specific probability distribution with a corresponding parameter fitting technique, forming a unique variant for analysis.

The performance of these variants is evaluated through multiple criteria. The metrics ω (reliability) and s (precision) are computed for each variant to enable quantitative comparison. Additionally, the analysis examines how closely each fitted probability curve remains within the 90% confidence interval boundaries (p_a, p_b) of the empirical probability curve across all probability points. This provides a measure of statistical consistency between the analytical and observed data. Complementing these quantitative assessments, the study incorporates expert evaluation of the hydrological plausibility of the peak discharge estimates obtained from each variant. This qualitative judgment considers whether the results align with physical expectations and regional hydrological characteristics.

Results of Extreme Flood Probability Assessment

Table 3 presents the evaluation of linear trend significance based on the Student’s t-test and correlation coefficient error across all analyzed observation series (refer to Table 1). The analysis reveals statistically significant trends in nearly one-third of the 19 examined time series. Specifically, five series exhibit downward trends indicating decreasing flood magnitudes, while one series demonstrates an upward trend suggesting increasing flood magnitudes. Importantly, all detected trends show consistent significance when assessed through both applied statistical criteria. The identified trends may stem from diverse underlying factors, including climatic influences such as shifting precipitation regimes or anthropogenic impacts like land-use modifications and water management practices. Given this complexity, a targeted follow-up study is recommended to elucidate the precise drivers behind these observed hydrological changes.

Figure 2 presents the development of probabilistic flood frequency models for the Kalzhyr River at Kalzhyr village (catchment area 434 km², 65-year observation period), which exhibits a statistically significant downward trend. The modeling approach involved an exhaustive evaluation of all possible combinations of theoretical distribution laws and parameter estimation methods described in Table 2, with final model selection based on both quantitative goodness-of-fit criteria and qualitative expert judgment. The model selection process was conducted in two distinct phases. The initial phase considered only those methods explicitly recommended by current regulatory guidelines, while the subsequent phase expanded the evaluation to include alternative computational approaches not covered by standard protocols. Given the presence of a statistically significant trend in the series, the analysis was performed separately for both the original observed data series (Fig. 2a) and a detrended series normalized to initial conditions (Fig. 2b). Model verification employed specialized normal probability paper, which transforms the cumulative normal distribution function into a linear plot, with discharge values plotted on a logarithmic scale. This visualization technique allows for immediate assessment of model performance by examining how closely the theoretical distribution curves align with the empirical data points and remain within the 5-95% confidence intervals across the entire probability range.

Appendix A presents a comprehensive comparison of optimal flood frequency models for annual maximum discharge series at 19 gauging with records exceeding 30 years (see Table 1). The table systematically organizes selection results for both standard methodology and alternative approaches, enabling direct comparison of their performance characteristics. For each station, the analysis provides parallel sets of parameters for the optimal standard methodology solution and the best-performing alternative approach. The tabulated parameters include key metrics for model evaluation: the observed maximum discharge (Q_max) as a reference value, the selected probability distribution type, the parameter estimation method, quantitative reliability (ω) and precision (s) estimates, and the calculated 1% exceedance discharge (Q_1%). These metrics collectively allow for assessment of both statistical adequacy and engineering safety implications, with particular attention to differences in extreme quantile estimation. For the six stations exhibiting statistically significant trends (previously identified in Table 3), Appendix A presents comparative results for both the original observed series and the detrended, normalized series.

Thus, Appendix A presents data from 25 individual calculations, each performed using both standard and alternative methodological approaches. The analysis reveals several key findings regarding the performance of different probability distributions and parameter estimation techniques. When applying the standard methodology, the Pearson Type III distribution demonstrated superior performance in 17 out of 25 cases, while the Kritsky-Menkel distribution proved optimal in the remaining 8 cases. The standard approach did not recommend any other probability distributions for these datasets. Regarding parameter estimation techniques, the method of moments and the graphical-analytical method each provided the best solution in 12 cases, whereas the maximum likelihood method yielded optimal results in only 1 out of 25 instances. An important engineering safety consideration emerges from the observation that in 20 out of 25 cases, the calculated 1% exceedance discharge (Q_1%)exceeds the observed maximum discharge (Q_max), which may be interpreted as an expert criterion for sufficient safety margin in engineering design.

For non-stationary series subjected to trend correction, the Q1% values exhibited an equal probability of either increasing or decreasing (3 cases each), despite the consistent increase in mean values following detrending procedures. In their original form, all non-stationary series were best described by the Pearson Type III distribution. However, after trend adjustment, half of these cases showed improved fit with the Kritsky-Menkel distribution, indicating the significant influence of trend treatment on distribution selection.

The alternative calculation approaches demonstrate substantially greater diversity in both probability distributions and parameter estimation methods, revealing noteworthy patterns in their application outcomes. Among the 25 cases analyzed, the three-parameter lognormal distribution emerged as optimal in 7 instances, followed by the generalized logistic distribution (6 cases), Pearson Type III distribution with L-moments parameterization (5 cases), Vinogradov’s C3 distribution and generalized extreme value distribution (3 cases each), and the two-parameter Gumbel distribution (1 case). Parameter estimation methods exhibited more consistent performance characteristics. From all available alternatives, either the L-moments method (15 cases) or the accuracy-based calibration (s-method, 10 cases) consistently provided optimal solutions. The latter approach, as previously noted, involves direct fitting of analytical curves to empirical data points by minimizing the mean squared differences in discharge values.

For the six non-stationary time series analyzed with both original and trend-corrected approaches, the distribution selection demonstrates specific patterns. In four cases, either the same distribution type was maintained or replaced with a functionally similar alternative (e.g., different parameterizations of gamma-type distributions). However, two cases exhibited more substantial changes - transitioning from Pearson Type III to Generalized Extreme Value distribution in one instance, and from three-parameter lognormal to Generalized Logistic distribution in another.

The parameter estimation methods showed greater stability during trend correction procedures, remaining unchanged in five out of six cases. Notably, the alternative calculation approaches produce more conservative engineering estimates. The computed 1% exceedance discharge (Q_1%) exceeds the observed maximum discharge (Q_max) in 22 of 25 cases (88%), with these exceedances being more pronounced than those obtained through standard regulatory methods. Furthermore, when applying alternative methods to detrended series, Q_1% values more frequently increased than decreased (4 cases versus 2), indicating an inherent tendency toward greater safety margins in the alternative framework.

The quantitative evaluation metrics (reliability ω and precision s) provide compelling evidence for the superiority of alternative approaches. In 20 of 25 cases (80%), these criteria unequivocally indicate better performance characteristics for probability models developed using alternative methodologies compared to standard regulatory solutions.

Conclusion

The results of maximum flow analysis in the study region generally correspond to previously established global patterns in flood frequency distributions (Gubareva 2011). The application of a comprehensive suite of stochastic hydrology tools for probabilistic estimation of extreme discharges in rivers of the East Kazakhstan region leads to several fundamental conclusions.

The analysis demonstrates the clear superiority of alternative computational approaches over the standard SP 33-101-2003 methodology. This conclusion is supported by three key factors:

— the alternative schemes employ heavy-tailed probability distributions that more accurately characterize extreme flood behavior, as documented in hydrological literature (Naidenov, 2004; Gubareva, 2010, 2011);

— the alternative schemes consistently achieve superior performance metrics, showing significantly better values for both reliability (ω) and precision (s) indicators;

— the alternative approaches systematically produce higher estimates for the 1% exceedance discharge (Q1%), thereby providing increased safety margins for flood protection infrastructure design.

The prevalence of different parameter estimation methods among the selected optimal models serves as indirect evidence supporting the greater adequacy of alternative approaches. Within the 25 optimal solutions obtained using the standard methodology, the approximate maximum likelihood method appears only once, despite its theoretical superiority as the most statistically efficient estimation technique (Gubareva, Gartsman 2010). This apparent contradiction can be explained by noting that the method’s advantages are strictly contingent upon correct distributional assumptions. The observed results therefore suggest potential inadequacies in the standard framework’s prescribed probability distributions for modeling extreme flood characteristics. In contrast, alternative methodologies demonstrate fundamentally different patterns. The L-moments method, which shows comparable precision and robustness to maximum likelihood estimation according to established research (Hosking & Wallis, 1997; Gubareva & Gartsman, 2010), appears in over half of the optimal alternative solutions. This striking methodological consistency strongly indicates that the alternative probability distributions better correspond to the actual statistical properties of maximum discharge series. The robust performance of L-moments estimation in this context provides compelling evidence for the theoretical soundness of the alternative distributional models, as the method’s statistical properties are known to be particularly sensitive to misspecification of the underlying probability distribution.

The second key finding of the regional analysis reveals the relatively limited impact of non-stationarity – at least in the form of downward trends – on the estimation of rare flood discharge quantiles. This conclusion emerges from a detailed examination of trend-adjusted series, where maximum discharge records were normalized to the highest moving average observed during the monitoring period. While this adjustment consistently produces mean values significantly higher than in the original series, the resulting changes in 1% exceedance discharge (Q1%) estimates remain within a modest range of -10% to +11% deviation from original values. This observed variation in Q_1% estimates due to trend correction falls well within the combined error envelope encompassing both measurement inaccuracies and computational uncertainties inherent in flood frequency analysis. Constructing the correct and adequate flood frequency curves based on relevant probability distributions may aide to future regional hydrological model development, which can be calibrated using signature measures based on flow duration curves (Gartsman, Solomatine, Gubareva, 2024).

The observed non-stationarity in river flow regimes across Eastern Kazakhstan represents a complex phenomenon requiring systematic investigation. Current evidence points to climate change as the primary driver, characterized by a well-documented rise in mean annual temperatures throughout the region. This warming trend has fundamentally altered the hydrological cycle through several interconnected mechanisms. A well-established redistribution of precipitation patterns has occurred, with meteorological records from RSE “Kazhydromet” and supporting studies documenting decreased summer rainfall alongside increased winter precipitation. This seasonal shift in moisture availability has significantly modified the hydrological behavior of regional rivers (Stambekov, Turulina, 2016, Salnikov et al., 2014). Increased aridity has led to more frequent forest fires, resulting in substantial reductions in forest cover – a critical natural regulator of both surface and subsurface flow regimes (Lebed’, Eserkepova, Suleimenov, 2020). While illegal logging activities represent an additional stressor, their contribution to overall forest loss appears secondary compared to climate-driven impacts. These transformations manifest most distinctly in Group I river basins, including the Oba, Turgysyn, Naryn, and Kalzhyr watersheds, where hydrological changes have been particularly pronounced.

Human economic activities represent a significant factor influencing alterations in hydrological regimes, particularly through agricultural expansion and water infrastructure development (Milanović Pešić 2024). In the Tarbagatai district, for instance, water from the Kandysu River is diverted through an extensive network of irrigation canals to support growing cropland areas. Under increasingly arid climatic conditions, agricultural water withdrawals have risen substantially, creating measurable impacts on river discharge patterns. The observed non-stationarity in flow records for the Turgysyn River may be directly attributed to the construction and commissioning of the Turgysyn Hydroelectric Power Station in 2021. Such hydraulic engineering projects fundamentally modify natural flow regimes through flow regulation, sediment trapping, and alteration of seasonal discharge patterns.

Appendix