A Stochastic Approach for the Analysis of Long Dry Spells with Different Threshold Values in Southern Italy

A non-homogeneous Poisson model was proposed to analyze the sequences of dry spells below prefixed thresholds as an upgrade of a stochastic procedure previously used to describe long periods of no rainfall. Its application concerned the daily precipitation series in a 60-year time span at four rain gauges (Calabria, southern Italy), aiming at testing the different behaviors of the dry spells below prefixed thresholds in two paired periods (1951–1980 and 1981–2010). A simulation analysis performed through a Monte Carlo approach assessed the statistical significance of the variation of the mean values of dry spells observed at an annual scale in the two 30-year periods. The results evidenced that the dry spells durations increased passing from the first 30-year period to the second one for all the thresholds analyzed. For instance, for the Cassano station, an increase of about 10% of the maximum dry spell duration was detected for a threshold of 5 mm. Moreover, the return periods evaluated for fixed long dry spells through the synthetic data of the period 1981–2010 were lower than the corresponding ones evaluated with the data generated for the previous 30-year period. Specifically, the difference between the two 30-year periods in terms of the return period of long dry spells occurrence increased with the growing thresholds. As an example, for the Cosenza rain gauge with a threshold of 1 mm, the return period for a dry spell length of 70 days decreased from 20 years (in 1951–1980) to about 10 years (in 1981–2010), while for a threshold of 5 mm, the return period for the dry spell lengths of 120 days decreases from 70 years to about 20 years. These results show a higher probability of the occurrence of long dry spells in the more recent period than in the past.


