插值多项式gydF4y2Ba

约翰D 'ErricogydF4y2Ba今天回来讨论多项式的插值。gydF4y2Ba

内容gydF4y2Ba

简介gydF4y2Ba

我用我的gydF4y2Ba上一篇博文gydF4y2Ba关于多项式回归模型。我们在多项式模型中看到了一些效用，并且增加多项式的阶数可以减少误差。我们能走多远?我们能把误差减小到零吗?gydF4y2Ba

在某些情况下，我可能会选择将多项式的次数增加到比数据点的数量小1。假设我的数据没有复制点，这是一个gydF4y2Ba插值多项式gydF4y2Ba这与我们的数据完全吻合，至少在我们计算的两倍精度范围内。直线可以通过任意两点，二次曲线可以通过三个点，三次曲线正好经过四个点，等等。如果你喜欢用统计自由度来思考，如果你有gydF4y2BangydF4y2Ba点，一个多项式的顺序gydF4y2BangydF4y2Ba1 (gydF4y2BangydF4y2Ba系数)准确地传递你的数据。gydF4y2Ba

我从最小二乘意义上的多项式建模开始，或者说是回归分析。但实际上，我想讨论的是插值，而不是多项式回归提供的近似。也许这是一种间接的方法，但我认为重要的事情之一是将插值与数学中使用的更通用的建模/曲线拟合工具区分开来。许多MATLAB用户通常不了解这种差异。gydF4y2Ba

我们从指数函数中抽取5个点开始。gydF4y2Ba

格式gydF4y2Ba短gydF4y2BaX = (-2:1:2)';Y = exp(x);情节(x, y,gydF4y2Ba“o”gydF4y2Ba网格)gydF4y2Ba在gydF4y2Ba

完全符合gydF4y2Ba

虽然这是事实gydF4y2BapolyfitgydF4y2Ba与命令完全吻合gydF4y2Ban - 1gydF4y2Ba多项式，直接法更有效。gydF4y2Ba

构造插值多项式使用gydF4y2Ba范德gydF4y2Ba配合使用gydF4y2Ba反斜杠gydF4y2Ba．gydF4y2Ba

P = vander(x)\y xev = linspace(-2,2,500)';yev = polyval(P,xev);情节(x, y,gydF4y2Ba“柯”gydF4y2Baxev yev,gydF4y2Ba的r -gydF4y2Ba网格)gydF4y2Ba在gydF4y2Ba标题(gydF4y2Ba' exp(x)上5点的四次多项式'gydF4y2Ba）gydF4y2Ba

P = 0.0492 0.2127 0.4939 0.9625 1.0000gydF4y2Ba

许多书教你如何使用gydF4y2Ba拉格朗日形式gydF4y2Ba插值。我在这里就不讲了因为我不建议一般使用多项式。gydF4y2Ba

插值应该在数据点上产生零残差gydF4y2Ba

验证残差在数据点上基本为零。gydF4y2Ba

格式gydF4y2Ba长gydF4y2BaggydF4y2Bapred = polyval(P,x) - ygydF4y2Ba

Pred = -1.11022302462516e-016 -1.11022302462516e-016 0 0 0gydF4y2Ba

数据是否有复制?gydF4y2Ba

假设我有两个自变量值相同的点，但因变量值不同?因为对于任何给定gydF4y2BaxgydF4y2Ba，单值函数仅为预测返回一个值，即thisgydF4y2Ba必须gydF4y2Ba制造麻烦。在这个例子中，我应该期望什么值插值时gydF4y2BaX = 2gydF4y2Ba？应该gydF4y2Baf (x)gydF4y2Ba是gydF4y2Ba3.gydF4y2Ba或gydF4y2Ba4gydF4y2Ba？gydF4y2Ba

