Next Article in Journal
A 6D Pose Estimation for Robotic Bin-Picking Using Point-Pair Features with Curvature (Cur-PPF)
Previous Article in Journal
A Population-Based Iterated Greedy Algorithm for Maximizing Sensor Network Lifetime
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A Modified Complex Variational Mode Decomposition Method for Analyzing Nonstationary Signals with the Low-Frequency Trend

College of Biomedical Engineering and Instrument Science, Zhejiang University, Hangzhou 310027, China
*
Author to whom correspondence should be addressed.
Sensors 2022, 22(5), 1801; https://doi.org/10.3390/s22051801
Submission received: 23 January 2022 / Revised: 15 February 2022 / Accepted: 21 February 2022 / Published: 24 February 2022
(This article belongs to the Section Intelligent Sensors)

Abstract

:
Complex variational mode decomposition (CVMD) has been proposed to extend the original variational mode decomposition (VMD) algorithm to analyze complex-valued data. Conventionally, CVMD divides complex-valued data into positive and negative frequency components using bandpass filters, which leads to difficulties in decomposing signals with the low-frequency trend. Moreover, both decomposition number parameters of positive and negative frequency components are required as prior knowledge in CVMD, which is difficult to satisfy in practice. This paper proposes a modified complex variational mode decomposition (MCVMD) method. First, the complex-valued data are upsampled through zero padding in the frequency domain. Second, the negative frequency component of upsampled data are shifted to be positive. Properties of analytical signals are used to get the real-valued data for standard variational mode decomposition and the complex-valued decomposition results after frequency shifting back. Compared with the conventional method, the MCVMD method gives a better decomposition of the low-frequency signal and requires less prior knowledge about the decomposition number. The equivalent filter bank structure is illustrated to analyze the behavior of MCVMD, and the MCVMD bi-directional Hilbert spectrum is provided to give the time–frequency representation. The effectiveness of the proposed algorithm is verified by both synthetic and real-world complex-valued signals.

1. Introduction

Nonstationary signal processing methods have attracted considerable interest over the last few decades. Typical methods such as short-time Fourier transform (STFT) [1], the continuous wavelet transform (CWT) [2], S transform (ST) [3], and Wigner–Ville distribution (WVD) [4] can provide a time–frequency or time–scale representation for analyzing nonstationary signals. However, these algorithms are based on specific bases. Their performance depends on the selection of basis and is restricted by the uncertainty principle [5]. In order to break the limitation of the uncertainty principle, the empirical mode decomposition (EMD) [6] was proposed by Huang et al. as a data-driven signal decomposition method. EMD can recursively decompose a nonstationary signal into data-dependent basis functions named the intrinsic mode function (IMF) [7]. Consequently, the instantaneous frequency characteristics of the nonstationary data can be directly obtained from the IMFs [7]. However, EMD lacks mathematical supports and has limitations such as mode mixing. Therefore, variational mode decomposition (VMD) [8] was proposed in 2014 as a newly developed method alternative to EMD, fully established on a mathematical framework [9,10]. It can adaptively decompose a multi-component signal into several sub-signals and is not restricted by the uncertainty principle. VMD has been proven to outperform EMD and many other time–frequency analysis techniques in many applications, such as detection, separation and fault diagnosis [8,11,12,13,14]. However, both EMD and VMD are designed for real-valued data, which do not satisfy the meet of complex-valued signal processing. Since complex-valued data are widely used nowadays [15,16,17,18], it is vital to extend existing real-valued signal decomposition approaches to complex-valued signals.
Several time–frequency decomposition techniques have been studied to deal with the nonstationary multivariate signals [19,20,21]. Bivariate EMD [22] proposed in 2007, multivariate EMD proposed [23] in 2010, and multivariate VMD [24] proposed in 2019 can be regarded as extensions of the mode decomposition methods to complex-valued data. However, this class of complex extension methods pays more attention to analyzing different channels of the data. Since complex-valued data are different from multivariate data, these approaches ignore the physical meaning of complex-valued data and lose the mutual information between the real and imaginary parts of the complex-valued data [25]. Therefore, another extended EMD method for complex-valued data named complex EMD (CEMD) [26] was developed by Toshihisa Tanaka and Danilo P. Mandic. The relationship between the Fourier spectrum and the analytical signal was fully used in CEMD. Similar ideas have been applied for extending standard VMD to the complex domain [25,27]. The complex extension of VMD (CVMD) divides the complex-valued signal into positive and negative frequency components. According to the properties of the analytical signal [28], two real-valued signals are generated from the positive and negative frequency components, respectively. Standard VMD is then applied to the two real-valued signals, and a set of sub-signals can be produced. Hilbert transform is finally used to get the complex-valued form of the sub-signals. However, this kind of approach encounters difficulties when the data contain low-frequency trends. Specifically, the low-frequency component may be around zero frequency and distributed over both positive and negative frequencies. CVMD separates it into two parts directly. Therefore, prior knowledge is required to rearrange the separated signal for reconstructing the low-frequency mode. Even if the reconstruction is successful with correct rearrangement, there is energy loss, which will be analyzed in Section 3. Another drawback in CVMD is related to the decomposition number. The decomposition number is needed as a priori parameter in VMD. Since CVMD separates the data into two parts, and the VMD algorithm should be applied twice, each of the parts’ decomposition number is required prior. Existing methods assume the number of sub-signals in the negative and the positive frequency domain is identical, while it is not always applicable in practice.
Complex-valued nonstationary data with low-frequency components are common in many areas, and there are great interests in analyzing them. For example, the low-frequency component in complex-valued time series represents a trend over time that can reveal a specific pattern [25,26]; in array processing, the low-frequency component in complex-valued array data represents the broadside direction of the incoming wave [29]; in wind forecasting, studying the low-frequency components in complex-valued wind data can help to forecast the wind profile [30]. However, there are no available methods to analyze complex-valued nonstationary data with low-frequency trends. On the one hand, VMD, one of the best novel nonstationary signal analysis methods, as well as its modified versions, is not suitable for complex-valued data. On the other hand, the existing complex extension of the VMD method, known as CVMD, has drawbacks when analyzing low-frequency signals. Moreover, no modified approaches focus on this problem to the authors’ knowledge. Therefore, proposing a modified CVMD method for analyzing complex-valued data with low-frequency trends is necessary and significant.
This paper proposes a modified complex variational mode decomposition method (MCVMD), giving a new and simple way to extend VMD to the complex domain, suitable for analyzing complex-valued nonstationary data with low-frequency trends. The rest of this paper is organized as follows. Firstly, Section 2 briefly recalls the related work about VMD and CVMD. Secondly, the proposed algorithm and its Hilbert spectrum and filter bank property are described in Section 3. Then, Section 4 provides the simulation results and analysis to verify the performance of the proposed method. At last, the conclusion is done in Section 5.