Introduction
The Mediterranean Basin is climatically affected by the interaction between mid-latitude and tropical processes, being located in a transition zone between the arid climate of North Africa and the temperate and rainy climate of central Europe. For this reason, it is considered a major hotspot of climate change, subject to strong warming and drying, with increasing consequences on drought occurrence and severity [1]. In fact, in recent decades, drought has been identified as a persistent problem in the Mediterranean Basin, inflicting serious damage, particularly in the agricultural sector. For example, many industries suffered from drought as a consequence of the severe events which affected some regions of the Mediterranean area in 2011, 2012, 2015 and 2017, including the production of fodder for livestock, vegetables, cereals, fruit growing, viticulture, beekeeping products and other livestock sectors [2,3]. Within this context, drought monitoring is becoming a critical factor in the mitigation of possible drought consequences and numerous drought analyses with different methodologies have been recently performed in the Mediterranean Basin [4][5][6][7], especially in Italy [8][9][10][11][12][13][14][15]. In particular, in order to monitor and assess drought, many indices have been proposed in literature (e.g., [16]), even though dry spells are also commonly applied with the aim to identify spatial patterns of drought risk [17,18]. Namely, a dry spell is defined as the number of consecutive days with a daily precipitation amount below a certain threshold, such as 0.1, 1.0, 5, 10, mm etc., preceded and followed by at least one day with rainfall exceeding the threshold [19,20].
The analysis of the probability distribution of dry spell characteristics (occurrence and duration) provides useful information for agricultural and environmental planning and management of drought-prone areas and it has always been a concern for researchers [21][22][23][24]. Several studies have been conducted worldwide to explore the variability of wet and dry spell characteristics. For example, Li et al. [25] provided a detailed framework to investigate the variability of wet and dry spells in Singapore in the period 1980-2013. Singh and Ranade [26] investigated the variation of wet and dry spell characteristics across India using gridded daily rainfall for the period 1951-2007. Daily precipitation records from 517 Russian stations  have been used to examine the relationships between continuous dry and wet day duration and surface air temperature for all four seasons [27]. Hussain and Lee [28] examined the number of wet days in Pakistan in the period 1950-2010 using precipitation data from 32 stations. In southern Africa, dry spells have been analyzed using gridded data by Usman and Reason [29]. Concerning the Mediterranean Basin, Raymond et al. [30] analyzed extreme dry spells during the wet season. Schmidli and Frei [31] investigated the trends of heavy precipitation and wet and dry spells in Switzerland during the 20 th century at a seasonal time scale based on daily precipitation data at 104 gauge stations. Pérez-Sánchez and Senent-Aparicio [32] analyzed meteorological droughts and dry spells in a semiarid region (Segura Basin, SE Spain). After an overview of the statistical distributions mostly used to model dry (wet) spell durations, Zolina et al. [33] analyzed the changes in the duration of wet and dry periods over Europe.
In the analysis of dry spells, one of the most common problems is the limited number of events in the historical series. In order to overcome such a problem, several authors derived analytically the probabilistic behavior of dry periods' characteristics by assuming a given stochastic structure. Within this context, the stochastic models can be subdivided into two different categories: driven data models, which reproduce the main characteristics of the available data series, and physically based models, which schematize the generating mechanism of atmospheric precipitation. The first category mainly refers to the Poisson model, while the second can be represented by Markov processes or alternating renewal processes. For example, Gupta and Duckstein [34], assuming the Poisson distribution for the number of dry events in a given time interval, analyzed the distribution of the maximum run length. Other authors used Markov processes to describe the persistence and the patterns of dry/wet periods [35][36][37][38][39][40]. Cancelliere and Salas [41] assumed a periodic Markov chain to estimate the occurrence probabilities of a dry event with a given length. Cindrić et al. [42] fitted two stochastic models to empirical dry spell distributions: the first-order Markov chain and the discrete autoregressive moving average model. Ochola and Kerkides [43] embedded first order Markov chains into a computer model to determine the critical wet and dry spells in the Kano Plains of Kenya. Sirangelo et al. [44,45] proposed a non-homogeneous Poisson model for the analysis of dry spells.
In this paper, following the data driven approach, a stochastic model for the analysis of long dry spells was proposed. In particular, this model allows the analysis of dry spells for different threshold values and can be considered as an upgrade of the one proposed by Sirangelo et al. [44] for the investigation of dry spells defined as the number of consecutive days without precipitation. The model applied in this study aims to evaluate the different behavior of the dry spells with five different thresholds in two 30-year periods, i.e., 1951-1980 and 1981-2010. The structure of this paper is the following: in Section 2, the stochastic model used to describe the dry spells with specified thresholds is presented as an upgrade of an existing model for dry spells with zero threshold. In Section 3, after a preliminary description of the study area, the model is applied to four series of daily rainfall characterized by a high quality of data, aiming at verifying the statistical variation of the properties of the dry spells in two consecutive 30-year periods. Finally, the results in terms of parameter calibration, model validation and occurrence probabilities of long dry spells are presented and discussed.

Materials and Methods
The probabilistic model employed to describe the properties of both the arrival and the duration processes of the dry spell for an assigned threshold value h s , is a generalization of the procedure introduced in a previous work [44]. Thus, the comprehensive framework of the model here proposed is only briefly described.
With reference to the sketch of Figure 1, by subdividing the temporal axes in a series of intervals of constant width ∆t > 0 (equal to one day), the number of consecutive intervals in which no rainfall events higher than the threshold value h s occur can be denoted as K , a discrete random variable (r.v.) with support {0, 1, 2, . . .}. Its probability mass function, p K (k ), indicates the probability that at the general instant i ∆t, a series of k intervals (with a total duration k ∆t) begins, before a total rainfall event h higher than h s stops the sequence. following: in Section 2, the stochastic model used to describe the dry spells with specified thresholds is presented as an upgrade of an existing model for dry spells with zero threshold. In Section 3, after a preliminary description of the study area, the model is applied to four series of daily rainfall characterized by a high quality of data, aiming at verifying the statistical variation of the properties of the dry spells in two consecutive 30-year periods. Finally, the results in terms of parameter calibration, model validation and occurrence probabilities of long dry spells are presented and discussed.

