# Fractal Analysis of the Long-Term Memory in Precipitation over Bénin (West Africa)

Authors: Medard Noukpo Agbazo; Gabin Koto N’Gobi; Eric Alamou; Basile Kounouhewa; Abel Afouda
Journal: Advances in Meteorology
 (2019)

---

## Abstract

This study analyzed the long-term memory (LTM) in precipitation over Bénin synoptic stations from 1951 to 2010 using the detrended fluctuation analysis (DFA) method. Results reveal the existence of positive long-term memory characteristic in rainfall field. DFA exponent values are different regarding the concerned synoptic stations, reflecting the effect of geographical position and climate on the LTM. These values were related to the type of climate. The best DFA1-4 method depends on the geographical position of the studied station. However, DFA2 is generally the best in terms of spatial average from DFA1 to DFA4. In Bénin synoptic stations, except the Parakou station, the long-term temporal correlations are systematically the source of multifractality in rainfall. Except Natitingou, the strength of long-term memory characteristic decreases each twenty years in the study period. Considering the fractal approach, our results show that the subperiod 1991–2010 is not really a transition period as shown before. Thus, the drought is prolonging until 2010. So, fractal theory reveals more Bénin climatic characteristics.


---

## Body

## 1. Introduction

According to several reports, such as those of the United Nations Framework Convention on Climate Change [1], the International Panel of Climate Change (IPCC) [2], several documents produced by the World Bank [3], and the Climate Change and Urban Vulnerability in Africa (CLUVA) project [4], Africa is considered as a continent particularly vulnerable to climate change. But, the real impact, especially on a local scale, is still poorly understood. Changes in climate have significant implications for societies, future generations, the economy, ecosystems, and agriculture [5]. Moreover, climate change has the self-memory characteristic [6]. This self-memory characteristic is widely known as long-range correlation (also called the long-term memory, LTM). So, a better understanding of long-term climate memory by improved methods is very important for predicting climate change in Africa. This would significantly reduce the vulnerability and prepare populations for climate change adaptation. Also, this understanding can help to understand more, a complex dynamic behavior in climate effects. Climate change is marked with variation in precipitations [7]. Long-term memory exists in precipitation records [8–11]. Precipitation memory characteristics have been widely studied in various regions of the world. For instance, Kantelhardt et al. [9] have studied, by the detrended fluctuation analysis (DFA) method, the LTM characteristics of precipitation records in Europe, Asia, and America. In Africa, Efstathiou and Varotsos [12] have applied the DFA method to explore the intrinsic properties of Sahel precipitation anomalies. They found that the fluctuations of the precipitation anomalies in small time intervals are positively correlated to those in longer time intervals in a power law fashion [12]. According to Taqqu et al. [13] and Kantelhardt et al. [9], the DFA method permits to estimate the Hurst exponent H. It is important to remind the distinction of the classic Hurst and the DFA exponent as announced in Taqqu et al. [13], Movahed et al. [14, 15], and more recently as clarified by Hall [16] and Varotsos and Efstathiou [17]. According to Hall [16] and Varotsos and Efstathiou [17], the relationship between the DFA exponent and the classic Hurst exponent H can be expressed as (i) for fractional Gaussian noise (fGn), H is equal to DFA exponent and (ii) for fractional Brownian motions (fBm), H is equal to DFA exponent 1. Moreover, for fGn processes (which is a stationary process), 0<DFA exponent<1 and for fBm processes (which is nonstationary process), 1<DFA exponent<2 [16]. The estimation of such exponent from a given dataset is an effective way to detect the LTM characteristics in it. However, Maraun et al. [18] and Varotsos and Efstathiou [17] suggested that two criteria are needed to ensure the existence of the long-range correlation: (i) the rejection of the simple exponential decay of the autocorrelation function and (ii) the establishment of the power law scaling. Although some claim that the DFA exponent or Hurst exponent H is not universal for meteorological and climatological data [19, 20], precipitation have characteristic values of DFA exponent or Hurst exponent H [21, 22]. One can be sure of the universality of the correlations in climatological time series but its DFA exponents can be related to local patterns.Bénin is one of the West Africa countries affected by the effects of climate change [3]. Long-term memory of precipitation is not yet studied in fractal framework in Bénin. In the literature, it has been shown that long-range correlations in weather and climate systems depend on geographical locations [11]. It is important to find out what is the long-term memory of precipitation at Bénin meteorological synoptic stations located in two main climatic areas (subequatorial climate and Sudanian climate).The purpose of this paper is to study, for the first time in Bénin, the long-term memory of precipitation. The paper is organized as follows: in Section2, the method of DFA and the dataset are described. Analysis of results and their interpretations are made in Section 3. Discussion of results is made in Section 4. Finally, we conclude the paper with a summary and outlook for further research in Section 5.


## 2. Materials and Methods

### 2.1. Materials

#### 2.1.1. Site Description

The study covers all the synoptic stations of Bénin, the geographical positions of which are shown in (Figure1). Bénin is characterized from the South to North by three climatic zones in which stations are located [23]: Cotonou, Bohicon, and Savè are located in the subequatorial region where March is the hottest month (∼26°C), while August is the coldest month (∼24°C). The daily and annual thermal amplitudes are, respectively, ∼10°C and ∼5°C. The relative humidity ranges between 70% and 95% because of the proximity to the Atlantic Ocean. The subequatorial climate has four seasons: a long rainy season (April to July) followed by a short dry season (August to September) and a short rainy season (October to November) followed by long dry season (December to March) in the year. However, the stations of Parakou, Kandi, and Natitingou are located in Sudanian region in the northern part of the country. The daily mean of air temperatures in Natitingou, Parakou, and Kandi are, respectively, ∼25°C, ∼27°C, and ∼35°C. Parakou and Savè are located in the transition area between the two kinds of climatic zone.Figure 1
Study sites location: (a) Bénin location in West Africa; (b) synoptic station locations in Bénin.


