INTRODUCTION

A suitable estimate of the expected extreme high flows of a river is a basic step for flood risk assessment and design, operation and management of hydraulic structures (^{PROSDOCIMI; KJELDSEN; SVENSSON, 2014}). Traditional methods developed for that purpose are based on the assumption of stationarity, which implies that the variable under analysis has a time invariant probability density function with fixed parameters (^{PETROW; MERZ, 2009}; READ; VOGEL, ^{2015}, ^{2016}; ^{SRAJ et al., 2016}; ^{VOGEL et al., 2015}; ^{VOGEL; YAINDL; WALTER, 2011}). However, the stationarity hypothesis might be rendered invalid due to several factors that influence streamflow, e.g., hydroclimatic changes (^{MILLY et al., 2008}), urbanization (^{VOGEL; YAINDL; WALTER, 2011}), agricultural management practices (^{FOUFOULA-GEORGIOU et al., 2015}) operation of hydropower plants (^{RÄSÄNEN et al., 2017}) and reservoirs (^{ZAJAC et al., 2017}).

The stationary approach is not appropriate to model hydrologic process in basins with large scale changes resulting from anthropogenic influences (^{KOUTSOYIANNIS, 2006}). ^{Milly et al. (2008)} suggested that “stationarity is dead” in face of hydroclimatic changes and that it is necessary to find ways to identify nonstationary probability models and use them in water resources risk assessment and planning (^{MILLY et al., 2015}). Nevertheless, ^{Montanari and Koutsoyiannis (2014)} advise caution in the use of nonstationary models and the expression “stationarity is dead”. It is considered that deterministic relationships are necessary to explain the evolution of a certain statistical process over time (^{KOUTSOYIANNIS, 2006}; ^{KOUTSOYIANNIS; MONTANARI, 2015}; ^{LINS; COHN, 2011}). The selection of a poor nonstationary model might turn out to be less efficient and less robust than a stationary one (^{MONTANARI; KOUTSOYIANNIS, 2014}).

Despite some criticism about the expression “stationarity is dead” (^{KOUTSOYIANNIS; MONTANARI, 2015}; ^{LINS; COHN, 2011}; ^{MONTANARI; KOUTSOYIANNIS, 2014}), ^{Milly et al. (2015)} recalls that the science suggests a substantial and growing antropoghenic climate change which cannot be readily assumed to be negligible over the decades long horizon of engineered water systems. In this context, the importance of large amount of historical data and the continuity of observations play an essential role and is a point of convergence between ^{Milly et al. (2008)} and ^{Montanari and Koutsoyiannis (2014)}.

Under nonstationary conditions, stationary models tend to underestimate maximum flows (^{SRAJ et al., 2016}). That fact is commonly ignored or a few times acknowledged through the simple use of multiplication factors applied to the results obtained by stationary models (^{PROSDOCIMI; KJELDSEN; MILLER, 2015}). Assuming the simplification of stationarity or only a multiplication factor in hydrological studies might be a risky practice. A change of concepts related to the estimation of the flood return periods and risks assumed in engineering projects is necessary and has led researchers to updated the standard statistical methods in order to account for possible nonstationarity in hydrological time series (^{READ; VOGEL, 2016}; ^{SRAJ et al., 2016}; ^{VOGEL et al., 2015}).

When hydrologic processes are nonstationary, the probability of exceedance associated with a certain event is changing over time, and thus the traditional formula in hydrology *T = 1/p* (where *T* is an average return period and p the probability of exceedance) no longer holds. Under nonstationary conditions, the exceedance probability *p* associated with a particular annual maximum flood discharge changes every year (^{READ; VOGEL, 2015}, ^{2016}). Thus, planning under nonstationary conditions is fundamentally different from planning under stationary conditions (^{VOGEL et al., 2015}).

