Multifractal Analysis of Standardized Precipitation Evapotranspiration Index in Serbia in the Context of Climate Change

1. Introduction

Drought is a slow-onset natural phenomenon that can occur anywhere in the world as a prolonged dry period relative to local climate conditions. It can affect agriculture and energy production [1,2], ecosystem function [3,4] human livelihood and health [5,6]. Knowledge of drought variability on multiple temporal and spatial scales is crucial for drought monitoring and development of mitigation methods.

Various standardized indices were developed to quantify water deficit and surplus with respect to normal conditions, and can be used to evaluate the severity, duration, and spatial extent of drought [7]. Recently, Yildirim et al. [8] performed bibliometric analysis of publication data on drought and found that the most frequently used indices for meteorological drought are Standardized Precipitation Index – SPI [9] based only on precipitation data, Palmer Drought Severity Index – PDSI [10], and Standardized Precipitation Evapotranspiration index – SPEI [11] which are based on precipitation and temperature data.

The changes in drought conditions have been intensively studied for different regions of world, including the analysis of historical data [12,13,14,15] and future projections [16,17,18]. Ionita and Nagavciuc [14] analyzed trends in drought indices in Europe for the period of 120 years (1901-2019). They found that SPEI12 showed significant negative trend in the Central Europe (CEU) and the Mediterranean region (MED) indicating that these regions are becoming dryer due to an increase in the potential evapotranspiration and mean air temperature. North Europe (NEU) is becoming wetter (SPI12 showed significant positive trend), for MED and CEU drying trend based on SPEI was confirmed at various timescales, while the drought evaluated by the SPI (that considers only precipitation deficit) shows the opposite or no changes. Spinoni et al. [17] investigated changes in drought occurrence, frequency, and severity in Europe until 2100 by using regional climate model (RCM). They found that under the moderate emission scenario (which for the end of 21^st century assumes CO₂ equivalent of about 650 ppm and global temperature increase of approximately 1.8 °C) severity and frequency of droughts will increase in the Mediterranean area, western Europe, and Northern Scandinavia, while the most severe emission scenario (CO₂ equivalent of about 1370 ppm and global temperature increase of approximately 4.0 °C) projects that whole European continent, with the exception of Iceland, will be affected by more frequent and severe extreme droughts, especially after 2070.

Serbia is a continental European country with southern and central part located in the Mediterranean region, and northern part in Central Europe. As other MED and CEU countries, over the last decades Serbia experienced increasing trend in drought frequency, severity and duration, causing a reduction in agricultural yield, especially in the northern province of Vojvodina, where more than half of arable land is used for maize, soybeans and sugar beets, and during the dry year the decrease of yield can reach more than 60% [19]. The concerns on evolution of drought conditions in Serbia resulted in increasing number of studies that analyze various aspects of this phenomenon on country level, such as spatial and temporal variability [20], drought regionalization [21] and trends [22]. In all these studies SPI was used to classify droughts categories. Recently Amiri and Gocic [23] compared the performance of five precipitation-based drought indices with SPI (as reference drought index) for drought analysis in Serbia on monthly, annual, and seasonal scale. Based on monthly precipitation data from 28 synoptic stations for the period 1946-2019, they found that rainfall anomaly index (RAI) has the most similar performance to the SPI index, while lowest similarity was found for China Z index (CZI) and modified China Z index (MCZI). Precipitation based indices such as SPI were extensively used to analyze drought spatial and temporal variability in different parts of the world [24,25,26,27]. However, the definition of SPI assumes that droughts are mainly controlled by the temporal variability in precipitation and do not consider other variables such as temperature, evapotranspiration, and soil water holding capacity, that can influence the development of drought conditions. Vicente-Serrano et al. [11] introduced standardized precipitation evapotranspiration index (SPEI) which is based on precipitation and temperature data, where temperature is used in calculation of potential evapotranspiration (PET). SPEI is simple for calculation and can be calculated for different timescales (same as SPI) but in addition it includes the effect of temperature and is considered more suitable for drought detection and monitoring, and for investigation of the consequences of global warming [28,29,30,31,32]. Bezdan et al. [33] developed SPEI based index for agricultural drought monitoring for different crops (AD-${SPEI}_{crop}$) in Vojvodina region in northern Serbia. The AD-${SPEI}_{crop}$ is obtained by replacing potential evapotranspiration ${ET}_0$ for potential crop evapotranspiration ${ET}_c$ (based on local agro-climatic conditions) in SPEI calculation. They found that AD-SPEI_crop showed a higher degree of correlation with yields of the 11 significant crops in the region than original SPEI index and can successfully detect dry and wet periods like other drought indices (SPI, SPEI and PDSI).