2. Related Work

2.1. Variational Mode Decomposition

VMD was first proposed by Dragomiretskiy and Zosso in 2014, a novel signal decomposition method for time–frequency analysis [8]. It can non-recursively decompose a real-valued multi-component signal f into K number of quasi-orthogonal band-limited sub-signals u k , also called band-limited intrinsic mode functions (BLIMFs). The VMD method is realized through the constrained variational problem:
min u k , ω k { k = 1 K t δ t + j π t * u k t e j ω k t 2 2 } , subject   to   k = 1 K u k = f ,
where { u k t } = u 1 , , u K are the BLIMFs to be solved, { ω k t } = ω 1 , , ω K are the center frequencies of the BLIMFs to be solved, and K is the decomposition number that is known as an a priori parameter. In Equation (1), the δ t + j π t * u k t part represents the transformation of the real-valued BLIMF u k t into its analytical signal. Then the analytical signals are shifted to the zero-frequency baseband, and finally, the effective bandwidth is obtained using the L2 norm of the time derivative. The quadratic penalty and Lagrangian multiplier are introduced to address Equation (1) by the augmented Lagrange L as follows:
L u k , ω k , λ = α k = 1 K t δ t + j π t * u k t e j ω k t 2 2 + f t k = 1 K u k t 2 2 + < λ t , f t k = 1 K u k t > ,
where α represents the balancing parameter of the data-fidelity constraint and λ t denotes the Lagrangian multipliers. Equation (2) is then solved through the alternate direction method of multipliers (ADMM) approach. Using the condition that the signal is real-valued, the BLIMFs in spectral domain can be solved by
u ^ k n + 1 ω = f ^ ω i k u ^ k ω + λ ^ ω / 2 1 + 2 α ω ω k 2 ,   ω 0 ,
where the ω k is calculated at the center of gravity of the corresponding mode’s power spectrum as
ω k n + 1 = 0 ω u ^ k ω 2 d ω 0 u ^ k ω 2 d ω .
The BLIMFs in time domain u k t are obtained from the real part of the inverse Fourier transform of the frequency spectra.
As the standard VMD has been a mature algorithm built as a function in MATLAB since 2020, the subsequent VMD procedure is based on the MATLAB function. Without a special note, all of the parameters, such as the way of initializing the center frequencies ω k and the parameter α , are set as default.

2.2. Complex Variational Mode Decomposition

Since the VMD is designed for real-valued data, the extension method of VMD for handling complex-valued data called CVMD was proposed. CVMD divides the complex-valued signal into positive and negative frequency components and uses the relationship between the Fourier spectrum and the analytical signal to get real signals for the standard VMD algorithm.
Let x t be the complex-valued data and X e j ω , where π ω < π is the Fourier transform of x t . The complex-valued data are separated into positive and negative components by using a band-pass filter as
H e j ω = 1 ,   0 ω < π 0 , π ω < 0 .
The following two real signals can be generated by
x + t = R { F 1 H e j ω X e j ω } ,
x t = R { F 1 H e j ω X * e j ω } ,
where * describes the complex conjugate operation, symbol R · defines the real part of the complex-valued data, and F 1 · denotes the inverse Fourier transform. According to the properties of the analytic signals, F 1 H e j ω X e j ω and F 1 H e j ω X * e j ω are analytic signals, and their real parts contain all of the information. The reconstruction of the original complex-valued signal is determined by x + t and x t as well as their Hilbert transform pairs, expressed as
x t = x + t + j x + t + x t + j x t * ,
where · denotes the Hilbert transform operator.
Since x + n and x n are real-valued signals, the corresponding BLIMFs can be obtained using standard VMD. The VMD is separately embraced to the two real-valued signals. Suppose the decomposition number of the positive and negative frequency planes is K + and K , respectively. The existing methods consider K + = K = K , so the decomposition is described as
x + t = i = 1 K x i t ,
x t = i = K 1 x i t ,
where x i t i = 1 K and x i t i = K 1 denote sets of BLIMFs corresponding to x + t and x t , respectively.
Based on the Equations (8)–(10), the CVMD for a complex-valued data x t is eventually expressed as
x t = i = K , i 0 K z i t ,
where z i t represents the i-th complex BLIMF, which can be written as
z i t = x i t + j x i t , i = 1 , , K x i t + j x i t * , i = K , , 1 .

3. Problem Statement

The VMD mathematical theory demonstrates that the algorithm is unsuitable for complex-valued signal analysis. CVMD extent VMD to complex-valued data, but here we will illustrate its limitations in decomposing the signal with a low-frequency trend. Take the complex-valued signal s 1 t = 0.5 e j 2 π 8 t + e j 2 π t + 0.6 t 2 + 0.8 e j 2 π 16 t + 0.2 t 2 for example. The blue solid line in Figure 1 shows the frequency spectrum of s 1 t . Three band-limited components can be seen distributed in different parts of the spectrum: the first component, x 1 = 0.5 e j 2 π 8 t , is mainly at negative frequency; the second component, named x 2 = e j 2 π t + 0.6 t 2 , is around zero frequency; and the third one, x 3 = 0.8 e j 2 π 16 t + 0.2 t 2 , is mainly at positive frequency. Note that the center frequency of the low-frequency signal x 2 is actually at the positive frequency, but there is energy spreading over negative frequencies.
In CVMD, the positive and negative frequency components of the signal are displayed on the positive and negative frequency planes, respectively. Both of the decomposition numbers of two planes are required. The decomposition number of the negative frequency plane is considered the same as that of the positive frequency plane. However, only the whole data’s decomposition number can be easily known in this example. It is unclear how many components are at the positive frequency plane or negative frequency plane without prior knowledge. Even if the correct decomposition number is known and the signals are rearranged well, the results by CVMD have energy loss in the components around zero frequency. The red dotted line in Figure 1 shows the reconstruction of the CVMD results. It can be seen clearly that signal x 2 is not recovered completely.
To further illustrate this problem, the behavior of CVMD as a filter bank is illustrated by applying the algorithm to the white Gaussian noise. Both the positive frequency and negative frequency components’ decomposition number is set as 6, for example. The equivalent filter bank of CVMD is shown in Figure 2, in which the parts near zero frequency, x 1 + and x 1 , are separated artificially by the filter. This separation goes against the data-driven concept of VMD. Although the low-frequency trend can be reconstructed by combining x 1 + and x 1 components, it will not be recovered ideally.

