Hilbert-Huang transform
| Norden E. Huang et al. (2008), Scholarpedia, 3(7):2544. | doi:10.4249/scholarpedia.2544 | revision #130614 [link to/cite this article] |
The Hilbert-Huang transform (HHT) is NASA's designated name for the combination of the empirical mode decomposition (EMD) and the Hilbert spectral analysis (HSA). It is an adaptive data analysis method designed specifically for analyzing data from nonlinear and nonstationary processes. The key part of the HHT is the EMD method with which any complicated data set can be decomposed into a finite and often small number of components, called intrinsic mode functions (IMF). An IMF is defined as any function having the same (or differing at most by one) numbers of zero-crossing and extrema, and also having symmetric envelopes defined by the local maxima and minima, respectively. The definition of an IMF guarantees a well-behaved Hilbert transform of the IMF. This decomposition method operating in the time domain is adaptive, and, therefore, highly efficient. Since the decomposition is based on the local characteristic time scale of the data, it is applicable to nonlinear and non-stationary processes. With the Hilbert transform, the IMF's yield instantaneous frequencies as functions of time that give sharp identifications of imbedded structures. The final presentation of the results is an energy-frequency-time distribution, designated as the Hilbert spectrum.
HHT is an empirical approach, and has been tested and validated exhaustively but only empirically. In almost all the cases studied, HHT gives results much sharper than any of the traditional analysis methods in time-frequency-energy representation. Additionally, it reveals true physical meanings in many of the data examined.
Contents |
Instantaneous frequency and the Hilbert transform
The development of HHT is motivated by the need to describe nonlinear and nonstationary distorted waves in details. Consider a very simple nonlinear system given by the non-dissipative Duffing Equation as, \[\tag{1} \frac{\partial ^2\,x}{\partial t^2}\,\,+\,x\,\,+\,\varepsilon \,x^3\,\,=\,\,\gamma \,\cos \,\omega t\,\,, \]
where \(\varepsilon\) is a parameter not necessarily small, and \(\gamma\) is the amplitude of a periodic forcing function with a frequency \(\omega\ .\) Equation (1) can be rewritten in a slightly different form as \[\tag{2} \frac{\partial ^2\,x}{\partial t^2}\,\,+\,x\,(\,1\,+\,\varepsilon \,x^2\,)\,\,=\,\,\gamma \,\cos \,\omega t\,\,, \]
Then the quantity within the parenthesis can be interpreted as a variable spring constant, or a variable pendulum length. Equation (2) represents a system with frequency changing from location to location, and from time to time, even within one oscillation cycle. This intra-frequency frequency variation is the hallmark of nonlinear systems. In the past, this intra-wave frequency variation is often ignored or dealt with using harmonics. Harmonics distortion is a mathematical artifact resulting from imposing a linear structure on a nonlinear system. They may have mathematic meaning, but make little physical sense. The physically meaningful representation is naturally the instantaneous frequency.
The instantaneous frequency can be computed through the Hilbert Transform, with which any real valued function \(x(t)\) of \(L^{p}\) class can be transformed into an analytic function by adding a complex part, \(y(t)\ ,\) given by \[\tag{3} y(t)\,\,=\,\,\frac{1}{\pi }\,\,P\,\,\int\limits_{-\infty }^\infty {\,\frac{x(\tau )}{t-\tau }} \,\,d\tau \,\,, \]
in which the \(P\) indicates the principal value of the singular integral. With the Hilbert transform, the analytic function is then, \[\tag{4} z(t)\,\,=x\,(t)\,\,+\,\,j\,\,y(t)\,\,=\,\,a\,(t)\,e^{j\,\theta (t)}\,\,, \]
where \[\tag{5} a(t)\,\,=\,\,\left( {x^2\,\,+\,\,y^2} \right)^{1/2}\,\,\,\,;\,\,\,\,\,\,\theta (t)\,\,=\,\,\tan ^{-1}\,\frac{y}{x}\,\,\,\,. \]
Here \(a\) is the instantaneous amplitude, and \(\theta\) is the phase function; and the instantaneous frequency is simply\[\tag{6} \omega \,\,=\,\,\,-\frac{d\theta }{dt}\,\,. \]
This represents the best local fit of an amplitude and phase varying trigonometric function to \(X(t)\ .\) But the Hilbert Transform can only produce physically meaningful results for `mono-component' signals (Huang et al. 1998).
The empirical mode decomposition (EMD)
In the past, no precise definition for the mono-component function was given until the introduction of the EMD and Intrinsic Mode Function (IMF). An IMF is defined as a function that satisfies the following requirements:
- In the whole data set, the number of extrema and the number of zero-crossings must either equal or differ at most by one.
- At any point, the mean value of the envelope defined by the local maxima and the envelope defined by the local minima is zero.
Therefore, an IMF represents a simple oscillatory mode as a counterpart to the simple harmonic function, but it is much more general: instead of constant amplitude and frequency in a simple harmonic component, an IMF can have variable amplitude and frequency along the time axis.
The EMD decomposes data in terms of IMFs. The decomposing process is as follows:
- Take the test data as given in Figure 1 (the blue line).
- Identify all the local extrema.
- Connect all the local maxima by a cubic spline line as the upper envelope (the upper green line).
- Repeat the procedure for the local minima to produce the lower envelope (the lower green line).
The upper and lower envelopes should cover all the data between them as shown in Figure 1. Their mean is designated as \(m_{1}\) (the red line) also shown in Figure 1. The difference between the data and \(m_{1}\) is the first component, \(h_{1}\) shown in Figure 2, i. e., \[\tag{7} X\,(t)\,\,-\,\,m_1 \,\,=\,\,h_1 \,\,. \]
Ideally, \(h_{1}\) should satisfy the definition of an IMF, for the construction of \(h_{1}\) described above should have made it symmetric and having all maxima positive and all minima negative. After the first round of sifting, a gentle hump on a slope may become a local maximum. New extrema generated in this way actually reveal the proper modes lost in the initial examination. In the subsequent sifting process, \(h_{1}\) can only be treated as a proto-IMF. In the next step, it is treated as the data, then
\[\tag{8} h_1 \,\,-\,\,m_{11} \,\,=\,\,h_{11} \,\,. \]
After repeated sifting, shown in Figure 3, up to k times, \(h_{1}\) becomes an IMF, that is \[\tag{9} h_{1\,(k-1)} \,\,-m_{1k} \,\,=\,\,h_{1k} \,\,; \]
then, it is designated as \[\tag{10} c_1 \,\,=\,\,h_{1k} \,\,, \]
the first IMF component from the data shown in Figure 4. Historically, two different criteria have been used: The first one was used in Huang et al. (1998). The stoppage criterion is determined by a Cauchy type of convergence test, i.e. \[\tag{11} SD_k \,=\,\frac{\sum\limits_{t=0}^T {\left| {h_{k-1} (t)\,-\,h_k (t)} \right|^2} }{\sum\limits_{t=0}^T {h_{k-1}^2 (t)} }\,\,\,, \]
which stops the sifting process when it becomes smaller than a pre-determined value. A second criterion is based on the agreement of the numbers of zero-crossings and extrema. Specifically, an \(S\)-number is pre-selected. The sifting process will stop only if for \(S\) consecutive times the numbers of zero-crossings and extrema stay the same, and are equal or at most differ by one.
Once a stoppage criterion is selected, the first IMF, \(c_{1}\ ,\) can be obtained. Overall, \(c_{1}\) should contain the finest scale or the shortest period component of the signal. We can, then, separate \(c_{1}\) from the rest of the data by \[\tag{12} X\,(t)\,\,-\,\,c_1 \,\,=\,\,r_1 \,\,. \]
Since the residue, \(r_{1}\ ,\) still contains longer period variations in the data, it is treated as the new data and subjected to the same sifting process as described above. This procedure can be repeated to all the subsequent \(r_{j}\)'s, and the result is \[\tag{13} \begin{array}{l} r_1 \,\,-\,\,c_2 \,\,=\,\,r_2 \,, \\ \,\,\,\,\,\,\,\,\,... \\ r_{n-1} \,\,-\,\,c_n \,\,=\,\,r_n \\ \end{array} \]
The sifting process stops finally when the residue, \(r_{n}\ ,\) becomes a monotonic function from which no more IMF can be extracted. By summing up Equations (12) and (13), it follows that \[\tag{14} X(t) = \sum\limits_{j=1}^n {c_j +r_n.} \]
Thus, a decomposition of the data into \(n\)-empirical modes is achieved. The components of the EMD are usually physically meaningful, for the characteristic scales are defined by the physical data. Flandrin et al. (2003) and Wu and Huang (2004) have shown that the EMD is equivalent to a dyadic filter bank.
Having obtained the intrinsic mode function components, the instantaneous frequency can be computed using the Hilbert Transform [see, equations (2) to (6)]. After performing the Hilbert transform to each IMF component, the original data can be expressed as the real part, RP, in the following form: \[\tag{15} X(t) = {\rm RP}\left\{\sum_{j=1}^n a_j(t) e^{i \int \omega_j(t) dt}\right\}. \]
Equation (15) gives both amplitude and frequency of each component as functions of time. The same data if expanded in Fourier representation would be \[\tag{16} X(t) = {\rm RP}\left\{\sum_{j=1}^\infty a_j e^{i \omega_j t}\right\}, \]
with both \(a_{j}\) and \(\omega_j\) constants. The contrast between Equations (15) and (16) is clear: the IMF, to a great degree, represents a generalized Fourier expansion. This frequency-time distribution of the amplitude is designated as the Hilbert amplitude spectrum, \(H(\omega, t)\ ,\) or simply Hilbert Spectrum. If amplitude squared is more preferred to represent energy density, then the squared values of amplitude can be substituted to produce the Hilbert energy spectrum just as well.
With the Hilbert Spectrum defined, the marginal spectrum, \(h(\omega)\ ,\) follows as \[\tag{17} h(\omega ) = \int_{0}^{T} H(\omega ,t) dt. \]
The marginal spectrum offers a measure of total amplitude (or energy) contribution from each frequency value. It represents the cumulated amplitude over the entire data span in a probabilistic sense.
Comparisons with other methods
A summary of comparison between Fourier, Wavelet and HHT analyses is given in the following table:
| Fourier | Wavelet | Hilbert | |
|---|---|---|---|
| Basis | a priori | a priori | adaptive |
| Frequency | convolution: global, uncertainty | convolution: regional, uncertainty | differentiation: local, certainty |
| Presentation | energy-frequency | energy-time-frequency | energy-time-frequency |
| Nonlinear | no | no | yes |
| Non-stationary | no | yes | yes |
| Feature extraction | no | discrete: no, continuous: yes | yes |
| Theoretical base | theory complete | theory complete | empirical |
After the invention of the basics of HHT, some related developments have been made, such as the Normalized Hilbert Transform, the confidence limit, and the statistical significance test of the IMF's. The details of these topics will be discussed in due course.
References
- Flandrin, P., Rilling, G. and Gonçalves, P., 2003: Empirical mode decomposition as a filterbank. IEEE Signal Proc Lett. 11 (2): 112-114.
- Huang, N. E., Long, S. R.and Shen, Z. 1996: The mechanism for frequency downshift in nonlinear wave evolution. Adv. Appl. Mech., 32, 59-111.
- Huang, et al. 1998: The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis. Proc. Roy. Soc. Lond., 454, 903-993.
- Huang, N. E., Z. Shen, R. S. Long, 1999: A New View of Nonlinear Water Waves -- The Hilbert Spectrum, Ann. Rev. Fluid Mech. 31, 417-457.
- Huang, N. E., M. L. Wu, S. R. Long, S. S. Shen, W. D. Qu, p. Gloersen, and K. L. Fan, 2003: A confidence limit for the Empirical Mode Decomposition and Hilbert Spectral Analysis. Proceedings Royal Society of London, A459, 2,317-2,345.
- Wu, Z. and N. E. Huang, 2004: A study of the characteristics of white noise using the empirical mode decomposition method, Proceedings Royal Society of London, A,460,1597-1611.
Internal references
- Takashi Kanamaru (2008) Duffing oscillator. Scholarpedia, 3(3):6327.
- Jeff Moehlis, Kresimir Josic, Eric T. Shea-Brown (2006) Periodic orbit. Scholarpedia, 1(7):1358.
Acknowledgements
N.E.H. was supported in part by a Chair at the National Central University endowed by Taiwan Semiconductor Manufacturing Company, Ltd., and by the National Research Council of Taiwan under grant NSC 95-2119-M-008-031-MY3. Z.W. was supported by the National Science Foundation (USA) under grants ATM-0342104 and ATM-0653123. S.R.L. was supported by the National Aeronautics and Space Administration, Earth Science Division, with special thanks to Dr. Eric Lindstrom of the Physical Oceanography Program. All the authors wish to express their sincere thanks to Andrew Whitford for his expert help in the final stages of the manuscript.
Recommended reading
- Norden E. Huang and Samuel S. P. Shen, editors (2005) Hilbert-Huang Transform and Its Applications. World Scientific
- Norden E. Huang and Nii O. Attoh-Okine, editors (2005) The Hilbert-Huang Transform in Engineering. Taylor and Francis