Despite high streamflow frequency analysis under nonstationarity being a topic widely explored in Europe (^{PROSDOCIMI; KJELDSEN; SVENSSON, 2014}; ^{SRAJ et al., 2016}; ^{VILLARINI et al., 2011}), North America (^{GADO; NGUYEN, 2016}; ^{SADRI; KAM; SHEFFIELD, 2016}; ^{VOGEL; YAINDL; WALTER, 2011}) and Asia (^{DU et al., 2015}; ^{XIONG et al., 2015}; ^{ZHANG et al., 2015}), it has been still timidly discussed in Brazil. That is troubling when taken into account the country’s extensive water network - which is directly related to power generation, flood control and water supply - the intense changes in land use that have been occurring in recent decades and the associated climate changes (^{SALAZAR et al., 2015}, ^{2016}). ^{Clarke (2007)} suggested that the assumptions of constant mean, variance and correlation structure in annual totals, means and extreme values of hydrological variables are invalid under some South American conditions.

^{Detzel et al. (2011)} analyzed inflow series of 146 hydroelectric plants located in Brazil and concluded that 51.4% of them have nonstationary characteristics. It is important to note that all stations from southern Brazil subsystem (30 in total) presented nonstationary characteristics. Nonstationary frequency analysis models associated to climate covariates for flood studies in Itajaí River Basin demonstrated that there is a relationship between climate covariates and floods (^{SILVA et al., 2015}; ^{SILVA; NAGHETTINI; PORTELA, 2016}). ^{Kruger, Kaviski and Muller (1998)} evaluated streamflow and precipitation series from upstream gauges of the Itaipu Dam and confirmed the presence of nonstationarity in fluviometric series. ^{Detzel and Mine (2014)} reviewed some proposed methods for nonstationarity studies and presented application of these in a fluviometric series corresponding to a gauge located in the Iguazu River. ^{Detzel, Fernandes and Mine (2016)} evaluated the effects of nonstationarity approach in the assessment of water availability in six Brazilian fluviometric gauges and found changes in the flow duration curves of the analyzed series.

Although not directly mentioning the term nonstationarity, there are some studies in Brazil that evaluated the influence of climate change on the fluviometric behavior using hydrologic models like MGB-IPH (^{ADAM et al., 2015}; ^{ADAM; COLLISCHONN, 2013}; ^{RIBEIRO JUNIOR; ZUFFO; SILVA, 2016}) and SWAT (^{ARROIO JUNIOR; MAUAD, 2015}; ^{VALÉRIO; FRAGOSO JUNIOR, 2015}). These researches have distinct conclusions, depending on the region and climate scenarios. In the same way, another group of studies focus in the identification of trends in streamflow series utilizing statistical techniques such as Mann-Kendall test for monotonic trend and Pettitt test for abrupt change detection (^{ROSIN; AMORIM; MORAIS, 2015}; ^{ULIANA et al., 2014}). The main difference between abrupt and gradual changes is that when a trend is detected, it is likely to continue in the future, while the presence of an abrupt change indicates a shift from one regime to another, and the status is likely to remain the same until a new shift occurs (^{VILLARINI et al., 2011}).

It is usually recommended that those analyses should be carried out from a regional perspective (^{HALL et al., 2014}; ^{PETROW; MERZ, 2009}; ^{VILLARINI et al., 2011}). When analyzing large regions and data sets, it is possible to reduce local noise and identify clearer spatial patterns in the observed changes. ^{Villarini et al. (2011)} found that abrupt changes, rather than monotonic trends, are responsible for violations of the stationary assumption in central Europe.

In the United States, there is a higher concentration in heavily urbanized regions of fluviometric stations with nonstationary characteristics (^{VOGEL; YAINDL; WALTER, 2011}). ^{Prosdocimi, Kjeldsen and Svensson (2014)} did not find a systematic spatial pattern for annual streamflows in United Kingdom; however, they found some scattered clusters of increasing and decreasing trends.

We could not find regional scale studies dealing the high streamflow frequency analysis under nonstationary conditions in Southern Brazil. In this context, some of the open questions that we have analyzed are:

i. Are there identifiable nonstationary signals in the flood series of Southern Brazil?

ii. Are there identifiable spatial patterns of nonstationarity in those series?

iii. How do the flood frequency curves change when considering the nonstationary approach?

iv. Is the presence of nonstationarity related to basin size?

The objective of this work was to evaluate the presence of nonstationarity conditions in maximum annual daily streamflows series of Southern Brazil. We used a nonstationary frequency analysis model developed by ^{Vogel, Yaindl and Walter (2011)} and we have estimated the return period associated with a high streamflow initially estimated as the 100 years flood in stationary conditions.

