# Bivariate EMD-Based Data Adaptive Approach to the Analysis of Climate Variability

Authors: Md. Khademul Islam Molla; Poly Rani Ghosh; Keikichi Hirose
Journal: Discrete Dynamics in Nature and Society
 (2011)

---

## Abstract

This paper presents a data adaptive approach for the analysis of climate variability using bivariate empirical mode decomposition (BEMD). The time series of climate factors: daily evaporation, maximum and minimum temperatures are taken into consideration in variability analysis. All climate data are collected from a specific area of Bihar in India. Fractional Gaussian noise (fGn) is used here as the reference signal. The climate signal and fGn (of same length) are combined to produce bivariate (complex) signal which is decomposed using BEMD into a finite number of sub-band signals named intrinsic mode functions (IMFs). Both of climate signal as well as fGn are decomposed together into IMFs. The instantaneous frequencies and Fourier spectrum of IMFs are observed to illustrate the property of BEMD. The lowest frequency oscillation of climate signal represents the annual cycle (AC) which is an important factor in analyzing climate change and variability. The energies of the fGn's IMFs are used to define the data adaptive threshold to separate AC. The IMFs of climate signal with energy exceeding such threshold are summed up to separate the AC. The interannual distance of climate signal is also illustrated for better understanding of climate change and variability.


---

## Body

## 1. Introduction

The climate variability and change (CVC) element emphasizes research to improve descriptions and understanding of past and current climate, as well as to advance national modeling capabilities to simulate climate and project how climate and related Earth systems may change in the future. The CVC refers to shifts in the mean state of the climate or in its variability, persisting for an extended period (decades or longer). Research under this element encompasses time scales ranging from short-term climate variations of a season or less to longer-term climate changes occurring over decades to centuries. The CVC element places a high priority on improving understanding and predictions of phenomena that may cause high impacts on society, the economy, and the environment. Climate variability refers to variations in the mean state of climate on all temporal and spatial scales beyond that of individual weather events. It may be due to natural changes or to persistent anthropogenic changes in the composition of the atmosphere or in land use. A lot of people have been observing drastic climate changes throughout certain portions of the world for the past few decades.There is a perception that extreme natural disasters such as floods, droughts, and heat waves, have become more frequent. This change in climate plays an important role in the Earth’s sustainability. The impact of hydrological change as a result of global warming caused by anthropogenic emissions of greenhouse gases is a fundamental concern [1–4]. In order to control water resources in the face of drought, flood and soil erosion, which frequently present a serious threat to human life and natural ecosystems, risk assessments are being required more frequently by policy makers [5]. Climate is currently changing in ways that mirror the effects of global warming [6]. There is also increasing demand for climate change information, particularly from policy makers for impact assessment studies [7]. Several linear statistical models have been applied to climate records, but the answers are not conclusive due to the high sensitivity of model results to model parameters [8–10], especially when stochastic processes are taken into account [11]. Various approaches have been employed to develop climate change scenarios at different scales [12]. In recent years, several weather events have caused large losses of life as well as a tremendous increase in economic losses from weather hazards. These life and property losses helped raise the alarm over the possibility that the recent increases were due to a shifting climate [13]. Geographical distribution of the change in annual mean surface air temperature rise is greater over continents than over oceans and greater in the Polar Regions than in tropical region.The features are consistent with previous studies. In the tropical Pacific, temperature rise is greater in the central-east part than in the western part. Such an anomaly pattern is observed in the El Nino years. The warming patterns among the scenarios are similar although the magnitudes are different [14]. Many critical impacts of climate are controlled by extreme events rather than mean values. Aspects of human activity and environment are usually well adjusted to mean climate conditions and show little sensitivity to moderate variations around these means. Systems that are particularly vulnerable are agricultural and forest ecosystems, coasts, and water resources. High temperatures can exacerbate drought conditions. Very high temperature also damage crops and reduce yields. However, low temperatures are important to some temperate horticulture and cereal crops in promoting yield and development [15].The surface runoff changes had regional differences and both the increase and decrease were suggested in summer. It might increase the risk of mismatch between water demand and water availability in the agricultural region. Under the global warming, both temperature and humidity were projected to increase in Asian region. The increases of annual mean surface runoff and its large fluctuations were projected in a lot of Asian regions. In some regions, the projected seasonal changes of hydrological cycles under the global warming potentially increase the risk of droughts and floods [16]. Coughlin and Tung [17] found that the atmosphere warms during the solar maximum almost everywhere over the globe. It should be pointed out that changes in the correlation values with latitude do not imply similar amplitude changes with latitude. However, the fact that the correlation with the solar flux is positive everywhere over the globe does imply that, on average, the temperatures increase during solar maxima at all latitude. This makes the phenomenon difficult to explain with dynamical mechanisms involving over turning meridional circulations in the troposphere. This is useful information to atmosphere dynamicity wishing to come up with a mechanism for the influence of solar cycle variations on the lower atmosphere. The El Nino and the southern oscillation phenomenon (ENSO) is the primary driver of interannual climate variability and have a large economic and social impact over the universe [18]. Fedorov and Philander [19] demonstrated that mean fluctuations of decadal timescale do contribute significantly to the later unusual ENSO events and suggested that global warming cannot be ruled out as a suspect [20]. El Nino-like differences can be found in the projected results by the AGCM under global warming [21]. As a result of global warming, the stronger temperature increase in the higher latitudes, the meridional gradients of air temperature became weak and waster caused by thermal wind were weakened. Monsoon westerly in the Arabian Sea, which is associated with Somali jet, was projected to be stronger and to bring more abundant water vapor to the south of India and the Bay of Bengal. Under global warming, increase of water vapor and higher air temperature over the land than over the ocean were projected. It enhanced a giant land-sea breeze, Asian monsoon circulation. As a result of the changes in the synoptic scale flow patterns and precipitation under the global warming, the increase of annual mean surface runoff was projected in a lot of Asian regions [22].Sampling is an important factor for modeling and analysis of climate data. Any climate signalx(t) at time instant t can be defined as x(t)⇒x(nT), where n is a nonnegative number representing the index of the discrete sample and T is the sampling period. For uniform sampling rate, T can be omitted and then x(t)⇒x(n), with n representing the discrete sample index. The data used in this study are collected with uniform sampling over total period of acquisition. The nonuniform sampling/discretization is always better to adapt with the signal variation. The basic concept of nonuniform sampling is to increase the sampling rate to acquire the signal with higher variation and the reverse when the lower variation of signal is observed in acquisition. The current research trend is to apply the nonuniform discretization to increase the efficiencies of different systems [23, 24]. The stabilization of linear time-invariant dynamic system is achieved by applying multirate discretization combined with fractional-order hold [23]. The stabilization of such system is not guaranteed without using multirate sampling. The controllability and observability are investigated for Caputo fractional differential linear system of any real order α in [24]. The developed model also supports those properties under nonuniform sampling. The authors have also proved that the choice of appropriate sampling instance is not restrictive as a result of the properties of associate Chebyshev’s systems. Since we are using the secondary data, the sampling rate is kept as it is although the nonuniform sampling has many advantages.A new nonlinear technique, empirical mode decomposition (EMD), has recently been pioneered by Huang et al. [25] for adaptively representing nonstationary time series data. Although it proved remarkably effective, the technique is faced with the difficulty of being essentially defined by an algorithm and therefore does not admitting an analytical formulation which would allow for a theoretical analysis and performance evaluation [25, 26]. The EMD is employed to analyze the climate change phenomena [27], where only decomposition efficiency and the spectral properties of the climate signal are observed. Three climate signals namely, daily evaporation, maximum, and minimum temperature are used in this study to analyze the climate variability. The data are collected from the area named Giridhi of Bihar in India. The details about the data are explained in Section 3.1. In this paper, recently developed bivariate EMD (BEMD) is applied to the climate signals to demonstrate its variability properly. The fGn is combined with one of the climate signal to make suitable for applying the BEMD. The data adaptive filtering approach is introduced to separate the annual cycle (AC) from the raw climate signal. The annual cycle is supposed to be unchanged over the years and sometimes changed due to the climate variability. Such type of variability is detected using BEMD and demonstrated by overcoming the limitations of the univariate EMD.