4. Modified Complex Variational Mode Decomposition

CVMD cannot reconstruct the low-frequency trend completely because it separates the data into positive and negative frequency components artificially and hurts the parts around zero frequency. A modified VMD method is proposed to overcome the defect based on upsampling and frequency shifting.

4.1. MCVMD Algorithm

For real-valued data, the positive- and negative frequency components are conjugate symmetric. For complex-valued data, its positive- and negative frequency components contain different information, which is not symmetric. The idea of the proposed MCVMD method is to convert complex-valued data into real-valued data without losing information so that standard VMD can be used. This idea is similar to CVMD and CEMD, but the operation is different. The well-known property of the analytic signal utilized in the proposed method is that the analytic signal of a real-valued signal only has positive frequency spectra but contains all of the information of the origin real-valued data. Alternatively, the real part of the analytic signal has symmetrical positive and negative frequencies and contains all of the origin complex-valued analytic signal information. According to the property, if the whole spectrum of the complex-valued data are shifted to the positive frequency, then the complex-valued data become an analytic signal. The real-valued data can be obtained by taking the real part of the analytic signal and contain all of the information of the original complex-valued data. Since the signal in practical applications is almost discrete data, the proposed method is designed for discrete signals. The steps of the algorithm are listed as follows:
Consider x n , where 0 n N 1 and n is an integer, as the discrete complex-valued data, and the sampling frequency rate is normalized as 1. Suppose the sampling rate satisfies the Nyquist sampling frequency. That is to say, the max signal frequency of the data do not exceed 1 / 2 , and the minimum signal frequency is not less than 1 / 2 .
First, the input data should be upsampled by the interpolation algorithm, among which the discrete sinc interpolation algorithm does not distort the signal as defined and is entirely reversible. However, it is often not used directly due to the boundary effects of sinc interpolation [31]. One of the simplest and most efficient methods for minimizing boundary effects is the mirror extension of the signal by half its length on each side. Accordingly, the mirror extension of the original input data are expressed as
x e = x e 1 , x e 2 , x e 3 ,
where x e 1 denotes the refection sequence of the first half of the original data as
x e 1 n = x N 2 n 1 ,   0 n N 2 1 ,
where x e 2 denotes the original sequence as
x e 2 n = x n ,   0 n N 1 ,
where x e 3 denotes the refection sequence of the last half of the original data as
x e 3 n = x N n 1 ,   0 n N 2 1 .
The extended signal can be rewritten as x e n , where 0 n 2 N 1 .
The discrete Fourier transform (DFT) is applied to the extended signal as
X e k = n = 0 2 N x e n e j 2 π 2 N k n   ,   0 k 2 N 1 ,
This represents the spectral samples. Then, 2N samples valued at zero are inserted after the first N spectral samples. A one-dimensional vector with length 4N can be generated as
X e z = X e 1 , Z , X e 2 ,
where X e 1 denotes the sequence X e k   ( 0 k N 1 ), Z denotes 2N samples valued at zero, and X e 2 denotes the sequence X e k   ( N k 2 N 1 ), which can be rewritten as X e z k ,   0 k 4 N 1 . The corresponding time–domain data can be obtained by
x e z n = 2 F 1 X e z k , 0 k 4 N 1 , 0 n 4 N 1 ,
where F 1 · denotes the inverse discrete Fourier transform.
The actual upsampled result is given by removing the mirror part as
x z n = x e z n + N , 0 n 2 N 1 .
It is known that zero padding in the frequency domain leads to an increased sampling rate in the time domain [32]. Since the data length is doubled by zero padding, the sampling rate of x z n now becomes 2. In addition, the Fourier spectrum of x z n now has definitions from the frequency of -1 to 1 and has values from the frequency of -1/2 to 1/2.
x z n is then shifted by multiplying exponential factors e j 2 π · n 4 as follows
x z s n = x z n e j 2 π · n 4 , 0 n 2 N 1 .
Since the frequency of x z n is within 1 / 2 to 1 / 2 , all of the components of x z s n are within 0 to 1 after frequency shifting, which means it has only positive frequency components. Thus, x z s n becomes an analytic signal. According to the property of the analytic signals, real parts of the analytic signal contain all of the information. A real-valued signal can be deduced as
x z s r n = R x z s n .
Since x z s r n is real valued, the BLIMFs can be obtained using standard VMD. Consider that the decomposition number is K . The mode decomposition procedure can be expressed as
x z s r n = i = 1 K u i n + r n ,
where u i n denotes the i-th BLIMF corresponding to x z s r n , and r n describes the residual.
After VMD, the Hilbert transform is applied to generate the corresponding complex and analytic i-th BLIMF as
y i n = u i n + j u i n .
Note the data are upsampled and modulated, so they are not the final result of MCVMD.
The final complex BLIMFs can be obtained by demodulation and downsampling. The i-th upsampled and modulated complex BLIMFs y i n are shifted back by multiplying exponential factors e j 2 π · n 4 as
z i n = y i n e j 2 π · n 4 ,
where z i n denotes the i-th upsampled complex BLIMFs, which has the double sample rate of the original data. The final i-th complex BLIMF q i n corresponding to the origin complex data x n can be obtained simply by downsampling.
q i n = z i 2 n + 1 , 0 n N 1 .