X = [0 1 2 2 3 4]';Y = [1 2 3 4 5 6]';情节(x, y,gydF4y2Ba“o”gydF4y2Ba）gydF4y2Ba

如果你看一下范德蒙矩阵的秩，你就能看到问题的线索。秩应该是6，但是复制的点导致了一个奇点。gydF4y2Ba

排名(范德(x))gydF4y2Ba

Ans = 5gydF4y2Ba

这个奇点的结果是一个相当无用的多项式。gydF4y2Ba

P = vander(x)\ygydF4y2Ba

警告:矩阵对工作精度是奇异的。P = NaN NaN NaN -Inf InfgydF4y2Ba

平均复制gydF4y2Ba

如果数据确实不是单值的，常用的方法是使用参数关系。(这是为将来可能的博客保存的东西。)gydF4y2Ba

然而，大多数时候，我只想做一些简单的事情，比如平均复制点。gydF4y2Ba集运商gydF4y2Ba(我贡献给File Exchange的一个工具)就可以做到这一点。如果您从未发现或查看过文件交换，请花些时间浏览一下那里的出色工具。gydF4y2Ba

高阶多项式不可靠gydF4y2Ba

有时多项式不太好用。了解你的数据，了解你的问题是很重要的。这里的数据来自一个单调递增的函数。多项式不喜欢单调，所以强迫高阶多项式插值这样的函数肯定会失败。gydF4y2Ba

X = (-5:1:5)';Y = erf(x);P = vander(x)\y xev = linspace(-5,5,1000)';yev = polyval(P,xev);情节(x, y,gydF4y2Ba“柯”gydF4y2Baxev yev,gydF4y2Ba的r -gydF4y2Ba网格)gydF4y2Ba在gydF4y2Ba标题(gydF4y2Ba在erf(x)上通过11个点进行多项式插值gydF4y2Ba）gydF4y2Ba

P = 3.49245965480804e-020 2.0991424886728e- 018 -0.00119783140998794 2.17562811018059e-017 -9.91270557701033e-017 -0.211403900317105 1.85037170770859e-016 1.03132935951897gydF4y2Ba

龙格现象gydF4y2Ba

龙格gydF4y2Ba发现了这个函数的例子，它会导致高阶插值多项式非常糟糕。gydF4y2Ba

X = (-5:1:5)';Y = 1 /(1+x.^2);P = vander(x)\y xev = linspace(-5,5,1000)';yev = polyval(P,xev);Ytrue = 1 /(1+xev.^2);情节(x, y,gydF4y2Ba“柯”gydF4y2Baxev yev,gydF4y2Bar -。gydF4y2Baxev ytrue,gydF4y2Ba“g -”gydF4y2Ba)传说(gydF4y2Ba的数据点gydF4y2Ba，gydF4y2Ba插值多项式的gydF4y2Ba，gydF4y2Ba'y = 1./(1+x.^2)'gydF4y2Ba网格)gydF4y2Ba在gydF4y2Ba标题(gydF4y2Ba“龙格例子”gydF4y2Ba）gydF4y2Ba

P = -2.26244343891403e-005 1.04238735515992e-019 0.00126696832579186 -4.47734652027112e-018 -0.0244117647058824 5.61044648135202e-017 0.19737556561086 -2.38551402729053e-016 -0.67420814479638 1.38777878078145e-016gydF4y2Ba

数值问题gydF4y2Ba

你还需要考虑高阶多项式插值的数值问题。你所创建的矩阵通常在数值上是奇异的。gydF4y2Ba

例如，矩阵的秩gydF4y2Ba范德(1:n)gydF4y2Ba应该是gydF4y2BangydF4y2Ba．然而，在有限精度算法中，这是失败的。gydF4y2Ba

Rank (vander(1:10) Rank (vander(1:20) Rank (vander(1:30))gydF4y2Ba

Ans = 10 Ans = 7 Ans = 6gydF4y2Ba

对数据进行居中和缩放可以在一定程度上提高高阶多项式的数值性能，但这也有局限性。你通常最好使用其他插值方法。gydF4y2Ba

你能做什么?你希望能够对数据进行插值。多项式看起来是个好地方，但它们也有自己的问题。高阶多项式插值通常存在一些问题，要么导致非单调插值，要么导致数值问题。你可以考虑其他函数族来构建你的插值函数，例如三角函数或贝塞尔函数，或正交多项式。这很大程度上取决于你想要建模的关系。了解你的问题!gydF4y2Ba

顺便说一下，对于那些对建模和曲线拟合的方法感兴趣的人，我在这里列出了一些想法gydF4y2Ba线程gydF4y2Ba这应该很有趣。gydF4y2Ba

我们仍然可以用多项式做很多事情。它们构成了许多建模工具的基础。gydF4y2Ba

一个很好的技巧是简单。为什么用一个高阶多项式来插值我们的数据?相反，将问题分解，使用更小的低阶多项式片段，拼接在一起。在我的下一篇博客中，我将开始讨论这些分段插值。在那之前请告诉我gydF4y2Ba在这里gydF4y2Ba当你发现插值多项式的使用，或者你发现它们无法解决的问题。gydF4y2Ba

获取MATLAB代码gydF4y2Ba

使用MATLAB®7.6发布gydF4y2Ba