Materials and Methods
The probabilistic model employed to describe the properties of both the arrival and the duration processes of the dry spell for an assigned threshold value s h , is a generalization of the procedure introduced in a previous work [44]. Thus, the comprehensive framework of the model here proposed is only briefly described.
With reference to the sketch of Figure 1, by subdividing the temporal axes in a series of intervals of constant width 0 > Δt (equal to one day), the number of consecutive intervals in which no rainfall events higher than the threshold value s h occur can be denoted as K′ , a discrete random variable . Its probability mass function, , indicates the probability that at the general instant , a series of k ′ intervals (with a total duration k ′ Δt) begins, before a total rainfall event h higher than s h stops the sequence.
where ( ) t λ is the intensity of the non-homogeneous Poissonian process describing the arrival of the rainfall events [44]. The function ( ) The probability mass function of the r.v.
′ ≥ , which describes the dry spells for assigned threshold value s h with duration greater than or equal to * k ′ , is To this aim, denoting with F H h; θ(t) , the cumulative distribution function of the total amount of each rainfall event, with parameters θ(t) depending on time, we can introduce the Poissonian reduced intensity: where λ(t) is the intensity of the non-homogeneous Poissonian process describing the arrival of the rainfall events [44]. The function for temporally homogeneous conditions, i.e., λ(t) = λ and F H h s ; θ(t) = F H , reproduces the known related to the maximum of the marked homogeneous Poissonian processes [46]. By using this notation, the probability mass function of the r.v. K , that is, p K (k ), can be represented as The probability mass function of the r.v. K |K ≥ k * , which describes the dry spells for assigned threshold value h s with duration greater than or equal to k * , is for which it is possible to express the expected value, µ K |K ≥k * , and the variance, σ 2 K |K ≥k * , as It is easy to verify that, by posing h s = 0, all the expressions reproduce those reported in Sirangelo et al. [44]. Moreover, so as to keep the full compatibility of the model here proposed with the non-homogeneous Poissonian process presented for the dry spells with zero threshold, the function λ * (t) is considered periodical, with period equal to the average duration of the solar year (D = 365.25 days), expressed through the truncated Fourier series: where n h is the number of harmonics with coefficients a * 0 , a * j and b * j (j = 0,1,2, . . . ,n h ). The estimation of these coefficients can be performed by means of the maximum likelihood method and using a numerical approach.
For the detection of the number of harmonics, the non-homogeneous Poisson process was first transformed into a homogeneous one with unitary intensity through a transformation in the time domain [47]. According to the parsimony principle, the number of harmonics to be chosen is the minimum one for which the homogeneity hypothesis of the arrival Poisson process in the modified time domain cannot be rejected. Finally, with the aim to verify the statistical significance of the difference between the dry spell's characteristics in the two sub-periods, a Monte Carlo procedure was applied. Specifically, synthetic data were generated for the first 30-year period in order to compare the statistics of the simulated world with the ones sampled for 1981-2010. Then, by applying again the Monte Carlo technique for both the sub-periods, the return periods of dry spell for specified lengths were assessed and compared for the different threshold values.

Study Area and Data Base
Calabria is the second southernmost region of Italy after Sicily, being located at the toe of the Italian peninsula. It is characterized by an elongated shape (about 248 km from north to south and a width ranging between 31 and 111 km) and covers a surface of 15,080 km 2 with a coastline of 738 km which reaches the Tyrrhenian Sea on the west side and the Ionian Sea on the southeast side of the region (Figure 2). Due to its position within the Mediterranean Sea and to its orography (42% of the regional area is over 500 m a.s.l. high) the climate of the region is typically Mediterranean but presents a high spatial variability of its climatic features and of extreme hydrological phenomena such as flood and drought [48]. In particular, rainfall deficits are not unusual in the region and thus, drought analysis is particularly important in Calabria where agriculture is one of the largest sectors of the tradable economy. In fact, several drought events hit the region in the 1980's [49,50] and many industries suffered from drought, including the production of fodder for livestock, vegetables, cereals, fruit growing, viticulture, beekeeping products, and other livestock sectors. In fact, the observation period for these rain gauges spans from 1916 to 2010 and presents the majority of the gaps in the series during the Second World War period. This daily database is almost complete after the year 1951, thus, only data from 1951 to 2010 were selected for the analysis, which was carried out with reference to two different 30-year periods (1951-1980 and 1981-2010) in order to investigate the temporal change in dry spells behavior. The main features of the four rain gauges are shown in Table 1.