#### 2.1.2. Data Records

Data were provided by the Agency for the Aerial Navigation’s Security in Africa and in Madagascar (ASECNA). Daily rainfall data from Bénin synoptic stations are used during a period from 1951 to 2010. Rainfall variation during the study period is presented in Figure2 for each synoptic station.Figure 2
Temporal variations of rainfall time series collected from 1951 to 2010: (a) Cotonou station, (b) Bohicon station, (c) Savè station, (d) Parakou station, (e) Natitingou station, and (f) Kandi station.
(a)
(b)
(c)
(d)
(e)
(f)


### 2.2. Methods

#### 2.2.1. Description of the MFDFA Method

The multifractal detrended fluctuation analysis (MFDFA) method is a generalization of the standard DFA method. It is generally preferred and widely used in the literature to detect the fractal scaling properties of themth-order moments and long-term memory in nonstationary and nonlinear time series as, meteorological and climatological series [6, 24, 25]. The modified MFDFA procedure consists of a sequence of steps and detailed information about computation can be found in [26]. The steps are essentially identical to the conventional DFA procedure. It can be briefly described as follow [6].Letxt∣t=1,2,…,N be an original time series of N equidistant measurements to which the procedure of the DFA method is applied. First, a news series named “profile” is determined as follows:(1)Yk=∑t=1kxt−x,where x is the mean value of xt and k=1,2,…,N.The profile is then divided intoN = int (N/s) equal-sized nonoverlapping windows with a length s. Since N is not the integral multiple of s in most cases, there might be a short part at the end of the profile that remains uncovered. To take full account of the series, the same procedure can be repeated starting from the end of the series. Hence, we can obtain 2Ns segments altogether. We then calculate the variance of each window as(2)F2ν,s=1s∑i=1sYv−1s+i−yvmi2,forv=1,2,…,Ns,1s∑i=1sYN−v−Nss+i−yvmi2,forv=Ns+1,…,2Ns.It should be noted that linear (i.e., polynomial orderm=1), quadratic (i.e., m=2), cubic (i.e., m=3), or higher-order polynomials yνmi can be used to fit the local trend, and MFDFA is noted as MFDFA1, MFDFA2, MFDFA3, MFDFA4, etc. In MFDFAm, possible m order trends are eliminated in the profile. By averaging over all windows, we obtain the fluctuation as(3)Fqs=12Ns∑ν=12NsF2ν,sq/21/q,forq≠0,exp14Ns∑ν=12NslnF2ν,s,forq=0.If the time series follows the power law, then we can obtain the scaling function:(4)Fqs∝shq,where hq is the generalized Hurst scaling function [9, 27]. In particular, h2 is the DFA exponent. It is very important to keep in mind the distinction between the classic Hurst and the DFA exponent as indicated earlier in the introduction. According to Hall [16] and Varotsos and Efstathiou [17], h2=H for fractional Gaussian noise (fGn) and h2=H+1 for fractional Brownian motions (fBm). For q=2, we have the standard DFA analysis and DFA1, DFA2, DFA3, and DFA4 are corresponding to MFDFA1, MFDFA2, MFDFA3, and MFDFA4 (i.e., for an m-order DFA process, an m-order polynomial is used to calculate the local trend in each segment). In this case, the scaling exponent h2 provides information about the average fluctuation of the series. The series can be categorized into one of the following three types depending on the h2 values. These are (i) 0 < h(2) < 0.5 for an antipersistent type long-range correlated process (negative long-term memory) where large values (compared to the average) are more likely followed by small values and vice versa, (ii) h(2) = 0.5 for an entirely random uncorrelated distribution, and (iii) 0.5 < h(2) < 1 for a persistent and long-range correlated process (positive long-term memory) where large values are more likely to be followed by large values and vice versa. In this study, h(2) is represented by h.


#### 2.2.2. The Best DFA Method

In order to choose the best DFA method, we have calculated the difference between theh values of two successive orders as follows:(5)Δhm,m+1=hDFAm−hDFAm+1,where hDFAm and hDFAm+1 are, respectively, the h values obtained by DFAm and DFAm+1m=1,2and3. The DFAm method minimizing Δhm,m+1 spatially averaged is then used as the best method.


#### 2.2.3. Origins of Multifractality

Generally, there are two types of source of multifractality in time series. One is due to different long-term temporal correlations for small and large fluctuations; the other is due to the fat-tailed distributions of variations [26, 28, 29]. The main methods to find the contributions of the two sources of multifractality are the shuffling procedure and the surrogating procedure, respectively [30]. The shuffling procedure destroys any temporal correlations in the data, but the distributions remain exactly the same. In order to quantify the influence of the broad probability density function, surrogate time series were generated from the original by randomizing their phases in Fourier space, so that the surrogate series are Gaussian. The shuffling procedure consists of generating a random permutation of the array elements of time series. In this paper, we use the common method of amplitude adjusted Fourier transform (AAFT) developed in [31–33] to obtain surrogates data. AAFT can be summarized according to Dong et al. [32] as follows: (1) a discrete Fourier transform of the original series is conducted, (2) the discrete Fourier transform of the data is multiplied by random phases, and (3) an inverse Fourier transform is performed to generate a phase randomized surrogate. For convenience, horigq, hshufq, and hsurrq denote the generalized Hurst scaling function obtained by the best MFDFAm corresponding to the best DFAm for the original, shuffled, and surrogated times series, respectively. In order to avoid a divergence of moments in the fat tails of the fluctuation distribution as mentioned by some authors [34, 35], we restrict the order q within the range −3 ≤ q ≤ 3 as in [34]. To better compare the impact of both types of multifractality in time, these following computations are done:(6)hCORRq=hshufq−horigq,hPDFq=hsurrq−horigq.The long-term temporal correlations have a slightly greater effect on multifractality than the broad probability density function, ifhCORRq>hPDFq; otherwise, hCORRq<hPDFq. However, nonzero values of hCORRqandhPDFq indicate that both of them influence the multifractality.


