16

私はこれこれを見てきました。

しかし、私には少し異なる問題があります。私のデータは、周期と振幅が不明な正弦曲線であり、非ガウス分布のノイズが追加されていることを知っています。

CのGSL非線形アルゴリズムを使用して適合させようとしていますが、適合は絶対にひどいものです。線形アルゴリズムを使用する必要がある非線形フィッティングアルゴリズムを(間違って)使用しているかどうか疑問に思っていますか?

特定のデータセットに線形アルゴリズムと非線形アルゴリズムのどちらが必要かを判断するにはどうすればよいですか?

編集:私の曲線は本当にうるさいので、周波数を把握するためにFFTを実行すると、誤検知や不適合が発生する可能性があります。もう少し頑丈なフィッティング方法を探しています。

約170点の曲線

ご覧のとおり、上のプロットには約170ポイント、下のプロットには約790ポイントがあります。

ここに画像の説明を入力してください

ノイズは明らかに非ガウスであり、データの振幅と比較して大きいです。ガウス分布のノイズでFFTを試しましたが、フィット感は素晴らしかったです。ここでは、かなりひどく失敗しています。

追加:最初の時系列データへのリンク。ファイルの各列は異なる時系列です。

4

4 に答える 4

7

データが正弦曲線であることがわかっている場合 (多数の複雑な指数関数として表すことができます)、Pisarkenko の調和分解を使用できます。http://en.wikipedia.org/wiki/Pisarenko_harmonic_decomposition

ただし、より多くのデータ ポイントにアクセスできる場合、私のアプローチでは引き続き DFT を使用します。

アップデート:

私はあなたのデータに Pisarenko の高調波分解 (PHD) を使用しました。信号が非常に短い (それぞれ 86 データポイントしかない) 場合でも、PHD アルゴリズムは、より多くのデータを利用できる場合、間違いなく可能性を秘めています。以下に含まれるのは、24 個の信号のうちの 2 つ (データの列 11 と 13) で、青色で示されています。赤色の正弦曲線は、PHD から推定された振幅/周波数値に対応しています。(位相シフトは不明であることに注意してください)

列 11 のデータのプロット 列 13 のデータのプロット

MATLAB (pisar.m) を使用して PHD を実行しました: http://www.mathworks.com/matlabcentral/fileexchange/74

% assume data is one single sine curve (in noise)
SIN_NUM = 1; 

for DATA_COLUMN = 1:24
    % obtain amplitude (A), and frequency (f = w/2*pi) estimate
    [A f]=pisar(data(:,DATA_COLUMN),SIN_NUM);

    % recreated signal from A, f estimate
    t = 0:length(data(:,DATA_COLUMN))-1;
    y = A*cos(2*pi*f*t);

    % plot original/recreated signal
    figure; plot(data(:,DATA_COLUMN)); hold on; plot(y,'r')
    title({'data column ',num2str(DATA_COLUMN)});

    disp(A)
    disp(f)
end

その結果、

1.9727     % amp. for  column 11
0.1323     % freq. for column 11
2.3231     % amp. for  column 13
0.1641     % freq. for column 13

博士号の証明:

また、振幅と周波数の値を知っている別のテストを行い、ノイズを追加して、PHD がノイズの多い信号から適切に値を推定できるかどうかを確認しました。信号は、周波数がそれぞれ 50 Hz、120 Hz、振幅がそれぞれ 0.7、1.0 の 2 つの追加された正弦曲線で構成されていました。下の図では、赤の曲線がオリジナルで、青がノイズを追加したものです。(図はトリミングされています)

PHD精度のテスト

Fs = 1000; % Sampling frequency
T = 1/Fs; % Sample time
L = 1000; % Length of signal
t = (0:L-1)*T; % Time vector

% Sum of a 50 Hz sinusoid and a 120 Hz sinusoid
x = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t);
y = x + 0.4*randn(size(t)); % Sinusoids plus noise

figure;
plot(Fs*t(1:100),y(1:100)); hold on; plot(Fs*t(1:100),x(1:100),'r')
title('Signal Corrupted with Zero-Mean Random Noise (Blue), Original (Red)')