4.2. MCVMD Hilbert Spectrum

The mode decomposition method is often used to get the Hilbert–Huang spectrum, which gives the time–frequency representation [33]. Specifically, mode decomposition is the first step in conventional Hilbert–Huang spectrum analysis. The Hilbert transform is then applied to the real-valued decomposed modes to get the analytic signal. Analytic signals are required to analyze each mode’s instantaneous phase and amplitude [34]. Here, the complex BLIMF from MCVMD is already analytical. Thus, the instantaneous phase and instantaneous amplitude can respectively be obtained directly by the following equations
i n = a r c t a n I q i n R q i n , 1 i K ,   0 n N 1 ,
A i n = I q i n 2 + R q i n 2 , 1 i K , 0 n N 1 ,
where symbol I · defines the imaginary part of a complex function, i n represents the instantaneous phase, and A i n corresponds to the instantaneous amplitude of the i-th complex BLIMF q i n . The instantaneous frequency can be derived as
f i n = 1 2 π d i n d n .
Accordingly, the Hilbert spectrum of each complex BLIMF can be obtained as
S i f , n = A i n .
The MCVMD Hilbert spectrum intuitively shows the instantaneous frequency and amplitude of each mode in the signal, which is not constrained by the uncertainty limitations. Unlike the VMD Hilbert spectrum that only contains positive frequency, the MCVMD Hilbert spectrum can represent the positive and negative frequencies. Thus, the time-varying characteristics of low-frequency trends varying from positive frequency to negative frequency or from negative to positive can also be seen.

4.3. MCVMD Equivalent Filter Bank

CVMD has been applied to white Gaussian noise to analyze its equivalent filter bank. In this subsection, the equivalent filter bank of MCVMD is also investigated. MCVMD is applied to Gaussian complex-valued data of 1024 randomly generated samples, and the decomposition number is set as 12. Thus, 12 complex BLIMFs and their spectra can be obtained. This process is repeated 2000 times independently, and then the resulting spectra are averaged. Figure 3 shows these averaged spectra for the 12 complex BLIMFs. It can be seen CVMD equivalent filter bank structures behave almost constant bandwidths. Compared to the equivalent filter bank of CVMD shown in Figure 2, the low-frequency components are no longer separated. Each resulting structure by MCVMD is natural and complete over the whole frequency range without artificial separation. The result exhibits the desired behavior of a filter bank, which confirms the proposed MCVMD being a natural extension of the standard VMD and outperforming the CVMD.

5. Results and Discussion

5.1. Synthetic Signal Analysis

5.1.1. Example 1

In this example, the performance of MCVMD in decomposing the signal with a low-frequency component is verified. Consider the signal s 1 t = 0.5 e j 2 π 8 t + e j 2 π t + 0.6 t 2 + 0.8 e j 2 π 16 t + 0.2 t 2 mentioned in Section 3. Suppose the sampling rate is 50 Hz and the number of samples is 150. The discrete signal can be written as s 1 n = 0.5 e j 2 π 8 · 0.2 n + e j 2 π 0.2 n + 0.6 · 0.2 n 2 + 0.8 e j 2 π 16 · 0.2 n + 0.2 · 0.2 n 2 , where n = 1 , 2 , , 150 . The mode decomposition methods aim to decompose the signal into an ensemble of sub-signals to analyze the nonstationary data. Applying the MCVMD to s n , three complex BLIMFs can be obtained. The solid blue line in Figure 4 shows the spectrum of the original signal, and the red dotted line represents the spectrum of the reconstruction signal obtained by summing up these three sub-signals. Compared with the corresponding CVMD result in Figure 1, MCVMD better reconstructs the original data, especially the second components around zero frequency.
The time–domain waveforms of the three BLIMFs obtained by CVMD and MCVMD are compared in Figure 5. Let x 1 ,     x 2 ,     a n d   x 3 describe the three constituent modes 0.5 e j 2 π 8 t , e j 2 π t + 0.6 t 2 , and 0.8 e j 2 π 16 t + 0.2 t 2 , respectively. Figure 5 shows the original modes (red solid line) and the decomposition results by CVMD (green dash-dot line) and MCVMD (blue dotted line). The real part and imaginary part of each component are displayed separately. It can be seen both CVMD and MCVMD give a proper decomposition of x 1 and x 3 . However, the performance of the two methods differs in reconstructing the second mode x 2 . MCVMD gives a better decomposition and reconstruction, especially in the imaginary part of x 2 . This phenomenon is consistent with the results of the spectrum display.
The root mean square error (RMSE) of each decomposed mode is shown in Table 1, in which the real and imaginary parts are considered separately. From Figure 5, both CVMD and MCVMD yield improper values at the beginning and end of the signal. Since this problem is because of the original VMD algorithm and is not the point here, the RMSE results in this paper are given by removing the first and last five samples. Compared with CVMD, the RMSE of the decomposed modes by MCVMD is generally smaller. As can be seen, the RMSE of the imaginary part of the low-frequency mode by MCVMD is almost one-tenth of that by CVMD.

5.1.2. Example 2

To further prove the proposed method’s ability to analyze the signal’s time–frequency characteristics, we give another more complex example of MCVMD developing the bi-directional Hilbert spectrum. The synthetic signal has two modes as follows
s 2 t = c o s 0.2 t · e j 2 π t + 0.5 t 2 + 0.2 t · e j 2 π 28 t s i n π t ,
where both amplitudes and frequencies are time varying. The signal length is 5 s, and the sampling rate is 100 Hz. The amplitude of the first mode varies as the cosine function, and the frequency varies linearly between -1 Hz and 4 Hz. The amplitude of the second mode increases linearly with time, and the frequency varies as the sine function. Hilbert transform is used to get the two modes’ instantaneous amplitudes and frequencies. Figure 6a shows the actual Hilbert spectrum of the two simulation sub-signals, which is used for comparison later.
The proposed MCVMD method is applied to s 2 t . After that, we use Equations (27)–(30) to get the corresponding Hilbert spectrum shown in Figure 6b. The MCVMD Hilbert spectrum presents the time–frequency characteristics and shows the two modes’ instantaneous amplitude and frequency as the actual Hilbert spectrum. The corresponding results by CVMD are shown in Figure 6c. As can be seen, the instantaneous amplitude of the low-frequency mode by CVMD behaves differently from the actual mode.
The RMSE of each mode’s instantaneous frequency and amplitude are presented in Table 2. It shows that both the instantaneous frequency and amplitude’s RMSE by MCVMD are smaller than by CVMD, especially in the low-frequency mode. It can be concluded that MCVMD not only gives a better decomposition of the multi-component signal, but also gives a better estimation for instantaneous frequency and instantaneous amplitude with the corresponding Hilbert spectrum.