#### 2.2.4. Assessment to Climate Change Effect by the DFA Method

To determine whether changes in the dynamics of rainfall time series can be assessed with the DFA method, the 60 years of data were divided into three separate datasets, the first is from 1951 to 1970; the second is from 1971 to 1990 and the third is from 1991 to 2010. In each synoptic station, the comparison ofhorig, hshuf, and hsurr of rainfall from these three periods are done. This division was made considering the classification mentioned by some studies in West Africa [36–38]: 1951–1970 is considered as the wet period, 1971–1990 is the drought period, and 1991–2010 is the end of the drought.


#### 2.2.5. Error Analysis of the DFA Exponents (h)

To show the efficiency of the linear regression used to fit the fluctuation function in a log-log representation for determining the DFA exponents (h), Pearson coefficient is computed. To find out whether the straight line fitting to the fluctuation function in a log-log representation is statistically significant, we determine the confidence interval (lower and upper bound of the 95% confidence interval) for Pearson coefficient. The error analysis is done from the obtained values of the mean squared error (MSE).


## 2.1. Materials

### 2.1.1. Site Description

The study covers all the synoptic stations of Bénin, the geographical positions of which are shown in (Figure1). Bénin is characterized from the South to North by three climatic zones in which stations are located [23]: Cotonou, Bohicon, and Savè are located in the subequatorial region where March is the hottest month (∼26°C), while August is the coldest month (∼24°C). The daily and annual thermal amplitudes are, respectively, ∼10°C and ∼5°C. The relative humidity ranges between 70% and 95% because of the proximity to the Atlantic Ocean. The subequatorial climate has four seasons: a long rainy season (April to July) followed by a short dry season (August to September) and a short rainy season (October to November) followed by long dry season (December to March) in the year. However, the stations of Parakou, Kandi, and Natitingou are located in Sudanian region in the northern part of the country. The daily mean of air temperatures in Natitingou, Parakou, and Kandi are, respectively, ∼25°C, ∼27°C, and ∼35°C. Parakou and Savè are located in the transition area between the two kinds of climatic zone.Figure 1
Study sites location: (a) Bénin location in West Africa; (b) synoptic station locations in Bénin.


### 2.1.2. Data Records

Data were provided by the Agency for the Aerial Navigation’s Security in Africa and in Madagascar (ASECNA). Daily rainfall data from Bénin synoptic stations are used during a period from 1951 to 2010. Rainfall variation during the study period is presented in Figure2 for each synoptic station.Figure 2
Temporal variations of rainfall time series collected from 1951 to 2010: (a) Cotonou station, (b) Bohicon station, (c) Savè station, (d) Parakou station, (e) Natitingou station, and (f) Kandi station.
(a)
(b)
(c)
(d)
(e)
(f)


## 2.1.1. Site Description

The study covers all the synoptic stations of Bénin, the geographical positions of which are shown in (Figure1). Bénin is characterized from the South to North by three climatic zones in which stations are located [23]: Cotonou, Bohicon, and Savè are located in the subequatorial region where March is the hottest month (∼26°C), while August is the coldest month (∼24°C). The daily and annual thermal amplitudes are, respectively, ∼10°C and ∼5°C. The relative humidity ranges between 70% and 95% because of the proximity to the Atlantic Ocean. The subequatorial climate has four seasons: a long rainy season (April to July) followed by a short dry season (August to September) and a short rainy season (October to November) followed by long dry season (December to March) in the year. However, the stations of Parakou, Kandi, and Natitingou are located in Sudanian region in the northern part of the country. The daily mean of air temperatures in Natitingou, Parakou, and Kandi are, respectively, ∼25°C, ∼27°C, and ∼35°C. Parakou and Savè are located in the transition area between the two kinds of climatic zone.Figure 1
Study sites location: (a) Bénin location in West Africa; (b) synoptic station locations in Bénin.


## 2.1.2. Data Records

Data were provided by the Agency for the Aerial Navigation’s Security in Africa and in Madagascar (ASECNA). Daily rainfall data from Bénin synoptic stations are used during a period from 1951 to 2010. Rainfall variation during the study period is presented in Figure2 for each synoptic station.Figure 2
Temporal variations of rainfall time series collected from 1951 to 2010: (a) Cotonou station, (b) Bohicon station, (c) Savè station, (d) Parakou station, (e) Natitingou station, and (f) Kandi station.
(a)
(b)
(c)
(d)
(e)
(f)


## 2.2. Methods

### 2.2.1. Description of the MFDFA Method

