在 python 上实现永无止境的矩阵公式以进行曲线拟合

Question

我正在尝试编写一个可以求解一般回归公式的程序： So I'm trying to implement this matrix equation，有没有办法做到这一点，比如让用户决定它可以有多大，而不用我做越来越多的 if 条件（所以只有一段代码可以折叠成用户希望的矩阵)? 代码：

    #Solving the general matrix for the coefficients
if 3 == n:
    a = np.array([[np.sum(np.multiply(FL[1],FL[1])),np.sum(np.multiply(FL[1],FL[2]))],
                    [np.sum(np.multiply(FL[1],FL[2])),np.sum(np.multiply(FL[2],FL[2]))]])
    b = np.array([np.sum(np.multiply(FL[0],FL[1])),np.sum(np.multiply(FL[0],FL[2]))])
    x = np.linalg.solve(a, b)
if 4 == n:
    a = np.array([[np.sum(np.multiply(FL[1],FL[1])),np.sum(np.multiply(FL[1],FL[2])),np.sum(np.multiply(FL[1],FL[3]))],
                    [np.sum(np.multiply(FL[1],FL[2])),np.sum(np.multiply(FL[2],FL[2])),np.sum(np.multiply(FL[2],FL[3]))],
                    [np.sum(np.multiply(FL[1],FL[3])),np.sum(np.multiply(FL[2],FL[3])),np.sum(np.multiply(FL[3],FL[3]))]])
    b = np.array([np.sum(np.multiply(FL[0],FL[1])),np.sum(np.multiply(FL[0],FL[2])),np.sum(np.multiply(FL[0],FL[3]))])
    x = np.linalg.solve(a, b)

1 这段代码中Phi_0对应FL[i=1]，FL[0]对应y.

Answer 1

您可以使算法独立于多项式的阶数。最简单的方法是使用 for 循环，尽管它们会很慢（因为它们不利用 NumPy 的矢量化）。
这是一个带有随机数据的可重现示例：

import numpy as np
# Order of polynomial
n = 5
# Random seed for reproducibility
np.random.seed(1)
# Input arrays
phi = np.random.random((100,n))
y = np.random.random(100)
# Output arrays
a = np.zeros((n,n))
b = np.zeros(n)
for i in range(n):
    b[i] = np.sum(y * phi[:,i])
    for j in range(i,n):
        # Exploit that matrix is diagonal
        a[i,j] = a[j,i] = np.sum(phi[:,i] * phi[:,j])
# Coefficients array
x = np.linalg.solve(a,b)

在 python 上实现永无止境的矩阵公式以进行曲线拟合

Implementing a never ending matrix formula on python for curve fitting

python

iteration

recursion

curve-fitting