5.1.3. Example 3

In this subsection, we discuss the prior requirements of the proposed method on the number of decompositions. The results obtained in Example 2 rely on the precise knowledge of the decomposition number K . In the MCVMD method, the correct decomposition number is K = 2 for the entire frequency range. In the CVMD method, the appropriate decomposition number is K + = 3 for the positive frequency range and K = 1 for the negative frequency range. Numerical experiments are presented to find out the influence of incorrect decomposition numbers.
Considering the same simulation scenario as Example 2, we set the decomposition number K from 1 to 6 in MCVMD. Figure 7a shows the second mode is lost for the case K = 1 , and when K > 1 , all of the modes can be obtained as in Figure 7b,c. Although one component may be divided into more than one mode when K > 2 as shown in Figure 7c, the separated modes can be summed up to recover the component without loss. This separation is not like the separation by the filter in CVMD, which is natural and data driven. The RMSEs of each decomposed mode versus decomposition number are summarized in Figure 8a, in which the RMSEs of the real and imaginary parts of each mode are averaged. The performance remains superior as long as the decomposition number exceeds two. Thus, the requirement of the decomposition number for MCVMD, in this case, is no less than two. This prior requirement is not difficult to meet in practical applications.
For CVMD, there are many more corresponding decomposition schemes to be considered for the same situation, as shown in Table 3. That is because the CVMD needs to consider the positive and negative frequencies separately. In Reference [26], K and K + are suggested to be the same number. Here, if we set K and K + both equal to 2, the CVMD results are obtained in Figure 9a, which can be seen the second mode is lost. The two modes can be obtained completely if we set K = 1 and K + = 3 , as shown in Figure 9b. Numerical experiments of the 15 schemes show that modes can be obtained completely only when K 1 and K + 3 . Precisely, only 3 of the 15 schemes succeed to recover the two modes, with the RMSE less than 0.11, as shown in Figure 8b. Thus, the prior requirement of CVMD is much more challenging to meet.

5.2. Real-World Data Analysis

5.2.1. Analysis of Float Drift Data

We apply the proposed method to a real-world complex-valued data set to demonstrate its ability further. The data were taken during the Eastern Basin experiments, which can be downloaded from the World Ocean Circulation Experiment Subsurface Float Data Assembly Center (WFDAC) http://wfdac.whoi.edu, accessed on 23 January 2022 [24]. It contains the position record of a subsurface oceanographic float deployed in the North Atlantic Ocean to track the trajectory of salty water flowing from the Mediterranean Sea [35]. The looping trajectories of the chosen float are shown in Figure 10a (blue solid line). Here, the east and north displacements can be regarded as the real and imaginary part of the complex-valued data. The application of CVMD for analyzing the float drift data are illustrated in Reference [26]. Each BLIMF associated with the real and imaginary parts indicates a different physical path during vortex evolution. However, to characterize low-frequency trends, CVMD must combine the positive frequency BLIMFs and negative frequency BLIMFs. It can be seen from Reference [26], even if the signal is rearranged well, the low-frequency physical pathway given by CVMD is quite different from the actual situation, which fails to describe the potential physical characteristic.
The complex-valued input data with 548 samples are shown in Figure 11 (top row), where the real part (blue line) and the imaginary part (red line) are displayed, respectively. The proposed MCVMD is applied to the data to analyze their principal modulated oscillations. The composition number is set as 4. Thus, four modes are shown as time plots (Figure 11) and 2D graphical plots (Figure 10). The 2D plot represents the drifting float’s modulated oscillations or rotating modes. The displacement east is the real part of the resulting BLIMF, and the displacement north is the imaginary part of the resulting BLIMF. The red solid line in Figure 10a shows the low-frequency drift trajectory of the buoy, and the other three decomposed modes in Figure 10b–d show its specific rotation characteristics. We also apply CVMD to the same dataset. The yellow dotted line in Figure 10 shows the decomposed rotating modes by CVMD. It can be seen that the low-frequency drift trajectory of the buoy in Figure 10a,e by CVMD does not fully show the true characteristics. Again, that is because CVMD uses band filters and loses the low-frequency information. This example demonstrates the ability of MCVMD to analyze physically meaningful real-world data while also verifying the successful extension of VMD to complex-valued data.

5.2.2. Analysis of Wind Data

In this subsection, we apply the proposed method to the real-world wind data, whose measurements consist of wind direction φ t and wind speed v t , which can be represented as a complex variable: x t = v t e j 2 π φ t / 360   [30]. If we apply VMD to speed and direction data, respectively, the signal components’ positive and negative frequencies cannot be distinguished [25,26]. Thus, the complex extension of VMD is required to analyze the complex-valued form of wind data. We apply the CVMD and MCVMD methods to the wind data recorded by the Automated Weather Observing System (AWOS), which can be downloaded from http://mesonet.agron.iastate.edu/ASOS/, accessed on 23 January 2022. Figure 12 shows the time series of the wind speed and wind direction. The decomposition number K is set as 8 in MCVMD, and the positive decomposition number K + and negative decomposition number K in CVMD are both set as 4. Figure 13 shows the complex-valued wind data decomposition results by MCVMD and CVMD, where the real part (blue line) and the imaginary part (red line) are displayed, respectively. From top to bottom corresponds to the mode with negative to positive frequency. In CVMD, the positive and negative frequency modes are obtained separately, in which the first four are the negative frequency, and the last four are positive. Summing up all of the decomposed modes, we can reconstruct the CVMD and MCVMD results. Figure 14 compares the spectrum of the original data, the reconstruction of CVMD results, and the reconstruction of MCVMD results. It can be seen that there is energy loss in the low-frequency component by CVMD, while MCVMD gives a better reconstruction. Similarly, that is because the CVMD uses band filters to separate the original data and go against the ‘data-driven’ idea of VMD. Thus, the low-frequency decomposed mode by CVMD cannot represent the wind’s real low-frequency trend ideally. MCVMD gives better decomposition and reconstruction of the signals with low-frequency components.

