memd: Multivariate empirical mode decomposition
In PierreMasselot/emdr: Tools for EMD-regression

Description Usage Arguments Details Value References Examples

Decomposes a multivariate series through empirical mode decomposition (EMD). Includes the possibility to perform ensemble EMD (EEMD) as well as noise-assisted multivariate EMD (NA-MEMD).

memd(
  x,
  tt = 1:NROW(x),
  ndirections = 64,
  stopping = c("absmean", "S"),
  tol = c(0.075, 0.075, 0.75),
  max.iter = 50,
  max.mimfs = NULL,
  l = 0,
  Ne = 1,
  wn.power = 0.02,
  keep.noise = FALSE,
  memd.stop = ifelse(p + l > 1, 3, 2)
)

`x`	A numeric matrix containing the multivariate signal. Each column corresponds to a variable.
`tt`	A numeric vector the same size of `row(x)` containing time indices for the data. Allows for irregularly sampled signals.
`ndirections`	Integer. The number of projections necessary to compute the multivariate enveloppe. Should be at least twice the number of variables (`ncol(x)`).
`stopping`	Character indicating the stopping criterion of the sifting process. When `stopping = "absmean"` (the default) the criterion of Rilling et al. (2003) based on the mean enveloppe is used. When `stopping = "S"`, the stopping criterion of Huang et al. (2003) based on the number of iteration is used. See details.
`tol`	A numeric vector givving the tolerance for the stopping criterion. A vector of length 3 when `stopping = "absmean"` and a single value when `stopping = "S"`. See details.
`max.iter`	Integer giving the maximum number of iterations for the sifting process.
`max.mimfs`	Integer giving the maximum number of IMFs to extract. If NULL (the default), IMFs are extracted until one extremum is left in the signal.
`l`	Integer giving the number of gaussian white noise variables to add for the noise-assisted MEMD.
`Ne`	Integer giving the ensemble number for EEMD.
`wn.power`	Numeric value > 0 giving the relative standard deviation of noise variables (see details).
`keep.noise`	Logical. If set to TRUE, the white noise channels created by the algorithm are returned as well. Note that it is useful only when `l > 0`.
`memd.stop`	Numeric value. The algorithm stops and consider that the IMF has been reached when the number of extrema of the remaining signal is below `memd.stop`. Note that in the multivariate case, the number of extrema of all the projections must below this value.

The EMD algorithm iteratively estimates IMF beginning with the highest frequency to the lowest frequency. Each time an IMF is estimated, it is retrieved from the signal and the next IMF is estimated. The algorithm continues until the remaining signal contains only one (or none) extremum. This remaining signal is then considered as the trend.

Each IMF is estimated through a sifting process, which consists in fitting the envelopes of the signal on its local extrema. The local mean is then computed as the mean of the envelopes and is retrieved from the signal to yield an IMF. If this IMF is not satisfying enough, the same process is repeated, until the stopping criterion is met. See references below for the full details of the algorithm.

The stopping criterion, given in stopping, is arguably the most important parameter of the algorithm. Two criteria are currently implemented. The one of Rilling et al. (2003) (stopping = "absmean") stops the sifting process when the relative mean envelope is below a predetermined threshold. In this case, the parameter tol must be a vector of length 3. tol[1] gives the threshold and tol[2] gives the maximum proportion of the signal allowed to be above this threshold. tol[3] indicates another threshold that none of the mean enveloppe is allowed to exceed. The second criteria implemented is the S of Huang et al. (2003) (stopping = "S"). It stops the sifting process when the difference between the number of local extrema and zero-crossings is at most 1, for tol steps. In this case, tol is a single integer value.

The function also includes the possibility to perform noise-assisted extensions to manage the mode-mixing issue. Ensemble EMD (EEMD) is performed by setting the number of ensembles Ne higher than 1. NA-MEMD is performed by setting the number of added white noise variables l higher than 0. When both are given, the NA-MEMD is performed. In all cases, the argument wn.power indicates the amplitude of added white noise.

At the end, IMFs with similar frequencies can be obtained. They can be summed using the function combine.mimf.

An object of class mimf, i.e an array with dimension nindividuals x nimfs x nvariables which also contains the attributes:

`x`	The original multivariate signal.
`tt`	The vector of time indices.
`nb.iter`	A vector containing the number of sifting iterations necessary to estimate each IMF. If `Ne > 1`, a `Ne` x nimfs matrix.
`call`	The function call.

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 of London. Series A: Mathematical, Physical and Engineering Sciences 454, 903-995.

Rehman, N., Mandic, D.P., 2010. Multivariate empirical mode decomposition. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Science 466, 1291-1302.

Rehman, N.U., Park, C., Huang, N.E., Mandic, D.P., 2013. EMD Via MEMD: Multivariate Noise-Aided Computation of Standard EMD. Advances in Adaptive Data Analysis 05, 1350007.

   library(dlnm)
   
   X <- chicagoNMMAPS[,c("temp", "rhum")]
   set.seed(3)
   
   # EMD
   imfs <- memd(X[,1])
   
   # EEMD
   imfs <- memd(X[,1], Ne = 100, wn.power = .05)
   
   # MEMD
   imfs <- memd(X) # Takes a couple of minutes
   
   # NA-MEMD
   imfs <- memd(X, l = 2, wn.power = .02) # Takes a couple of minutes
   
   # Plot resulting (M)IMFs
   plot(imfs)