如何使用曲线拟合的标准偏差误差来绘制误差条

Question

我想知道我是否理解标准差的工作原理。

曲线拟合函数返回函数的参数及其协方差，用于查找标准差。

我可以使用标准差来绘制带有误差线的拟合曲线吗？

这是一个例子。自从我学习以来，我尝试让一切尽可能简单。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

x = np.arange(0,10)
y = np.array([5,4.8,4.7,4.5,4.4,4.5,4.8,4.9,5,5.1])

def function(x,a,b,c):
    return a*x**2+ b*x+c

plt.scatter(x,y)

parameters, pcov = curve_fit(function, x, y)

perr = np.sqrt(np.diag(pcov))

print(perr)

plt.scatter(x, function(x, parameters[0], parameters[1], parameters[2]))

Answer 1

您所要求的一定是可能的，因为 MATLAB 中有一个函数可以执行此操作（请参阅polyconf）。我不知道如何在 Python 中做到这一点。之前有人问过这个问题here，但我认为这个答案并不令人满意。 2022 年的一条评论说“我认为 sklearn 上没有任何内置功能来计算置信区间”。

希望有人能为您提供答案，但您可能需要考虑使用高斯过程回归。这可能是实现您目标的另一种方式。

GP 背后的数学知识非常复杂，但如果你能克服这个障碍，它们实际上是非常有用和实用的工具。

由于您选择了一条您认为是真实函数的良好模型的曲线，因此您可以将其用作高斯过程的固定均值函数（请参阅 Rasmussen & C. K. I. Williams，2006 年的第 2.7 节）。正如他们在书中所说：

简单地将通常的零均值 GP 应用于观测值与固定均值函数之间的差异

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import ConstantKernel, RBF

x = np.arange(0, 10)
y = np.array([5, 4.8, 4.7, 4.5, 4.4, 4.5, 4.8, 4.9, 5, 5.1])

def function(x, a, b, c):
    return a * x**2+ b * x+c

# Fit polynomial (your code)
parameters, pcov = curve_fit(function, x, y)
m_func = function(x, parameters[0], parameters[1], parameters[2])

# Get the differences between measurements and the fitted curve
# We will use a GP to model these
y_gp = y - m_func

# Gaussian process parameters - adjust these to reflect 
# your beliefs about the true system
noise_std = 0.1  # measurement noise level
length_scale = 6  # the smoothness of the true function
kernel = RBF(length_scale, length_scale_bounds="fixed")
gpr = GaussianProcessRegressor(kernel=kernel, alpha=noise_std**2)
gpr.fit(x.reshape(-1, 1), y_gp)

# Plot the data points
plt.scatter(x, y, label="Data")

# Plot the mean function
X_plot = np.linspace(x[0], x[-1], 101)
m_func = function(X_plot, parameters[0], parameters[1], parameters[2])
plt.plot(X_plot, m_func, linestyle='--', color='k', label='Mean function')

# Compute GP predictions
y_mean_pred, y_std_pred = gpr.predict(X_plot.reshape(-1, 1), return_std=True)

# Plot mean prediction and the confidence interval
plt.plot(X_plot, y_mean_pred + m_func, label="Mean prediction")
plt.fill_between(
    X_plot,
    y_mean_pred - 1.96 * y_std_pred + m_func,
    y_mean_pred + 1.96 * y_std_pred + m_func,
    alpha=0.5,
    label=r"95% confidence interval",
)
plt.xlabel("$x$")
plt.show()