6. Conclusions

This paper presents a modified complex variational mode decomposition (MCVMD) method for decomposing a complex-valued signal into an ensemble of complex band-limited intrinsic mode functions (BLIMFs). The proposed approach extends standard variational mode decomposition (VMD) to complex-valued data. In contrast to the existing complex variational mode decomposition (CVMD) method, we refrain from separating the data into positive and negative frequency planes, applying VMD twice, and rearranging the positive and negative frequency parts. Instead, we upsample and modulate the input data to have only positive frequency components to apply the standard VMD. The properties of analytical signals are used to get the final complex-valued BLIMFs. Compared with CVMD, MCVMD has the following advantages. (1) The low-frequency components can be decomposed adaptively like other frequency parts, without artificial separation. (2) The MCVMD gives better decomposition and reconstruction of complex-valued signals with lower RMSE, especially of the low-frequency component. (3) VMD is only used once, resulting in a minor computation and less prior requirement for the decomposition number. (4) It is unnecessary to rearrange the decomposition results to reconstruct the low-frequency mode. The numerical experiment illustrates that the proposed MCVMD method behaves as an ideal filter bank structure over the whole positive and negative frequency plane. Applications in synthetic and real-world complex-valued signals demonstrate the effectiveness of the proposed MCVMD method.

Author Contributions

