IMF decomposition: TS2IMF (EMD / EEMD / VMD)

In gmsp: Ground Motion Signal Processing

knitr::opts_chunk$set(collapse = TRUE, comment = "#>", eval = FALSE)

# Minimal executable example
library(gmsp)
library(data.table)
t_vec <- seq(0, 5, by = 0.01)
dt_sig <- data.table(t = t_vec, s = sin(2 * pi * 2 * t_vec) + 0.3 * sin(2 * pi * 8 * t_vec))
imf_res <- TS2IMF(dt_sig, method = "vmd", K = 4, output = "IMF")
print(imf_res)

TS2IMF() decomposes a time series into Intrinsic Mode Functions (IMFs) and a residue. A common goal is to separate oscillatory components by scale / band and then reconstruct a filtered signal by removing or keeping selected IMFs.

Canonical references: EMD (Huang et al., 1998), EEMD (Wu & Huang, 2009), VMD (Dragomiretskiy & Zosso, 2014).

What is an IMF?

In the EMD framework, a component $c(t)$ is an IMF if it satisfies the classic definition:

1. the number of zero crossings and extrema differ at most by 1, and
2. the local mean of the upper and lower envelopes is approximately zero.

This enables a meaningful instantaneous frequency (via the Hilbert transform) for non-stationary signals.

Input and grouping

• .x: one signal as a data.table with canonical columns t and s.

Note: internally TS2IMF() regularizes the time grid and rebuilds a time vector ts from 0. The workflow requires t / s as the working column names.

TS2IMF() is a worker and does not detect groups. For a canonical TSL, group outside the helper:

TSL[ID == "AT",
    TS2IMF(.SD, method = "vmd", K = 6, output = "TSL"),
    by = .(RecordID, OCID),
    .SDcols = c("t", "s")]

Available engines

Parameter method selects one of three engines.

VMD (method = "vmd", default)

Calls VMDecomp::vmd(...). Returns a matrix $U \in \mathbb{R}^{N \times K}$ ($K$ modes). The residue is defined as

$$r[n] = s[n] - \sum_{k=1}^{K} U_k[n].$$

EMD (method = "emd")

Calls EMD::emd(...) with boundary and stop.rule controlling sifting. The output AUX$imf and AUX$residue are taken as-is.

EEMD (method = "eemd")

Calls hht::EEMD(...) followed by hht::EEMDCompile(...). Runs an ensemble with noise (noise.amp, noise.type, trials). Uses a tempdir() as a working folder and then compiles the averaged IMFs.

Output

If output = NULL, TS2IMF() returns a list with three elements.

| Element | Shape | Description | |---|---|---| | TSW | wide | t, signal, IMF1..IMFk, residue | | TSL | long | melt of TSW with columns t, IMF, s | | IMF | metrics | per-IMF center frequency / bandwidth / relative energy |

If output is "TSW", "TSL" or "IMF", only that component is returned.

Per-IMF metrics

For each IMF the code estimates the center frequency Fo[Hz], the bandwidth BW[Hz], and the relative energy E[k]/Eo [%]. The procedure is:

1. apply a Hann window $w[n]$,
2. compute FFT with zero-padding to a power of 2,
3. compute a power spectrum $P(f)$,
4. apply a −10 dB threshold relative to the maximum and compute Fo / BW as moments of $P(f)$ over "active" bins.

Reconstruction filters

Inside TS2IMF() a local helper .applyIMFFilter() adds a reconstructed signal signal.R when any of the following arguments is set.

• imf.remove — if character, removes those columns (IMF1, IMF2, ...); if numeric, positive values remove IMFs counted from the first mode, negative values remove IMFs counted from the last mode, and both selections are unioned. For example, c(1L, -1L, -2L) removes IMF1 plus the last two IMFs. Zero, non-finite, and out-of-range numeric values are ignored.

• remove.Fo = c(f1, f2) — removes IMFs whose frequency interval $[F_o - BW,\,F_o + BW]$ lies entirely within $[f_1, f_2]$.

• keep.Fo = c(f1, f2) — keeps only IMFs whose frequency interval $[F_o - BW,\,F_o + BW]$ lies entirely within $[f_1, f_2]$.

• keep.Residue — if TRUE, signal.R adds the residue on top of the kept IMFs.

Composing with AT2TS / VT2TS / DT2TS

In practice the typical pattern is:

1. Start from a TSL (e.g. output of AT2TS(..., output = "TSL")).
2. Filter a specific ID + OCID (e.g. ID == "DT" & OCID == "H1").
3. Run TS2IMF() on that slice, or inside a data.table grouped call.
4. Use RES$TSL[IMF == "signal.R"] as the reconstructed (filtered) signal.

AT_H1 <- TSL[ID == "AT" & OCID == "H1", .(t, s)]
RES   <- TS2IMF(AT_H1, method = "vmd", K = 6,
                keep.Fo = c(2, 15), keep.Residue = TRUE)

# Reconstructed signal (only IMFs whose [Fo-BW, Fo+BW] fits in [2, 15] Hz):
RES$TSL[IMF == "signal.R"]

References

• 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. (1998). The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis. Proceedings of the Royal Society A, 454(1971), 903–995.
• Wu, Z., & Huang, N. E. (2009). Ensemble Empirical Mode Decomposition: A Noise-Assisted Data Analysis Method. Advances in Adaptive Data Analysis, 1(1), 1–41.
• Dragomiretskiy, K., & Zosso, D. (2014). Variational Mode Decomposition. IEEE Transactions on Signal Processing, 62(3), 531–544.

gmsp documentation built on July 18, 2026, 5:07 p.m.