The multifractal detrended fluctuation analysis (MFDFA) method is a generalization of the standard DFA method. It is generally preferred and widely used in the literature to detect the fractal scaling properties of themth-order moments and long-term memory in nonstationary and nonlinear time series as, meteorological and climatological series [6, 24, 25]. The modified MFDFA procedure consists of a sequence of steps and detailed information about computation can be found in [26]. The steps are essentially identical to the conventional DFA procedure. It can be briefly described as follow [6].Letxt∣t=1,2,…,N be an original time series of N equidistant measurements to which the procedure of the DFA method is applied. First, a news series named “profile” is determined as follows:(1)Yk=∑t=1kxt−x,where x is the mean value of xt and k=1,2,…,N.The profile is then divided intoN = int (N/s) equal-sized nonoverlapping windows with a length s. Since N is not the integral multiple of s in most cases, there might be a short part at the end of the profile that remains uncovered. To take full account of the series, the same procedure can be repeated starting from the end of the series. Hence, we can obtain 2Ns segments altogether. We then calculate the variance of each window as(2)F2ν,s=1s∑i=1sYv−1s+i−yvmi2,forv=1,2,…,Ns,1s∑i=1sYN−v−Nss+i−yvmi2,forv=Ns+1,…,2Ns.It should be noted that linear (i.e., polynomial orderm=1), quadratic (i.e., m=2), cubic (i.e., m=3), or higher-order polynomials yνmi can be used to fit the local trend, and MFDFA is noted as MFDFA1, MFDFA2, MFDFA3, MFDFA4, etc. In MFDFAm, possible m order trends are eliminated in the profile. By averaging over all windows, we obtain the fluctuation as(3)Fqs=12Ns∑ν=12NsF2ν,sq/21/q,forq≠0,exp14Ns∑ν=12NslnF2ν,s,forq=0.If the time series follows the power law, then we can obtain the scaling function:(4)Fqs∝shq,where hq is the generalized Hurst scaling function [9, 27]. In particular, h2 is the DFA exponent. It is very important to keep in mind the distinction between the classic Hurst and the DFA exponent as indicated earlier in the introduction. According to Hall [16] and Varotsos and Efstathiou [17], h2=H for fractional Gaussian noise (fGn) and h2=H+1 for fractional Brownian motions (fBm). For q=2, we have the standard DFA analysis and DFA1, DFA2, DFA3, and DFA4 are corresponding to MFDFA1, MFDFA2, MFDFA3, and MFDFA4 (i.e., for an m-order DFA process, an m-order polynomial is used to calculate the local trend in each segment). In this case, the scaling exponent h2 provides information about the average fluctuation of the series. The series can be categorized into one of the following three types depending on the h2 values. These are (i) 0 < h(2) < 0.5 for an antipersistent type long-range correlated process (negative long-term memory) where large values (compared to the average) are more likely followed by small values and vice versa, (ii) h(2) = 0.5 for an entirely random uncorrelated distribution, and (iii) 0.5 < h(2) < 1 for a persistent and long-range correlated process (positive long-term memory) where large values are more likely to be followed by large values and vice versa. In this study, h(2) is represented by h.


### 2.2.2. The Best DFA Method

In order to choose the best DFA method, we have calculated the difference between theh values of two successive orders as follows:(5)Δhm,m+1=hDFAm−hDFAm+1,where hDFAm and hDFAm+1 are, respectively, the h values obtained by DFAm and DFAm+1m=1,2and3. The DFAm method minimizing Δhm,m+1 spatially averaged is then used as the best method.


### 2.2.3. Origins of Multifractality

Generally, there are two types of source of multifractality in time series. One is due to different long-term temporal correlations for small and large fluctuations; the other is due to the fat-tailed distributions of variations [26, 28, 29]. The main methods to find the contributions of the two sources of multifractality are the shuffling procedure and the surrogating procedure, respectively [30]. The shuffling procedure destroys any temporal correlations in the data, but the distributions remain exactly the same. In order to quantify the influence of the broad probability density function, surrogate time series were generated from the original by randomizing their phases in Fourier space, so that the surrogate series are Gaussian. The shuffling procedure consists of generating a random permutation of the array elements of time series. In this paper, we use the common method of amplitude adjusted Fourier transform (AAFT) developed in [31–33] to obtain surrogates data. AAFT can be summarized according to Dong et al. [32] as follows: (1) a discrete Fourier transform of the original series is conducted, (2) the discrete Fourier transform of the data is multiplied by random phases, and (3) an inverse Fourier transform is performed to generate a phase randomized surrogate. For convenience, horigq, hshufq, and hsurrq denote the generalized Hurst scaling function obtained by the best MFDFAm corresponding to the best DFAm for the original, shuffled, and surrogated times series, respectively. In order to avoid a divergence of moments in the fat tails of the fluctuation distribution as mentioned by some authors [34, 35], we restrict the order q within the range −3 ≤ q ≤ 3 as in [34]. To better compare the impact of both types of multifractality in time, these following computations are done:(6)hCORRq=hshufq−horigq,hPDFq=hsurrq−horigq.The long-term temporal correlations have a slightly greater effect on multifractality than the broad probability density function, ifhCORRq>hPDFq; otherwise, hCORRq<hPDFq. However, nonzero values of hCORRqandhPDFq indicate that both of them influence the multifractality.


### 2.2.4. Assessment to Climate Change Effect by the DFA Method

To determine whether changes in the dynamics of rainfall time series can be assessed with the DFA method, the 60 years of data were divided into three separate datasets, the first is from 1951 to 1970; the second is from 1971 to 1990 and the third is from 1991 to 2010. In each synoptic station, the comparison ofhorig, hshuf, and hsurr of rainfall from these three periods are done. This division was made considering the classification mentioned by some studies in West Africa [36–38]: 1951–1970 is considered as the wet period, 1971–1990 is the drought period, and 1991–2010 is the end of the drought.


### 2.2.5. Error Analysis of the DFA Exponents (h)

To show the efficiency of the linear regression used to fit the fluctuation function in a log-log representation for determining the DFA exponents (h), Pearson coefficient is computed. To find out whether the straight line fitting to the fluctuation function in a log-log representation is statistically significant, we determine the confidence interval (lower and upper bound of the 95% confidence interval) for Pearson coefficient. The error analysis is done from the obtained values of the mean squared error (MSE).


## 2.2.1. Description of the MFDFA Method