MATERIAL AND METHODS

Study area

Southern Brazil is located between the latitudes 22°S and 34°S with a total area of 576,774 km^{2} (Figure 1). The population density of Southern Brazil is 48.58 inhabitants km^{-2}, the second highest value among the Regions of Brazil (^{IBGE, 2010}). The mean annual precipitation varies between 1200 and 1900 mm. The climate is characterized by high contrasts in precipitation and temperature regime (^{GRIMM, 2009}). The northern part of the region has a more accentuated seasonality, with a larger difference in total winter and summer precipitations. In the south, predominates a characteristic regime of the midlatitudes, with precipitation distributed uniformly during the year. Topographic effects are present in some areas. The two main biomes of the region are Atlantic Forest and Pampas (^{IBGE, 2004}).

The installed hydroelectric power generation capacity in the Southern region corresponds to 26.9% of the national total – 24.622 MW. This value does not consider yet the 7.000 MW installed in the Itaipu Dam that belong to Paraguay, but much of it is sold to Brazil (^{EPE, 2016}).

Fluviometric data selection

We obtained streamflow data from the Hidroweb Portal of the Brazilian National Water Agency (ANA) – a total of 765 gauge stations – and from the Brazilian National Grid Operator (ONS) website – a total of 38 gauge stations (Figure 1), the ONS series refer to natural flows, which means that all effects of the installation and operation of the reservoirs, as well as anthropic actions in the water courses are removed (^{DETZEL et al., 2011}). All historical series refer to the maximum period available for each fluviometric gauge station. Only those series that had at least 30 years of data were used, as recommended by World Meteorological Organization (^{WMO, 1994}) for statistical studies. We used consisted data as much as possible; however, non-consisted data recently added to the series were also used, as well as series that had only this type of data.

The maximum daily value discharge for each year was selected in order to deal with the analysis of annual maxima series (AMS). In addition to the minimum record length, we checked the AMS for the presence of failure. We adopted the methodology proposed by ^{Papalexiou and Koutsoyiannis (2013)}, in which it is evaluated the percentage of missing values registered in the years corresponding to the 40% lower values of maximum annual daily streamflow series. If in any of these years the failure percentage is equal or superior to 30%, the series was discarded.

Frequency analysis – stationary case

In order to explore trends in annual maximum flow frequency analysis and their consequences on flood frequency analysis is necessary to assume a probabilistic model (^{COSTA; FERNANDES, 2015}; ^{VOGEL; YAINDL; WALTER, 2011}). In this work, the two parameter lognormal probability distribution function – LN was employed. The lognormal distribution is the one recommended for use in Brazil (^{COSTA; FERNANDES 2015}). ^{Vogel, Yaindl and Walter (2011)} also used the lognormal distribution and showed that when it is combined with a log-linear trend model, it produces a simple flood frequency nonstationary model.

In this work a frequency analysis model developed by ^{Vogel, Yaindl and Walter (2011)} was used, which was applied by the same authors in United States streamflow series and with some modifications by ^{Prosdocimi, Kjeldsen and Svensson (2014)} in United Kingdom fluviometric series.

Based on the assumption that a maximum annual daily streamflow series *Q _{t}* (composed by

*Q*values, AMS approach) follows a lognormal distribution probability with two parameters, the relative quantile to annual maximum flow

_{1}, Q_{2},…,Q_{n}*p*or return period

*T=1/p*can be given by (stationary case):

where
*Q*, respectively.
*p*. This equation is said to be stationary because it assumes that the moments of *y* (where *y = ln Q*) given by

Frequency analysis – nonstationary case

^{Vogel, Yaindl and Walter (2011)} report that a simple exponential model of *Q* (annual maximum flow) versus *t* (time) provides an excellent approximation between the two variables. Besides that, the use of an ordinary least squares (OLS) method is an easily way to fit a log linear trend model that describe the relationship between the two variables. Considering the model,

**ε _{t}** the errors model.