## 2. BEMD for Climate Variability Analysis

Climate signals in the atmosphere are often lost among the noise and other data collection and assimilation problems. Although climate signals are mainly recognized by their time scales, previous analysis, such as the use of empirical orthogonal functions (EOFs), relies on spatial decompositions to remove excess noise. In some cases this is useful, but in general, there is no a priori reason to expect temporal signals to be defined in terms of spatial patterns ordered by their variances (as required by EOF analysis). In fact, time series analysis is more appropriate for the decomposition of climate signals. The traditional time series analysis tools usually rely on Fourier transforms in one way or another. However, Fourier transforms lead to inconclusive interpretations due mainly to the global nature (in the time domain) of the transforms.Even wavelet analysis, developed to deal with non stationary and local frequency changes, produces confusing and sometimes contradicting results when applied to climate signals [28, 29]. The authors use similar types of wavelet analysis on global, surface climate data and come up with three different results. Using wavelet analysis, it is sometimes difficult to distinctly define local frequency changes because the spectra are created by stepping through various predetermined frequencies producing an often blurred result. Wavelets have the additional problem of shift variance. If the starting point is varied, like by dropping the initial point, wavelet analysis can produce completely different results. In comparison, the EMD method makes no assumption about linearity or stationary and the intrinsic mode functions (IMFs) are usually easy to interpret and relevant to the physical system being studied [25, 26] for further discussion of the method and comparison to other time series analysis techniques.Empirical mode decomposition (EMD) can be a useful time series analysis tool. One example where EMD is found to be particularly useful is in analyzing climate records of the atmosphere beyond annual time scales. Global climate phenomena are often separated in temporal, rather than spatial, scales. Therefore, time series analysis is more appropriate for the initial decomposition of climate data than spatial methods that have previously been used to find climate components. The nonstationary and nonlinear nature of climate signals make the EMD appropriate to be used in analyzing such signals. One difficulty encountered when using this method is the sensitivity to end point treatments. The envelopes are calculated using a cubic spline; however, splines are notoriously sensitive to end points. It is important to make sure that the end effects do not propagate into the interior solution. The BEMD is employed here for robust and data adaptive analysis of climate variability.


### 2.1. Bivariate EMD

The empirical mode decomposition (EMD) is a signal processing decomposition technique that decomposes the signal into waveforms modulated in both amplitude and frequency by extracting all of the oscillatory modes embedded in the signal [25]. The decomposition is an intuitive and adaptive signal-dependent decomposition and does not require any conditions about the stationary and linearity of the signal. The waveforms extracted by EMD are named Intrinsic Mode Functions (IMFs). Each IMF is symmetric, is assumed to yield a meaningful local frequency, and different IMFs do not exhibit the same frequency at the same time. The traditional EMD is only suitable for univariate (real-valued) signals.The Complex Empirical Mode Decomposition (Complex-EMD) is an extension of the basic EMD suitable for dealing with complex signals [30]. The motivation to extend EMD is that a large number of signal processing applications have complex signals. In addition, this extension is applied on both the real and imaginary parts simultaneously because complex signals have a mutual dependence between the real and imaginary parts. Thus, if the decomposition is done separately, the mutual dependency will be lost.The Bivariate Empirical Mode Decomposition (BEMD) is more generalized extension of the EMD to complex signals. The main difference between the BEMD and the Complex-EMD is that the latter uses the basic EMD to decompose complex signals, whereas the BEMD adapts the rationale underlying the EMD to a bivariate framework [31, 32]. In BEMD, two variables are decomposed simultaneously based on their rotating properties. The algorithm of the BEMD, as proposed in [31], is as follows.(1)
For1<m<M,
(a)
projectx(n) on direction ϕm:pϕm(n)=Re(e-jϕmx(n)),(b)
extract the maxima ofpϕm(n):(nim,pim),(c)
interpolate the set of points(nim,ejϕmpim) to obtain the partial envelope curve in direction ϕm named  eϕm(n),(2)
compute the mean of all tangents:e(n)=2/M∑meϕm(n),(3)
subtract the mean to obtainc(n)=x(n)-e(n),(4)
test ifc(n) is an IMF:
(a)
if yes, repeat the procedure from the step (2.1) on the residual signal,(b)
if not, replacex(n) with c(n) and repeat the procedure from step (2.1).The Bivariate-EMD can now be expressed as(2.1)x(n)=∑k=1Kck(n)+rK(n),
where ck(n) denotes the kth extracted complex empirical mode and rK(n) the final residual.Herex(n) can be modeled as x(n)=a+jb, where a represents one of the weather signal and the b is the fractional Gaussian noise (fGn). The BEMD is employed to decompose the climate signals to filter the annual cycles. The reference signal fGn together with the daily evaporation signals is decomposed using BEMD as shown in Figure 1. The normalized fGn and the climate signals are considered as real and imaginary parts of the complex signal, respectively.BEMD of climate signal (daily evaporation) associated with the fGn. The left and right columns represent the IMFs corresponding to climate signals (daily evaporation) and fGn, respectively.
(a)(b)It is well known that the EMD of fGn is acting as dyadic filter bank [33] and hence it is considered as the reference signal to determine the IMFs representing the actual climate signals. The spectra of the individual IMF of fGn as well as the climate signal are shown in Figure 2. The IMFs are interpreted as the basis vectors representing the data [25]. The Fourier spectra of the IMF components are identical in shape and cover the same area on a semilogarithmic period scale. It can effectively be applied to the data without harmonics, whereas any harmonic analysis (Fourier, wavelet, and so on) would end up in the same context with a much less compact and physically less meaningful decomposition. The focus in this study will be mainly on the most recent climate records.The Fourier spectra of the IMFs of (a) fractional Gaussian noise (fGn) and (b) climate signal (daily evaporation).
(a)(b)


### 2.2. Instantaneous Frequency