As a complement to these classical approaches, in this work we investigated the complexity of dry/wet conditions in Serbia through the multifractal analysis of the Standardized precipitation evapotranspiration index (SPEI). Multifractality of climate phenomena is well known and has been identified for climate variables such as rainfall [34,35], temperature [36,37,38], wind [39,40], and humidity [41]. While multifractality of rainfall temporal series was extensively studied [34,35], much less is known about multifractal properties of rainfall-based drought indices. Recently Adarsh et al. [42] analyzed SPI time series at three aggregation timescales (3, 6, 12 months) from 30 meteorological subdivisions of India, for the period 1871-2016. By using multifractal detrended fluctuation analysis - MFDFA [43], they found that SPI time series exhibit multifractal properties that vary spatially throughout studied area, with stronger persistence and higher multifractality degree for larger time scales. Ogunho [44] applied MFDFA on SPI series (at 1, 6, 12, and 24-month temporal scale) for the period 1980-2010 over different locations across Nigeria. He found that all SPI series show multifractal properties that vary spatially but are unrelated to climatic delineation, and that the strength of multifractality (quantified by the width of the multifractal spectrum) increases with increasing temporal scale. Da Silva et al. [45] used MFDFA to analyze SPI (for 1, 3, 6, and 12-month accumulation periods) at 133 locations across the state of Pernambuco, Brazil, with climate varying from tropical humid in coastal area to semiarid in deep inland region. The results show that all SPI series display multifractal dynamics, with persistent long-term correlations, and multifractality that becomes stronger as accumulation time increases. For all temporal scales the SPI series display stronger persistence in the coastal region with tropical humid climate than in in dry semiarid inland region. For SPI-12, which describes long term dry/wet conditions, multifractality is strongest in the semiarid region. Recently Zhan et al. [46] investigated spatiotemporal characteristics of multiscale drought in the Yellow River Basin (YRB) by analyzing multifractal properties of SPEI (using MFDFA) at different timescales (1, 3, 6, and 12-month), derived for 80 meteorological stations from 1961 to 2017. They found that all SPEI series demonstrate multifractality with persistent long-term correlations and dominance of small fluctuations. With the increase of temporal scale drought persistence increases, being strongest in the mid-temperate semi-arid zone. The strength of multifractality measured by the width of multifractal spectrum also increased with the aggregation timescale.

In this work we investigate the complexity of temporal fluctuations of dry/wet conditions in Serbia by applying MFDFA method on multiscale SPEI (1, 3, 6 and 12-month), derived from the high-resolution gridded climate dataset. The time span is from 1961 to 2020 which permits to investigate the impact of climate change. We calculate multifractal parameters that describe following properties of SPEI temporal series: persistence (quantified by the position of the maximum of the multifractal spectrum), degree of multifractality (quantified by the spectrum width) and the contribution of small/large fluctuations (quantified by the spectrum asymmetry). By comparing the spatial distribution of these parameters over Serbia for two periods of 30 years (1961-1990 and 1991-2020) we investigate: i) the relation of multifractal parameters with accumulation timescale of SPEI; ii) the relation of multifractal parameters with terrain topology, and iii) the changes from first to second sub-period.

2. Materials and Methods

2.1. Study Area and Dataset