Intensity Function Modelled through Truncated Fourier Series
The estimation of the model parameters was performed through the maximum likelihood method. For both the 30-year periods and the different threshold values, the number of harmonics for which the Poisson process can be considered as a homogeneous one was equal to 1 for all the rain gauges. The maximum likelihood estimates of the Fourier coefficients , and for each In this study, high-quality daily precipitation data were selected from the database of the Multi-Risk Functional Centre of the Regional Agency for Environmental Protection of the Calabria Region. In particular, four homogeneous and almost complete rainfall series were collected: Cosenza, Cassano allo Jonio, San Pietro in Guarano and Campotenese rain gauges ( Figure 2).
In fact, the observation period for these rain gauges spans from 1916 to 2010 and presents the majority of the gaps in the series during the Second World War period. This daily database is almost complete after the year 1951, thus, only data from 1951 to 2010 were selected for the analysis, which was carried out with reference to two different 30-year periods (1951-1980 and 1981-2010) in order to investigate the temporal change in dry spells behavior. The main features of the four rain gauges are shown in Table 1.

Intensity Function Modelled through Truncated Fourier Series
The estimation of the model parameters was performed through the maximum likelihood method. For both the 30-year periods and the different threshold values, the number of harmonics for which the Poisson process can be considered as a homogeneous one was equal to 1 for all the rain gauges.
The maximum likelihood estimates of the Fourier coefficients a 0 , a 1 and b 1 for each different combination of threshold, period and rain gauge are shown in Table 2. In particular, five threshold values were selected (h 1 = 1 mm, . . . , h 5 = 20 mm). By means of these values evaluated from the data series, the expected value µ K |K ≥k * of the r.v. K |K ≥ k * for the different days of the year can be estimated for both the periods and for various thresholds.
As an example referred to the Cosenza and Cassano rain gauges, Figure 3 shows, for the periods 1951-1980 and 1981-2010, the behaviors of the observed values of the r.v. K |K ≥ k * calculated as 30-year means for every day of the year, m K |K ≥k * , and the corresponding expected values µ K |K ≥k * for two different thresholds (2 mm and 5 mm, respectively). The values are represented in a logarithmic scale with the aim to better graphically show the differences between the observed and estimated values. As Figure 3 shows, the maximum values of the expected values were among the 170th and 190th Julian day (19 June and 10 July); for the Cassano station, the maximum values of the means of observed values for both the sub-periods were anticipated (among the 120 th and 150 th day, that are days in the month of May). Figure 3 also presents evidence that there was a low increase of the maximum values of the expected values for both the stations, passing from the first sub-period to the second one. As an example, for the Cassano station, fixing a 5-mm threshold, the maximum value for the 1951-1980 sub-period was equal to 20.7 days (about 1.31 on the y-axis) increases to 22.3 days (about 1.35 on the y-axis) for the successive sub-period.         The differences between these values for a prefixed threshold slightly decreased for increasing thresholds. In fact, the decreasing behavior passing from the first 30-year period to the second one is noticeable and quite similar for all the rain gauge. Specifically, also considering the values of a * 0 previously obtained for the zero threshold [44], this parameter, averaged on all the considered threshold values, shows a decrease of about 11% for Cosenza and 13% for the other three rain gauges (Cassano allo Ionio, San Pietro in Guarano and Campotenese). The highest decreasing rates occurred for the maximum threshold here considered (20 mm), and this happened in particular for the Cosenza and the Cassano allo Ionio rain gauges (rates of about 18% and 20%, respectively).
In order to evaluate the performance of the developed theoretical model relative to the corresponding empirical results, a comparison between simulated and observed values of the maximum dry spell length was performed for the different thresholds and for each 30-year period. The results obtained for the Cosenza rain gauge ( Figure 5) clearly show that the observed values fell within the range obtained through the Monte Carlo procedure, thus evidencing a good performance of the stochastic model in reproducing the sample data, also in the case of extreme values.
Water 2019, 11, 2026 8 of 13 The differences between these values for a prefixed threshold slightly decreased for increasing thresholds. In fact, the decreasing behavior passing from the first 30-year period to the second one is noticeable and quite similar for all the rain gauge. Specifically, also considering the values of * previously obtained for the zero threshold [44], this parameter, averaged on all the considered threshold values, shows a decrease of about 11% for Cosenza and 13% for the other three rain gauges (Cassano allo Ionio, San Pietro in Guarano and Campotenese). The highest decreasing rates occurred for the maximum threshold here considered (20 mm), and this happened in particular for the Cosenza and the Cassano allo Ionio rain gauges (rates of about 18% and 20%, respectively).
In order to evaluate the performance of the developed theoretical model relative to the corresponding empirical results, a comparison between simulated and observed values of the maximum dry spell length was performed for the different thresholds and for each 30-year period. The results obtained for the Cosenza rain gauge ( Figure 5) clearly show that the observed values fell within the range obtained through the Monte Carlo procedure, thus evidencing a good performance of the stochastic model in reproducing the sample data, also in the case of extreme values.