Instantaneous frequency (IF) represents the signal’s frequency at an instance. The concept of IF is physically meaningful only when applied to monocomponent signals. In order to apply the concept of IF to arbitrary signals, it is necessary to decompose the signal into a series of monocomponent contributions. In the recent approaches [25, 26], EMD technique decomposes a time domain signal into a series of mono-component contributions designated as intrinsic mode functions (IMFs). Then the IF derived for each component provides the meaningful physical information. The IF is defined as the rate of change of the phase angle at the instant of the analytic version of the signal. Then the complex (analytic) version of the kth IMF ck(n) is defined as(2.2)zk(n)=ck(n)+jℏ[ck(n)]=ak(n)ejθk(n)
wherej=√(-1), ℏ[·] represents the Hilbert transform; ak(n) and θk(n) are instantaneous amplitude and phase, respectively, of the kth IMF. The analytic signal is advantageous in determining the instantaneous quantities such as energy, phase, and frequency. The discrete-time IF of kth IMF is computed by taking the derivative of discrete phase vector θk(n) and defined as
(2.3)ωk(n)=|dθ̃k(n)dn|,
where θ̃k(n) represents the unwrapped version of phase θk(n) at sample index n. More precisely, the derivative in (2.3) is evaluated to be approximated |θ̃k′(n)|  by the absolute value of forward difference that is, ωk(n)=|θ̃k(n+1)-θ̃k(n)|. Phase unwrapping is necessary in IF computation. During the calculation of discrete phase using analytic signal method given in (2.2), the phase is returned in a form that suffers from 2π phase jump yielding an inherently wrapped output. Such wrapped phase is not useful until the discontinuities are removed—and this is practically achieved by employing a phase unwrapping algorithm (Matlab unwrapping function is used here).The IF of individual IMF given in Figure1 is shown in Figure 3. It should be noted that the IF obtained by direct derivative introduces the abrupt fluctuations and hence the smoothing is required. Here, the moving average smoothing filter of three sample length is used, to produce a slowly varying estimate of frequency.Instantaneous frequencies of the IMFs shown in Figure1.
(a)(b)


### 2.3. Bandpass Filtering

One of the useful properties of EMD as well as BEMD is that the first IMF contains the signal component with the highest oscillation frequency and theKth (K represents the total number of IMF) IMF comprises the lowest frequency components of the analyzing signal. Thus the average frequency ordering is made from highest to lowest (descending) with the ascending order of indices of the IMFs. Another way to explain EMD is that the first IMF extracts out the highest frequency oscillation that remains in the signal. Thus locally, each IMF contains lower frequency oscillation than the one extracted just before. The IMFs are considered as the basis function with ordered frequency components representing the original signal. Such basis components are used to design time-space filtering [26]. Traditionally, filtering is carried out in frequency space only. Using IMFs, however, we can devise a time-space filtering. For example, the result of high, band, and low pass filters by using EMD of any signal s(n) having K IMF components can be simply expressed as(2.4)sh(n)=∑k=1y1ck(n),sb(n)=∑k=y1+1y2ck(n),sl(n)=∑k=y2+1Kck(n)+rK(n),
where 1<y1<y2<K; sh(n), sb(n), and sl(n) represent high-pass, band-pass, and low-pass signals, respectively; rK(n) is the final residue after extracting all the IMFs. The decomposition ends up with that residual signal. The final residue is defined as a signal with no oscillation or the signal with constant amplitude close to zero. It represents the mean trend of the signal. The main advantage of this time-space filtering is that the results preserve full meaning, nonlinearity, and nonstationary in physical space. In this paper, the low frequency annual cycle is to be separated by applying the bandpass filtering. The rest is to select the criteria to select the values of y1 and y2. The linear model of IMF’s energy of fGn and confidence limits are used to index of the IMFs for the band pass filtering. The three subband signals of daily evaporation data derived from the BEMD shown in Figure 1 with (2.4) are illustrated in Figure 4; the values of y1 and y2 are assigned as 4 and 8, respectively.Figure 4
Three subband signals obtained by applying BEMD on daily evaporation data.


## 2.1. Bivariate EMD

The empirical mode decomposition (EMD) is a signal processing decomposition technique that decomposes the signal into waveforms modulated in both amplitude and frequency by extracting all of the oscillatory modes embedded in the signal [25]. The decomposition is an intuitive and adaptive signal-dependent decomposition and does not require any conditions about the stationary and linearity of the signal. The waveforms extracted by EMD are named Intrinsic Mode Functions (IMFs). Each IMF is symmetric, is assumed to yield a meaningful local frequency, and different IMFs do not exhibit the same frequency at the same time. The traditional EMD is only suitable for univariate (real-valued) signals.The Complex Empirical Mode Decomposition (Complex-EMD) is an extension of the basic EMD suitable for dealing with complex signals [30]. The motivation to extend EMD is that a large number of signal processing applications have complex signals. In addition, this extension is applied on both the real and imaginary parts simultaneously because complex signals have a mutual dependence between the real and imaginary parts. Thus, if the decomposition is done separately, the mutual dependency will be lost.The Bivariate Empirical Mode Decomposition (BEMD) is more generalized extension of the EMD to complex signals. The main difference between the BEMD and the Complex-EMD is that the latter uses the basic EMD to decompose complex signals, whereas the BEMD adapts the rationale underlying the EMD to a bivariate framework [31, 32]. In BEMD, two variables are decomposed simultaneously based on their rotating properties. The algorithm of the BEMD, as proposed in [31], is as follows.(1)
For1<m<M,
(a)
projectx(n) on direction ϕm:pϕm(n)=Re(e-jϕmx(n)),(b)
extract the maxima ofpϕm(n):(nim,pim),(c)
interpolate the set of points(nim,ejϕmpim) to obtain the partial envelope curve in direction ϕm named  eϕm(n),(2)
compute the mean of all tangents:e(n)=2/M∑meϕm(n),(3)
subtract the mean to obtainc(n)=x(n)-e(n),(4)
test ifc(n) is an IMF:
(a)
if yes, repeat the procedure from the step (2.1) on the residual signal,(b)
if not, replacex(n) with c(n) and repeat the procedure from step (2.1).The Bivariate-EMD can now be expressed as(2.1)x(n)=∑k=1Kck(n)+rK(n),
where ck(n) denotes the kth extracted complex empirical mode and rK(n) the final residual.Herex(n) can be modeled as x(n)=a+jb, where a represents one of the weather signal and the b is the fractional Gaussian noise (fGn). The BEMD is employed to decompose the climate signals to filter the annual cycles. The reference signal fGn together with the daily evaporation signals is decomposed using BEMD as shown in Figure 1. The normalized fGn and the climate signals are considered as real and imaginary parts of the complex signal, respectively.BEMD of climate signal (daily evaporation) associated with the fGn. The left and right columns represent the IMFs corresponding to climate signals (daily evaporation) and fGn, respectively.
(a)(b)It is well known that the EMD of fGn is acting as dyadic filter bank [33] and hence it is considered as the reference signal to determine the IMFs representing the actual climate signals. The spectra of the individual IMF of fGn as well as the climate signal are shown in Figure 2. The IMFs are interpreted as the basis vectors representing the data [25]. The Fourier spectra of the IMF components are identical in shape and cover the same area on a semilogarithmic period scale. It can effectively be applied to the data without harmonics, whereas any harmonic analysis (Fourier, wavelet, and so on) would end up in the same context with a much less compact and physically less meaningful decomposition. The focus in this study will be mainly on the most recent climate records.The Fourier spectra of the IMFs of (a) fractional Gaussian noise (fGn) and (b) climate signal (daily evaporation).
(a)(b)


## 2.2. Instantaneous Frequency

