重拾基礎 - 線性回歸

05-24

重拾基礎 - 線性回歸

背景

好像該找實習了，建了一個倉庫：

Ceruleanacg/Descent?

github.com

專門寫寫預習筆記，預習演算法實現什麼的，嗯。

所以第一篇就是線性回歸了。

問題設定

已知有 $N$ 個 $x, y$ 對構成數據集 $X, Y$ ，他們在坐標軸上的分布如下圖：

%matplotlib inlineimport matplotlib.pyplot as pltimport numpy as np# 生成100對x, ydata_count = 100x_data = np.linspace(-20, 20, data_count)y_data = np.multiply(2, x_data) + 3 + np.random.normal(loc=0, scale=8.0, size=(data_count,))plt.figure(figsize=(16, 6))plt.scatter(x_data, y_data, s=10, color=g)plt.xlabel(x)plt.ylabel(y)plt.show()

我們希望找到一個函數：

$h(x) = wx+b$

這個函數會儘可能的擬合數據集 $X, Y$ ，為了做到這點，我們希望這個函數 $h(x)$ 在 $X$ 上每一個取值 $x_i$ 的函數值 $h(x_i)$ 與 $Y$ 上每一個對應的 $y_i$ 的平方差儘可能小。進一步，我們可以寫出平方損失函數：

$loss(w, b) = frac{1}{N}sum^{N}_{i=0}(wx_i+b-y_i)^2$

我們希望找到一組 $w, b$ ，能使得 $loss(w, b)$ 最小，我們畫出損失函數的圖像如下所示：

from mpl_toolkits.mplot3d import Axes3Dw_sample = np.linspace(-10, 10, data_count).reshape((-1, 1))b_sample = np.linspace(-10, 10, data_count).reshape((-1, 1))x_data = x_data.reshape((-1, 1))y_data = y_data.reshape((-1, 1))loss = np.square(np.dot(w_sample, x_data.T) + b_sample - y_data) / data_countfigure = plt.figure(figsize=(16, 6))axes = Axes3D(figure)axes.set_xlabel(w)axes.set_ylabel(b)axes.plot_surface(w_sample.T, b_sample, loss, cmap=rainbow)<mpl_toolkits.mplot3d.art3d.Poly3DCollection at 0x10e83f588>

梯度下降

我們通過梯度下降找到目標 $w, b$ ，我們先隨機初始化一對 $w_0, b_0$ ，並在坐標軸上表示它：

w = np.random.rand()b = np.random.rand()y_predict = w * x_data + bplt.figure(figsize=(16, 6))plt.scatter(x_data, y_data, s=10, color=g)plt.plot(x_data, y_predict)plt.xlabel(x)plt.ylabel(y)plt.show()

函數的負梯度方向是函數值下降最快的方向，為此，我們首先對 $w, b$ 求其偏微分：

$egin{aligned} frac{partial loss(w, b)}{partial w} &= frac{2}{N}sum^{N}_{i=0}(wx_i+b-y_i)cdot x_i, \ frac{partial loss(w, b)}{partial b} &= frac{2}{N}sum^{N}_{i=0}(wx_i+b-y_i) end{aligned}$

然後，我們通過下式在每次迭代中更新 $w, b$ ：

$egin{aligned} w_{t+1} &= w_t - eta frac{partial l(w_t, b_t)}{partial w_t} \ b_{t+1} &= b_t - eta frac{partial l(w_t, b_t)}{partial b_t} end{aligned}$

其中， $eta$ 是學習率，我們迭代3000次，繪製出損失函數值的變化和迭代結束後 $w, b$ 的方程：

w_cache, b_cache, l_cache, = [], [], []# 迭代3000次for iteration in range(3000): y_predict = w * x_data + b diff = y_predict - y_data error = np.sum(np.square(diff)) / data_count grad_w = np.mean(diff * x_data) grad_b = np.mean(diff) w -= 0.003 * grad_w b -= 0.003 * grad_b w_cache.append(w) b_cache.append(b) l_cache.append(error)y_predict = w * x_data + bplt.figure(figsize=(16, 6))plt.scatter(x_data, y_data, s=10, color=g)plt.plot(x_data, y_predict)plt.title(y=2x+3)plt.xlabel(x)plt.ylabel(y)plt.show()

事情就這麼成了。

後續

爭取一天一篇。