The multifractal detrended fluctuation analysis (MFDFA) method is a generalization of the standard DFA method. It is generally preferred and widely used in the literature to detect the fractal scaling properties of themth-order moments and long-term memory in nonstationary and nonlinear time series as, meteorological and climatological series [6, 24, 25]. The modified MFDFA procedure consists of a sequence of steps and detailed information about computation can be found in [26]. The steps are essentially identical to the conventional DFA procedure. It can be briefly described as follow [6].Letxt∣t=1,2,…,N be an original time series of N equidistant measurements to which the procedure of the DFA method is applied. First, a news series named “profile” is determined as follows:(1)Yk=∑t=1kxt−x,where x is the mean value of xt and k=1,2,…,N.The profile is then divided intoN = int (N/s) equal-sized nonoverlapping windows with a length s. Since N is not the integral multiple of s in most cases, there might be a short part at the end of the profile that remains uncovered. To take full account of the series, the same procedure can be repeated starting from the end of the series. Hence, we can obtain 2Ns segments altogether. We then calculate the variance of each window as(2)F2ν,s=1s∑i=1sYv−1s+i−yvmi2,forv=1,2,…,Ns,1s∑i=1sYN−v−Nss+i−yvmi2,forv=Ns+1,…,2Ns.It should be noted that linear (i.e., polynomial orderm=1), quadratic (i.e., m=2), cubic (i.e., m=3), or higher-order polynomials yνmi can be used to fit the local trend, and MFDFA is noted as MFDFA1, MFDFA2, MFDFA3, MFDFA4, etc. In MFDFAm, possible m order trends are eliminated in the profile. By averaging over all windows, we obtain the fluctuation as(3)Fqs=12Ns∑ν=12NsF2ν,sq/21/q,forq≠0,exp14Ns∑ν=12NslnF2ν,s,forq=0.If the time series follows the power law, then we can obtain the scaling function:(4)Fqs∝shq,where hq is the generalized Hurst scaling function [9, 27]. In particular, h2 is the DFA exponent. It is very important to keep in mind the distinction between the classic Hurst and the DFA exponent as indicated earlier in the introduction. According to Hall [16] and Varotsos and Efstathiou [17], h2=H for fractional Gaussian noise (fGn) and h2=H+1 for fractional Brownian motions (fBm). For q=2, we have the standard DFA analysis and DFA1, DFA2, DFA3, and DFA4 are corresponding to MFDFA1, MFDFA2, MFDFA3, and MFDFA4 (i.e., for an m-order DFA process, an m-order polynomial is used to calculate the local trend in each segment). In this case, the scaling exponent h2 provides information about the average fluctuation of the series. The series can be categorized into one of the following three types depending on the h2 values. These are (i) 0 < h(2) < 0.5 for an antipersistent type long-range correlated process (negative long-term memory) where large values (compared to the average) are more likely followed by small values and vice versa, (ii) h(2) = 0.5 for an entirely random uncorrelated distribution, and (iii) 0.5 < h(2) < 1 for a persistent and long-range correlated process (positive long-term memory) where large values are more likely to be followed by large values and vice versa. In this study, h(2) is represented by h.


## 2.2.2. The Best DFA Method

In order to choose the best DFA method, we have calculated the difference between theh values of two successive orders as follows:(5)Δhm,m+1=hDFAm−hDFAm+1,where hDFAm and hDFAm+1 are, respectively, the h values obtained by DFAm and DFAm+1m=1,2and3. The DFAm method minimizing Δhm,m+1 spatially averaged is then used as the best method.


## 2.2.3. Origins of Multifractality

Generally, there are two types of source of multifractality in time series. One is due to different long-term temporal correlations for small and large fluctuations; the other is due to the fat-tailed distributions of variations [26, 28, 29]. The main methods to find the contributions of the two sources of multifractality are the shuffling procedure and the surrogating procedure, respectively [30]. The shuffling procedure destroys any temporal correlations in the data, but the distributions remain exactly the same. In order to quantify the influence of the broad probability density function, surrogate time series were generated from the original by randomizing their phases in Fourier space, so that the surrogate series are Gaussian. The shuffling procedure consists of generating a random permutation of the array elements of time series. In this paper, we use the common method of amplitude adjusted Fourier transform (AAFT) developed in [31–33] to obtain surrogates data. AAFT can be summarized according to Dong et al. [32] as follows: (1) a discrete Fourier transform of the original series is conducted, (2) the discrete Fourier transform of the data is multiplied by random phases, and (3) an inverse Fourier transform is performed to generate a phase randomized surrogate. For convenience, horigq, hshufq, and hsurrq denote the generalized Hurst scaling function obtained by the best MFDFAm corresponding to the best DFAm for the original, shuffled, and surrogated times series, respectively. In order to avoid a divergence of moments in the fat tails of the fluctuation distribution as mentioned by some authors [34, 35], we restrict the order q within the range −3 ≤ q ≤ 3 as in [34]. To better compare the impact of both types of multifractality in time, these following computations are done:(6)hCORRq=hshufq−horigq,hPDFq=hsurrq−horigq.The long-term temporal correlations have a slightly greater effect on multifractality than the broad probability density function, ifhCORRq>hPDFq; otherwise, hCORRq<hPDFq. However, nonzero values of hCORRqandhPDFq indicate that both of them influence the multifractality.


## 2.2.4. Assessment to Climate Change Effect by the DFA Method

To determine whether changes in the dynamics of rainfall time series can be assessed with the DFA method, the 60 years of data were divided into three separate datasets, the first is from 1951 to 1970; the second is from 1971 to 1990 and the third is from 1991 to 2010. In each synoptic station, the comparison ofhorig, hshuf, and hsurr of rainfall from these three periods are done. This division was made considering the classification mentioned by some studies in West Africa [36–38]: 1951–1970 is considered as the wet period, 1971–1990 is the drought period, and 1991–2010 is the end of the drought.


## 2.2.5. Error Analysis of the DFA Exponents (h)