Instantaneous frequency (IF) represents the signal’s frequency at an instance. The concept of IF is physically meaningful only when applied to monocomponent signals. In order to apply the concept of IF to arbitrary signals, it is necessary to decompose the signal into a series of monocomponent contributions. In the recent approaches [25, 26], EMD technique decomposes a time domain signal into a series of mono-component contributions designated as intrinsic mode functions (IMFs). Then the IF derived for each component provides the meaningful physical information. The IF is defined as the rate of change of the phase angle at the instant of the analytic version of the signal. Then the complex (analytic) version of the kth IMF ck(n) is defined as(2.2)zk(n)=ck(n)+jℏ[ck(n)]=ak(n)ejθk(n)
wherej=√(-1), ℏ[·] represents the Hilbert transform; ak(n) and θk(n) are instantaneous amplitude and phase, respectively, of the kth IMF. The analytic signal is advantageous in determining the instantaneous quantities such as energy, phase, and frequency. The discrete-time IF of kth IMF is computed by taking the derivative of discrete phase vector θk(n) and defined as
(2.3)ωk(n)=|dθ̃k(n)dn|,
where θ̃k(n) represents the unwrapped version of phase θk(n) at sample index n. More precisely, the derivative in (2.3) is evaluated to be approximated |θ̃k′(n)|  by the absolute value of forward difference that is, ωk(n)=|θ̃k(n+1)-θ̃k(n)|. Phase unwrapping is necessary in IF computation. During the calculation of discrete phase using analytic signal method given in (2.2), the phase is returned in a form that suffers from 2π phase jump yielding an inherently wrapped output. Such wrapped phase is not useful until the discontinuities are removed—and this is practically achieved by employing a phase unwrapping algorithm (Matlab unwrapping function is used here).The IF of individual IMF given in Figure1 is shown in Figure 3. It should be noted that the IF obtained by direct derivative introduces the abrupt fluctuations and hence the smoothing is required. Here, the moving average smoothing filter of three sample length is used, to produce a slowly varying estimate of frequency.Instantaneous frequencies of the IMFs shown in Figure1.
(a)(b)


## 2.3. Bandpass Filtering

One of the useful properties of EMD as well as BEMD is that the first IMF contains the signal component with the highest oscillation frequency and theKth (K represents the total number of IMF) IMF comprises the lowest frequency components of the analyzing signal. Thus the average frequency ordering is made from highest to lowest (descending) with the ascending order of indices of the IMFs. Another way to explain EMD is that the first IMF extracts out the highest frequency oscillation that remains in the signal. Thus locally, each IMF contains lower frequency oscillation than the one extracted just before. The IMFs are considered as the basis function with ordered frequency components representing the original signal. Such basis components are used to design time-space filtering [26]. Traditionally, filtering is carried out in frequency space only. Using IMFs, however, we can devise a time-space filtering. For example, the result of high, band, and low pass filters by using EMD of any signal s(n) having K IMF components can be simply expressed as(2.4)sh(n)=∑k=1y1ck(n),sb(n)=∑k=y1+1y2ck(n),sl(n)=∑k=y2+1Kck(n)+rK(n),
where 1<y1<y2<K; sh(n), sb(n), and sl(n) represent high-pass, band-pass, and low-pass signals, respectively; rK(n) is the final residue after extracting all the IMFs. The decomposition ends up with that residual signal. The final residue is defined as a signal with no oscillation or the signal with constant amplitude close to zero. It represents the mean trend of the signal. The main advantage of this time-space filtering is that the results preserve full meaning, nonlinearity, and nonstationary in physical space. In this paper, the low frequency annual cycle is to be separated by applying the bandpass filtering. The rest is to select the criteria to select the values of y1 and y2. The linear model of IMF’s energy of fGn and confidence limits are used to index of the IMFs for the band pass filtering. The three subband signals of daily evaporation data derived from the BEMD shown in Figure 1 with (2.4) are illustrated in Figure 4; the values of y1 and y2 are assigned as 4 and 8, respectively.Figure 4
Three subband signals obtained by applying BEMD on daily evaporation data.


## 3. Analysis and Results

BEMD  is  an  adaptive  method  to  decompose  any  complex/bivariate signal  into  a  finite set  of complex intrinsic mode  function  (IMF)  components,  which  become  the  basis functions representing  the data.  As the bases are fully data adaptive, it usually offers a physically meaningful representation of the underlying processes. Also because of  the adaptive  nature of  the  bases,  there  is  no need  for  harmonics consideration;  therefore,  BEMD  is  ideally  suited  for  analyzing  data  from  nonstationary and nonlinear processes.  Even with these nice properties, BEMD still cannot resolve signal from noise in the most complicated cases, when the processes are nonlinear and  the  noise  also  has  the  same  time  scale  as  the  signal;  their  separation  becomes difficult.  Nevertheless, BEMD offers a totally different approach to data decomposition, and we apply it to the study of the characteristics of white noise.  We will show that with this approach, we can offer some measure of the information content of climate signals.


### 3.1. Data for CVC Analysis

The Indian Statistical Institute (ISI), Kolkata, India has a large database of weather parameters. The daily evaporation, maximum temperature and minimum temperature, data are collected from that database. The study area is the place near Giridhi of Bihar, India and the data are acquired for the period of 1989 to 2004. All the weather parameters are measured with instruments (without any touch of human factors) specified by Indian Meteorological Department, Government. of India. Climatologically, the area under study is located in the tropical Indian monsoon region. The climate of the area before the monsoon is characterized by a hot summer seasons; this is called premonsoon season. However, in early March, the area also experiences the impact of western disturbances. The data are collected from a single weather station and hence it is very local for a specific area. The (uniform) sampling period of the data is one day that is, the weather parameters are acquired once a day. No global data is used in this study.


### 3.2. Effectiveness of BEMD in CVC Analysis

The climate variables are always treated as discrete time series. The nonlinear trend, nonstationary nature, irregular frequencies, and amplitude modulation which are inherently present in climate signals. Such types of properties of climate signal requires a data adaptive signal analysis tool which is suitable to analyze the time series of nonlinear and nonstationary systems. The annual cycle (AC) of climate signal is effectively separated using BEMD-based time domain filtering. Fourier transform represents any stationary signal as a combination of predefined basis functions (also called harmonics), that is, the analyzing signal is fitted with a set of harmonics. Sometimes it makes a distortion of the original signal to properly fit with the bases. The wavelet transform is a slightly improved version to analyze nonstationary signals. It is also troublesome to select the basis wavelet depending on the analyzing signal. To overcome the mentioned limitations, BEMD is applied here to implement data adaptive time domain filtering. It does not require any predefined basis vector. It derives the bases directly from the data and hence it is called data adaptive approach.


### 3.3. Significance of Mean Periods of IMFs

Before examining any results, it is necessary to list the properties of an IMF as follows:  an IMF is any function having symmetric envelops defined by the local maxima and minima separately, and also having the same number of zero-crossing and extrema.  Based on this definition, we can determine the mean period of the function by counting the number of the peaks (local maxima) of the function. The mean periodτm of any IMF is calculated as(3.1)τm=1L-1∑l=1Lρl′,
where L is the total number of maxima, ρl′ is the sample length between lth and (l+1)th maxima, ρ′=∂ρ/∂n, and ρ is the vector containing the sample indices of maxima. It is noted that each sample represents the weather data of a specific day, that is, t indicates a day.  The mean average periods (in terms of the number of days) of the IMFs of climate signal (daily evaporation) and fGn are shown in Figure 5.Figure 5
The mean periods of the IMFs of daily evaporation and fGn.An interesting property of EMD is that the averaged period of any IMF component almost exactly doubles that of the previous one, suggesting that the EMD is a dyadic filter.  This finding is consistent with the recent result by Flandrin et al. [33]. Such property of EMD on climate signals shown in Tables 1, 2, and 3 justified the evidence of multiband analysis. The second row of the table is the total number of local maxima, the peaks of individual IMF over the entire data points. The third row is the averaged period measured in terms of days.Table 1
Period statistics of daily evaporation data.
IMFs1234567891011No. of maxima1992114865636620911763291583Mean period (day)2.804.908.5015.3026.9047.9189.23189.12364.75773.341579.50Table 2
Period statistics of daily maximum temperature data.
IMFs1234567891011No. of maxima1980110460633218210447251473Mean period (day)2.815.139.2016.9430.7253.23118.67225.81400.70767.511363.1Table 3
Period statistics of daily minimum temperature data.
IMFs1234567891011No. of maxima2045116562834419810351221374Mean period (day)2.744.828.9016.2228.2354.41108.85238.0446.26816.321409.3


