matlab - 如何绘制多项式拟合参数中的误差 - MATLAB

Question

我有一个数据列表，我试图将其拟合到多项式中，并且我也在尝试绘制参数的 95% 置信带(在 Matlab 中)。如果我的数据是 x 和 y

f=fit(x,y,'poly2')
plot(f,x,y)
ci=confint(f,0.95);
a_ci=ci(1,:);
b_ci=ci(2,:);

我不知道在那之后如何获得我的数据周围的最小和最大范围。有谁知道这是怎么做到的吗？

Answer 1

我可以看到您安装了曲线拟合工具箱，这很好，因为您需要它才能使以下代码正常工作。

示例数据的基本拟合

让我们定义一些示例数据和一个可能的拟合函数。 (我也可以在这里使用 poly2，但我想让它更通用一些。)

xdata = (0:0.1:1)';               % column vector!
noise = 0.1*randn(size(xdata));
ydata = xdata.^2 + noise;
f = fittype('a*x.^2 + b'); 
fit1 = fit(xdata, ydata, f, 'StartPoint', [1,1])
plot(fit1, xdata, ydata)

旁注:plot() 不是我们常用的绘图函数，而是一个 method of the cfit-object适合1.

拟合参数的置信区间

我们的拟合使用数据来确定基础模型 f(x)=ax2+b 的系数 a、b。您已经这样做了，但为了完整起见，这里是您如何读出任何置信区间的系数不确定性的方法。系数按字母顺序排列，这就是为什么我可以将 ci(1,:) 用于 a，等等。

names = coeffnames(fit1)   % check the coefficient order!
ci = confint(fit1, 0.95);  % 2 sigma interval
a_ci = ci(1,:)
b_ci = ci(2,:)

默认情况下，Matlab 使用 2σ (0.95) 置信区间。有些人(物理学家)更喜欢引用 1σ (0.68) 区间。

置信带和预测带

画图是个好习惯confidence bands or prediction bands围绕数据——尤其是当系数相关时!但是你应该花点时间考虑一下你想要绘制两个中的哪一个:

• 预测带:如果我采用新的测量值，我预计它位于何处？在 Matlab 术语中，这称为“观察带”。
• 置信区间:我认为真正的值(value)在哪里？在 Matlab 术语中，这称为“功能带”。

与系数的置信区间一样，Matlab 默认使用 2σ 波段，我们中的物理学家将其切换为 1σ 区间。从本质上讲，预测带更大，因为它是模型误差(置信带!)和测量误差的组合。

还有一个我不完全理解的命运。 Matlab 和维基百科 make that distinction .

• Pointwise:单个测量值/真实值的预测/置信区间有多大？在我能想到的几乎所有情况下，这就是您作为物理学家想要问的问题。
• 同时:如果您希望一组所有新测量值/所有预测点都位于给定置信度的区间内，您必须将预测/置信区间设为多大？

我个人认为，“同步乐队”不是乐队!对于 n 个点的测量，它应该是 n 个单独的误差条!

预测/置信度区别和逐点/同时区别为您提供了围绕图的“the”波段的总共四个选项。 Matlab 使 2σ 逐点预测带 易于访问，但您似乎感兴趣的是 2σ 逐点置信带。这是一个bit more cumbersome绘制，因为您必须指定要评估预测带的虚拟数据:

x_dummy = linspace(min(xdata), max(xdata), 100);
figure(1); clf(1);
hold all
plot(xdata,ydata,'.')
plot(fit1)    % by default, evaluates the fit over the currnet XLim
% use "functional" (confidence!) band; use "simultaneous"=off
conf1 = predint(fit1,x_dummy,0.95,'functional','off');
plot(x_dummy, conf1, 'r--')
hold off

请注意，x=0 处的置信带等于拟合系数 b 的置信区间!

外推

如果您想外推到数据范围未涵盖的 x 值，您可以评估更大范围的拟合和预测/置信带:

x_range = [0, 2];
x_dummy = linspace(x_range(1), x_range(2), 100);
figure(1); clf(1);
hold all
plot(xdata,ydata,'.')
xlim(x_range)
plot(fit1)
conf1 = predint(fit1,x_dummy,0.68,'functional','off');
plot(x_dummy, conf1, 'r--')
hold off