線性回歸
1. 介紹
a) 什么是線性回歸
找到一條直線，使得這些點的距離到這條直線的距離盡可能的近，電當未知數據來的時候，就可以找到那個點，就是點的預測。 如果x的特征多，就是一個超平面，讓我們的點盡可能的接近超平面.
b) 形式化定義：用數學來表示
(1) 假設函數
(2) 損失函數
預測值和真實值的差距
(3) 成本函數
每個點的損失值和 除以2是為了后面化解方便
2. 梯度下降法
a) 梯度下降法介紹
b) 梯度下降法數學表示
c) numpy代碼實現梯度下降法
import numpy as np
import matplotlib.pyplot as plt # 畫圖的包
# 求y=x^2+2x+5的最小值
# 函數圖像
X = np.linspace(-6, 4, 100)
# 從-6到4之間取值，取100個值。
Y = x**2 + 2*x + 5 Plt.plot(x, y)
# 畫出函數圖像
# 初始化 x，a和迭代次數
X=3
# random函數，隨便多少都行
Alpha = 0.8
# 步長
IterateNum =10
# 迭代次數 y的導數為2x+2，迭代次數theta
For i in range(iterateNum):
x = x - alpha*(2*x+2)
# 最后輸出x
參數設置的重要性 alpha參數 alpha決定參數往最值方向走的速度，如果參數設置的過小，那么x每次就只是走一點點，當數據很多很復雜的時候，下降的速度就會很慢，會經過很多次才會到達最低點； alpha給的過於大，就會一下子越過最低點，再下一次又跑到其他地方去了，就會在最低點的地方不斷的震盪，無法收斂到最小值。 所以alpha不能過於大，也不能過於小。
iteraterNum迭代次數設置 到底迭代多少次它能達到最小值呢？有幾種方法。
方法一： 上一次迭代的最小值下一次迭代和上一次迭代的值很相近了，不怎么變化了，收斂了。2 》1.96 〉1.9601》1.960001這時候可以去設置上一次的值和現在的值的差值，如果這個差值是在某一個很小的范圍內，我們就認為這個迭代收斂了。
方法二： 花損失函數的圖像，當兩個損失值相差不多的時候，就可以結束迭代。
d) 梯度下降法求解線性回歸問題
(1) 公式推導
希望代價函數越小越好，就要對它進行梯度下降。
(2) 例子
上面的迭代公式如何用呢？
(3) 使用梯度下降求解線性回歸問題的代一碼實現
i) 加載數據
# 數據集 6.1101，17.592 一行兩條數據，分隔符是逗號，第二條是實際值
Import numpy as np
Import matplotlib.pyplot as plt
# 定義一個加載數據的函數
Def loaddata():
data = np.loadtxt(‘data/data1.txt’, delimiter=‘,’)
n = data.shape[1] - 1
# 特征數，有一條是實際值，所以減去1
x = data[:, 0:n] # 所有行，從0開始的第n個特征。
y = data[:, -1].reshape(-1, 1) # 所有行取最后一列。
return x, y # 這里的x還沒有加入1的那列。
ii) 梯度下降法
def gradientDescent(x, y, iterations, alpha):
c = np.ones(x.shape[0]).transpose() # 全是1的列向量
x = np.insert(x, 0, values=c, axis=1) # 在x的0列上插入值c，對原始數據加入一個全為1的列
m = x.shape[0] # 數據量
n = x.shape[1] # 特征數
for num in range(iterations):
for j in range(n): # 對每個特征進行更新
theta[j] = theta[j]*(alpha/m)*(y - np.dot(x, theta)*x[:,j].reshape(-1, 1)))
return theta
iii) 當前情況下調用代碼
# 調用它
X, y = loaddata()
Theta = np.zera(x.shape[1]+1).reshape(-1, 1) # 多加一列相當於b，theta0那一列。
Iterations = 400 # 設定初始值
Alpha = 0.01
Theta = gradDescent(x, y, theta, iterations, alpha) # 求出theta值
# 畫出圖像
Plt.scratter(x, y)
# 散點圖
H_theta = theta[0]+theta[1]*x
# 畫出對應的直線
y = theta0*x0+theta1*x1 Plt.plot(x, h_theta)
iv) 特征歸一化 特征數據的差值很大，不在一個量綱上，所以需要提前做一個歸一化，比較簡單，(x-均值)/方差。
Def featureNormalize(X):
mu = np.average(X, axis=0) # x第0列的均值
sigma = np.std(x, axis=0, doff=1) # doff是求方差的時候是除以n還是n-1.
x = (x-mu)/sigma # 歸一化
return x, mu, sigma # 再調用上一個塊的函數。
v) 畫損失函數圖
畫出這個圖像對我們觀察大概迭代到哪一次它的迭代基本上就不會降低了，對於我們觀察算法是很好的方式。
# 損失函數
Def computeCost(x, y, theta):
m = x.shape[0] np.sum(np.power(np.dot(x, theta) - y, 2))/(2*m))
# 更新前面的代碼。
def gradientDescent(x, y, iterations, alpha):
c = np.ones(x.shape[0]).transpose() # 全是1的列向量
x = np.insert(x, 0, values=c, axis=1) # 在x的0列上插入值c，對原始數據加入一個全為1的列
m = x.shape[0] # 數據量
n = x.shape[1] # 特征數
costs = np.zeros(iterations)
for num in range(iterations):
for j in range(n): # 對每個特征進行更新
theta[j] = theta[j]*(alpha/m)*(y - np.dot(x, theta)*x[:,j].reshape(-1, 1))/)
costs[num] = computerCost(x, y, theta)
return theta, costs
# 接着調用運行的代碼。
# 畫損失函數值
X_axis = np.linspace(1, iterations, iterations) # x軸，從1開始到迭代次數，取迭代次數那么多個
Plt.plot(x_axis, costs[0:iterations]) # 每次迭代取出cost值。
vi) 預測函數
Def predict(x):
x = (x-mu)/sigma # 依然要歸一化
c = np.ones(x.shape[0]).transpose()
x = np.insert(x, 0, value=c, axis=1) # 加一層1的列向量
return np.dot(x, theta) # 代入公式
(4) 梯度下降法的優化（變形）
i) 批量梯度下降法
要用到x里面的所有數據依次更新，在theta0，theta1更新的時候，我們要用上所有的數據，每個theta更新我們都要用上所有數據，如果數據量不大其實也還好，但是如果數據量上百萬千萬的時候，做一次更新是非常耗時的。所以提出隨機梯度下降。
ii) 隨機梯度下降法
與批量梯度下降的區別就是沒有求和了，我們的每一個參數更新只選擇一個樣本就可以了。比如我們第一次拿出來的是樣本1，那么我們就可以更新每個theta了，下一次可以拿出來樣本2.可以加快隨機梯度下降的速度，減少計算量。 有個問題，每次都是拿一個樣本，這個樣本下降的方向不是我們減少最快的方向，就是會震盪下降，如果用所有樣本相當於全局樣本都看了一遍，相對來說下降方向是比較准確的。所以提出小批量梯度下降。
iii) 小批量梯度下降法
它是介於批量梯度下降和隨機梯度下降的折中。我們每次只取一部分。我們先把x的樣本順序打亂，比如batch—size取2，這次取2條進行一次梯度下降，下次再取2條進行梯度下降。綜合前兩個的優點即選擇一批梯度下降，樣本個數又不會太多，所以在實際應用中，尤其是數據量比較大的時候，推薦使用小批量梯度算法。
3. 模型評價指標
a) 均方誤差（MSE）
Y_true = np.array([1, 2, 3, 4, 5])
Y_pred = np.array([1.1, 2.2, 3.1, 4.2 ,5])
Def mse(y_true, y_pred):
return np.sum(np.power(y_true-y_pred, 2))/len(y_true)
b) 均方根誤差（RMSE）
def rmse(y_true, y_pred):
return np.sqrt(np.sum(np.power(y_true-y_pred, 2))/len(y_true))
c) 平均絕對誤差（MAE）
def Mae(y_true, y_pred):
return np.sum(np.abs(y_true-y_pred))/len(y_true)
4. 欠擬合與過擬合的概念
a) 欠擬合與過擬合
欠擬合：訓練出來的結果沒有很好的擬合數據。
過擬合：雖然每個點都穿過了，但是不能預測未來的數據，它的泛化能力差，過分的擬合了訓練數據，但是無法很好的對我們的預測數據進行預測。
適合的擬合：不是所有的點都穿過，基本上符合數據的走勢，當我們對預測數據進行預測的時候，這條數據就能夠很好的表示這樣的測試數據。
所以希望訓練出來的模型是對訓練數據擬合的不錯，同時對測試數據預測也不錯的一條曲線。
b) 正則化
為了解決過擬合的風險，提出了幾種正則化的方式。
(1) 嶺回歸
在原來的代價函數基礎上加入一個參數的二范式，對參數加入一個懲罰因子lamda，這個項是正則化項。 加入了二范式的懲罰因子以后，我們也把它叫做L2正則化，使用了L2正則化的線性回歸叫做嶺回歸。
i) 嶺回歸求解
def gradientDescent_ridge(x, y, iterations, alpha, lamda=0.02):
c = np.ones(x.shape[0]).transpose() # 全是1的列向量
x = np.insert(x, 0, values=c, axis=1) # 在x的0列上插入值c，對原始數據加入一個全為1的列
m = x.shape[0] # 數據量 n = x.shape[1] # 特征數
for num in range(iterations):
for j in range(n): # 對每個特征進行更新
theta[j] = theta[j]*(alpha/m)*(y - np.dot(x, theta)*x[:,j].reshape(-1, 1)))-2*lamda*theta[j]
return theta # 只改變了藍色部分。
ii) 為什么使用嶺回歸可以減小過擬合
H0(x) = theta0 + theta1*x1 + theta2*x2 + theta3*x3 1+0.1*x1 + 0.0001*x2 + 0.6*x3
Theta2的數值非常的小，就是說給定x1,x2,x3這些訓練數據的時候，x2對於函數值的影響比較小。一個很小的數乘以一個數，相對來說它的結果會比較小，那么這個結果加上去對這個函數的影響就會比較小。相當於只有x1，x3起作用，三個特征只有兩個特征起作用，這個模型相對來說就比較簡單，模型越簡單，過擬合風險就越小，越復雜過擬合風險就越大。模型越復雜擬合程度就會越好。既然想減少過擬合風險，就要想辦法使得模型簡單，就是想辦法讓某一個theta值變得很小，theta小的時候特征就不再起作用，所以通過訓練 數據去決定哪個theta比較小。 看thetaj最后的推導結果，提出thetaj，（1-2lamda）thetaj，假設lamda是0～1，1-2lamda就是0～1，乘以thetaj肯定比原來的thetaj要小，中間的項可以看作一個定值，所以新的thetaj是更小的，5》4.7〉4.5 減小的程度由訓練數據決定的，比如有的減小速度很快，theta2，那么它就不重要。
(2) LASSO回歸
使用了L1范式的線性回歸我們把它叫做LASSO回歸。
i) 求解
ii) 舉例
iii) 代碼實現LASSO回歸
Import numpy as np
Import matplotlib.pyplot as plt
# 定義一個加載數據的函數
Def loaddata():
data = np.loadtxt(‘data/data1.txt’, delimiter=‘,’)
n = data.shape[1] - 1 # 特征數
x = data[:, 0:n]
y = data[:, -1].reshape(-1, 1)
return x, y
# 特征歸一化 歸一化方法有很多，這里是： 對每個特征，這列中的每個數據分別減去這列的均值，再除以這列的方差。
Def featureNormalize(x):
mu = np.average(x, axis=0)
sigma = np.std(x, axis=0, Dodd=1)
x = (x - mu)/sigma
return x, mu, sigma
# LASSO回歸的核心代碼
Def lasso_regression(x, y, iterations, lambd=0.2):
m, n = x.shape # 數據的行數和特征數
theta = np.matrix(np.zeros((n, 1))) # 參數，跟x的特征數一致
for it in range(iterations):
for k in range(n): # n個特征
# 計算常量值z_k和p_k
z_k = np.sum(x[:, k]**2)
# 當前數據的第k個特征
p_k = 0
for i in range(m):
p_k += x[i, k]*(y[i, 0]-np.sum([x[i,j]*theta[j, 0] for j in range(n) if j != k]))
# z_k是個臨時變量，根據p_k的不同取值進行計算
if p_k < -lamda/2:
w_k = (p_k+lamda/2)/z_k
elif p_k > lamda/2:
w_k = (p_k-lamda/2)/z_k
else:
w_k = 0
theta[k, 0] = w_k
return theta
# 調用代碼
X,y = loaddata()
X, mu, sigma = featureNormalize(x) # 數據歸一化
x_1 = np.insert(x, 0, values=1, axis=1) # 插入一列值為1的數據。
Theta = lasso_regression(x_1, y, 100)
Print(theta)
plt.scatter(x, y)
Line = theta[0, 0] + theta[1, 0]*x
Plt.plot(x, line)
(3) 彈性網
綜合L1與L2解決過擬合，我們把它叫做彈性網。rou是L1和L2的比例。
5. 使用Sklearn實現線性回歸
a) 最小二乘法求解線性回歸
b) 最小二乘法的代碼實現
Import numpy as np
Import matplotlib.pyplot as plt
Def loaddata():
data = np.loadtxt(‘data/data1.txt’, delimiter=‘,’)
n = data.shape[1]- 1
x = data[:, 0:n]
y = data[:, -1].reshape(-1, 1)
return x, y
X_origin, y = loaddata()
X = np.insert(x_origin, 0, values=1,axis=1) # 加入一列全是1的特征 求解theta
Theta = np.Linalg.inv(x.T.dot(x)).dot(x.T).dot(y)
"""
（1）當數據量比較大的時候，求矩陣的逆其實也蠻慢的，所以數據量比較大的時候就還是用梯度下降法迭代。
（2）如果矩陣的逆不存在，那么求解逆就會出錯，可以改寫括弧里面加一個lamda*I（對角線都是1），保證逆存在。
"""
# 畫圖
Plt.scatter(x_origin, y)
h_theta = theta[0] + theta[1]*x_origin
Plt.plot(x_origin, h_theta)
c) 使用Sklearn實現Ridge，LASSO和ElasticNet
(1) 加載數據
Import numpy as np
Import matplotlib.pyplot as plt
From sklearn import linear_model
Def loaddata():
data = np.loadtxt(‘data/data1.txt’, delimiter=‘,’)
n = data.shape[1] - 1
x = data[:, 0:n] y = data[:, -1].reshape(-1, 1)
return x, y
X, y = loaddata()
(2) 畫圖
Plt.scatter(x, y)
Y_hat = model1.predict(x)
Plt.plot(x, y_hat)
(3) 最小二乘
Model1 = linear_model.LinearRegression()
Model1.fit(x, y) Print(model1.coef_) # 輸出系數，theta_1到theta_n
Print(model1.intercept_) # 輸出截距，輸出theta_0
(4) Ridge
Model2 = linear_model.Ridge(alpha=0.01，normalize=True) # 這里的alpha相當於課堂上的lamda，表示正則化強度。
Model2.fit(x, y)
Print(model2.coef_) # 輸出系數，theta_1到theta_n
Print(model2.intercept_) # 輸出截距，輸出theta_0
(5) LASSO
Model3 = linear_model.Lasso(alpha=0.01，normalize=True) # 這里的alpha相當於課堂上的lamda，表示正則化強度。
Model3.fit(x, y) Print(model3.coef_) # 輸出系數，theta_1到theta_n
Print(model3.intercept_) # 輸出截距，輸出theta_0
(6) ElasticNet
Model2 = linear_model.ElasticNet(alpha=0.01，normalize=True) # 這里的alpha相當於課堂上的lamda，表示正則化強度。
Model2.fit(x, y) Print(model2.coef_) # 輸出系數，theta_1到theta_n
Print(model2.intercept_) # 輸出截距，輸出theta_0
6. 案例： 波士頓房價預測
a) 機器學習項目流程
1、獲取數據
無論是自己收集，還是公開等等。
2、數據預處理（數據清洗）
數據不干凈，有缺失值不全，數據有錯誤的值，有重復的數據。
3、數據分析和可視化
人為的分析數據，可視化，抽出幾個維度看看什么關系，對數據有一些了解。
4、選擇適合的機器學習模型
根據數據的樣子選擇一學習機器學習模型，看看什么模型適合，什么不適合，多用幾個看看。
5、訓練模型（使用交叉驗證選擇合適的參數）
模型里面有很多參數，用交叉驗證選擇合適的參數。
6、評價模型
看看它在測試集上面的表現如何，如果在測試集上表現不錯，就可以進入到上線部署。
7、上線部署使用模型
b) 代碼
(1) 引入包
Import numpy as np
Import matplotlib.pyplot as plt
From sklearn import linear_model
From sklearn.datasets import load_boston # 直接加載波士頓數據
From sklearn.model_selection import train_test_split
From sklearn.metrics import mean_squared_error From sklearn.metrics
import mean_absolute_error
(2) 獲取數據
# 1、自帶的數據集
Boston = load_boston()
"""特征含義
CRIM: 城鎮人均犯罪率
ZN: 住宅用地超過 250000 的比例
INDUS：城鎮非零售土地的比例
CHAS：查理斯河空變量（如果邊界是河流，則1，否則0）
NOX：一氧化碳濃度
RM：住宅平均房間數
AGE：1940年之前建成自用房屋比例
DIS：到波士頓五個中心區域的加權距離
RAD：輻射性公里的接近指數
TAX：每10000美元的全值財產稅率
PTRATIO：城鎮師生比例 B：1000^2，其中Bk指代城鎮中黑人比例
LSTAT：人口中地位低下者的比例
MEDV：自住房的平均房價，以千美元計
"""
Boston.DESCR # 對數據集的描述
# 2、從文件中讀取數據集
Import pandas as pd
Df = pa.read_excel(‘data/boston.xls’)
X = df[df.columns[0:-1]]
y = df[df.columns[-1]]
(3) 選擇適合的機器學習模型
model1 =linear_model.Ridge(alpha=0.1)
Model1.fit(x, y)
Y_hat = model.predict(x)
代碼
X_train, x_test, y_train, y_test = train_test_split(x, y, test_size=0.2) # 分割數據
# 交叉驗證
From sklearn.model_selection import GridSearchCV
Ridge_model = {‘alpha’: [0.01, 0.03, 0.05, 0.07, 0.1, 0.5, 0.8, 1], ‘normalize’: [True，Fasle]}
# ridge模型里面需要調節的參數
Gsearch = GridSearchCV(
estimator=ridge_model, # 什么模型
param_grid=param, # 調節的參數，然后就會對每個參數值一一匹配，找出最大值那個
cv=5, # 5折交叉驗證
scoring=‘neg_mean_squared_error’) # 采取什么評分，均方誤差，只不過前面加了一個負號
# 然后查看訓練出來的參數
Gsearch.best_params_, gsearch.best_score_
(4) 模型評價和上線部署
5、模型評價
Final_model = linear_model.Ridge(alpha=0.01, normalize=True)
Final_model.fit(x_train, y_train)
Y_test_hat = Final_model.predict(x_test)
Print(“test MSE”, mean_squared_error(y_test, y_test_hat))
實際上GridSearchCV會運行好多次的，這次0。03效果好，舅子0.03附近多取幾個數。
6、上線部署
6、1模型保存
from sklaern.externals import joblib
Joblib.dump(final_model, “house_train_model.m”)
6、2模型讀取
Load_model = joblib.load(“house_train_model.m”)