Test on the Statistical Significance of the Variation of Dry Spells between the Two Periods
Aiming at verifying whether the variation of the properties of dry spells with assigned thresholds between the two 30-year periods 1951-1980 and 1981-2010 is statistically significant, as assessed for the zero threshold [44], the following equation is introduced:

Test on the Statistical Significance of the Variation of Dry Spells between the Two Periods
Aiming at verifying whether the variation of the properties of dry spells with assigned thresholds between the two 30-year periods 1951-1980 and 1981-2010 is statistically significant, as assessed for the zero threshold [44], the following equation is introduced: which expresses the expected value of the means µ K |K ≥k * with thresholds varying in the range 0 ≤ h s ≤ h s,max for the 30-year period 1951-1980. Through a Monte Carlo simulation, the distribution function F ∆ (δ) is evaluated and its percentile δ 1−α , with 1 − α = 0.95, is compared to the sample value d estimated for the 30-year period 1981-2010: in which h 0 = 0, h 1 = 1 mm, . . . , h 5 = 20 mm, and ω s are the weighting coefficients, which allow to obtain the best approximation of the integral.
The hypothesis H 0 , for which the variation among the statistical properties of the dry spells for assigned thresholds between the two 30-year periods is no statistically significant, should be rejected with a significance level α = 0.05 if d > δ 0.95 . The results of the test, synthetized for all the four rain gauges here employed in Table 3, confirm that the variation already assessed as significant for h s = 0 is also the same in the thresholds range 0 ≤ h s ≤ 20 mm, resulting in every case d > δ 0.95 , that is F ∆ (d) > 0.95. Therefore, it can be assessed that for all the rain gauges, the observed differences between 1951-1980 and 1981-2010 are statistically significant.

Occurrence Probability of Long Dry Spells
By means of Equation (4) and using the synthetic series generated through Monte Carlo techniques, it is possible to quantify, for each day of the year, the occurrence probability of long dry spells with durations greater than the fixed threshold values k * . Moreover, the difference of occurrence of long dry spells can be highlighted by the evaluation of the return periods (T). Being k max the maximum value of k over a time interval of length 1 year, the return period T of k max , in years, can be estimated by means of where F K max (k max ) is the cumulative distribution function of k max [51]. Thus, through Equation  The differences for the two periods are clearly shown. As an example, with a threshold of 1 mm and for a dry spell length of 70 days, the return period decreased from 20 years (in 1951-1980) to about 10 years (in 1981-2010). Similarly, the return periods decreased from about 70 years to about 20 years for dry spell lengths of 120 days and a threshold of 5 mm. By increasing the threshold values, for even higher values of dry spell length, the return periods evaluated with synthetic data of 1981-2010 were much less than half the corresponding return periods evaluated with the data generated for 1951-1980.

