Multivariate Empirical Mode Decomposition (MEMD)

• MEMD is a data-driven, multiscale decomposition technique that generalizes empirical mode decomposition to multichannel nonstationary signals.
• It employs directional projection, envelope interpolation, and iterative sifting to extract scale-aligned intrinsic mode functions across all channels.
• MEMD is applied in fields like neuroscience and flow analysis to improve feature extraction and classification accuracy over traditional univariate approaches.

Multivariate Empirical Mode Decomposition (MEMD) is a fully data-driven, multiscale decomposition technique that generalizes Empirical Mode Decomposition (EMD) from univariate to multichannel (vector-valued) nonstationary signals. MEMD produces a set of joint, scale-aligned Intrinsic Mode Functions (IMFs), thereby enabling adaptive, mode-synchronous time-frequency analysis, especially for applications where oscillatory modes manifest coherently across multiple measurement channels. The algorithm’s key operations are based on directional projection, envelope interpolation, and iterative sifting, supporting multidimensional instantaneous frequency analysis and facilitating robust feature extraction in multivariate settings (Islam et al., 2022, Islam et al., 2022, Eriksen et al., 2022, You, 2021).

1. Mathematical Foundations and Algorithmic Structure

MEMD is defined for a $d$-dimensional signal $$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$. The decomposition takes the form:

$$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$

where $\mathbf{c}_j(t) \in \mathbb{R}^d$ are the multivariate IMFs and $\mathbf{r}(t)$ is the multivariate residue (You, 2021, Islam et al., 2022, Islam et al., 2022, Eriksen et al., 2022).

The canonical MEMD algorithm proceeds as follows:

1. Direction Vector Generation: Sample $K$ unit vectors $\{\mathbf{v}_k\}\subset S^{d-1}$ quasi-uniformly over the $(d-1)$-sphere (often using Hammersley or Halton sequences) (You, 2021, Eriksen et al., 2022).
2. Projection: For each direction, the signal is projected: $p_k(t) = \mathbf{v}_k^\top\mathbf{X}(t)$.
3. Extremum Identification: For each projected $p_k(t)$, find times of local maxima $$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$0.
4. Envelope Construction: At each $$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$1, record $$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$2 and, for each channel, interpolate multivariate envelopes (typically with cubic splines) (Islam et al., 2022, Islam et al., 2022, Eriksen et al., 2022).
5. Computation of Local Mean: Compute the local mean by averaging all $$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$3 envelopes:

$$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$4

1. Sifting: The “detail” is $$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$5. This is iteratively sifted until the multivariate IMF criterion is satisfied, i.e., the local mean is near zero and extrema/zero-crossing counts differ by at most one in all projections (You, 2021, Islam et al., 2022, Islam et al., 2022).
2. Deflation: The extracted IMF is subtracted, and sifting continues on the residue until it is monotonic or contains no more oscillatory energy.

This process produces a set of $$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$6 joint IMFs that are aligned in scale and structure across all channels (Eriksen et al., 2022, You, 2021).

2. Intrinsic Mode Function Definition in the Multivariate Context

A multivariate IMF in the MEMD framework is defined by two key properties (Islam et al., 2022, Islam et al., 2022):

• The multivariate local mean, formed by projection-averaged envelopes, is approximately zero in every channel:

$$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$7

• For every projected direction, the number of extrema and zero-crossings differ at most by one.
• Stopping criteria involve normalized standard deviation thresholds or a fixed count of sifting iterations.

These constraints ensure that each IMF corresponds to a physically meaningful oscillatory mode, suitable for the application of Hilbert spectral analysis (Islam et al., 2022).

3. Mode Alignment and Comparison with Univariate EMD

A critical advantage of MEMD over univariate EMD is its ability to enforce mode alignment—i.e., oscillatory modes occupying the same spectral band appear at the same IMF index across all channels. Univariate EMD, when applied separately per channel, often results in scale-misalignment and loss of joint channel relationships (Eriksen et al., 2022, You, 2021, Islam et al., 2022).

MEMD maintains common time indices for extrema across all channels per direction, leading to robust extraction of coherent modes—essential in ambient oscillation analysis, neurophysiology, and multichannel dynamics. In contrast, mode-mixing and loss of phase relationship are endemic to univariate approaches, especially under noise.

4. Key Algorithmic Parameters and Implementation Considerations

Parameter	Typical Range	Role
Number of directions $$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$8	$$\mathbf{X}(t) = [x_1(t), x_2(t), \dots, x_d(t)]^\top$$9–$$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$0	Spheres sampling granularity; higher $$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$1 improves mode separation at higher computational cost (Eriksen et al., 2022, Gul et al., 2020)
Envelope Interpolation	Cubic spline (per channel)	Ensures smooth, physically plausible envelopes; linear can introduce errors, quartic may overfit (Eriksen et al., 2022, You, 2021)
Sifting SD threshold $$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$2	$$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$3–$$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$4	Controls IMF purity and sifting convergence (Islam et al., 2022, Islam et al., 2022, Eriksen et al., 2022)
Max sifting iterations	$$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$5–$$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$6	Limits over-sifting and resource usage