To show the efficiency of the linear regression used to fit the fluctuation function in a log-log representation for determining the DFA exponents (h), Pearson coefficient is computed. To find out whether the straight line fitting to the fluctuation function in a log-log representation is statistically significant, we determine the confidence interval (lower and upper bound of the 95% confidence interval) for Pearson coefficient. The error analysis is done from the obtained values of the mean squared error (MSE).


## 3. Results

### 3.1. Long-Term Memory Characteristic Study

In the present study, we calculate the results for all the six synoptic stations in the range of 8–1024 days. Figure3 shows South-North distribution of DFA exponent (ℎ) obtained by the DFA1-4 method. The main results of this figure are the following: (a) h values depend on the applied DFAm method; (b) except Cotonou, whatever the DFAm method, h values present spatial variation, there is a South-North gradient, (c) over all synoptic stations, h values are systematically between 0.5 and 1, this means that there exists positive long-term memory characteristic in the rainfall field, and (d) among synoptic stations, the highest and lowest values of h are, respectively, obtained at Cotonou and Bohicon stations, indicating the existence of strong and weak long-term memory characteristics in Cotonou and Bohicon.Figure 3
South-North distribution of DFA exponents (ℎ) obtained by theDFA1-4 method. Cotonou (Co), Bohicon (Bo), Savè (Sa), Parakou (Pa), Natitingou (Na), and Kandi (Ka). For a DFAm, an m-order polynomial is used to calculate the local trend in each segment.


### 3.2. The Best DFA Method

Figure4 presents the South-North distribution of Δhm,m+1 obtained for (m = 1, 2, and 3). It is observed that among DFA1-4 methods, none is systematically the best over all synoptic stations. This result indicates that different kind of local trend exist in rainfall time series over Bénin. In addition, the DFA1-4 method affects significantly and differently their multifractal characteristics. Thus, the best DFA method depends on the geographical position of the studied station. The comparison of the spatial average of Δhm,m+1 (shown in Figure 4) reveals that DFA2 minimizes Δhm,m+1 spatially averaged among DFA1-DFA4. Therefore, DFA2, is adopted in the rest of this study.Figure 4
Spatial distribution ofΔhm,m+1 obtained by the DFA1-4 method (Δhm,m+1=hDFAm−hDFAm+1). For a DFAm, an m-order polynomial is used to calculate the local trend in each segment.


### 3.3. Origins of Multifractality

Figure5 shows the absolute difference of generalized Hurst scaling function obtained by MFDFA2 corresponding to DFA2 for original and shuffled data horigq−hshufq = hCORRq and original and surrogate data horigq−hsurrq = hPDFq as a function of q (q = −3 : 0.5 : 3) for the studied stations. The results reveal that in Bénin, excepted Parakou, at the synoptic station, the long-term temporal correlations are systematically the source of multifractality in rainfall.Figure 5
Absolute difference of Hurst exponentshCORRq and hPDFq as a function of (q) (where q) is the order of the fluctuation function), q = −3 : 0.5 : 3 at each synoptic station: (a) Cotonou station, (b) Bohicon station, (c) Savè station, (d) Parakou station, (e) Natitingou station, and (f) Kandi station.
(a)
(b)
(c)
(d)
(e)
(f)


### 3.4. Temporal Evolution of Long-Term Memory Characteristics

Figure6 presents South-North distribution of DFA exponents (ℎ) obtained over 1951–1970, 1971–1990, and 1991–2010 by the DFA2 method. For all synoptic stations except Natitingou, h1951−1970<h1971−1990<h1991−2010. This means that except Natitingou, the strength of long-term memory characteristic decreases each twenty years of the study period. The reason behind the exception done at Natitingou could be linked to the fact that Natitingou is the rainiest region of Bénin because of the presence of mountains (∼800 m of high). However, DFA exponents are between 0.5 and 1 on each period. Thus, there exists positive long-term memory characteristics in the rainfall field at all stations on these periods.Figure 6
Spatiotemporal variation of DFA exponents (ℎ) obtained by theDFA2 method.


### 3.5. Error Analysis of DFA Exponents (h)

Figure7 presents the latitudinal variation of Pearson coefficient for Hurst exponent for each subperiod. The results show that whatever the studied subperiod, the Hurst exponent is estimated by a Pearson coefficient greater than 97%, indicating the robustness of the linear regression used to fit the fluctuation function in a log-log representation.Figure 7
Spatiotemporal variation of Pearson coefficient for DFA exponent: (a) from 1951–2010, (b) from 1951–1970, (c) from 1971–1990, and (d) from 1991–2010. Red and black represent the upper bound and the lower bound of the 95% confidence interval. Green represents the Pearson coefficient.
(a)
(b)
(c)
(d)The upper and the lower bound of the 95% confidence interval have the same sign meaning that the computation is statistically significant.The error in the computation above for each subperiod is shown in Figure8. One can remark that all the values of the MSE are lower than 8% showing a robustness of the computation process.Figure 8
Spatiotemporal variation of the mean squared error (MSE).


### 3.6. Assessment to Climate Change Effect by the DFA Method