### 3.4. Analysis of Annual Cycles

The annual cycle (AC) is an important factor to take into consideration, while the climate variability is to be analyzed [34]. Unlike more traditional methods, where the twelve-month climatology is subtracted from the data and only the anomalies are studied, this method decomposes the total climate signal, not just the anomalies. In fact, the removal of climatology, which is a linear operation, can actually degrade the nonlinear analysis. The annual cycle is an important component and is retrieved here by partial reconstruction of the IMFs in the decomposition. Usually the modes contain a mixture of frequencies and these mixed modes are much more difficult to interpret. Since we are not interested in intraseasonal variations, a minimal amount of presmoothing is performed so that single or multiple modes collectively contain the annual cycle. The IMFs will generally be ordered from high to low frequency and the last IMF often contains a trend. It is important to determine the significance of the IMFs. We should not expect all modes to be significant.


#### 3.4.1. Trend Extraction

The effect of annual cycle is considered as the low frequency trend of climate signals. The trends of the recorded climate signals are detected using the energy distribution of the signal over the individual IMF. The analyzing signals(n) consisting of a slowly varying trend superimposed to a fluctuating process x(n); the trend is captured by IMFs of large indices including the final residue [35]. Detrending s(n), which corresponds to estimate x(n), may therefore amount to capturing the partial reconstruction
(3.2)x̃C(n)=∑k=1Cck(n),
where C is the larger IMF index prior contamination by the trend. Each of the IMF {ck(n);k=1,2,…,C} being close to zero-mean. A rule of thumb for choosing C is to observe the evolution of the empirical mean of x̃c(n) as a function of a test order c, and to identify for which c=C it departs significantly from zero. An example of this approach is illustrated in Figure 6, where a noisy voiced speech signal (0 dB SNR) is taken into consideration.Trend detection noisy speech signal. (a) Standardized empirical mean of the fine-to-coarse EMD reconstruction, evidencingC=2 as the change point; (b) top: original signal, middle: estimated trend obtained from the partial reconstruction with IMFs 3 to 7 + the residual, bottom: noise (detrended) signal obtained from the partial reconstruction with IMFs 1 to 2.
(a)(b)The detrending method described in [35] is only suitable to extract low frequency trend signals contaminated by fGn. The IMFs’ energies of the analyzing signal are compared with fGn-based energy model. It is not assumed that the climate signals are contaminated with Gaussian noise. We are using the advantage of BEMD in which the fGn is only used as the reference signal to compare the frequency and corresponding energy as well at each IMF level. It is already statistically proven that the EMD on fGn noise acts as dyadic filterbank. The fGn is also used for statistical validation of the frequency and corresponding energy of the IMFs of climate signal. As EMD does not have any basis function to parameterize the frequency response also, we can use the IMF of fGn as like the basis function. Then we get proper frequency tracing from fGn part of BEMD (including climate signal). The IMF-dependent energy of fGn and its confidence limits are used as reference to detect the upper frequency limit of the climate signals.


#### 3.4.2. Separation of AC Using BEMD

The energy of the IMF of climate exceeds the confidence limit as it represents the trend of the signal. Some reprocessing on the raw climate signal is performed to make peaks of the annual cycle prominent. Any climate signals(n) is preprocessed as(3.3)ŝ(t)=I[ξ[s(t)]],
where ξ[·] is the function of computing the Hilbert envelope of the climate signal s(t) and  ℑ[·]  is the function of smoothing filter of the Hilbert envelope. Figure 7 shows the original climate signal s(t), its Hilbert envelope, and the smoothed Hilbert envelope.Preprocessing of the raw climate signal. The original climate signal (a), the Hilbert envelope (b), and the smoothed version of the Hilbert envelope (c).
(a)(b)(c)The BEMD is applied to the complex signalx(n)=ŝevp(n)+jηfGn(n), where ηfGn(n) is the fractional Gaussian noise and ŝevp(n) is the climate signal for daily evaporation. Both signals are normalized. After decomposition, the IMF vector has two parts: the real part corresponds to the IMFs of fGn and the imaginary parts represent the IMFs of the ŝevp(n) The IMFs’ energies of ŝevp(n) are compared with the confidence limits of that of the fGn signal. The IMF with energy exceeding the upper confidence limit is tracked as the starting of the trend that is effects of annual cycle to the analyzing climate signal. In Figure 8(a), illustrating the IMFs’ energies with confidence limits, the kth (here, k=8) IMF of ŝevp(n) is selected as the starting point of reconstructing the annual cycle (AC) effects. The AC effect is separated by adding up the IMFs of ŝevp(n) signal (i.e., the real parts of the IMF vector derived by BEMD) starting from kth one up to the residue. The separation results of AC effects from the preprocessed climate signal (daily evaporation) are shown in Figure 8(b). Similarly, the starting IMF selection and the separation of ACs from preprocessed daily maximum temperature and daily minimum temperature signals are illustrated in Figures 9 and 10, respectively. The proposed AC separation algorithm can be summarized as follows.(i)
Use climate signal  sayŝevp(n) and fractional Gaussian noise (fGn) as the real and imaginary components, respectively, of a complex signal x(n).(ii)
Decomposex(n) using BEMD. Each IMF has real and imaginary components which correspond to the IMFs of ŝevp(n) and fGn, respectively.(iii)
Compute the energy of each of imaginary IMFs and also its 99% confidence limits.(iv)
Compute the energy of each of real IMFs. Find which IMF (i.e., its energy) exceeds the upper confidence limits of the energies of fGn, that is, the staring index (k) of IMFs responsible to produce the trend.(v)
The annual cycle effect is separated by summing up the real IMFs starting fromkth one up to the residue.(a) The energy distribution of fGn over the IMFs, its upper and lower bounds with 99% confidence interval. The energies of individual IMF of daily evaporation signal indicate that the 8th one is the starting IMF to represent the annual cycle. (b) The separation of annual cycle of daily evaporation. The preprocessed climate signal (top), extracted annual cycle (middle), and the residual signal (bottom).
(a)(b)(a) The energy distribution of fGn over the IMFs, its upper and lower bounds with 99% confidence interval. The energies of individual IMF of daily max. temperature signal indicate that the 9th one is the starting IMF to represent the annual cycle. (b) The separation of annual cycle of daily maximum temperature. The preprocessed climate signal (top), extracted annual cycle (middle), and the residual signal (bottom).
(a)(b)(a) The energy distribution of fGn over the IMFs, its upper and lower bounds with 99% confidence interval. The energies of individual IMF of daily minimum temperature signal indicate that the 8th one is the starting IMF to represent the annual cycle. (b) The separation of annual cycle of daily maximum temperature. The preprocessed climate signal (top), extracted annual cycle (middle), and the residual signal (bottom).
(a)(b)It is difficult to implement the-model based approach described in [35] to separate AC from recorded climate signals using univariate EMD (UEMD). If the reference signal (fGn) and the climate signal are decomposed separately using UEMD, it is difficult to synchronize IMFs corresponding to fGn and climate signals in terms of energy and spectral characteristics. Hence the newly developed BEMD is effectively applied here for separating the annual cycle from the raw climate signals.There is another issue of analyzing the annual cycle (AC) of different climate signals as a function of climate variability and change. When there is a climate change, there is a  good chance of nonstationary of the interannual distances of the climate signals. The interannual distances of the three climate signals (daily evaporation, maximum, and minimum temperature) are shown in Figure11. The cross-correlation among the annual cycles of the three climate signals is shown in Table 4. It is found that the ACs are not fully correlated and hence there are some changes in climate.Table 4
Correlation Coefficients among the annual cycles (ACs) of three climate signals.
AC of evaporationAC of max temperatureAC of min temperatureAC of evaporation10.742080.34837AC of max temperature0.7420810.75252AC of min temperature0.348370.752521Figure 11
The interannual distances computed from the extracted annual cycles of  of three climate signals.It is expected that the interannual distance is close to 360 days, whereas such distances for the mentioned three climate signals are varied up to two months (as illustrated in Figure11). So there is much amount of climate change with those climate signals.


