python实现最小二乘拟合/平方逼近问题

本文介绍python实现最小二乘拟合/平方逼近问题的方法

最优平方逼近问题的定义为：

可以使用正规方程组求解：

解方程即得$(c_0,c_1,...,c_n)$，即为拟合式的系数。

当拟合式为多项式时：

正规方程组可以转化为：

用python语言描述正规方程组的求解过程如下：

def mox(n):
    sum = 0
    for i in x:
        sum += pow(i, n)
    return sum

def moy(n):
    sum = 0
    for i, j in zip(x, y):
        sum += pow(i, n) * j
    return sum

def p(i):
    p = 0
    for idx, k in enumerate(c):
        p += k * pow(i, idx)
    return p

def compute_error(c):
    error = 0
    for i, j in zip(x, y):
        error += pow(p(i)-j, 2)
    return math.sqrt(error)

def compute(n):
    global c
    A = np.zeros((n, n))
    for i in range(n):
        for j in range(n):
            A[i, j] = A[j, i] = mox(i+j)
    b = np.zeros(n)
    for i in range(n):
        b[i] = moy(i)
    print('A = ', A)
    print('b = ', b)
    c = np.matmul(np.linalg.inv(A), b)
    print('c = ', c)
    error = compute_error(c)
    print('error = ', error)

使用一道例题来检验一下：

输入例题：

 x = [0.25*i for i in range(5)]
 y = [1.0000, 1.2840, 1.6487, 2.1170, 2.7183]
 compute(n=3)

输出结果：

A =  [[5.        2.5       1.875    ]
 [2.5       1.875     1.5625   ]
 [1.875     1.5625    1.3828125]]
b =  [8.768     5.4514    4.4015375]
c =  [1.00513714 0.86418286 0.84365714]
error =  0.016556949339433698

拟合曲线和数据点的可视化结果为：

和手算的结果对比可发现他们是一致的：