The comparison of DFA exponents for original, shuffled, and surrogated data obtained over 1951–1970, 1971–1990, 1991–2010, and 1951–2010 is shown in Figure9. For a good understanding of the results analysis, it is important to keep in mind that the red and black curves provide information on climatic and geographical influence, respectively. The analysis of results (shown by red curve) reveal the following: (1) In 1951–1970, considered as humid period in literature, the minimum and maximum values of DFA exponents are obtained at Savè and Natitingou stations. These findings could be explained by the fact that Savè is in the transitional zone between the Sudanian and subequatorial climates while Natitingou is in a mountainous area. Natitingou is therefore more watered, which is not surprising since Natitingou is known as the most watered area of Bénin. These results associated with the increasing of DFA exponents obtained and presented in blue curves at all stations confirm that 1951–1970 is a humid period. (2) In 1971–1990, considered as drought period in the literature, the transition zone is much more affected by drought, and this is noted by the greater values of h observed at Savè. The result confirms that a great drought has occurred during this subperiod. (3) In 1991–2010, considered as the end of the drought, it is observed that only the Natitingou station was actually watered, and this can be facilitated by the presence of mountains. However, no real recovery has been made in the other stations. Thus, in this period the, drought is still present in several regions of Bénin, and it is strongly severe in Savè and Parakou. These findings are confirmed by the randomly variation of DFA exponents over synoptic stations presented in blue curves, contrary to what is obtained by some authors in the literature for West Africa in general [36–38]. In the fractal approach, our results show that the subperiod 1991–2010 is not really a transition period, and the drought is prolonged until 2010. This result is not too surprising since there is some unanimity to date the onset of the drought episode around the 1970s, with the support of breakthrough statistical tests like that of Pettitt [39], there is a debate about its end [36].Figure 9
h
orig (black line), hshuf (red line), and hsurr (blue line) for original, shuffled, and surrogated data obtained over 1951–1970 (a), 1971–1990 (b), 1991–2010 (c), and 1951–2010 (d).
(a)
(b)
(c)
(d)


## 3.1. Long-Term Memory Characteristic Study

In the present study, we calculate the results for all the six synoptic stations in the range of 8–1024 days. Figure3 shows South-North distribution of DFA exponent (ℎ) obtained by the DFA1-4 method. The main results of this figure are the following: (a) h values depend on the applied DFAm method; (b) except Cotonou, whatever the DFAm method, h values present spatial variation, there is a South-North gradient, (c) over all synoptic stations, h values are systematically between 0.5 and 1, this means that there exists positive long-term memory characteristic in the rainfall field, and (d) among synoptic stations, the highest and lowest values of h are, respectively, obtained at Cotonou and Bohicon stations, indicating the existence of strong and weak long-term memory characteristics in Cotonou and Bohicon.Figure 3
South-North distribution of DFA exponents (ℎ) obtained by theDFA1-4 method. Cotonou (Co), Bohicon (Bo), Savè (Sa), Parakou (Pa), Natitingou (Na), and Kandi (Ka). For a DFAm, an m-order polynomial is used to calculate the local trend in each segment.


## 3.2. The Best DFA Method

Figure4 presents the South-North distribution of Δhm,m+1 obtained for (m = 1, 2, and 3). It is observed that among DFA1-4 methods, none is systematically the best over all synoptic stations. This result indicates that different kind of local trend exist in rainfall time series over Bénin. In addition, the DFA1-4 method affects significantly and differently their multifractal characteristics. Thus, the best DFA method depends on the geographical position of the studied station. The comparison of the spatial average of Δhm,m+1 (shown in Figure 4) reveals that DFA2 minimizes Δhm,m+1 spatially averaged among DFA1-DFA4. Therefore, DFA2, is adopted in the rest of this study.Figure 4
Spatial distribution ofΔhm,m+1 obtained by the DFA1-4 method (Δhm,m+1=hDFAm−hDFAm+1). For a DFAm, an m-order polynomial is used to calculate the local trend in each segment.


## 3.3. Origins of Multifractality

Figure5 shows the absolute difference of generalized Hurst scaling function obtained by MFDFA2 corresponding to DFA2 for original and shuffled data horigq−hshufq = hCORRq and original and surrogate data horigq−hsurrq = hPDFq as a function of q (q = −3 : 0.5 : 3) for the studied stations. The results reveal that in Bénin, excepted Parakou, at the synoptic station, the long-term temporal correlations are systematically the source of multifractality in rainfall.Figure 5
Absolute difference of Hurst exponentshCORRq and hPDFq as a function of (q) (where q) is the order of the fluctuation function), q = −3 : 0.5 : 3 at each synoptic station: (a) Cotonou station, (b) Bohicon station, (c) Savè station, (d) Parakou station, (e) Natitingou station, and (f) Kandi station.
(a)
(b)
(c)
(d)
(e)
(f)


## 3.4. Temporal Evolution of Long-Term Memory Characteristics

Figure6 presents South-North distribution of DFA exponents (ℎ) obtained over 1951–1970, 1971–1990, and 1991–2010 by the DFA2 method. For all synoptic stations except Natitingou, h1951−1970<h1971−1990<h1991−2010. This means that except Natitingou, the strength of long-term memory characteristic decreases each twenty years of the study period. The reason behind the exception done at Natitingou could be linked to the fact that Natitingou is the rainiest region of Bénin because of the presence of mountains (∼800 m of high). However, DFA exponents are between 0.5 and 1 on each period. Thus, there exists positive long-term memory characteristics in the rainfall field at all stations on these periods.Figure 6
Spatiotemporal variation of DFA exponents (ℎ) obtained by theDFA2 method.


## 3.5. Error Analysis of DFA Exponents (h)

Figure7 presents the latitudinal variation of Pearson coefficient for Hurst exponent for each subperiod. The results show that whatever the studied subperiod, the Hurst exponent is estimated by a Pearson coefficient greater than 97%, indicating the robustness of the linear regression used to fit the fluctuation function in a log-log representation.Figure 7
Spatiotemporal variation of Pearson coefficient for DFA exponent: (a) from 1951–2010, (b) from 1951–1970, (c) from 1971–1990, and (d) from 1991–2010. Red and black represent the upper bound and the lower bound of the 95% confidence interval. Green represents the Pearson coefficient.
(a)
(b)
(c)
(d)The upper and the lower bound of the 95% confidence interval have the same sign meaning that the computation is statistically significant.The error in the computation above for each subperiod is shown in Figure8. One can remark that all the values of the MSE are lower than 8% showing a robustness of the computation process.Figure 8
Spatiotemporal variation of the mean squared error (MSE).