According to ^{Vogel, Yaindl and Walter (2011)}, if the estimated model parameters are significantly different from 0, the regression model can provide an estimate of the average natural logarithms of *Q* (defined as

An estimate of the intercept term of the trend model can be calculated by ordinary least squares (OLS) by

where

The trend model can be derived by combining the Equation 3 with Equation 4

where

Substituting the nonstationary trend model (7) for

where
*p*.

In this work, with the purpose of verifying if a series can be considered nonstationary, only those series wich present the slope trend model significantly bigger than zero were used, calculated by one-sided Student’s t-test. As required for Student’s t-test, it was verified whether the model residuals are normally distributed by Anderson-Darling test. Independence and Homocedasticity of the residuals were evaluated by Durbin-Watson and White’s tests, respectively. Level of significance p<0.05 was used for all these tests.

Recurrence reduction

^{Vogel, Yaindl and Walter (2011)} define the term “recurrence reduction” as the average time (
_{f} associated with the flood with an average recurrence interval of *T _{0}* in some reference year t

_{0}. According to the authors, the recurrence reduction can be given by

where the function

This work assumed, as in ^{Prosdocimi, Kjeldsen and Svensson (2014)} and ^{Vogel, Yaindl and Walter (2011)}, a recurrence reduction associated with a flood which currently has a T_{0} = 100 years in year t_{0} and the planning horizon was a decade, so

Abrupt changes – Pettitt Test

As ^{Gebremicael et al. (2017)}, ^{Ma et al. (2008)}, ^{Nka et al. (2015)}, ^{Ribeiro Junior, Zuffo and Silva (2016)}, ^{Villarini et al. (2011)}, ^{Zhang et al. (2014)}, ^{Wagesho, Goel and Jain (2012)}, Pettitt test (^{PETTITT, 1979}) was used to detect abrupt changes in time series and a probable date for the abrupt change. According to ^{Villarini et al. (2011)} Pettitt test is a rank-based test which allows for the testing of whether abrupt changes in the mean of the variable of interest at an unknown point in time. It is based on the Mann–Whitney sample test, and allows testing whether two samples (X_{1}, …, X_{M} and X_{M+1}, …, X_{N}) come from the same population.

Pettitt statistic test (

where

The Equation 10 is applied to 1≤t<T, where T is the size of the series. The null hypothesis (H_{0}) of the Pettitt test admits absence of abrupt changes in the series, while de alternative hypothesis (H_{1}) admits a change point. Its statistic

In this work was adopted p<0.05 for Pettitt test. It was verified the period that occurs abrupt changes

RESULTS AND DISCUSSION

Preliminary data selection and trend analysis

There were 157 fluviometric series that attended all the minimum 30-year length and the missing value percentage criteria, Figure 2 shows the coverage period for each of these streamflow series, while the location of gauges are in Figure 3, which also shows the gauges whose series had the linear trend model slope coefficient significantly bigger than zero. Figure 4 represents the absolute frequency histogram corresponding to the basin area of the selected series.

We can observe that those selected gauges are distributed uniformly across the study region, with a few gauges in the south of Rio Grande do Sul State. The series that presented the linear trend model slope coefficient significantly bigger than zero are located mainly in the Uruguay, Iguazu and Paranapanema basins.

The fact that the slope coefficient of the linear trend model is significantly bigger than zero represents an important information when evaluating the presence of nonstationarity in historic fluviometric series. According to ^{Prosdocimi, Kjeldsen and Svensson (2014)} the violation of the normality premises, fact which occurred to some series evaluated in this work, can be related to the presence of the extremely high values in the series, for example.

As an example we selected the streamflow series corresponding to Itaipu (Parana River) and Canoas II (Paranapanema River) hydroelectric power plants affluence, whose series are presented in Figure 5. There is an increasing trend of the maximum annual streamflow corresponding to Itaipu station and Canoas II station, both in terms of mean and variance. In fact, there is a similar pattern in the two series, which are gauges located in the Paraná River basin. Additionally, an abrupt change in the maximum streamflows is observed starting in the 70’s-80’s; which is discussed in the next section.

Frequency analysis and recurrence reduction

There were 47 series identified as nonstationary with slope coefficient of the log linear trend model significantly different from zero and the model residuals obey the assumption of normality, homoscedasticity and independence (Figure 6). Only one of the analyzed stations presented negative trend in the maximum annual flow. The recurrence reduction (for a planning horizon equal to ten years) associated with a current 100-year were calculated according to the model of ^{Vogel, Yaindl and Walter (2011)}. The recurrence reduction term was previously described and can be understood as a new return period, associated to nonstationary condition.

When taken into account the nonstationary frequency analysis model, there are significant changes in return period initially estimated as 100 years for the stationary condition (Figure 6). In extreme cases, the return period for the nonstationary condition becomes less than 50 years. Similar results were obtained by proponents of the model in fluviometric series of the United States, their conclusion is that high streamflows initially associated to T_{0} = 100 year will be much more frequent in the future.

^{Sraj et al. (2016)} also found a great difference between maximum streamflow values estimated by stationary and nonstationary models for series from Europe. They used the Extreme Value Generalized probabilistic distribution with their respective parameters estimated as function of time. They found that a high streamflow initially estimated as 100 years recurrence interval in stationary conditions might become 20 years in the non-stationary conditions.

As an example of the analysis, we show the high streamflow frequency curves for the stationary and nonstationary conditions corresponding to fluviometric series located upstream to Taquaruçu (Paranapanema River), Salto Osório (Iguazu River) and Itá (Rio Uruguay) hydroelectric power plants (Figure 7). Considering the nonstationary condition, there is an increase in estimated streamflow for all nonstationary series. For example, a high streamflow associated with a return period equal to 100 years in Salto Osorio changes from 14.800 m^{3}/s to 16.400 m^{3}/s approximately, corresponding to 11% increase.

Basin area vs nonstationarity

Figure 8 show the area distribution of all basins whose series were evaluated and from those considered nonstationary. We can observe that the relative frequency of the higher areas is great for the group of series considered as nonstationary (All criteria and p-value criteria) in relation to all basins evaluated (All).

The possible drivers of changes in fluviometric regimes are atmosphere (e.g. changes in precipitation), catchment (e.g. land use change) and river (e.g. construction of dams) (^{BLÖSCHL et al., 2015}; ^{MERZ et al., 2012}; ^{VIGLIONE et al., 2016}). ^{Blöschl et al. (2015)} discuss that changes in synoptic precipitation characteristics are relevant for large basins (hundred to hundred thousands of square kilometers), while changes in convective precipitations are relevant in small basins (hundreds of square kilometers or less). The effect of changes in land use on runoff and consequently on the flood frequency and magnitude decreases with the size of the basin (^{BLÖSCHL et al., 2007}, ^{2015}; ^{PETROW; MERZ, 2009}; ^{VIGLIONE et al., 2016}).

In this work, despite the limited availability of fluviometric series from small basins, it is observed that the presence of nonstationarity is predominant in basins with area greater than one thousand square kilometers. Therefore, the considerations of the previous paragraphs lead to the hypothesis that the nonstationarity of most of the series evaluated are predominantly associated to climatic factors. That is also the widely accepted hypothesis in the literature (^{DOYLE; BARROS, 2011}; ^{GENTA; PEREZ-IRIBARREN; MECHOSO, 1998}; ^{SAURRAL; BARROS; LETTENMAIER, 2008}). The use of the term “predominantly” is emphasize, since different drivers act in parallel in a basin, and there are interactions between them (^{MERZ et al., 2012}; ^{VOGEL; YAINDL; WALTER, 2011}).

Other studies related to fluviometric regimes changes in Southern Brazil and their possible causes focused on trend identification, especially in relation to mean streamflow values (^{DETZEL; MINE, 2014}; ^{DOYLE; BARROS, 2011}; ^{GENTA; PEREZ-IRIBARREN; MECHOSO, 1998}; ^{KRUGER; KAVISKI; MULLER, 1998}; ^{SAURRAL; BARROS; LETTENMAIER, 2008}).

^{Detzel and Mine (2014)} evaluated a fluviometric series from upper Iguazu basin (Porto Amazonas) and found positive trend in mean and high streamflows. The authors consider that this fact could be interpreted result of land use change that occurred in the last 30 years in the basin, mainly substitution of forest for agricultural crops. However, they emphasize that precipitation and streamflow series are similar, which raises questions about the real influence of land change on the fluviometric regime.

^{Doyle and Barros (2011)} studied series from lower Iguazu and Uruguay basin (Salto Caxias and Salto Grande gauges, respectively) and found an significant increase in the mean streamflow between 1960-1979 and 1980-1999 periods for both gauges (approximately 40% for Salto Caxias). The authors hypothesized that this behavior is due to an increase in the precipitation and, of lesser importance, due to land change. This is the same conclusion obtained from ^{Saurral, Barros and Lettenmaier (2008)} for Uruguay basin.

^{Kruger, Kaviski and Muller (1998)} studied trends in mean streamflow of the rivers affluent to Itaipu Dam (which includes some series evaluated in this work, mainly in the Paranapanema basin) and concluded that the significant increase in this component can be explained by increase in precipitation and infiltration (infiltration facilitated by adequate soil practices) and decrease in evapotranspiration (due to removal of native forest). According to the authors, the streamflow series corresponding to hydroletric power plants located in the Paranapanema River can not be considered stationary since the ratio between the annual mean streamflow before and after the year 1970 is 1.45 in this basin.

Abrupt changes

We identified abrupt changes in 54 out of the total series we analyzed using the Pettitt test (Figure 9). Most of those are located in Iguazu, Uruguay and Paranapanema basins. Similar results were found for the gauges that attended the trend model premises developed by ^{Vogel, Yaindl and Walter (2011)}.

Many of the series that attended the trend model premises also presented significant abrupt change (25 of the total 46). There were 21 stations in which we found a significant trend model and there was no abrupt changes found. Therefore, when adding all this subsets, the total of 75 series can be considered as nonstationary in Southern Brazil.

Abrupt changes are concentrated around the 70’s (Figure 10). According to ^{Carvalho et al., (2014)}, there was an increase in the total maximum precipitation average in the southern Brazil in the 70’s. That phenomenon can be explained by abrupt changes in the El Niño-Southern Oscillation (ENSO) behavior that occurred in this period. ^{Doyle and Barros (2011)} observed a relationship between the positive phase of the Pacific Decadal Oscillation (PDO) and an intense El Niño. The last positive phase of the PDO started in the 70’s and remained until at least to the end of the century.

Fluviometric regimes can be altered as a consequence of climatic phenomena, like ENSO and PDO (^{DAI et al., 2009}; ^{ALVES; SOUZA FILHO; SILVEIRA, 2013}). El Niño tends to increase the streamflow in many rivers around the world, including the Paraná and Uruguay rivers (^{DAI et al., 2009}). ^{Alves, Souza Filho and Silveira (2013)}, when analyzing average and maximum streamflow series from all parts of Brazil found that positive changes agreed with El Niño and PDO occurrence for most of the maximum series evaluated.

In a regional perspective, ^{Detzel and Mine (2014)} identified in a series from Iguazu River that the two largest streamflow events occurred in years that had an intense manifestation of the El Niño, in 1983 and 1998. ^{Silva et al. (2015)} point out that the increase in ENSO activity after the 70’s has led to the occurrence of higher maximum streamflow values in Itajai basin.

CONCLUSIONS

There are many series with nonstationary characteristics in Southern Brazil, mainly of fluviometric gauges located in Iguazu, Uruguay and Paranapanema basins. For many of those series, there is a significant difference between the return periods estimated using the stationary and nonstatioionary models. This fact must be taken into account in flood frequency analysis for hydraulic structures projects and water management. Data indicated that series from large basins are more likely to present nonstationarity.

The quantification of the weight of each possible variable inducing the nonstationarity (climate or land change, eg.) is a complex question to answer, since several drivers act in parallel in a basin. Climatic factors are pointed out as the predominant cause for the fluviometric regime alterations in all South America. We emphasize that these factors are not necessarily related to anthropogenic activities.

In relation to abrupt changes, we identified that they occurred mainly in the 70’s, the same period in which many researchers point out that a significant change in the precipitation regime occurred in Southern Brazil.

The questions about the adoption of nonstationary models are relevant. However, the applicability of this approach in hydraulic structures projects and water management must be further evaluated. Nonstationarity is still little discussed in relation to the stationary traditional methods adopted in the high streamflow frequency analysis. In order to determine which type of model is most appropriate for flood frequency analysis, further evaluation and comparision of nonstationary and stationary models and the continuity of hydrological monitoring across several basins play a key role in determining the most appropriate methods for flood frequency analysis.