ドキュメンテーション

最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

emd

構文

[imf,residual] = emd(X)
[imf,residual,info] = emd(X)
[___] = emd(___,Name,Value)
emd(___)

説明

[imf,residual] = emd(X) は、X の経験的モード分解に対応する固有モード関数 imf と残差信号 residual を返します。emd を使用して、複雑な信号をヒルベルト スペクトル解析に必要な有限数の固有モード関数に分解して単純化します。

[imf,residual,info] = emd(X) は、診断の目的のために IMF の追加情報 info と残差信号を返します。

[___] = emd(___,Name,Value) では、1 つ以上の Name,Value のペア引数で指定された追加オプションを使用して、emd を推定します。

emd(___) では、同じ Figure 内で元の信号、IMF、残差信号をサブプロットとしてプロットします。

すべて折りたたむ

この例では、周波数にはっきりした変化が含まれる正弦波で構成される非定常連続信号について考えます。ジャックハンマーの振動、または花火の音は、非定常連続信号の例です。

サンプリング周波数 fs で非定常信号のデータを読み込み、混合正弦波信号を可視化します。

load('sinusoidalSignalExampleData.mat','X','fs');
t = (0:length(X)-1)/fs;
plot(t,X);
xlabel('Time(s)');

混合信号に異なる振幅と周波数値をもつ正弦波が含まれることを観察します。

ヒルベルト スペクトル プロットを作成するには、信号の IMF が必要です。経験的モード分解を実行して、固有モード関数と信号の残差を計算します。信号が滑らかではないため、'pchip' を Interpolation メソッドとして指定します。

[imf,residual,info] = emd(X,'Interpolation','pchip');
Current IMF  |  #Sift Iter  |  Relative Tol  |  Stop Criterion Hit  
      1      |        2     |     0.026352   |  SiftMaxRelativeTolerance
      2      |        2     |    0.0039573   |  SiftMaxRelativeTolerance
      3      |        1     |     0.024838   |  SiftMaxRelativeTolerance
      4      |        2     |      0.05929   |  SiftMaxRelativeTolerance
      5      |        2     |      0.11317   |  SiftMaxRelativeTolerance
      6      |        2     |      0.12599   |  SiftMaxRelativeTolerance
      7      |        2     |      0.13802   |  SiftMaxRelativeTolerance
      8      |        3     |      0.15937   |  SiftMaxRelativeTolerance
      9      |        2     |      0.15923   |  SiftMaxRelativeTolerance
The decomposition stopped because the number of extrema of the residual signal is less than 'MaxNumExtrema'.

コマンド ウィンドウ内に生成されたテーブルは、生成された各 IMF のシフト反復の数、相対許容誤差、およびシフト停止基準を示します。この情報は info にも含まれます。Display0 と指定することによってテーブルを非表示にできます。

経験的モード分解を使用して取得した imf 成分を使用してヒルベルト スペクトル プロットを作成します。

hht(imf,fs);

周波数対時間のプロットは、IMF の各点における瞬間エネルギーを示す垂直方向のカラー バーがあるスパース プロットです。このプロットは、元の混合信号から分解された各成分の瞬時周波数スペクトルを表します。プロットから、1 秒での周波数に明瞭な変化がある 3 つの IMF が観察されます。

この例では、周波数にはっきりした変化が含まれる正弦波で構成される非定常連続信号について考えます。ジャックハンマーの振動、または花火の音は、非定常連続信号の例です。

非定常信号のデータを読み込み、混合正弦波信号を可視化します。

load('sinusoidalSignalExampleData.mat','X','fs');
plot(X);

混合信号に異なる振幅と周波数値をもつ正弦波が含まれることを観察します。

経験的モード分解を実行して、固有モード関数と信号の残差をプロットします。信号が滑らかではないため、'pchip' を Interpolation メソッドとして指定します。

emd(X,'Interpolation','pchip');
Current IMF  |  #Sift Iter  |  Relative Tol  |  Stop Criterion Hit  
      1      |        2     |     0.026352   |  SiftMaxRelativeTolerance
      2      |        2     |    0.0039573   |  SiftMaxRelativeTolerance
      3      |        1     |     0.024838   |  SiftMaxRelativeTolerance
      4      |        2     |      0.05929   |  SiftMaxRelativeTolerance
      5      |        2     |      0.11317   |  SiftMaxRelativeTolerance
      6      |        2     |      0.12599   |  SiftMaxRelativeTolerance
      7      |        2     |      0.13802   |  SiftMaxRelativeTolerance
      8      |        3     |      0.15937   |  SiftMaxRelativeTolerance
      9      |        2     |      0.15923   |  SiftMaxRelativeTolerance