## 3.6. Assessment to Climate Change Effect by the DFA Method

The comparison of DFA exponents for original, shuffled, and surrogated data obtained over 1951–1970, 1971–1990, 1991–2010, and 1951–2010 is shown in Figure9. For a good understanding of the results analysis, it is important to keep in mind that the red and black curves provide information on climatic and geographical influence, respectively. The analysis of results (shown by red curve) reveal the following: (1) In 1951–1970, considered as humid period in literature, the minimum and maximum values of DFA exponents are obtained at Savè and Natitingou stations. These findings could be explained by the fact that Savè is in the transitional zone between the Sudanian and subequatorial climates while Natitingou is in a mountainous area. Natitingou is therefore more watered, which is not surprising since Natitingou is known as the most watered area of Bénin. These results associated with the increasing of DFA exponents obtained and presented in blue curves at all stations confirm that 1951–1970 is a humid period. (2) In 1971–1990, considered as drought period in the literature, the transition zone is much more affected by drought, and this is noted by the greater values of h observed at Savè. The result confirms that a great drought has occurred during this subperiod. (3) In 1991–2010, considered as the end of the drought, it is observed that only the Natitingou station was actually watered, and this can be facilitated by the presence of mountains. However, no real recovery has been made in the other stations. Thus, in this period the, drought is still present in several regions of Bénin, and it is strongly severe in Savè and Parakou. These findings are confirmed by the randomly variation of DFA exponents over synoptic stations presented in blue curves, contrary to what is obtained by some authors in the literature for West Africa in general [36–38]. In the fractal approach, our results show that the subperiod 1991–2010 is not really a transition period, and the drought is prolonged until 2010. This result is not too surprising since there is some unanimity to date the onset of the drought episode around the 1970s, with the support of breakthrough statistical tests like that of Pettitt [39], there is a debate about its end [36].Figure 9
h
orig (black line), hshuf (red line), and hsurr (blue line) for original, shuffled, and surrogated data obtained over 1951–1970 (a), 1971–1990 (b), 1991–2010 (c), and 1951–2010 (d).
(a)
(b)
(c)
(d)


## 4. Discussion

The obvious reasons behind the South-North gradient observed in Figure3 except Cotonou may be explained by the fact that Bohicon is in the transition region between Sudanian and subequatorial region and Kandi is near to the Sahelian region. The highest values of h, obtained at Cotonou may be due firstly to the presence of Atlantic Ocean, which can affect Cotonou rainfall by their sea surface temperature and secondly to the anthropic effect created by the presence of many industries around.The comparison of the spatial average ofΔhm,m+1 (shown in Figure 4) shows that meanΔh2,3<meanΔh3,4<meanΔh1,2; thus, DFA2 minimizes Δhm,m+1 spatially averaged among DFA1-DFA4. In Bénin, at the synoptic scale, DFA2 is the best method in terms of spatial average. This result confirms those of [6, 24]. These authors have found the second-order DFA (DFA2) method as the best and the most commonly used order, respectively, for outgoing longwave and temperature records.It is clearly observed in Figure5 that, for all synoptic stations, for almost all values of q, with the exception of Parakou station, hCORRq is systematically larger than hPDFq. This means that the effect of the broad probability density function on multifractality is more important in Parakou than that of long-term temporal correlations. Thus, in Bénin, except Parakou, at the synoptic station the long-term temporal correlations are systematically the source of multifractality in rainfall. Therefore, the source of multifractality in rainfall depends on the geographical position and the type of climate of the studied station. However, nonzero values of hCORRq and hPDFq indicate that both of them influence the multifractality. In specific cases, the long-range correlations have a slightly greater effect on rainfall multifractality than the broad probability density function, as hPDFq<0>hPDFq>0.


## 5. Conclusions

Traditional approaches cannot detect natural variability due to the dynamic behavior of climate variables. In this paper, the DFA method is used to analyze the long-term memory characteristic of rainfall over Bénin. The results reveal that over all synoptic stations, DFA exponent (h) is systematically between 0.5 and 1, meaning that there is a positive long-term memory characteristic in the rainfall field. However, among synoptic stations, highest and lowest values of h are, respectively, obtained at Cotonou and Bohicon stations, indicating strong and weak long-term memory characteristics in Cotonou and Bohicon. It is observed that among DFA1-4 methods, none is systematically the best over all synoptic stations. This result indicates that different kinds of local trend exist in rainfall time series over Bénin. Thus, the best DFA method depends on the geographical position of the studied station. The comparison of the spatial average of Δh between two successive orders of DFA reveals that DFA2 minimizes more Δh among DFA1-DFA4. Therefore, DFA2 is considered as the best method and is used to calculate LTM in this study. For all stations, for almost all values of q, with the exception of Parakou station,hCORRq is larger than hPDFq. This means that the effect of the broad probability density function on multifractality is more important than that of long-term temporal correlations, only for Parakou. Thus, in Bénin, except Parakou, at the synoptic station, the long-term temporal correlations are systematically the source of multifractality in rainfall. However, nonzero values of hCORRq and hPDFq indicate that both broad probability density function and long-term temporal correlations influence the multifractality. Except Natitingou, the strength of long-term memory characteristic decreases each twenty years of the study period. In the fractal approach, our results show that the subperiod 1991–2010 is not really a transition period, and the drought is prolonged until 2010. Fractal theory is adapted to characterize Bénin climate. It should be good to confirm the existence of long-range correlations by studying the autocorrelation function (ACF) and the variability of the local slopes for large time scales as suggested in Maraun et al. [18] and clarified by Varotsos and Efstathiou [17]. These aspects of the analysis will be the objectives of future studies.



---

*Source: 1353195-2019-01-09.xml*