vmd
Decomposition by variational modes.
| Library |
|
Syntax
Function call
-
___ = vmd(x,Name,Value)— decomposes the signal into variational modes with additional parameters specified by one or more arguments of the typeName,Value.
-
vmd(___;out=:plot)— Plots the original signal, IMF, and residual signal as separate graphs in one figure.
Arguments
Input arguments
# x is a time domain signal with uniform sampling
+
vector
Details
A time domain signal with uniform sampling, defined as a vector of real values.
Name-value input arguments
Specify optional argument pairs as Name,Value, where Name — the name of the argument, and Value — the appropriate value. Name-value arguments should be placed after other arguments, but the order of the pairs does not matter.
Use commas to separate the name and value, and Name put it in quotation marks.
# AbsoluteTolerance — absolute error of convergence of modes
+
5e−6 (by default) | positive real scalar
Details
The absolute error of mode convergence, given as a positive real scalar. Argument AbsoluteTolerance It is one of the criteria for stopping optimization, that is, optimization stops when the mean quadratic absolute improvement in the convergence of the IMF in two consecutive iterations becomes less than the value AbsoluteTolerance.
# RelativeTolerance — relative error of mode convergence
+
AbsoluteTolerance*1e3 (by default) | positive real scalar
Details
The relative error of mode convergence, given as a positive real scalar. Argument RelativeTolerance It is one of the criteria for stopping optimization, that is, optimization stops when the average relative improvement in IMF convergence in two consecutive iterations becomes less than the value RelativeTolerance.
The optimization process stops when the conditions are met simultaneously. AbsoluteTolerance and RelativeTolerance.
|
# MaxIterations — maximum number of optimization iterations
+
500 (by default) | a positive integer scalar
Details
The maximum number of optimization iterations, set as a positive integer scalar. Argument MaxIterations It is one of the criteria for stopping optimization, that is, optimization stops when the number of iterations exceeds the value MaxIterations.
Argument MaxIterations It can only be specified with positive integers.
# NumIMFs — number of extracted IMFs
+
5 (by default) | a positive integer scalar
Details
The number of extracted IMFs, set as a positive integer scalar
# InitialIMFs — initial IMF
+
the zero matrix (by default) | the real matrix
Details
The initial IMFs, given as a real matrix. The rows correspond to time counts, and the columns correspond to modes.
# PenaltyFactor — penalty factor
+
1000 (by default) | positive real scalar
Details
The penalty coefficient, set as a positive real scalar. This argument determines the accuracy of the recovery. Use lower penalty factor values to achieve more rigorous data accuracy.
# InitialLM — initial Lagrange multiplier
+
the complex vector of zeros (by default) | the complex vector
Details
The initial Lagrange multiplier, defined as a complex vector. The range of the initial Lagrange multiplier in the frequency domain is [0, 0.5] cycle/countdown. The multiplier ensures that recovery limits are met. The length of the multiplier depends on the size of the input signal.
# LMUpdateRate — Lagrange multiplier update rate
+
0.01 (by default) | the real scalar
Details
The update rate of the Lagrange multiplier at each iteration, given as a positive real scalar. Higher speed leads to faster convergence, but increases the likelihood of the optimization process getting stuck in the local optimum.
# InitializeMethod is a method for initializing central frequencies
+
"peaks" (by default) | "random" | "grid"
Details
The central frequency initialization method, defined as "peaks", "random" or "grid".
-
"peaks"— initialization of the center frequencies as the peak values of the signal in the frequency domain (by default). IfNumIMFsmore than the number of peaks in the frequency domain, the functionvmdinitializes the remaining central frequencies randomly; -
"random"— initialization of the central frequencies as random numbers evenly distributed in the interval[0, 0.5]cycle/countdown; -
"grid"— initialization of the central frequencies as a uniformly sampled grid in the interval[0, 0.5]cycle/countdown.
# Display — toggles the progress display in the command window
+
false or 0 (default) | true or 1
Details
Switching the progress display in the command window, set as true (or 1) or false (or 0). If you specify true The function will display the average absolute and relative improvement of modes and center frequencies every 20 iterations, as well as show information about the final stop.
Set for the argument Display meaning 1 to display the data, or 0 to hide the data.
# out — type of output data
+
:none (by default) | :plot
Details
Type of output data:
-
:none— the function returns data; -
:plot— the function returns a graph.
Output arguments
# imf — decomposition function into internal modes
+
the matrix
Details
Decomposition functions into Intrinsic Mode Function (IMF), returned as a matrix. Each one imf It is a signal with amplitude and frequency modulation, having positive and slowly changing envelopes. Each mode has an instantaneous frequency that does not decrease, changes slowly, and is centered around a central value. Use imf to apply the Hilbert—Huang transform for spectral signal analysis.
Argument imf it is returned as a matrix, each column of which is imf if x — vector.
# info — additional diagnostic information
+
structure
Details
Additional diagnostic information returned as a structure with the following fields:
-
ExitFlag— completion flag. Meaning0indicates that the algorithm stops when the maximum number of iterations is reached. Meaning1indicates that the algorithm stops when the absolute and relative tolerances are reached.; -
CentralFrequencies— central IMF frequencies; -
NumIterations— total number of iterations; -
AbsoluteImprovement— RMS absolute improvement in IMF convergence between the last two iterations; -
RelativeImprovement— average relative improvement in IMF convergence between the last two iterations; -
LagrangeMultiplier— the Lagrange multiplier in the frequency domain at the last iteration.
Examples
Decomposition by variational modes of a piecewise signal
Details
We will generate a piecewise composite signal consisting of a quadratic trend, a chirp, and a cosine with a sharp transition between two constant frequencies at .
The signal is sampled at a frequency 1 kHz during 1 seconds. We will plot graphs for each individual component and for the composite signal.
fs = 1e3
t = 0:1/fs:1-1/fs
t_half = length(t) ÷ 2
t_first_half = t[1:t_half]
t_second_half = t[t_half+1:end]
quadratic_component = 6 .* t.^2
chirp_component = cos.(4π .* t .+ 10π .* t.^2)
cos_first_half = cos.(60π .* t_first_half)
cos_second_half = cos.(100π .* (t_second_half .- 10π))
cos_component = vcat(cos_first_half, cos_second_half)
x = quadratic_component .+ chirp_component .+ cos_component
p1 = plot(t, vcat(zeros(t_half), cos_second_half),
xlabel="Time (s)",
ylabel="Cosine",
legend=false)
p2 = plot(t, vcat(cos_first_half, zeros(length(t_second_half))),
xlabel="Time (s)",
ylabel="Cosine",
legend=false)
p3 = plot(t, chirp_component,
xlabel="Time (s)",
ylabel="Chirp",
legend=false)
p4 = plot(t, quadratic_component,
xlabel="Time (s)",
ylabel="Quadratic trend",
legend=false)
p5 = plot(t, x,
xlabel="Time (s)",
ylabel="Signal",
legend=false,
linewidth=2)
l = @layout [a{0.333w} b{0.333w} c{0.333w}
d{0.333w} e{0.667w}]
plot(p1, p2, p3, p4, p5,
layout=l,
size=(800, 600))
Let’s perform the expansion in variational modes to calculate the four decomposition functions into internal modes. Four different signal components are restored.
import EngeeDSP.Functions: vmd
imf, res = vmd(x, "NumIMFs", 4)
p1 = plot(t, imf[:, 1],
ylabel = "IMF1",
showlegend = false,
grid = true,
xlabel = "Time (s)")
p2 = plot(t, imf[:, 2],
ylabel = "IMF2",
showlegend = false,
grid = true,
xlabel = "")
p3 = plot(t, imf[:, 3],
ylabel = "IMF3",
showlegend = false,
grid = true,
xlabel = "Time (s)")
p4 = plot(t, imf[:, 4],
ylabel = "IMF4",
showlegend = false,
grid = true,
xlabel = "Time (s)")We will restore the signal by adding the mode functions and the remainder. Let’s build a graph and compare the original and restored signals.
sig = sum(imf, dims=2) + res
p5 = plot(t, sig,
linewidth = 1.5,
label = "Reconstructed signal",
xlabel = "Time (s)",
ylabel = "Signal",
legendfontsize = 7,
legend = (0.45, 0.5))
plot!(p5, t, x,
linestyle = :dash,
linewidth = 2,
label = "Original signal")
plot(p1, p2, p3, p4, p5,
layout = l,
size = (800, 600))
Let’s calculate the norm of the difference between the original and restored signals.
import EngeeDSP.Functions: norm
norm(x-sig, Inf)0.0Additional Info
Decomposition functions into internal modes
Details
Function vmd decomposes the signal by a small number narrow-band internal modes (IMF):
IMFs have the following characteristics:
-
Every fashion it is a signal with amplitude and frequency modulation of the form
where — the fashion phase, eh — its envelope.
-
The modes have positive and slowly changing envelopes.
-
Each mode has an instantaneous frequency , which does not decrease, changes slowly, and is centered around a central value .
The variational mode decomposition method simultaneously calculates all modes of oscillation and their central frequencies. The process is to find a set of values and which minimize the variational problem with constraints.
Optimization
Details
To calculate and the procedure finds the optimum of the extended Lagrangian
where is the scalar product and the second norm . The regularization member includes the following steps:
-
Using the Hilbert transform to calculate the analytical signal associated with each mode, where denotes a convolution. As a result, each mode has a purely positive spectrum.
-
Demodulating an analytical signal to the baseband by multiplying it by a complex exponential.
-
Estimation of bandwidth by calculating the square of the second gradient norm of the demodulated analytical signal.
Members and ensure compliance with the restriction by imposing a quadratic penalty and including a Lagrange multiplier. Argument PenaltyFactor measures the relative contribution of a member compared to members and .
The algorithm solves the optimization problem using the alternating direction multiplier method described in [1].
Algorithms
Function vmd calculates the IMF in the frequency domain by reconstructing using the functions . To eliminate edge effects, the algorithm expands the signal, mirroring half of its length on each side.
The Lagrange multiplier, introduced in Optimization, has a Fourier transform . The length of the Lagrange multiplier vector is equal to the length of the extended signal.
Unless otherwise specified in the argument InitialIMFs, IMFs are initialized to zero. Initialize the argument CentralFrequencies using one of the methods specified in the argument InitializeMethod. Function vmd iteratively updates mods until one of the following conditions is met:
-
and performed jointly, while and they are set using arguments
RelativeToleranceandAbsoluteToleranceaccordingly; -
The algorithm exceeds the maximum number of iterations specified in the argument.
MaxIterations.
On At the 3rd iteration, the algorithm performs the following steps:
-
He’s going through the signal mode (specified using the argument
NumIMFs) to calculate:-
Frequency domain signals for each mode using
where — the Fourier transform -th mode, calculated on - th iteration.
-
-th central frequency using
-
-
Updates the Lagrange multiplier using , where — the update rate of the Lagrange multiplier, set using the argument
LMUpdateRate.
Literature
-
Boyd, Stephen, Neal Parikh, Eric Chu, Borja Peleato, and Jonathan Eckstein. «Distributed Optimization and Statistical Learning via the Alternating Direction Method of Multipliers.» Foundations and Trends® in Machine Learning. Vol 3, Number 1, 2011, pp. 1–122.
-
Dragomiretskiy, Konstantin, and Dominique Zosso. «Variational Mode Decomposition.» IEEE® Transactions on Signal Processing. Vol. 62, Number 3, 2014, pp. 531–534.
-
Moody, George B., and Roger G. Mark. «The impact of the MIT-BIH Arrhythmia Database.» IEEE Engineering in Medicine and Biology Magazine. Vol. 20, No. 3, May–June 2001, pp. 45–50.