## 3.1. Data for CVC Analysis

The Indian Statistical Institute (ISI), Kolkata, India has a large database of weather parameters. The daily evaporation, maximum temperature and minimum temperature, data are collected from that database. The study area is the place near Giridhi of Bihar, India and the data are acquired for the period of 1989 to 2004. All the weather parameters are measured with instruments (without any touch of human factors) specified by Indian Meteorological Department, Government. of India. Climatologically, the area under study is located in the tropical Indian monsoon region. The climate of the area before the monsoon is characterized by a hot summer seasons; this is called premonsoon season. However, in early March, the area also experiences the impact of western disturbances. The data are collected from a single weather station and hence it is very local for a specific area. The (uniform) sampling period of the data is one day that is, the weather parameters are acquired once a day. No global data is used in this study.


## 3.2. Effectiveness of BEMD in CVC Analysis

The climate variables are always treated as discrete time series. The nonlinear trend, nonstationary nature, irregular frequencies, and amplitude modulation which are inherently present in climate signals. Such types of properties of climate signal requires a data adaptive signal analysis tool which is suitable to analyze the time series of nonlinear and nonstationary systems. The annual cycle (AC) of climate signal is effectively separated using BEMD-based time domain filtering. Fourier transform represents any stationary signal as a combination of predefined basis functions (also called harmonics), that is, the analyzing signal is fitted with a set of harmonics. Sometimes it makes a distortion of the original signal to properly fit with the bases. The wavelet transform is a slightly improved version to analyze nonstationary signals. It is also troublesome to select the basis wavelet depending on the analyzing signal. To overcome the mentioned limitations, BEMD is applied here to implement data adaptive time domain filtering. It does not require any predefined basis vector. It derives the bases directly from the data and hence it is called data adaptive approach.


## 3.3. Significance of Mean Periods of IMFs

Before examining any results, it is necessary to list the properties of an IMF as follows:  an IMF is any function having symmetric envelops defined by the local maxima and minima separately, and also having the same number of zero-crossing and extrema.  Based on this definition, we can determine the mean period of the function by counting the number of the peaks (local maxima) of the function. The mean periodτm of any IMF is calculated as(3.1)τm=1L-1∑l=1Lρl′,
where L is the total number of maxima, ρl′ is the sample length between lth and (l+1)th maxima, ρ′=∂ρ/∂n, and ρ is the vector containing the sample indices of maxima. It is noted that each sample represents the weather data of a specific day, that is, t indicates a day.  The mean average periods (in terms of the number of days) of the IMFs of climate signal (daily evaporation) and fGn are shown in Figure 5.Figure 5
The mean periods of the IMFs of daily evaporation and fGn.An interesting property of EMD is that the averaged period of any IMF component almost exactly doubles that of the previous one, suggesting that the EMD is a dyadic filter.  This finding is consistent with the recent result by Flandrin et al. [33]. Such property of EMD on climate signals shown in Tables 1, 2, and 3 justified the evidence of multiband analysis. The second row of the table is the total number of local maxima, the peaks of individual IMF over the entire data points. The third row is the averaged period measured in terms of days.Table 1
Period statistics of daily evaporation data.
IMFs1234567891011No. of maxima1992114865636620911763291583Mean period (day)2.804.908.5015.3026.9047.9189.23189.12364.75773.341579.50Table 2
Period statistics of daily maximum temperature data.
IMFs1234567891011No. of maxima1980110460633218210447251473Mean period (day)2.815.139.2016.9430.7253.23118.67225.81400.70767.511363.1Table 3
Period statistics of daily minimum temperature data.
IMFs1234567891011No. of maxima2045116562834419810351221374Mean period (day)2.744.828.9016.2228.2354.41108.85238.0446.26816.321409.3


## 3.4. Analysis of Annual Cycles

The annual cycle (AC) is an important factor to take into consideration, while the climate variability is to be analyzed [34]. Unlike more traditional methods, where the twelve-month climatology is subtracted from the data and only the anomalies are studied, this method decomposes the total climate signal, not just the anomalies. In fact, the removal of climatology, which is a linear operation, can actually degrade the nonlinear analysis. The annual cycle is an important component and is retrieved here by partial reconstruction of the IMFs in the decomposition. Usually the modes contain a mixture of frequencies and these mixed modes are much more difficult to interpret. Since we are not interested in intraseasonal variations, a minimal amount of presmoothing is performed so that single or multiple modes collectively contain the annual cycle. The IMFs will generally be ordered from high to low frequency and the last IMF often contains a trend. It is important to determine the significance of the IMFs. We should not expect all modes to be significant.


### 3.4.1. Trend Extraction

The effect of annual cycle is considered as the low frequency trend of climate signals. The trends of the recorded climate signals are detected using the energy distribution of the signal over the individual IMF. The analyzing signals(n) consisting of a slowly varying trend superimposed to a fluctuating process x(n); the trend is captured by IMFs of large indices including the final residue [35]. Detrending s(n), which corresponds to estimate x(n), may therefore amount to capturing the partial reconstruction
(3.2)x̃C(n)=∑k=1Cck(n),
where C is the larger IMF index prior contamination by the trend. Each of the IMF {ck(n);k=1,2,…,C} being close to zero-mean. A rule of thumb for choosing C is to observe the evolution of the empirical mean of x̃c(n) as a function of a test order c, and to identify for which c=C it departs significantly from zero. An example of this approach is illustrated in Figure 6, where a noisy voiced speech signal (0 dB SNR) is taken into consideration.Trend detection noisy speech signal. (a) Standardized empirical mean of the fine-to-coarse EMD reconstruction, evidencingC=2 as the change point; (b) top: original signal, middle: estimated trend obtained from the partial reconstruction with IMFs 3 to 7 + the residual, bottom: noise (detrended) signal obtained from the partial reconstruction with IMFs 1 to 2.
(a)(b)The detrending method described in [35] is only suitable to extract low frequency trend signals contaminated by fGn. The IMFs’ energies of the analyzing signal are compared with fGn-based energy model. It is not assumed that the climate signals are contaminated with Gaussian noise. We are using the advantage of BEMD in which the fGn is only used as the reference signal to compare the frequency and corresponding energy as well at each IMF level. It is already statistically proven that the EMD on fGn noise acts as dyadic filterbank. The fGn is also used for statistical validation of the frequency and corresponding energy of the IMFs of climate signal. As EMD does not have any basis function to parameterize the frequency response also, we can use the IMF of fGn as like the basis function. Then we get proper frequency tracing from fGn part of BEMD (including climate signal). The IMF-dependent energy of fGn and its confidence limits are used as reference to detect the upper frequency limit of the climate signals.