Conceptualization, Q.M. and X.S.; Data curation, B.W.; Formal analysis, Q.M. and K.S.; Investigation, Q.S.; Methodology, Q.M. and Q.S.; Software, B.W.; Validation, K.S.; Writing—original draft, Q.M.; Writing—review and editing, X.S. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The data that support the findings of this study are available from the corresponding author upon reasonable request.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Allen, J.B.; Rabiner, L.R. A unified approach to short-time fourier analysis and synthesis. Proc. IEEE 1977, 65, 1558–1564. [Google Scholar] [CrossRef]
  2. Daubechies, I. Orthonormal bases of compactly supported wavelets. Commun. Pure Appl. Math. 1988, 41, 909–996. [Google Scholar] [CrossRef] [Green Version]
  3. Stockwell, R.G.; Mansinha, L.; Lowe, R.P. Localization of the complex spectrum: The S transform. IEEE Trans. Signal Process. 1996, 44, 998–1001. [Google Scholar] [CrossRef]
  4. Cohen, L. Time-frequency distributions-a review. Proc. IEEE 1989, 77, 941–981. [Google Scholar] [CrossRef] [Green Version]
  5. Hlawatsch, F.; Boudreaux-Bartels, G.F. Linear and quadratic time-frequency signal representations. IEEE Signal Process. Mag. 1992, 9, 21–67. [Google Scholar] [CrossRef]
  6. Huang, N.E.; Shen, Z.; Long, S.R.; Wu, M.C.; Shih, H.H.; Zheng, Q.; Yen, N.-C.; Tung, C.C.; Liu, H.H. The empirical mode decomposition and the hilbert spectrum for nonlinear and non-stationary time series analysis. Proc. R. Soc. Lond. Ser. A 1998, 454, 903–998. [Google Scholar] [CrossRef]
  7. Huang, N.E.; Shen, S.S.P. Hilbert-Huang Transform and Its Applications; Interdisciplinary Mathematical Sciences; World Scientific: Singapore, 2005; pp. 2–15. ISBN 978-981-256-376-7. [Google Scholar]
  8. Huang, N.E.; Wu, Z.; Long, S.R.; Arnold, K.C.; Chen, X.; Blank, K. On instantaneous frequency. Adv. Adapt. Data Anal. 2009, 1, 177–229. [Google Scholar] [CrossRef]
  9. Dragomiretskiy, K.; Zosso, D. Variational mode decomposition. IEEE Trans. Signal Process. 2014, 62, 531–544. [Google Scholar] [CrossRef]
  10. Nazari, M.; Sakhaei, S.M. Successive variational mode decomposition. Signal Process. 2020, 174, 107610. [Google Scholar] [CrossRef]
  11. Wang, Y.; Markert, R. Filter bank property of variational mode decomposition and its applications. Signal Process. 2016, 120, 509–521. [Google Scholar] [CrossRef]
  12. Mohanty, S.; Gupta, K.K.; Raju, K.S. Hurst based vibro-acoustic feature extraction of bearing using EMD and VMD. Measurement 2018, 117, 200–220. [Google Scholar] [CrossRef]
  13. Liu, Y.; Yang, G.; Li, M.; Yin, H. Variational mode decomposition denoising combined the detrended fluctuation analysis. Signal Process. 2016, 125, 349–364. [Google Scholar] [CrossRef]
  14. Lu, Q.; Shen, X.; Wang, X.; Li, M.; Li, J.; Zhang, M. Fault diagnosis of rolling bearing based on improved VMD and KNN. Math. Probl. Eng. 2021, 2021, e2530315. [Google Scholar] [CrossRef]
  15. Civera, M.; Surace, C. A comparative analysis of signal decomposition techniques for structural health monitoring on an experimental benchmark. Sensors 2021, 21, 1825. [Google Scholar] [CrossRef] [PubMed]
  16. Krim, H.; Viberg, M. Two decades of array signal processing research: The parametric approach. IEEE Signal Process. Mag. 1996, 13, 67–94. [Google Scholar] [CrossRef]
  17. Adali, T.; Schreier, P.J.; Scharf, L.L. Complex-valued signal processing: The proper way to deal with impropriety. IEEE Trans. Signal Process. 2011, 59, 5101–5125. [Google Scholar] [CrossRef]
  18. Magarey, J.; Kingsbury, N. Motion estimation using a complex-valued wavelet transform. IEEE Trans. Signal Process. 1998, 46, 1069–1084. [Google Scholar] [CrossRef]
  19. Yu, S.; Ma, J. Complex variational mode decomposition for slop-preserving denoising. IEEE Trans. Geosci. Remote Sens. 2018, 56, 586–597. [Google Scholar] [CrossRef]
  20. Stanković, L.; Mandić, D.; Daković, M.; Brajović, M. Time-frequency decomposition of multivariate multicomponent signals. Signal Process. 2018, 142, 468–479. [Google Scholar] [CrossRef]
  21. Sony, S. Towards Multiclass Damage Detection and Localization Using Limited Vibration Measurements. Ph.D. Thesis, Western University, London, ON, Canada, 2021. [Google Scholar]
  22. Huang, G.; Su, Y.; Kareem, A.; Liao, H. Time-frequency analysis of nonstationary process based on multivariate empirical mode decomposition. J. Eng. Mech. 2016, 142, 04015065. [Google Scholar] [CrossRef]
  23. Rilling, G.; Flandrin, P.; Goncalves, P.; Lilly, J.M. Bivariate empirical mode decomposition. IEEE Signal Process. Lett. 2007, 14, 936–939. [Google Scholar] [CrossRef] [Green Version]
  24. Rehman, N.; Mandic, D.P. Multivariate empirical mode decomposition. Proc. R. Soc. Math. Phys. Eng. Sci. 2010, 466, 1291–1302. [Google Scholar] [CrossRef]
  25. ur Rehman, N.; Aftab, H. Multivariate variational mode decomposition. IEEE Trans. Signal Process. 2019, 67, 6039–6052. [Google Scholar] [CrossRef] [Green Version]
  26. Wang, Y.; Liu, F.; Jiang, Z.; He, S.; Mo, Q. Complex variational mode decomposition for signal processing applications. Mech. Syst. Signal Process. 2017, 86, 75–85. [Google Scholar] [CrossRef]
  27. Tanaka, T.; Mandic, D.P. Complex empirical mode decomposition. IEEE Signal Process. Lett. 2007, 14, 101–104. [Google Scholar] [CrossRef]
  28. Xia, S.; Yang, J.; Cai, W.; Zhang, C.; Hua, L.; Zhou, Z. Adaptive complex variational mode decomposition for micro-motion signal processing applications. Sensors 2021, 21, 1637. [Google Scholar] [CrossRef]
  29. Doroslovački, M.I. On nontrivial analytic signals with positive instantaneous frequency. Signal Process. 2003, 83, 655–658. [Google Scholar] [CrossRef]
  30. Xu, W.; Stewart, W.K. Multiresolution-signal direction-of-arrival estimation: A wavelets approach. IEEE Signal Process. Lett. 2000, 7, 66–68. [Google Scholar] [CrossRef]
  31. Goh, S.L.; Chen, M.; Popović, D.H.; Aihara, K.; Obradovic, D.; Mandic, D.P. Complex-valued forecasting of wind profile. Renew. Energy 2006, 31, 1733–1750. [Google Scholar] [CrossRef]
  32. Yaroslavsky, L. Fast Transform Methods in Digital Signal Processing: Theory, Applications, Efficient Algorithms; Digital Signal Processing in Experimental Research; Bentham eBooks: Oak Park, IL, USA, 2011; ISBN 978-1-60805-230-1. [Google Scholar]
  33. Streamlining Digital Signal Processing: A Tricks of the Trade Guidebook|IEEE EBooks|IEEE Xplore. Available online: https://ieeexplore.ieee.org/book/6241055 (accessed on 23 January 2022).
  34. Huang, N.E.; Wu, M.L.; Qu, W.; Long, S.R.; Shen, S.S.P. Applications of Hilbert-Huang Transform to Non-Stationary Financial Time Series Analysis. Appl. Stoch. Models Bus. Ind. 2010, 19, 245–268. [Google Scholar] [CrossRef]
  35. Richardson, P.L.; Price, J.F.; Walsh, D.; Armi, L.; Schröder, M. Tracking three meddies with SOFAR floats. J. Phys. Oceanogr. 1989, 19, 371–383. [Google Scholar] [CrossRef]
