sgolay函数如何在Matlab R2013a中工作? [英] How does the sgolay function work in Matlab R2013a?

查看:282
本文介绍了sgolay函数如何在Matlab R2013a中工作?的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我对Matlab中的 sgolay 函数有疑问R2013a.我的数据库有165个光谱,带有2884个变量,我想取它们的一阶和二阶导数.我如何定义KFsgolay的输入?

I have a question about the sgolay function in Matlab R2013a. My database has 165 spectra with 2884 variables and I would like to take the first and second derivatives of them. How might I define the inputs K and F to sgolay?

下面是一个示例:

sgolay用于平滑嘈杂的正弦曲线,并将生成的一阶和二阶导数与使用diff计算的一阶和二阶导数进行比较.请注意,使用diff如何放大噪声并产生无用的结果.

sgolay is used to smooth a noisy sinusoid and compare the resulting first and second derivatives to the first and second derivatives computed using diff. Notice how using diff amplifies the noise and generates useless results.

K = 4;                 % Order of polynomial fit
F = 21;                % Window length
[b,g] = sgolay(K,F);   % Calculate S-G coefficients

dx = .2;
xLim = 200;
x = 0:dx:xLim-1;

y = 5*sin(0.4*pi*x)+randn(size(x));  % Sinusoid with noise

HalfWin  = ((F+1)/2) -1;
for n = (F+1)/2:996-(F+1)/2,
  % Zero-th derivative (smoothing only)
  SG0(n) =   dot(g(:,1), y(n - HalfWin: n + HalfWin));

  % 1st differential
  SG1(n) =   dot(g(:,2), y(n - HalfWin: n + HalfWin));

  % 2nd differential
  SG2(n) = 2*dot(g(:,3)', y(n - HalfWin: n + HalfWin))';
end

SG1 = SG1/dx;         % Turn differential into derivative
SG2 = SG2/(dx*dx);    % and into 2nd derivative

% Scale the "diff" results
DiffD1 = (diff(y(1:length(SG0)+1)))/ dx;    
DiffD2 = (diff(diff(y(1:length(SG0)+2)))) / (dx*dx);

subplot(3,1,1);
plot([y(1:length(SG0))', SG0'])
legend('Noisy Sinusoid','S-G Smoothed sinusoid')

subplot(3, 1, 2);
plot([DiffD1',SG1'])
legend('Diff-generated 1st-derivative', 'S-G Smoothed 1st-derivative')

subplot(3, 1, 3);
plot([DiffD2',SG2'])
legend('Diff-generated 2nd-derivative', 'S-G Smoothed 2nd-derivative')

推荐答案

在本质上嘈杂的过程中进行派生.因此,如果您的数据中确实已经有一些噪声,那么当您采用高阶导数时,噪声会被放大. Savitzky-Golay 是将平滑和微分合并为一个操作的一种非常有用的方法.这是一种通用方法,它可以将导数计算为任意顺序.但是,需要权衡取舍.对于具有特定结构的数据,还存在其他特殊方法.

Taking derivatives in an inherently noisy process. Thus, if you already have some noise in your data, indeed, it will be magnified as you take higher order derivatives. Savitzky-Golay is a very useful way of combining smoothing and differentiation into one operation. It's a general method and it computes derivatives to an arbitrary order. There are trade-offs, though. Other special methods exist for data with a certain structure.

关于您的申请,我没有任何具体答案.在很大程度上取决于数据的性质(采样率,噪声比等).如果使用过多的平滑处理,则会涂抹数据或产生锯齿.如果通过使用高阶多项式系数K对数据进行过拟合,也会发生同样的事情.在演示代码中,您还应该绘制sin函数的解析导数.然后使用不同数量的输入噪声和平滑滤波器进行播放.如果您可以近似真实数据的各个方面,则具有已知确切答案的此类工具可能会有所帮助.在实践中,为了产生不太嘈杂的导数,我尝试使用尽可能少的平滑处理.通常,这意味着三阶多项式(K = 3)和窗口大小F尽可能小.

In terms of your application, I don't have any concrete answers. Much depends on the nature of the data (sampling rate, noise ratio, etc.). If you use too much smoothing, you'll smear your data or produce aliasing. Same thing if you over-fit the data by using high order polynomial coefficients, K. In your demo code you should also plot the analytical derivatives of the sin function. Then play with different amounts of input noise and smoothing filters. Such a tool with known exact answers may be helpful if you can approximate aspects of your real data. In practice, I try to use as little smoothing as possible in order to produce derivatives that aren't too noisy. Often this means a third-order polynomial (K = 3) and a window size, F, as small as possible.

是的,很多人建议您用眼睛来调整它们参数.但是,最近也有一些关于自动选择系数的研究: 在选择上最佳Savitzky-Golay滤波器 (2013). Savitzky-Golay也有其他选择,例如,本文基于正规化,但是您可能需要在Matlab中自己实现.

So yes, many suggest that you use your eyes to tune these parameters. However, there has also been some very recent research on choosing the coefficients automatically: On the Selection of Optimum Savitzky-Golay Filters (2013). There are also alternatives to Savitzky-Golay, e.g., this paper based on regularization, but you may need to implement them yourself in Matlab.

顺便说一句,前一段时间我写了一个sgolay的替代品.像您一样,我只需要第二个输出,即微分滤波器G,这就是它的全部计算结果.此功能也更快(大约2到4倍):

By the way, a while back I wrote a little replacement for sgolay. Like you, I only needed the second output, the differentiation filters, G, so that's all it calculates. This function is also faster (by about 2–4 times):

function G=sgolayfilt(k,f)
%SGOLAYFILT  Savitzky-Golay differentiation filters
s = vander(0.5*(1-f):0.5*(f-1));
S = s(:,f:-1:f-k);
[~,R] = qr(S,0);
G = S/R/R';

具有输入验证功能的此功能的完整版本为可在我的GitHub上找到.

A full version of this function with input validation is available on my GitHub.

这篇关于sgolay函数如何在Matlab R2013a中工作?的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