Serbia is a continental country in the southeast part of Europe on the Balkan peninsula. The geographical location of Serbia shown in Figure 1 is between latitudes 41°50′ and 46°10′ N which belongs to the temperate climate zone. The northern part of the country with flat and low elevation terrain is entirely located within Pannonian plane, the elevation increases towards central and southern part which are covered by hills and mountains intersected by river valleys. Its climate is shaped by the influence of Alps (Dinaric Mountain), Mediterranean Sea, Pannonian Plains, Morava River valley, Carpathian, Rhodope and Balkan Mountains. There are three climate types: moderate continental in the northern part, continental climate in the central part and modified Mediterranean climate in the southern and southwestern parts [47]. In the lowlands mean annual temperature is between 11 and 12 °C and mean annual precipitation is between 500 and 700 mm. In the mountains, the mean annual temperature is below 8 °C, while mean annual precipitation is above 1000 mm [48].

This study uses daily mean temperature and daily accumulated precipitation data from the E-OBS gridded dataset, version 27.0e with a horizontal resolution of 0.1 degrees (approximately 12 km x 12 km) in regular latitude/longitude coordinates during the period of 1950–2020 [49]. The dataset was created by interpolating the most comprehensive collection of station data available across Europe [50], including 31 stations with temperature observations, and 29 stations with precipitation observations from Serbia.

2.2. Standardized Precipitation Evapotranspiration Index

Standardized precipitation evapotranspiration index (SPEI) was introduced by Vincente Serrano et al. [11] as a multiscalar drought index that is based on precipitation (P) and potential evapotranspiration (PET). It is temperature sensitive and thus more suitable for drought studies, especially those related to observed climate change, future climate scenarios and global warming. The procedure for calculation of SPEI is described in [11]. It is based on the calculation of SPI (albeit with a different probability distribution), which uses monthly precipitation as an input variable. In SPEI input variable is the monthly difference between precipitation and PET

$$D_j=P_j-{PET}_j$$

(1)

which quantifies water deficit or surplus (climate water balance) for month j. In this work PET is calculated using Thornthwaite method [51] which requires only monthly-mean temperature and geographical position of a certain point. The accumulated difference $X_{i,j}^k$ in a month j and year i on time scale k (k=3,6,12 month, etc.) is calculated as

$$X_{i,j}^k=\left\{\begin{array}{cc}\begin{array}{c}{\displaystyle \sum _{l=13-k+j}^{12}}{D}_{i-1,j}+{\displaystyle \sum _{l=1}^j}{D}_{i,j} \\ \end{array}& ,j<k\\ {\displaystyle \sum _{l=j-k+1}^j}{D}_{i,j}& ,j\ge k\end{array}\right. .$$

(2)

The SPEI was obtained by standardizing D series utilizing the three-parameter Gamma distribution, with the use of the python library SPEI [52]. According to [11], in three-parameter distributions $x$ can assume negative values, which are common in $D$ series. We analyzed SPEI on accumulation time scales of 1 month (SPEI-1), 3 months (SPEI-3), 6 months (SPEI-6) and 12 months (SPEI-12), in order to access meteorological, agricultural and hydrological droughts. The drought classification based on SPEI which is the same as for SPI [9] is presented in Table 1.

2.3. Multifractal Detrended Fluctuation Analysis

Kantelhardt et al. [43] introduced the method Multifractal detrended fluctuation analysis (MFDFA) which was shown to produce reliable results for multifractal time series that require a hierarchy of scaling exponents for full description of long-term correlations [53]. Over the last two decades the MFDFA was widely used for data analysis in various fields, such as physiology [54], hydrology [55], geophysics [56], engineering [57], ecology [58] and finances [59].

The MFDFA algorithm is performed through the following steps [43]. It starts with the integration of original time series $x\left(i\right), i=1,\dots ,N$, to produce $X\left(k\right)={\sum }_{i=1}^k\left[x\left(i\right)-〈x〉\right], k=1,\dots ,N$, where $〈x〉$ is the average of $x\left(i\right).$ Next, the integrated series $X\left(k\right)$ is divided into $⌊N/n⌋$ (here $⌊·⌋$ stands for the floor function) non-overlapping segments of length $n$ starting from $X\left(1\right)$, and another $⌊N/n⌋$ segments of size $n$ starting from $X\left(N-n⌊N/n⌋+1\right)$, so that the last segment ends at $X\left(N\right)$ and thus all values of the series $X\left(k\right)$ are included in segments division. In each segment $\nu =1,\dots ,N_n$, where $N_n=2⌊N/n⌋$ is the total number of segments, the local trend$X_{n,\nu }\left(k\right)$ (linear or higher order polynomial least square fit) is estimated and subtracted from$X\left(k\right)$ to calculate the detrended variance

$$F^2\left(n,\nu \right)=\frac{1}{n}\sum _{k=\left(\nu -1\right)n+1}^{\nu n}{\left[X\left(k\right)-X_{n,\nu }\left(k\right)\right]}^2 .$$

(3)

The average over all segments is used to obtain the $q$th order fluctuation function

$$F_q\left(n\right)={\left\{\frac{1}{N_n}\sum _{\nu =1}^{N_n}{\left[F^2\left(n,\nu \right)\right]}^{q/2}\right\}}^{1/q} ,$$

(4)

where, $q$ can take on any real value except zero. For the special case q = 0 the fluctuation function is given by

$$F_0\left(n\right)=exp\left\{\frac{1}{N_n}\sum _{\nu =1}^{N_n}ln\left[F\right(n,\nu \left)\right]\right\} .$$

(5)

Parameter $q$ serves as a “magnifying glass”, where positive values of $q$ emphasize large fluctuations, and negative values enhance small fluctuations.

This procedure is repeated for various box sizes, providing the relationship between fluctuation function $F_q\left(n\right)$ and box size$n$. If the long-term correlations are present, the fluctuation function increases with $n$ as a power law $F_q\left(n\right) ~ n^{h\left(q\right)}$. The scaling exponent $h\left(q\right)$ is obtained as the slope of linear regression of $\log F_q\left(n\right)$ versus $\log n$. It is called the generalized Hurst exponent (as for stationary time series $h\left(2\right)$ is identical to the classical Hurst exponent $H)$. The exponents $h\left(q\right)$ describe the scaling behavior of subsets of series with large fluctuations ($q>0$) and subsets with small fluctuations ($q<0$). For multifractal series subsets with small and large fluctuations of the series scale differently and $h\left(q\right)$ is a decreasing function of $q$, while for monofractal time series, $h\left(q\right)$ is constant. The series that generates larger range of $h\left(q\right)$ values is considered to display a higher degree of multifractality.

Multifractality of time series can also be described by using the properties of multifractal spectrum $f\left(\alpha \right)$ which is obtained through the Legendre transform

$$\alpha \left(q\right)=\frac{d\tau \left(q\right)}{dq} ,$$

(6)

$$f\left(\alpha \left(q\right)\right)=q\alpha \left(q\right)-\tau \left(q\right) ,$$

(7)

where $\tau \left(q\right)=qh\left(q\right)-1 .$ A monofractal process in this approach is represented by a single point in the $f\left(\alpha \right)$ plane, while a multifractal process generates a single humped function [43].

Multifractal spectrum $f\left(\alpha \right)$ contains the information about complexity of the underlying stochastic process, described by the following complexity parameters: the position of maximum ${\alpha }_0$, the width of the spectrum $W={\alpha }_{max}-{\alpha }_{min}$, and the skew parameter$r=({\alpha }_{max}+{\alpha }_{min}{-2\alpha }_0)/ {(\alpha }_{max}-{\alpha }_{min})$ where $-1\le r\le 1$ $(r=0$ for symmetric shapes, $r>0$ for right-skewed shapes, and $r<0$ for left-skewed shapes) [60,61,62]. If ${\alpha }_0>0.5,$the underlying process is overall persistent (larger value of ${\alpha }_0$ indicates stronger persistency), and if ${\alpha }_0<0.5$ the process is overall antipersistent (smaller value of ${\alpha }_0$ indicates stronger antipersistency). The width $W$ of the spectrum measures the degree of multifractality of the process (the larger $W$ value indicates the wider the range of the fractal exponents, and stronger multifractality). The skew parameter $r$ indicates which fractal exponents are dominant: if $f\left(\alpha \right)$ spectrum is right-skewed $\left(r>0\right)$ the process is characterized by “fine structure” (large scaling exponents describing small fluctuations), and if the $f\left(\alpha \right)$ spectrum is left-skewed $\left(r<1\right)$ the process is dominated by the scaling of large fluctuations (small scaling exponents). The signals that generate multifractal spectrum with a wide range $W$ of fractal exponents (higher degree of multifractality), and a right-skewed shape ($r>0$) may be considered more complex than those with opposite characteristics.

3. Results and Discussion

We applied MFDFA method on SPEI-1, SPEI-3, SPEI-6 and SPEI-12 temporal series for all the 1278 grid cells for the two periods (1961-1990 and 1991-2020), which totals 10224 MFDFA runs. The multifractal spectra for SPEI-1 for the two periods, for a sample grid point at latitude 43.15 and longitude 22.45, corresponding to the city of Pirot, are shown in Figure 2 together with the position of this grid point on the map of Serbia. The multifractal parameters ${\alpha }_0$, $W$ and $r$ are then extracted from these spectra, for all the grid points, and all the accumulation periods.

The spatial distribution of multifractal parameters over the Serbian territory for four SPEI accumulation time scales and for two analyzed periods, interpolated across the map of Serbia, is shown in Figure 3, Figure 4 and Figure 5, where we observed the following patterns. The values of ${\alpha }_0$ (position of maximum of $f\left(\alpha \right)$ spectrum) shown in Figure 3 indicate that: i) SPEI-1, SPEI-3 and SPEI-6 series are persistent (${\alpha }_0>0.5$, persistence becomes stronger as accumulation period increases), while SPEI-12 showed $1.0<{\alpha }_0<1.5$ indicating that increments of SPEI-12 are anti-persistent ($0< {\alpha }_0-1<0.5)$; ii) in the first period SPEI-1, SPEI-3 and SPEI-6 series are more persistent in northern lowland part of the country indicated by higher ${\alpha }_0$ values, while SPEI-12 series are more anti-persistent in central and southern part with mountainous terrain, indicated by lower ${\alpha }_0$ values; iii) in the second period the values of ${\alpha }_0$ increased for all accumulation time scales, indicating that short term and medium term dry/wet conditions described by SPEI-1, SPEI-3 and SPEI-6, which are related to short and medium term soil moisture, become more persistent, which is relevant for agricultural planning, while SPEI-12 series which are related to streamflow, reservoir levels, and groundwater levels become less anti-persistent (approaching random regime).