The decomposition stopped because the number of extrema of the residual signal is less than 'MaxNumExtrema'.

元の信号、最初の 3 つの IMF、残差信号を含む対話型プロットが生成されます。コマンド ウィンドウ内に生成されたテーブルは、生成された各 IMF のシフト反復の数、相対許容誤差、およびシフト停止基準を示します。Display0 と指定することによってテーブルを非表示にできます。

プロット内の空白を右クリックし、[IMF セレクター] ウィンドウを開きます。[IMF セレクター] を使用して、生成された IMF、元の信号、および残差を選択的に表示します。

リストから表示する IMF を選択します。プロットに元の信号と残差を表示するかどうかを選択します。

これで選択した IMF がプロットに表示されます。

プロットを使用して、残差と共に元の信号から分解された個々の成分を可視化します。残差は IMF の総数に対して計算され、[IMF セレクター] ウィンドウで選択された IMF によって変化はしないことに注意してください。

入力引数

すべて折りたたむ

一様にサンプリングされた時間領域。ベクトルまたは単一データ列 timetable として指定します。

名前と値のペアの引数

オプションの引数 Name,Value のコンマ区切りペアを指定します。Name は引数名で、Value は対応する値です。Name は引用符で囲まなければなりません。Name1,Value1,...,NameN,ValueN のように、複数の名前と値のペアの引数を、任意の順番で指定できます。

例: 'MaxNumIMF',5

コーシー型収束基準。'SiftRelativeTolerance' と正のスカラーで構成されるコンマ区切りのペアとして指定します。SiftRelativeTolerance はシフティング停止基準の 1 つです。つまり、現在の相対許容誤差が SiftRelativeTolerance 未満の場合にシフティングが停止します。

シフティング反復の最大数。'SiftMaxIterations' と正のスカラー整数で構成されるコンマ区切りのペアとして指定します。SiftMaxIterations はシフティング停止基準の 1 つです。つまり、現在の反復数が SiftMaxIterations より大きい場合にシフティングが停止します。

SiftMaxIterations は 0 または正の整数のみを使用して指定できます。

抽出された IMF の最大数。'MaxNumIMF' と正のスカラー整数で構成されるコンマ区切りのペアとして指定します。MaxNumIMF は分解停止基準の 1 つです。つまり、生成された IMF の数が MaxNumIMF に等しい場合に分解が停止します。

MaxNumIMF は 0 または正の整数のみを使用して指定できます。

残差信号の極値の最大数。'MaxNumExtrema' と正のスカラー整数で構成されるコンマ区切りのペアとして指定します。MaxNumExtrema は分解停止基準の 1 つです。つまり、極値の数が MaxNumExtrema 未満の場合に分解が停止します。

MaxNumExtrema は 0 または正の整数のみを使用して指定できます。

信号対残差エネルギー比。'MaxEnergyRatio' とスカラーで構成されるコンマ区切りのペアとして指定します。MaxEnergyRatio はシフティングの開始における信号のエネルギーと平均包絡線エネルギーとの比です。MaxEnergyRatio は分解停止基準の 1 つです。つまり、現在のエネルギー比が MaxEnergyRatio より大きい場合に分解が停止します。

包絡線構築のための内挿法。'Interpolation' と 'spline' または 'pchip' のいずれかで構成されるコンマ区切りのペアとして指定します。

以下のように Interpolation を指定します。

  • X が滑らかな信号の場合は、'spline'

  • X が滑らかでない信号の場合は、'pchip'

'spline' 内挿法は 3 次スプラインを使用します。一方で 'pchip' は区分的 3 次エルミート内挿多項式法を使用します。

コマンド ウィンドウの情報表示の切り替えです。'Display' と 1 または 0 で構成されるコンマ区切りのペアで指定します。コマンド ウィンドウ内に生成されたテーブルは、生成された各 IMF のシフト反復の数、相対許容誤差、およびシフト停止基準を示します。テーブルを表示するには Display に 1 を、非表示にするには 0 を指定します。

出力引数

すべて折りたたむ

