将多项式拟合到函数的最大值

我想将一个多项式拟合成有噪声的数据,这样近似多项式始终> =原始数据。 例如:

x = linspace (-2, 6);
y = (x-2).^2 + 1 + 2 * randn (size (x));

function ret = delta (P, x, y)
  yP = polyval (P, x);
  d = yP - y;
  d (d < 0) *= 1000;
  ret = sumsq (d);
endfunction

P0 = polyfit (x, y, 2);
f = @(P) delta (P, x, y);
[P, FVAL] = sqp (P0, f)

xi = linspace (min(x), max(x), 100);
yi = polyval (P, xi);

plot(x, y, xi, yi);
grid on

原始(蓝色)和约。 (绿色)数据

有更好的方法/方法也适用于高阶多项式吗?

easies的方法是只使用polyfit,然后计算max(y-yi)并将其作为偏移量添加,但这不是最优的...

编辑:我想使用GNU OCtave,但添加了“matlab”作为标签,因为语言和功能是相似的。

编辑:根据tvo的答案和实际数据:

x = [10 20 30 40 50 60 80 100];
y = [0.2372, 0.1312, 0.0936, 0.0805, 0.0614, 0.0512, 0.0554, 0.1407];

function ret = delta (P, x, y)
  ret = sumsq (polyval (P, x) - y);
endfunction

f = @(P) delta (P, x, y);
h = @(P) polyval(P, x) - y;

P0 = polyfit (x, y, 3);
[P] = sqp (P0, f, [], h)

xi  = linspace (min(x), max(x));
yi  = polyval (P0, xi);
yio = polyval (P, xi);

plot(x, y, xi, yi, ";initial guess;", xi, yio, ";optimized;");
grid on

tvos结果图

但正如你所看到的那样,经过优化和评估的poly已经指向<不允许的原始点。


你的方法看起来很好,我没有看到它不能用于实际的高阶多项式。 请解释为什么如果你认为它不能使用。

你正在使用Octave的'sqp'求解器。 文档在这里:http://www.gnu.org/software/octave/doc/v4.0.1/Nonlinear-Programming.html

您可能希望避免将错误乘以任意数字(在您的示例中为1000),否则将为负数。 对于不同的数据集,这可能会失败,尤其是如果它们更大,即更多的数据点。

您可以尝试使用Octave'sqp'提供的非线性不等式约束选项,即h(x)> = 0(请参阅doc)。

作为目标函数phi,您可以使用平方规范错误(如您的示例中所示),并为每个数据点添加格式h(x)> = 0的约束。 请注意,'x'是您想要拟合的多项式系数,h(x)是在特定数据点处评估的多项式。

例如:

phi = @(P) delta_mod (P, x, y); % mod: Don't increase the importance of negative residuals!!
h = @(P) polyval(P, x1) - y1;

Psol = sqp(P0, phi, [], h);

请注意,约束函数'h'确保多项式位于(x1,y1)之上,并且目标函数'phi'将尽量保持它尽可能接近。 您可以扩展'h'来为集合中的每个数据点包含一个约束(请参阅doc)。

链接地址: http://www.djcxy.com/p/91417.html

上一篇: Fit polynomial to the maximums of a function

下一篇: Difference between sphinxcontrib.napoleon and numpy.numpydoc