线性/非线性拟合到正弦曲线

我看过这个和这个 。

但我的问题略有不同。 我知道我的数据是一条正弦曲线,具有未知周期和未知幅度,具有加性非高斯分布噪声。

我试图在C中使用GSL非线性算法来拟合它,但是拟合非常糟糕。 我想知道我是否(错误地)使用非线性拟合算法,我应该使用线性算法?

如何判断特定数据集是否需要线性算法或非线性算法?

编辑:我的曲线非常嘈杂,因此使用FFT来计算频率可能会导致误报和错误拟合。 我正在寻找一种更健壮的贴合方式。

曲线约170点

如您所见,上面的情节大约有170个点,下面的情节大约有790个点。

在此处输入图像描述

噪声明显是非高斯的,并且与数据的幅度相比较大。 我已经尝试过高斯分布式噪声的FFT,而且我很适合。 在这里,它的失败非常严重。

增加:链接到第一次系列数据。 文件中的每一列都是不同的时间序列。

如果您知道您的数据是正弦曲线(可以表示为多个复指数),那么您可以使用Pisarkenko的谐波分解; http://en.wikipedia.org/wiki/Pisarenko_harmonic_decomposition

但是,如果您可以访问更多数据点,我的方法仍然是使用DFT。

更新:

我在你的数据上使用了Pisarenko的谐波分解(PHD),即使你的信号非常短(每个只有86个数据点),如果有更多的数据,PHD算法肯定有潜力。 下面包括24个信号中的两个(数据的第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的validation:

我还做了另一个测试,我知道振幅和频率的值,然后增加噪声,看看PHD是否可以从噪声信号中正确估计值。 该信号由两条相加的正弦曲线组成,频率分别为50Hz,120Hz,振幅分别为0.7,1.0。 在下图中,红色曲线是原始曲线,蓝色曲线是增加的噪声。 (图被裁剪)

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估计放大器/频率值为:

 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) 

对于相当多的噪音并不坏,只知道信号所包含的波数。

回复@Alex:

是的,这是一个很好的算法,我在我的DSP研究期间遇到过它,我认为它运作得很好,但重要的是要注意Pisarenko的Harm.Dec。 将任何信号建模为N> 0个正弦曲线,N从开始指定,并使用该值忽略噪声。 因此,根据定义,只有当您大致了解数据的人体正弦曲线时,它才有用。 如果您不知道N的值,并且您需要运行该算法以获得一千个不同的值,那么肯定会建议使用不同的方法。 也就是说,此后评估是直截了当的,因为它返回N个幅度和频率值。

多信号分类(MUSIC)是Pisarenko离开的另一种算法。 http://en.wikipedia.org/wiki/Multiple_signal_classification

Kitchi:你能提供一些样本数据吗? 您必须使用多长时间的典型信号? (就样本数量和正弦波周期数而言)信噪比(dB)是多少?

  1. 在您知道algo将起作用之前,我建议您使用python / numpy / scipy(或matlab / octave,或R / S或Mathematica ……)进行原型设计,无论您最喜欢哪种原型语言/工具集,除了C.将节省大量时间,您将使用更丰富的工具。

  2. 你确定噪音会严重影响FFT吗? 这不一定是一个好的假设,特别是如果噪声相对“白”,并且分析窗口很长。 如果正弦波的频率非常稳定,您可以进行大量的FFT并从比信号强的噪声数量级提取信号。 尝试大约几百到几百万个预期正弦波周期的东西。

  3. 曲线拟合正弦波不能很好地工作。 我猜周期性会产生很多局部最小值,而相移变量也会使问题显着非线性。 您可以看到其他遇到同一问题的人提出的问题。 除非你预先线性化问题,否则你最好尝试几乎任何其他东西,而不是非线性最小二乘拟合,这让我…

  4. 自相关对于这种事物来说非常好。 尝试立即计算整个信号的自相关(如果源频率稳定,则数据越多越好)。 正弦波周期作为自相关中的高峰应该非常明显,并且您可能会获得比FFT更准确的频率估计(除非您使用极大的FFT)。 此外,您可以从第一个高自相关峰值的高度计算平均幅度。

编辑 :经过进一步研究,有更多的技术可能比FFT更适合您的问题。 Pisarenko的谐波分解(下面由Fredrik Rubin首先提出)是一个; 另一种是最小二乘谱分析 (LSSA),它与您原来的问题思想非常相似。 LSSA有许多变体,例如Lomb-Scargle,基础追踪等,它们以各种方式处理我上面描述的拟合问题。 但是我想如果你在大型FFT中看不到任何信号,那么其他方法都不可能找到任何东西:)

PS有关无法正确拟合正弦波的其他问题,请参阅:

  • curve_fit甚至在正弦波上失败

  • SciPy leastsq适合正弦波失败

  • 我知道scipy curve_fit可以做得更好

这实际上是一个谱估计问题。 您正在尝试估算一个“线谱”,在那里您可以知道正弦波的数量(在您的情况下,一个)。 像MUSIC或ESPRIT这样的方法应该能够解决问题。

作为参考, Stoica的书将派上用场。 本书的第4章是线谱的参数化方法,其中包含用于查找所需信号的幅度,相位和频率的算法。 本书还附带了在MATLAB中实现的算法 ,它们也很容易自己实现。

如果你正在对sin进行回归,你可以使用FFT进行傅立叶变换。

编辑

尝试用filter去除噪音。 如果您有像传感器这样的物理来源,请在传感器上放置低通滤波器。 FFT是相对较差的滤波器。

编辑2 – 这种测量是完全错误的

可能是,你做错了测量。 根据Nyquist-Shannon采样定理,您的采样频率太低,或者输入频率太高。 这导致WRONG解决方案,因为如果您采样5kHz采样,例如5kHz采样,您将根据此定理测量2kHz。

我确信您无法通过此类测量来确定正确的输入频率。