As seen in Figure 4 the width $W$ of multifractal spectrum also varies spatially: i) for SPEI-1, SPEI-3 and SPEI-6 the width of the spectrum increases (indicating stronger multifractality) with accumulation time, and becomes smaller for SPEI-12 in the first period; ii) in the first period the clear north to south positive gradient in $W$ is observed for SPEI-1 and SPEI-3 indicating stronger multifractality in short and medium dry/wet conditions in mountainous regions; for SPEI-6 stronger multifractality is observed in southwestern Dinaric Mountain region, in eastern Carpathian Mountains, and for SPEI-12 in border southern region that encompass Rhodope Mountains, Balkan Mountains and southern part of Dinaric Mountains; iii) in the second period the width of the spectrum decreased for SPEI-1, SPEI-3 and SPEI-6 indicating that multifractality of the short-term and medium-term dry/wet conditions become weaker; for SPEI-12 multifractality becomes stronger ($W$ increased) in northeast and central regions and weaker ($W$ decreased) for northwest and southern region. For all SPEI series (all accumulation time scales and both sub-periods) the skew parameter $r>0$ indicating that small fluctuations are dominant in multifractality of the process leading to right-skewed $f\left(\alpha \right)$ spectrum. The spatial distribution of skewness parameter $r$ shown in Figure 5 reveals similar pattern for all accumulation time scales for the first sub-period: lower values in northern and southern than in central region of Serbia, indicating that in central part the dominance of small fluctuation is stronger. In the second period the value of $r$increased in northern and southern part and decreased in central part for SPEI-1, SPEI-3, SPEI-6 while for SPEI-12 the increase of $r$ is observed in southern and south-western part while in other regions the value of $r$ decreased.

