上篇部落格寫完之後，終於發現自己線性迴歸入門！

然後洗澡的時候就在想一個問題，線性迴歸會了，寫線性擬合是完全沒問題的，但是np庫的多項式擬合到底是怎麼做出來的呢？突然靈光一閃多項式擬合？多變數的線性迴歸？好像發現了什麼？重新理清一下思路。什麼是多項式擬合？對的，這個問題以前沒有好好思考過，以前的直觀感覺是，不就是給一些x，y的點，讓一個多項式去擬合，使得這個多項式的曲線看起來大致符合那些點。等等，好像忽略了什麼，什麼叫大致符合？對的，什麼叫大致符合？答案已在心中了！就是讓這些點距離這個擬合曲線的和最小！慢著，這其實就是線性迴歸的損失函式的定義嘛！最小化，損失函式！多項式多個引數怎麼確定？多變數線性迴歸！每個一個次方項當成獨立的一個變數！比如a*x^2+b*x+c這個標準的二次項，x^2當成x1，x當成x2，1當成x3！好了，整個思路有了，嘗試一下實現吧！

import numpy as np

def my_fit(x, y, power):
    size = len(x)
    c = np.ones((1, size))
    x_arr = np.array(x)
    x_a=x_arr.copy()
    x_a.resize((1,size))
    x_mat = np.append(x_a, c, axis=0)
    y_arr = np.array(y)
    y_arr.resize((1, size))
    y_mat = np.mat(y_arr)
    for i in range(2, power+1):
        temp_x = x_arr**i
        temp_x.resize((1, size))
        x_mat = np.append(temp_x, x_mat, axis=0)
    x_mat=np.mat(x_mat)
    w = (x_mat * x_mat.T).I * x_mat * y_mat.T
    w0=w.T
    w0.resize(w0.size)
    return w0

def f(x):
    return 2*x ** 2 +x+3

x = np.linspace(-5, 5)
y = f(x) + 0.5*np.random.randn(len(x))  # 加入噪音
y_fit = np.polyfit(x, y, 2)  # 二次多項式擬合
print(y_fit)
w=my_fit(x,y,2)
print(w)
 

[1.99205043 1.00302882 3.08335708]
[1.99205043 1.00302882 3.08335708]

竟然和np庫的結果完全一樣！驚了！雖說實現程式碼不算優美，但是基本思路沒啥區別！懶得調參用梯度下降了，直接解出來。

很好，還有一個疑問，scipy的曲線擬合怎麼做的呢？一樣嘛！當成多變數線性迴歸就好！

code：

import numpy as np
from scipy.optimize import curve_fit


def f_fit(x, a, b, c):
    return a*np.sin(x)+b*x+c


def f_test(x):
    return 2*np.sin(x)+3*x+1


def my_curve_fit1(fit_fun, x, y, p_num):  #p_num引數的個數
    size = len(x)
    if p_num <= 0 or p_num > size:
        print('no parameter to fit')
        return
    test_list = [0]*p_num
    test_list[0] = 1
    x_arr = np.array(fit_fun(x, *test_list))
    x_arr.resize((1, size))
    for i in range(1, p_num):
        test_list[i-1] = 0
        test_list[i] = 1
        temp_x = np.array(fit_fun(x, *test_list))
        temp_x.resize((1, size))
        x_arr = np.append(x_arr, temp_x, axis=0)
    x_mat = np.mat(x_arr)
    y_arr = np.array(y)
    y_arr.resize((1, size))
    y_mat = np.mat(y_arr)
    w = (x_mat * x_mat.T).I * x_mat * y_mat.T
    w0 = w.T
    w0.resize(w0.size)
    return w0

#為了強行接近庫函式的介面,無需引數的個數，但實現有點撈
#很容易知道引數個數多於樣本x的個數無法擬合
def my_curve_fit2(fit_fun, x, y):
    test_list = [0]
    size = len(x)
    p_num = 1  #引數個數預設為1
    for i in range(0, size):
        try:
            fit_fun(x, *test_list)
            break
        except:
            p_num += 1
            test_list.append(0)
    if p_num > size:
        print('can not fit')
        return
    test_list[0] = 1
    x_arr = np.array(fit_fun(x, *test_list))
    x_arr.resize((1, size))
    for i in range(1, p_num):
        test_list[i - 1] = 0
        test_list[i] = 1
        temp_x = np.array(fit_fun(x, *test_list))
        temp_x.resize((1, size))
        x_arr = np.append(x_arr, temp_x, axis=0)
    x_mat = np.mat(x_arr)
    y_arr = np.array(y)
    y_arr.resize((1, size))
    y_mat = np.mat(y_arr)
    w = (x_mat * x_mat.T).I * x_mat * y_mat.T
    w0 = w.T
    w0.resize(w0.size)
    return w0





x = np.linspace(-2*np.pi, 2*np.pi)
y = f_test(x)+0.3*np.random.randn(len(x))   #加入噪音
p_fit, prov = curve_fit(f_fit, x, y)  #曲線擬合
my_fit1 = my_curve_fit1(f_fit, x, y, 3)
my_fit2 = my_curve_fit2(f_fit, x, y)
print('sicpy庫的曲線擬合')
print('a,b,c', p_fit)
print('手寫曲線擬合1')
print('a,b,c', my_fit1)
print('手寫曲線擬合2')
print('a,b,c', my_fit2)
 

sicpy庫的曲線擬合
a,b,c [1.99131057 3.00473361 1.05173153]
手寫曲線擬合1
a,b,c [1.99131057 3.00473361 1.05173153]
手寫曲線擬合2
a,b,c [1.99131057 3.00473361 1.05173153]

結果竟然也就是一樣！厲害了，竟然一樣。好了，實現的小竅門用到了列表引數傳遞，想單獨提出來的那個項，就令該項的引數為1，其他引數為0，然後就把每個項都單獨提出來。一用多變數線性迴歸，搞掂！然後為了強行接近scipy庫的曲線擬合的介面，不知道有幾個引數擬合，怎麼辦呢？好像沒有什麼好的方法，只好用try來試了，反正傳的引數個數不對就會報錯的，就一個一個引數往上加咯。

以上。