FPGA-based implementations of MEMD leverage directional projection units, fixed-point arithmetic, hardware accelerators for cubic-spline interpolation, and pipelined dataflow architectures to enable real-time throughput (e.g., $$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$7 for four channels) (Gul et al., 2020). On software platforms, scaling with channel count and time points is linear in $$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$8. Order-statistics-based filters (used in FA-MVEMD) further improve computational speed for high-dimensional applications (Souza et al., 2023).

5. Hilbert Spectral Analysis and Feature Extraction

After MEMD, channel-wise IMFs are subjected to the Hilbert transform to derive analytic signals:

$$\mathbf{X}(t) = \sum_{j=1}^M \mathbf{c}_j(t) + \mathbf{r}(t)$$9

where $\mathbf{c}_j(t) \in \mathbb{R}^d$0 is the Hilbert transform. The instantaneous amplitude $\mathbf{c}_j(t) \in \mathbb{R}^d$1, phase $\mathbf{c}_j(t) \in \mathbb{R}^d$2, and frequency $\mathbf{c}_j(t) \in \mathbb{R}^d$3 provide a joint, channel-synchronous time–frequency representation (Islam et al., 2022).

Time-marginalization yields the Marginal Hilbert Spectrum (MHS), a channel-wise spectral density:

$\mathbf{c}_j(t) \in \mathbb{R}^d$4

These features are used extensively in neuroscience (e.g., EEG-based emotion and mental state recognition) and physical flow analysis (Islam et al., 2022, Islam et al., 2022, Souza et al., 2023).

After IMF extraction, multivariate non-linear features such as Hjorth parameters, coefficient of variation, fluctuation index, entropy measures, and fractal dimension are computed per channel from the high-oscillation IMFs. These features demonstrate significantly superior classification accuracy relative to DFT or DWT-derived features (e.g., $\mathbf{c}_j(t) \in \mathbb{R}^d$5 with MEMD features vs. $\mathbf{c}_j(t) \in \mathbb{R}^d$6 for DFT in EEG-based mental state detection) (Islam et al., 2022).

6. Variants, Extensions, and Comparative Performance

Extensions encompass noise-assisted MEMD (NA-MEMD), which adds artificial noise channels to improve mode alignment in high-noise regimes (Eriksen et al., 2022, Souza et al., 2023); vector-valued EMD (VEMD), using back-projected smooth envelopes for improved accuracy in high dimensions (Huang et al., 2015); spatiotemporal IMF decomposition (STIMD) for coupled spatiotemporal structure (Hirsh et al., 2018); and fast adaptive MEMD (FA-MVEMD) with order-statistics filtering for multidimensional flow data (Souza et al., 2023).

Experimental comparisons with other mode decomposition methods highlight that MEMD provides:

• Fully data-driven, non-parametric decomposition for nonstationary, nonlinear signals.
• Mode alignment (for SNR $\mathbf{c}_j(t) \in \mathbb{R}^d$7 dB), but loss of alignment or mode mixing in lower SNR regimes; MVMD (Multivariate VMD) is more noise robust (Eriksen et al., 2022).
• Higher computational cost than VMD and especially MVMD; FA-MVEMD mitigates this with fast algorithms and noise assistance (Souza et al., 2023).
• Superior spatial and temporal support in flow analysis versus SVD/truncated DMD, with modes better localized and more physically interpretable (Souza et al., 2023, Hirsh et al., 2018).

7. Limitations, Practical Challenges, and Open Questions

MEMD is limited by empirical, non-rigorous convergence criteria, lack of uniqueness guarantees, and sensitivity to algorithmic parameters. Computational demands are significant for large $\mathbf{c}_j(t) \in \mathbb{R}^d$8. Mode mixing remains an issue, especially for signals with intermittent characteristics or under high noise, although NA-MEMD and post-sifting combination of adjacent IMFs can mitigate this to some extent (Eriksen et al., 2022, Souza et al., 2023).

Key challenges include scalability to many channels, theoretical analysis of convergence and uniqueness, robustness of alignment in pathological signals, and optimized envelope construction—spurring research into variants and accelerated implementations (Huang et al., 2015, Gul et al., 2020).

---

References

(Hirsh et al., 2018, Islam et al., 2022, Islam et al., 2022, Eriksen et al., 2022, You, 2021, Huang et al., 2015, Souza et al., 2023, Gul et al., 2020)