For northern part which is principal agricultural region the width of the spectrum decreased while the skewness parameter increased in the second period for SPEI-1, SPEI-3 and SPEI-6 indicating that the short-term and medium-term dry/wet conditions that are related to soil moisture evolved toward weaker multifractality with the increase of the right side of the spectrum (stronger influence of small fluctuations), thus the decrease of the spectrum width is the result of the decrease of the left side of the spectrum and loss of the influence of large fluctuations. This indicates that agricultural dry/wet conditions become more stable and more persistent (parameter ${\alpha }_0$ increased) which can affect agricultural yield. SPEI-12, which is tied to streamflow, reservoirs, and groundwater level response to longer-term precipitation anomalies, presents changes in the second period towards larger width of multifractal spectrum and lower skewness parameter in eastern and central regions, indicating that the width of multifractal spectrum in these regions increased because of the stronger influence of large fluctuations. In the southern region the width of multifractal spectrum became smaller and skewness parameter became larger, indicating that the influence of small fluctuations increased. These results are relevant for the electricity production as most of Serbian hydroelectric plants and accumulation lakes are located in eastern, central, and southern part of country. These results are summarized in Table 2 which presents descriptive statistics for multifractal parameters for two considered periods. The mean values calculated from all the grid cells show clear patterns for SPEI-1, SPEI-3 and SPEI-6: i) the values of ${\alpha }_0$ are greater than 0.5, increase with time scale, and become greater in the second period; ii) the width of the spectrum W also increases with the time scale and become smaller in the second period; iii) the skew parameter $r>0,$decreases with the time scale and increases in the second period. For SPEI-12, ${\alpha }_0>1$ and increases in the second period, W also increases in the second period while $r>0$ and decreases in the second period. These results indicate that considering the entire Serbian territory, in the second period short term and medium-term dry/wet conditions evolve towards the regime with stronger persistence, weaker multifractality and stronger influence of small fluctuations, while long term dry/wet conditions shift towards more random regime, stronger multifractality and reduced dominance of small fluctuations.

Considering that for all analyzed series the $f\left(\alpha \right)$ spectrum is right skewed, the complexity of SPEI series in two sub-periods can be addressed by comparing the degree of multifractality (spectrum width $W$). For SPEI-1, SPEI-3 and SPEI-6 the width of the $f\left(\alpha \right)$ spectrum decreased in the second sub-period in most of the country, indicating the decrease of complexity in temporal dynamics of short and medium dry/wet conditions. The complexity of SPEI-12 series increased ($W$ increased) in the northern and decreased ($W$ decreased) in the southern part of Sebia. The joint information of changes in ${\alpha }_0$ and $W$ indicates that SPEI-1, SPEI-3 and SPEI-6 time series become more persistent and less complex in the second sub-period and thus easier for development of predictive models. SPEI-12 series become less anti-persistent (approaching the random regime) and more complex in northern part of the country, indicating that long term dry/wet conditions become more difficult for modeling and prediction. Recently we investigated the changes in multifractality of air temperature across Serbia for the same sub-periods and using the same gridded dataset [63], and we found that in the second sub-period the persistence of temperature series increased across the whole country, the width of the spectrum increased in the northern and decreased in the central and southern parts of country, and that in both periods the multifractality is dominated by small fluctuations (right skewed $f\left(\alpha \right)$ spectrum), indicating the increase of complexity in the northern part of the country. The SPEI time series also showed increase in persistence in the second period for short term and medium-term dry/wet conditions (SPEI-1, SPEI-3 and SPEI-6), while long term dry/wet conditions described by SPEI-12 showed similar changes in the degree of multifractality as temperature series: an increase in the multifractality degree (increase of $W$) in the northern part, and decrease in multifractality degree (decrease of $W$) in the southern part of country. Ionita and Nagavciuc [16] found that negative trends in drought (based on SPEI-12) in the Central Europe (CEU) and the Mediterranean region (MED) where Serbia is located are due to an increase in the potential evapotranspiration and mean air temperature. Our results support these findings in the sense that changes in multifractal properties from the first to the second sub-period of SPEI temporal series coincide with changes of multifractal properties of temperature series as reported in [63]. The multifractal properties of drought indices were investigated only recently. All results obtained for SPI [42,44,45] and SPEI [46] are qualitatively in agreement with what was observed for Serbia. Temporal series of drought indices show multifractal properties with persistent long-term correlations and the dominance of small fluctuation. As accumulation temporal scale increases, the drought indices become more persistent and display stronger multifractality.

