在Python中获得Hermite级数对数据的最小二乘拟合

在科学计算中，最小二乘拟合是一种常用的方法，用于找到最适合一组数据的函数。以下将介绍如何在Python中使用Hermite级数进行最小二乘拟合。

Hermite级数简介

Hermite级数是一种基于Hermite多项式的级数，用于在某个区间内逼近连续函数。如果我们想要在区间 $[-1,1]$ 内逼近一个一次函数 $f(x)$ ，可以使用形如下面这个级数：

$f(x) = c_0 + c_1 H_1(x) + c_2 H_2(x)$

其中， $H_n(x)$ 是Hermite多项式，它们可以使用递归形式计算得出：

$H_0(x) = 1$
$H_1(x) = 2x$
$H_n(x) = 2xH_{n-1}(x) – 2(n-1)H_{n-2}(x)$

$c_0$ 、 $c_1$ 、 $c_2$ 等系数可以通过最小二乘法得到。更高阶的Hermite级数可以使用类似的方法来定义和拟合，但是计算起来会更复杂。

最小二乘拟合

在这里，我们假设我们已经有了一些数据点 ${(x_i,y_i)}$ ，我们想要找到一个适合这些点的Hermite级数 $f(x)$ 。可以将这个问题表示为以下的最小二乘拟合问题：

$\min_{c_0,c_1,c_2} \sum_{i=1}^n (y_i – f(x_i))^2$

这个问题的解可以通过矩阵运算求得：

$\begin{bmatrix} 1&H_1(x_1)&H_2(x_1)\\1&H_1(x_2)&H_2(x_2)\\\vdots&\vdots&\vdots\\1&H_1(x_n)&H_2(x_n)\\\end{bmatrix} \begin{bmatrix} c_0\\c_1\\c_2 \end{bmatrix} = \begin{bmatrix} y_1\\y_2\\\vdots\\y_n \end{bmatrix}$

即：

$A\mathbf{c} = \mathbf{y}$

其中， $A$ 是矩阵， $\mathbf{c}$ 是我们要求的系数向量， $\mathbf{y}$ 是原始数据中的 $y$ 值组成的向量。

可以使用NumPy中的polyfit函数来计算拟合系数。在这里，我们设定要拟合的数据为：

x = np.linspace(-1, 1, 10)
y = np.sin(x*np.pi/2) + np.random.normal(size=10)*0.1

使用Hermite级数进行最小二乘拟合的代码如下所示：

def hermite(x, n):
    if n == 0:
        return np.ones_like(x)

    if n == 1:
        return 2*x

    return 2*x*hermite(x, n-1) - 2*(n-1)*hermite(x, n-2)

def hermite_fit(x, y, degree):
    X = np.zeros((len(x), degree+1))

    for i in range(degree+1):
        X[:,i] = hermite(x, i)

    c = np.linalg.lstsq(X, y, rcond=None)[0]

    return c

c = hermite_fit(x, y, 2)
xx = np.linspace(-1, 1, 1000)
yy = np.zeros_like(xx)

for i in range(len(c)):
    yy += c[i]*hermite(xx, i)

plt.plot(xx, yy, label='fit')
plt.plot(x, y, 'o', label='data')
plt.legend()
plt.show()

在这段代码中，首先定义了一个计算Hermite多项式 $H_n(x)$ 的递归函数hermite。然后，我们定义一个新的函数hermite_fit，它计算拟合的系数向量 $\mathbf{c}$ 并返回。在函数中，我们首先创建一个系数矩阵 $A$ ，矩阵的每一列都是递归计算得到的Hermite多项式。然后，使用np.linalg.lstsq函数求解拟合系数向量。最后，使用拟合系数向量生成拟合曲线。

在这个例子中，拟合了一个三次Hermite级数，生成的拟合曲线如下所示：