### 3.4.2. Separation of AC Using BEMD

The energy of the IMF of climate exceeds the confidence limit as it represents the trend of the signal. Some reprocessing on the raw climate signal is performed to make peaks of the annual cycle prominent. Any climate signals(n) is preprocessed as(3.3)ŝ(t)=I[ξ[s(t)]],
where ξ[·] is the function of computing the Hilbert envelope of the climate signal s(t) and  ℑ[·]  is the function of smoothing filter of the Hilbert envelope. Figure 7 shows the original climate signal s(t), its Hilbert envelope, and the smoothed Hilbert envelope.Preprocessing of the raw climate signal. The original climate signal (a), the Hilbert envelope (b), and the smoothed version of the Hilbert envelope (c).
(a)(b)(c)The BEMD is applied to the complex signalx(n)=ŝevp(n)+jηfGn(n), where ηfGn(n) is the fractional Gaussian noise and ŝevp(n) is the climate signal for daily evaporation. Both signals are normalized. After decomposition, the IMF vector has two parts: the real part corresponds to the IMFs of fGn and the imaginary parts represent the IMFs of the ŝevp(n) The IMFs’ energies of ŝevp(n) are compared with the confidence limits of that of the fGn signal. The IMF with energy exceeding the upper confidence limit is tracked as the starting of the trend that is effects of annual cycle to the analyzing climate signal. In Figure 8(a), illustrating the IMFs’ energies with confidence limits, the kth (here, k=8) IMF of ŝevp(n) is selected as the starting point of reconstructing the annual cycle (AC) effects. The AC effect is separated by adding up the IMFs of ŝevp(n) signal (i.e., the real parts of the IMF vector derived by BEMD) starting from kth one up to the residue. The separation results of AC effects from the preprocessed climate signal (daily evaporation) are shown in Figure 8(b). Similarly, the starting IMF selection and the separation of ACs from preprocessed daily maximum temperature and daily minimum temperature signals are illustrated in Figures 9 and 10, respectively. The proposed AC separation algorithm can be summarized as follows.(i)
Use climate signal  sayŝevp(n) and fractional Gaussian noise (fGn) as the real and imaginary components, respectively, of a complex signal x(n).(ii)
Decomposex(n) using BEMD. Each IMF has real and imaginary components which correspond to the IMFs of ŝevp(n) and fGn, respectively.(iii)
Compute the energy of each of imaginary IMFs and also its 99% confidence limits.(iv)
Compute the energy of each of real IMFs. Find which IMF (i.e., its energy) exceeds the upper confidence limits of the energies of fGn, that is, the staring index (k) of IMFs responsible to produce the trend.(v)
The annual cycle effect is separated by summing up the real IMFs starting fromkth one up to the residue.(a) The energy distribution of fGn over the IMFs, its upper and lower bounds with 99% confidence interval. The energies of individual IMF of daily evaporation signal indicate that the 8th one is the starting IMF to represent the annual cycle. (b) The separation of annual cycle of daily evaporation. The preprocessed climate signal (top), extracted annual cycle (middle), and the residual signal (bottom).
(a)(b)(a) The energy distribution of fGn over the IMFs, its upper and lower bounds with 99% confidence interval. The energies of individual IMF of daily max. temperature signal indicate that the 9th one is the starting IMF to represent the annual cycle. (b) The separation of annual cycle of daily maximum temperature. The preprocessed climate signal (top), extracted annual cycle (middle), and the residual signal (bottom).
(a)(b)(a) The energy distribution of fGn over the IMFs, its upper and lower bounds with 99% confidence interval. The energies of individual IMF of daily minimum temperature signal indicate that the 8th one is the starting IMF to represent the annual cycle. (b) The separation of annual cycle of daily maximum temperature. The preprocessed climate signal (top), extracted annual cycle (middle), and the residual signal (bottom).
(a)(b)It is difficult to implement the-model based approach described in [35] to separate AC from recorded climate signals using univariate EMD (UEMD). If the reference signal (fGn) and the climate signal are decomposed separately using UEMD, it is difficult to synchronize IMFs corresponding to fGn and climate signals in terms of energy and spectral characteristics. Hence the newly developed BEMD is effectively applied here for separating the annual cycle from the raw climate signals.There is another issue of analyzing the annual cycle (AC) of different climate signals as a function of climate variability and change. When there is a climate change, there is a  good chance of nonstationary of the interannual distances of the climate signals. The interannual distances of the three climate signals (daily evaporation, maximum, and minimum temperature) are shown in Figure11. The cross-correlation among the annual cycles of the three climate signals is shown in Table 4. It is found that the ACs are not fully correlated and hence there are some changes in climate.Table 4
Correlation Coefficients among the annual cycles (ACs) of three climate signals.
AC of evaporationAC of max temperatureAC of min temperatureAC of evaporation10.742080.34837AC of max temperature0.7420810.75252AC of min temperature0.348370.752521Figure 11
The interannual distances computed from the extracted annual cycles of  of three climate signals.It is expected that the interannual distance is close to 360 days, whereas such distances for the mentioned three climate signals are varied up to two months (as illustrated in Figure11). So there is much amount of climate change with those climate signals.


## 3.4.1. Trend Extraction

The effect of annual cycle is considered as the low frequency trend of climate signals. The trends of the recorded climate signals are detected using the energy distribution of the signal over the individual IMF. The analyzing signals(n) consisting of a slowly varying trend superimposed to a fluctuating process x(n); the trend is captured by IMFs of large indices including the final residue [35]. Detrending s(n), which corresponds to estimate x(n), may therefore amount to capturing the partial reconstruction
(3.2)x̃C(n)=∑k=1Cck(n),
where C is the larger IMF index prior contamination by the trend. Each of the IMF {ck(n);k=1,2,…,C} being close to zero-mean. A rule of thumb for choosing C is to observe the evolution of the empirical mean of x̃c(n) as a function of a test order c, and to identify for which c=C it departs significantly from zero. An example of this approach is illustrated in Figure 6, where a noisy voiced speech signal (0 dB SNR) is taken into consideration.Trend detection noisy speech signal. (a) Standardized empirical mean of the fine-to-coarse EMD reconstruction, evidencingC=2 as the change point; (b) top: original signal, middle: estimated trend obtained from the partial reconstruction with IMFs 3 to 7 + the residual, bottom: noise (detrended) signal obtained from the partial reconstruction with IMFs 1 to 2.
(a)(b)The detrending method described in [35] is only suitable to extract low frequency trend signals contaminated by fGn. The IMFs’ energies of the analyzing signal are compared with fGn-based energy model. It is not assumed that the climate signals are contaminated with Gaussian noise. We are using the advantage of BEMD in which the fGn is only used as the reference signal to compare the frequency and corresponding energy as well at each IMF level. It is already statistically proven that the EMD on fGn noise acts as dyadic filterbank. The fGn is also used for statistical validation of the frequency and corresponding energy of the IMFs of climate signal. As EMD does not have any basis function to parameterize the frequency response also, we can use the IMF of fGn as like the basis function. Then we get proper frequency tracing from fGn part of BEMD (including climate signal). The IMF-dependent energy of fGn and its confidence limits are used as reference to detect the upper frequency limit of the climate signals.


## 3.4.2. Separation of AC Using BEMD

The energy of the IMF of climate exceeds the confidence limit as it represents the trend of the signal. Some reprocessing on the raw climate signal is performed to make peaks of the annual cycle prominent. Any climate signals(n) is preprocessed as(3.3)ŝ(t)=I[ξ[s(t)]],
where ξ[·] is the function of computing the Hilbert envelope of the climate signal s(t) and  ℑ[·]  is the function of smoothing filter of the Hilbert envelope. Figure 7 shows the original climate signal s(t), its Hilbert envelope, and the smoothed Hilbert envelope.Preprocessing of the raw climate signal. The original climate signal (a), the Hilbert envelope (b), and the smoothed version of the Hilbert envelope (c).
(a)(b)(c)The BEMD is applied to the complex signalx(n)=ŝevp(n)+jηfGn(n), where ηfGn(n) is the fractional Gaussian noise and ŝevp(n) is the climate signal for daily evaporation. Both signals are normalized. After decomposition, the IMF vector has two parts: the real part corresponds to the IMFs of fGn and the imaginary parts represent the IMFs of the ŝevp(n) The IMFs’ energies of ŝevp(n) are compared with the confidence limits of that of the fGn signal. The IMF with energy exceeding the upper confidence limit is tracked as the starting of the trend that is effects of annual cycle to the analyzing climate signal. In Figure 8(a), illustrating the IMFs’ energies with confidence limits, the kth (here, k=8) IMF of ŝevp(n) is selected as the starting point of reconstructing the annual cycle (AC) effects. The AC effect is separated by adding up the IMFs of ŝevp(n) signal (i.e., the real parts of the IMF vector derived by BEMD) starting from kth one up to the residue. The separation results of AC effects from the preprocessed climate signal (daily evaporation) are shown in Figure 8(b). Similarly, the starting IMF selection and the separation of ACs from preprocessed daily maximum temperature and daily minimum temperature signals are illustrated in Figures 9 and 10, respectively. The proposed AC separation algorithm can be summarized as follows.(i)
Use climate signal  sayŝevp(n) and fractional Gaussian noise (fGn) as the real and imaginary components, respectively, of a complex signal x(n).(ii)
Decomposex(n) using BEMD. Each IMF has real and imaginary components which correspond to the IMFs of ŝevp(n) and fGn, respectively.(iii)
Compute the energy of each of imaginary IMFs and also its 99% confidence limits.(iv)
Compute the energy of each of real IMFs. Find which IMF (i.e., its energy) exceeds the upper confidence limits of the energies of fGn, that is, the staring index (k) of IMFs responsible to produce the trend.(v)
The annual cycle effect is separated by summing up the real IMFs starting fromkth one up to the residue.(a) The energy distribution of fGn over the IMFs, its upper and lower bounds with 99% confidence interval. The energies of individual IMF of daily evaporation signal indicate that the 8th one is the starting IMF to represent the annual cycle. (b) The separation of annual cycle of daily evaporation. The preprocessed climate signal (top), extracted annual cycle (middle), and the residual signal (bottom).
(a)(b)(a) The energy distribution of fGn over the IMFs, its upper and lower bounds with 99% confidence interval. The energies of individual IMF of daily max. temperature signal indicate that the 9th one is the starting IMF to represent the annual cycle. (b) The separation of annual cycle of daily maximum temperature. The preprocessed climate signal (top), extracted annual cycle (middle), and the residual signal (bottom).
(a)(b)(a) The energy distribution of fGn over the IMFs, its upper and lower bounds with 99% confidence interval. The energies of individual IMF of daily minimum temperature signal indicate that the 8th one is the starting IMF to represent the annual cycle. (b) The separation of annual cycle of daily maximum temperature. The preprocessed climate signal (top), extracted annual cycle (middle), and the residual signal (bottom).
(a)(b)It is difficult to implement the-model based approach described in [35] to separate AC from recorded climate signals using univariate EMD (UEMD). If the reference signal (fGn) and the climate signal are decomposed separately using UEMD, it is difficult to synchronize IMFs corresponding to fGn and climate signals in terms of energy and spectral characteristics. Hence the newly developed BEMD is effectively applied here for separating the annual cycle from the raw climate signals.There is another issue of analyzing the annual cycle (AC) of different climate signals as a function of climate variability and change. When there is a climate change, there is a  good chance of nonstationary of the interannual distances of the climate signals. The interannual distances of the three climate signals (daily evaporation, maximum, and minimum temperature) are shown in Figure11. The cross-correlation among the annual cycles of the three climate signals is shown in Table 4. It is found that the ACs are not fully correlated and hence there are some changes in climate.Table 4
Correlation Coefficients among the annual cycles (ACs) of three climate signals.
AC of evaporationAC of max temperatureAC of min temperatureAC of evaporation10.742080.34837AC of max temperature0.7420810.75252AC of min temperature0.348370.752521Figure 11
The interannual distances computed from the extracted annual cycles of  of three climate signals.It is expected that the interannual distance is close to 360 days, whereas such distances for the mentioned three climate signals are varied up to two months (as illustrated in Figure11). So there is much amount of climate change with those climate signals.


## 4. Discussion and Conclusions

The EMD method is highly data adaptive and efficient for nonlinear and nonstationary time series. Other decompositions, for example, Fourier-based filtering and wavelet transform are very much model dependent. There are some assumptions on data to be adapted with the required model and hence there occurs the change in original properties of the analyzing data. Such types of loss or gain of climate data affect the climate analysis, greatly whereas, the EMD-based filtering, being fully data adaptive, does not cause any loss of original data.The main superiority of this method is to apply the EMD method yielding IMFs based on local properties of the signal and the instantaneous frequencies for complicated data sets. The use of BEMD makes easy to represent the nonstationary and nonlinear climate signals without considering as a collection of harmonics (as in Fourier transformation). The EMD is a new approach to many researchers in climate research. This study plays a vital role for analysis of the properties of nonlinear and nonstationary daily maximum and minimum temperature and evaporation time series data. This study focuses on the determination of climate change and variability based on three climate signals namely, daily maximum temperature, minimum temperature, and evaporation using EMD data analyzing method.The proposed BEMD is the extension of ordinary EMD to generalize as bivariate decomposition. Using EMD, it is required to decompose the fGn and the climate signals separately and hence different number of IMFs could be produced. The frequency correspondence is difficult between such types of heterogeneous sets of IMFs. In BEMD, the climate signal and the reference signals (i.e., fGn) are considered as the real and imaginary components and decomposed simultaneously. There is no chance to produce different numbers of IMFs for one climate signal and the corresponding fGn. Hence it performs better and more appropriate than the ordinary EMD for the analysis of climate signals.The IMFs partition the time series as a function of time-scale (frequency) in a statistically significant way. The residual series show that the data is overall fitted though a slight underprediction of extreme values which occurred due to small underlying trends caused climate change. Some further statistical research would be needed to address these problems. The IMFs, each carrying its own time scales, could be used in statistical prediction of future climate scenarios. The annual cycle extraction in a data adaptive way is the most significant achievement of this research. The climate variability is crystal clear from the scenario of the interannual distance. The deviation of the interannual distance from a year (360 days) illustrates the climate variability over the years.



---

*Source: 935034-2011-07-25.xml*