Adarsh et al. [42] used MFDFA to analyze SPI time series (at 3-,6-, and 12-months aggregation timescales) from 30 meteorological subdivisions of India, for the period 1871-2016. They found that SPI time series exhibit multifractal properties that vary spatially throughout the studied area, with stronger persistence and higher multifractality degree for larger time scales, which agrees with our results. By comparing two periods separated by the Pacific climatic shift in 1976/1977, they found that in the second period SPI series become more persistent with stronger multifractality for about 70% of subdivision for SPI-3 and SPI-6 and 80% of subdivision for SPI-12. Our results also indicated that multifractal properties of wet/dry conditions could be affected by climate change. In the recent climate period (1991-2020) short- and medium-term conditions (SPEI-, SPEI -3 and SPEI -6) become more persistent and with stronger multifractality, while long- term conditions (SPEI-12) become more random and with weaker multifractality. Climate change affects both rainfall and temperature which are included in SPEI definition, while SPI includes only monthly rainfall. Thus, the comparison of our results with other studies will be clearer when new studies addressing SPEI index changes with similar climate sub periods are published.

4. Conclusions

In this work we investigated the changes in multifractal properties of time series of Standardized precipitation evapotranspiration index (SPEI) in Serbia between two standard climate periods (1961-1990 and 1991-2020). The SPEI is calculated for 1, 3, 6, and 12-month accumulation time scales in order to assess agricultural and hydrological dry/wet conditions. By applying the method Multifractal detrended fluctuation analysis (MFDFA) we compared multifractal spectrum parameters for different accumulation time scales among the two considered periods. These parameters describe specific properties of SPEI time series: long term correlations (position of the spectrum maximum ${\alpha }_0$), degree of multifractality (the spectrum width W) and the contribution of small and large fluctuations (the spectrum asymmetry $r$). We found that, for both periods and for all accumulation time scales SPEI series display multifractal properties with dominance of small fluctuations. Short and medium dry/wet conditions (SPEI-1, SPEI-3 and SPEI-6) showed persistent correlations that become stronger with increase of accumulation time scale, while long term dry/wet conditions (SPEI-12) are antipersistent. In the second period, in the northern agricultural region the multifractal properties of SPEI-1, SPEI-3 and SPEI-6 changed towards stronger persistency, weaker multifractality and stronger influence of small fluctuations, indicating that these series that describe short- and medium-term soil moisture become less complex and easier for development of predictive models, which is relevant for agricultural planning and production. SPEI-12 series on the other hand become less antipersistent (demonstrating an approach towards random regime) in the second period, and more complex (displaying stronger multifractality) in the eastern and central regions, indicating that these series that describe long term dry/wet conditions (related to streamflow, reservoir levels, and groundwater levels) become more difficult for modeling and prediction. These changes are relevant for the electricity production as most of Serbian hydroelectric plants and accumulation lakes are located in the eastern and central parts of country. The results of this study should be taken into account for agricultural planning, better planning of use of water resources for renewable energy, and to develop more efficient measures for the protection and mitigation of the impact from water related natural hazards (drought and floods).