Discussion and Conclusions
A reliable analysis of the theoretical distributions of dry spells is important for water resources management. In particular, in agriculture, prolonged dry spells during cropping seasons can impact the yield production. For these reasons, knowledge of the length of dry spells and their occurrence probability can aid the planning for supplementary risk aversion strategies. With the aim to perform an analysis of dry spell durations with observed data, high quality records are required in order to overcome completeness problems, which can lead to underestimating the exact number of dry days. Since the early 20th century, several probability models for the distribution of dry spell have been developed based on single site or multisite data (e.g., [52]). In particular, in order to overcome the problem of the limited number of events in the historical series, the probabilistic behavior of dry periods characteristics has been often derived analytically, assuming a given stochastic structure of The differences for the two periods are clearly shown. As an example, with a threshold of 1 mm and for a dry spell length of 70 days, the return period decreased from 20 years (in 1951-1980) to about 10 years (in 1981-2010). Similarly, the return periods decreased from about 70 years to about 20 years for dry spell lengths of 120 days and a threshold of 5 mm. By increasing the threshold values, for even higher values of dry spell length, the return periods evaluated with synthetic data of 1981-2010 were much less than half the corresponding return periods evaluated with the data generated for 1951-1980.

Discussion and Conclusions
A reliable analysis of the theoretical distributions of dry spells is important for water resources management. In particular, in agriculture, prolonged dry spells during cropping seasons can impact the yield production. For these reasons, knowledge of the length of dry spells and their occurrence probability can aid the planning for supplementary risk aversion strategies. With the aim to perform an analysis of dry spell durations with observed data, high quality records are required in order to overcome completeness problems, which can lead to underestimating the exact number of dry days. Since the early 20th century, several probability models for the distribution of dry spell have been developed based on single site or multisite data (e.g., [52]). In particular, in order to overcome the problem of the limited number of events in the historical series, the probabilistic behavior of dry periods characteristics has been often derived analytically, assuming a given stochastic structure of the underlying hydrological series. With this goal, in this paper, a non-homogeneous Poisson model was applied for the analysis of the duration of long dry spells with different threshold values. Specifically, the analysis focused on the different behaviors of the dry spell durations observed in two different 30-year periods (1951-1980 and 1981-2010). Within this aim, a Monte Carlo procedure was developed to generate two 30-year time spans, each reproducing the main characteristics of the two sub-periods. First, for each threshold value, the mean expected values of the dry spells evaluated in the two sub-periods were compared, and the statistical significance was been verified. As a result, passing from the 1951-1980 to the 1981-2010 sub-period, an increase of the maximum duration was detected for all the stations, thus evidencing a longer duration of recent dry spells than in the past. Moreover, the analysis of the behavior of the Fourier parameters evidenced that the differences between the two sub-periods tended to decrease with the increase of the threshold values. Finally, in order to better point out the different behaviors in the two sub-periods, the return periods of high values of dry spell length were estimated. The results show that for all the thresholds, the return periods evaluated for fixed long dry spells with the synthetic data of 1981-2010 d lower results than the corresponding ones evaluated with the data generated for 1951-1980. Specifically, with the increase of both the threshold values and the dry spell lengths, the return periods of the latest sub-period were much less than half the corresponding return periods evaluated for 1951-1980.
Overall, the stochastic approach proposed in this paper proved to be highly reliable for explaining the seasonal variability of dry spell series in southern Italy. In fact, although not straightforward to apply, the model was powerful in the estimation of parameters. Moreover, as opposed to other data-driven models, which reproduce only some characteristics of the available data series, this model allows the reproduction of all the main features of the dry spell evaluated for each day of the year. This investigation could be applied for assessing change in probabilities of long dry spells, thus allowing a better drought prediction for agricultural and environmental issues.

Conflicts of Interest:
The authors declare no conflict of interest.