## 2 MEMD Method  
The EMD adaptively decomposes the original signal into a finite set of signals, called intrinsic mode functions, abbreviated IMFs [9]. Each IMF represents different vibrational modes embedded in the data, and also has a localized frequency content by preventing frequency spreading because of asymmetric waveforms. Obtaining the IMFs with a frequency range from high to low frequency is an iterative sifting procedure, as follows:  
1. Compute the mean envelope \(m_{n}=(env_{n}^{up}+env_{n}^{down})/2\), where \(env_{n}^{up}\) and \(env_{n}^{down}\) are the upper and lower envelopes that pass through the maxima and minima of the input signal \(s_{n}\), respectively;
2. Set the residue \(r_{n}=s_{n}-m_{n}\) as the new input and repeat step 1. This two-step iterative process is repeated until the stopping criterion is satisfied. The resulting residue of this iteration can be regarded as the first IMF with the lowest frequency content. In standard EMD, the local mean is computed by taking an average of the upper and lower envelopes, and IMFs in turn are obtained between the local maxima and minima. For multicomponent signals, the local maximum and minimum values are not directly defined by the EMD. [PERSON] and [PERSON] [15] proposed to generate multiple \(n\)-dimensional envelopes by taking signal projections along different directions in \(n\)-dimensional spaces. Cubic spline interpolations in different directions are adopted to form multiple local envelops, and then these values are averaged. Later, the MEMD proposed was used to process the nonlinear \(n\)-dimensional time-series signals. Let the sequence \(\{V(t)\}=\{v_{1}(t),v_{2}(t),\ldots,v_{n}(t)\}\) represent a multivariate signal with \(n\) components, and \(X^{\theta_{k}}=\left\{x_{1}^{k},x_{2}^{k},\ldots,x_{n}^{k}\right\}\) denote a set of direction vectors along the directions given by angles\(\theta^{k}=\left\{\theta_{1}^{k},\theta_{2}^{k},\ldots,\theta_{(n-1)}^{k}\right\}\) on an \((n-1)\)-dimensional sphere. The main steps of the MEMD method are as follows,  
1. Take advantage of Hammersley sequence sampling method, choose a suitable pointset for sampling on an \((n-1)\)-dimensional sphere;
2. Calculate a projection set, denoted by \(\left\{p^{\theta_{k}}(t)\right\}\), of the input signal \(\left\{V(t)\right\}\) along the direction vector \(X^{\theta_{k}}\), for all \(k\), which is the whole set of direction vectors, giving \(\left\{p^{\theta_{k}}(t)\right\}\) as the set of projections;
3. Find the time instants \(\left\{t_{i}^{\theta^{k}}\right\}\) corresponding to the maxima of the set of projected signals \(\left\{p^{\theta_{k}}(t)\right\}_{k=1}^{K}\);
4. Interpolate \(\left[t_{i}^{\theta^{k}},V\left(t_{i}^{\theta^{k}}\right)\right]\) to perform multiple spline interpolation and obtain multivariate envelope curves \(\left\{e^{\theta_{k}}(t)\right\}_{k=1}^{K}\);
5. For a set of \(K\) direction vectors, the mean value \(m(t)\) of the envelope curves is calculated as \(m(t)=\frac{1}{K}\sum\limits_{k=1}^{K}e^{\theta_{k}}(t)\);
6. Extract the 'detail' \(d(t)\) using \(d(t)=v(t)-m(t)\). If the 'detail' \(d(t)\) fulfills the stopping criterion for a multivariate IMF, apply the above procedure to \(v(t)-d(t)\), otherwise apply the above procedure to \(d(t)\). Through a series of MEMD processes, like the EMD decomposition algorithm, the \(n\)-dimensional multi-signal \(\left\{V(t)\right\}=\left\{v_{1}(t),v_{2}(t),\ldots,v_{n}(t)\right\}\) is decomposed into a series of IMF components and a residue \(r(t)=V(t)-\sum\limits_{i=1}^{q}d_{i}(t)\), where \(q\) denotes the number of IMF components. The ability of MEMD which aligns to the intrinsic mode functions is demonstrated in Figure 1. We tested the MEMD method on synthetic data using a three-component signal s = [X, Y, Z] whose components are mixtures of a 16 Hz sinusoid common to all data channels, a 64 Hz and 4 Hz tone in the X component, a 32 Hz tone in the Y component and a 4 Hz tone in the Z components.