[A, f] = pisar(y',2); 
disp(A)
disp(f/Fs)

PHD は、amp/freq 値を次のように推定しました。

0.7493    % amp wave 1  (actual 0.7)
0.9257    % amp wave 2  (actual 1.0)
58.5      % freq wave 1 (actual 50)
123.8     % freq wave 2 (actual 120)

かなりのノイズがあり、信号が構成されている波の数だけを知っているのは悪くありません.

@アレックスに返信:

ええ、それは素晴らしいアルゴリズムです。私は DSP の研究中にこれに出会い、非常にうまく機能すると思いましたが、ピサレンコの Harm.Dec. は任意の信号を N > 0 の正弦波としてモデル化し、N は最初から指定され、その値を使用してノイズを無視します。したがって、定義上、データがどのように正弦波で構成されているかを大まかに知っている場合にのみ役立ちます。N の値の手がかりがなく、1000 の異なる値に対してアルゴリズムを実行する必要がある場合は、別のアプローチが間違いなく推奨されます。とはいえ、N 個の振幅値と周波数値が返されるため、その後の評価は簡単です。

多重信号分類 (MUSIC) は、ピサレンコが中断したところから続く別のアルゴリズムです。http://en.wikipedia.org/wiki/Multiple_signal_classification

于 2013-01-24T17:11:25.670 に答える
4

Kitchi: サンプルデータを提供していただけますか? 処理しなければならない典型的な信号の長さはどれくらいですか? (サンプル数と正弦波周期数に関して) 信号対雑音比 (dB) は?

  1. どのアルゴリズムが機能するかを知る前に、python/numpy/scipy (または matlab/octave、R/S、または Mathematica など) でプロトタイプを作成することをお勧めします。C 以外の、お気に入りのプロトタイプ作成言語/ツールセットは何でもかまいません。大幅な時間の節約になり、より豊富なツールを使用できるようになります。

  2. ノイズが FFT に悪影響を与えると確信していますか? 特にノイズが比較的「ホワイト」で、分析ウィンドウが長い場合、これは必ずしも適切な仮定ではありません。正弦波の周波数が非常に安定している場合は、巨大な FFT を実行して、信号よりも数桁強いノイズから信号を引き出すことができます。予想される正弦波の数百から数百万サイクルのオーダーで何かを試してください。

  3. カーブ フィッティング サイン波はうまく機能しません。周期性によって多くの極小値が作成され、位相シフト変数によって問題が大幅に非線形になると思います。以下にリンクされている同じ問題に遭遇した他の人々からのいくつかの質問を見ることができます. 問題を事前に線形化しない限り、非線形最小二乗法ではなく、他のほとんどすべてを試したほうがよいでしょう...

  4. 自己相関は、この種のものに優れています。一度に信号全体の自己相関を計算してみてください (ソース周波数が安定している場合は、データが多いほど良い)。正弦波の周期は、自己相関の高いピークとして非常に明白である必要があり、おそらく FFT よりも正確な周波数の推定値が得られます (非常に大きな FFT を使用しない限り)。また、最初の高い自己相関ピークの高さから平均振幅を計算できます。

編集:さらに調査すると、FFTよりも問題に適した手法が他にもあります。ピサレンコの調和分解 (以下の Fredrik Rubin によって最初に提案された) はその 1 つです。もう1つは、元の質問のアイデアと非常によく似た最小二乗スペクトル分析(LSSA)です。LSSA にはさまざまなバリエーションがあり、たとえば Lomb-Scargle や基底追跡などがあります。これらは、上で説明したフィッティングの問題をさまざまな方法で処理します。ただし、大きなFFTで信号がまったく見えない場合は、他の方法でも何も見つからない可能性が高いと思います:)

PS 正弦波にうまく適合できないことに関連するその他の質問については、以下を参照してください。

于 2013-01-21T05:23:10.967 に答える
2

sin に対する回帰を行っている場合は、FFT を使用してフーリエ変換を適用できます。

編集

フィルタでノイズを除去してみてください。センサーのような物理的なソースがある場合は、センサーにローパス フィルターを付けます。FFTは比較的悪いフィルターです。

EDIT2 - この測定は単純に間違っています

間違った測定をしている可能性があります。ナイキスト・シャノンのサンプリング定理によれば、サンプリング周波数が低すぎるか、入力周波数が高すぎます。たとえば、5kHz サンプリングで 3kHz をサンプリングしている場合、この定理に従って 2kHz を測定するため、間違ったソリューションになります。

このような測定では、正しい入力周波数を判断できないと思います。

于 2013-01-07T14:21:56.363 に答える
2

これは実際にはスペクトル推定の問題です。あなたが持っている正弦波の数(あなたの場合は1つ)を知っている「線スペクトル」を推定しようとしています。MUSICESPRITなどの方法で問題を解決できるはずです。

参考までに、ストイカの本が重宝します。この本の第 4 章は、目的の信号の振幅、位相、および周波数を見つけるためのアルゴリズムを含む線スペクトルのパラメトリック法です。この本には、MATLAB で実装されたアルゴリズムも付属しており、独自に実装するのも簡単です。

于 2013-01-27T13:59:58.760 に答える