固有モード関数。行列または timetable として返されます。imf は、包絡線がゼロに関して対称である同数の極値とゼロクロッシングをもつ任意の関数です。imf を使用してヒルベルト・ファン変換を適用し信号のスペクトル解析を実行します。

imf は以下として返されます。

  • X がベクトルの場合、各列が imf の行列

  • X が単一データ列の timetable の場合、timetable

信号の残差。列ベクトルまたは単一データ列の timetable として返されます。residualemd によって分解されていない元の信号 X の部分を表します。

residual は以下として返されます。

  • X がベクトルの場合、列ベクトル

  • X が単一データ列の timetable の場合、単一データ列の timetable

診断の詳細。以下のフィールドを含む構造体として返します。

  • NumIMF - 信号から抽出された IMF の数

  • NumExtrema - 各 IMF における極値の数

  • NumZeroCrossing - 各 IMF のゼロクロッシングの数

  • NumSifting - 各 IMF の実行されたシフティングの数

  • MeanEnvelopeEnergy - 各 IMF 計算において取得された上側および下側の包絡線の平均のエネルギー

  • RelativeTolerance - 各 IMF の得られた相対許容誤差

アルゴリズム

Empirical Mode Decomposition

emd は、シフティング プロセスを使用して、信号 X(t)k の数の固有モード関数 (IMF) と残差 rk(t) に分解します。[1][2]にリストされているシフティング プロセスの簡単な概要は以下のとおりです。

  1. 信号 X(t) の局所的最大値と局所的極小値を検索して上側包絡線 s+(t) と、下側包絡線 s-(t) を構築します。

  2. ith 反復の平均包絡線、mk,i(t) を計算します。

    mk,i(t)=12[s+(t)+s(t)]

  3. 最初の反復については ck(t) = X(t) を使用して、残差信号から平均包絡線を引きます。

    ck(t)=ck(t)mk,i(t)

    ck(t) が IMF の基準に一致しない場合は、手順 4 および 5 がスキップされます。この手順を ck(t) の新しい値で手順 1 において再度繰り返します。

  4. ck(t) が IMF の基準に一致する場合、新しい残差が計算されます。残差信号を更新するには、kth IMF を前の残差信号から減算します。

    rk(t)=rk1(t)ck(t)

  5. 次に新しい信号 rk(t) として取得した残差を使用して手順 1 から開始し、ck(t) を固有モード関数として保存します。

固有モード関数 N の場合、元の信号は次のように表されます。

X(t)=i=1Nci(t)+rN(t)

シフティング プロセスの詳細は、[1]および[2]を参照してください。

SiftRelativeTolerance

SiftRelativeTolerance[4]で提案されているコーシー型の停止基準です。現在の相対許容誤差が SiftRelativeTolerance 未満の場合にシフティングが停止します。現在の相対許容誤差は、以下のように定義されます。

RelativeTolerancec(t)previousc(t)current2c(t)current2.

MaxEnergyRatio

エネルギー比は、シフティングの開始における信号のエネルギーと平均包絡線エネルギーの比です。[3]現在のエネルギー比が MaxEnergyRatio より大きい場合に分解が停止します。k IMF では、EnergyRatio

EnergyRatio10log10(X(t)2rk(t)2).

と定義されます。

参照

[1] Norden E. Huang, Zheng Shen, Steven R. Long, Manli C. Wu, Hsing H. Shih, Quanan Zheng, Nai-Chyuan Yen, Chi Chao Tung, Henry H. LiuProc. R. Soc. Lond. A. "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 1998. 454. 903-995. 10.1098/rspa.1998.0193.

[2] Rilling, G & Flandrin, Patrick & Gonçalves, Paulo. "On empirical mode decomposition and its algorithms." IEEE-EURASIP workshop on nonlinear signal and image processing 2003. NSIP-03. Grado, Italy. 8-11.

[3] Rato, R.T. & Ortigueira, Manuel & Batista, Arnaldo. "On the HHT, its problems, and some solutions." Mechanical Systems and Signal Processing 2008. 22. 1374-1394. 10.1016/j.ymssp.2007.11.028.

[4] Wang, Gang & Chen, Xianyao & Qiao, Fang-Li & Wu, Zhaohua & Huang, Norden. "On Intrinsic Mode Function." Advances in Adaptive Data Analysis 2010. 2. 277-293. 10.1142/S1793536910000549.

拡張機能

C/C++ コード生成
MATLAB® Coder™ を使用して C および C++ コードを生成します。

参考

R2018a で導入