Figure 1. Spectrum of the original signal (blue solid line) and the reconstruction signal by CVMD (red dash-dot line).
Figure 1. Spectrum of the original signal (blue solid line) and the reconstruction signal by CVMD (red dash-dot line).
Sensors 22 01801 g001
Figure 2. CVMD equivalent filter banks.
Figure 2. CVMD equivalent filter banks.
Sensors 22 01801 g002
Figure 3. MCVMD equivalent filter banks.
Figure 3. MCVMD equivalent filter banks.
Sensors 22 01801 g003
Figure 4. Spectrum of the original signal (blue solid line) and reconstruction signal by MCVMD (red dash-dot line).
Figure 4. Spectrum of the original signal (blue solid line) and reconstruction signal by MCVMD (red dash-dot line).
Sensors 22 01801 g004
Figure 5. The original signal (red solid line), CVMD decomposed modes (green dash-dot line), and MCVMD decomposed modes (blue dotted line). From top to bottom: the real part of x 1 , the imaginary part of x 1 , the real part of x 2 , the imaginary part of x 2 , the real part of x 3 , and the imaginary part of x 3 .
Figure 5. The original signal (red solid line), CVMD decomposed modes (green dash-dot line), and MCVMD decomposed modes (blue dotted line). From top to bottom: the real part of x 1 , the imaginary part of x 1 , the real part of x 2 , the imaginary part of x 2 , the real part of x 3 , and the imaginary part of x 3 .
Sensors 22 01801 g005
Figure 6. The Hilbert spectrum: (a) the actual Hilbert spectrum, (b) the MCVMD Hilbert spectrum, and (c) the CVMD Hilbert spectrum.
Figure 6. The Hilbert spectrum: (a) the actual Hilbert spectrum, (b) the MCVMD Hilbert spectrum, and (c) the CVMD Hilbert spectrum.
Sensors 22 01801 g006
Figure 7. Spectrum of the original signal (blue solid line) and the decomposed modes (dotted line): (a) decomposition number K = 1 , (b) decomposition number K = 2 , and (c) decomposition number K = 3 .
Figure 7. Spectrum of the original signal (blue solid line) and the decomposed modes (dotted line): (a) decomposition number K = 1 , (b) decomposition number K = 2 , and (c) decomposition number K = 3 .
Sensors 22 01801 g007
Figure 8. The RMSE of mode 1 (blue line) and mode 2 (red line) versus the decomposition number: (a) MCVMD results and (b) CVMD results.
Figure 8. The RMSE of mode 1 (blue line) and mode 2 (red line) versus the decomposition number: (a) MCVMD results and (b) CVMD results.
Sensors 22 01801 g008
Figure 9. Spectrum of the original signal (blue solid line) and the decomposition modes (dotted line): (a) decomposition number K = 2 and K + = 2 , and (b) decomposition number K = 1 and K + = 3 .
Figure 9. Spectrum of the original signal (blue solid line) and the decomposition modes (dotted line): (a) decomposition number K = 2 and K + = 2 , and (b) decomposition number K = 1 and K + = 3 .
Sensors 22 01801 g009
Figure 10. The trajectories of the original float data (blue solid line), its MCVMD results (red dotted line), and CVMD results (yellow dotted line): (a) original trace, MCVMD decomposed mode 1 and CVMD decomposed mode 1; (b) MCVMD decomposed mode 2,;(c) MCVMD decomposed mode 3; (d) MCVMD decomposed mode 4; (e) CVMD decomposed mode 2; (f) CVMD decomposed mode 3; and (g) CVMD decomposed mode 4.
Figure 10. The trajectories of the original float data (blue solid line), its MCVMD results (red dotted line), and CVMD results (yellow dotted line): (a) original trace, MCVMD decomposed mode 1 and CVMD decomposed mode 1; (b) MCVMD decomposed mode 2,;(c) MCVMD decomposed mode 3; (d) MCVMD decomposed mode 4; (e) CVMD decomposed mode 2; (f) CVMD decomposed mode 3; and (g) CVMD decomposed mode 4.
Sensors 22 01801 g010
Figure 11. The time–domain waveform of the original drifting float signal and its MCVMD decomposition. From top to bottom: original signal, decomposed mode 1, decomposed mode 2, decomposed mode 3, and decomposed mode 4.
Figure 11. The time–domain waveform of the original drifting float signal and its MCVMD decomposition. From top to bottom: original signal, decomposed mode 1, decomposed mode 2, decomposed mode 3, and decomposed mode 4.
Sensors 22 01801 g011
Figure 12. Time series of wind speed and direction.
Figure 12. Time series of wind speed and direction.
Sensors 22 01801 g012
Figure 13. The decomposition results of the complex-valued wind data: (a) MCVMD results and (b) CVMD results.
Figure 13. The decomposition results of the complex-valued wind data: (a) MCVMD results and (b) CVMD results.
Sensors 22 01801 g013
Figure 14. Spectrum of the original wind data (blue solid line), the reconstruction signal of MCVMD results (red dotted line), and the reconstruction signal of CVMD results (green dash-dot line).
Figure 14. Spectrum of the original wind data (blue solid line), the reconstruction signal of MCVMD results (red dotted line), and the reconstruction signal of CVMD results (green dash-dot line).
Sensors 22 01801 g014
Table 1. The RMSE of the decomposed modes by CVMD and MCVMD.
Table 1. The RMSE of the decomposed modes by CVMD and MCVMD.
x 1 R   1   x 1 I   x 2 R   x 2 I   x 3 R   x 3 I  
CVMD0.0200.0200.0420.3340.0120.015
MCVMD0.0100.0130.0270.0330.0050.006
1  x 1 R denotes the real part of mode x 1 , and x 1 I denotes the imaginary part.
Table 2. The RMSE of the instantaneous frequency and amplitude by CVMD and MCVMD.
Table 2. The RMSE of the instantaneous frequency and amplitude by CVMD and MCVMD.
Mode 1-F 1Mode 1-AMode 2-FMode 2-A
CVMD0.05383.39540.12020.3566
MCVMD0.01490.93340.01040.0780
1 Mode 1-F and Mode 1-A denote the instantaneous frequency and amplitude of mode 1, respectively.
Table 3. Decomposition schemes for CVMD.
Table 3. Decomposition schemes for CVMD.
K = 0   K = 1   K = 2   K = 3  
K + = 0 -(1,0)(2,0)(3,0)
K + = 1 (0,1)(1,1)(2,1)(3,1)
K + = 2 (0,2)(1,2)(2,2)(3,2)
K + = 3 (0,3)(1,3)(2,3)(3,3)
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Miao, Q.; Shu, Q.; Wu, B.; Sun, X.; Song, K. A Modified Complex Variational Mode Decomposition Method for Analyzing Nonstationary Signals with the Low-Frequency Trend. Sensors 2022, 22, 1801. https://doi.org/10.3390/s22051801

AMA Style

Miao Q, Shu Q, Wu B, Sun X, Song K. A Modified Complex Variational Mode Decomposition Method for Analyzing Nonstationary Signals with the Low-Frequency Trend. Sensors. 2022; 22(5):1801. https://doi.org/10.3390/s22051801

Chicago/Turabian Style

Miao, Qiuyan, Qingxin Shu, Bin Wu, Xinglin Sun, and Kaichen Song. 2022. "A Modified Complex Variational Mode Decomposition Method for Analyzing Nonstationary Signals with the Low-Frequency Trend" Sensors 22, no. 5: 1801. https://doi.org/10.3390/s22051801

APA Style

Miao, Q., Shu, Q., Wu, B., Sun, X., & Song, K. (2022). A Modified Complex Variational Mode Decomposition Method for Analyzing Nonstationary Signals with the Low-Frequency Trend. Sensors, 22(5), 1801. https://doi.org/10.3390/s22051801

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop