您好,登錄后才能下訂單哦!
這篇文章給大家介紹python中怎么推導(dǎo)線性回歸模型,內(nèi)容非常詳細(xì),感興趣的小伙伴們可以參考借鑒,希望對大家能有所幫助。
首先,先看一張圖:
圖是我們在初中學(xué)習(xí)過的直角坐標(biāo)系二維平面,上面遍布著一些點(diǎn)。從整體趨勢看,y隨x的增大而增大。如果曾經(jīng)你和我一樣,數(shù)學(xué)每次考試都是90的話,那么接下來,我相信你會情不自禁地做一件事:
沒錯(cuò),我們會以(0,0)和(10,10)為兩點(diǎn),畫出一條貫穿其中的線,從視覺上,這條紅線正好把所有點(diǎn)一分為二,其對應(yīng)的數(shù)學(xué)表達(dá)式為:
y=x
而這就是我們線性回歸所要做的事:找到一組數(shù)學(xué)表達(dá)式(圖中的紅線),用來反映數(shù)據(jù)(圖中的點(diǎn))的變化規(guī)律。
目標(biāo)有了,問題也來了:
貫穿圖中密密麻麻點(diǎn)的線有無數(shù)條,為什么不是y=2x,y=x+1,偏偏是y=x呢?
我們又是通過何種方法去找到這條線呢?
先解決第一個(gè)問題,上天書:
這個(gè)式子就是第一個(gè)問題的解,沒見過的符號太多,看不懂是吧?那么我來翻譯一下:
先求出(每個(gè)點(diǎn)的Y值-以每個(gè)點(diǎn)的X值通過函數(shù)求出的Y值)的平方
求和;
乘以1/2
再通俗點(diǎn):
把每個(gè)點(diǎn)的實(shí)際y值與它通過某個(gè)函數(shù)求出的y值的差的平方加起來,再乘以1/2。
而文章開篇中的均方差損失,MSE,平方損失函數(shù),二次代價(jià)函數(shù)其實(shí)都指的是它。這個(gè)式子其實(shí)計(jì)算的是真實(shí)值和用函數(shù)預(yù)測的值之間的誤差之和。那么第一個(gè)問題就迎刃而解了:哪一個(gè)表達(dá)式所求出的誤差和最小,就是我們要找的那條“紅線”。
我們繼續(xù)解決第二個(gè)問題,先上圖:
這個(gè)問題還要簡單,我們只要從斜率為0的那條“紅線”(y=0*X)開始畫線,然后一點(diǎn)點(diǎn)增大斜率,每條線求一個(gè)誤差值,找出其中誤差值最小的那條線,就大功告成了。而中間有著巨大計(jì)算量的遍歷過程,我們可以通過python,瞬間完成。
————————
線性回歸的Python實(shí)現(xiàn)
————————
重點(diǎn):梯度下降法!
導(dǎo)入一些包,待用:
import numpy as np
import pandas as pd
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import MinMaxScaler
import matplotlib.pyplot as plt
import seaborn as sns
sns.set_context('notebook')
sns.set_style('white')
導(dǎo)入案例數(shù)據(jù):
model_data = pd.read_csv('model_data.csv',engine='python')
model_data.head()
數(shù)據(jù)是一份上海的房價(jià)數(shù)據(jù),我們要把房屋價(jià)格作為因變量y,房屋面積,房間數(shù)附近餐飲POI數(shù)量,評論,距離市中心距離等作為自變量,擬合一個(gè)線性回歸模型,用于預(yù)測房價(jià)。
根據(jù)要求提取自變量和因變量:
feat_cols = model_data.columns.tolist()[1:]
print(feat_cols)
X = model_data[feat_cols].values
y = model_data['價(jià)格'].values
構(gòu)建損失函數(shù):
def Cost_Function(X,y,theta):
'''
需要傳入的參數(shù)為
X:自變量
y:應(yīng)變量
theta:權(quán)重
使用均方誤差(MSE),作為損失函數(shù)
'''
m = y.size #求出y的個(gè)數(shù)(一共多少條數(shù)據(jù))
t = X.dot(theta) #權(quán)重和變量點(diǎn)乘,計(jì)算出使用當(dāng)前權(quán)重時(shí)的預(yù)測值
c = 0 #定義損失值
c = 1/(2*m) * (np.sum(np.square(t-y)))
#預(yù)測值與實(shí)際值的差值,平方后除以數(shù)據(jù)的條數(shù),計(jì)算出均方誤差。最后乘以1/2(無實(shí)際意義,方便以后計(jì)算)
return c #返回?fù)p失值
*θ為每個(gè)變量前的權(quán)重,什么是權(quán)重?比如y=2x,2就是自變量x的權(quán)重
求損失值我們就用先前說到的損失函數(shù)。如果你夠仔細(xì),可能會有一個(gè)問題,我們的損失函數(shù)前需要乘以一個(gè)1/2,似乎沒有特別的意義。恭喜你很機(jī)智,1/2的確沒有任何意義,只是為了接下來方便求導(dǎo)。
構(gòu)建梯度下降法:
def GradientDescent(X, y, feat_cols, alpha=0.3, num_iters=10000):
'''
需傳入?yún)?shù)為
X:自變量
y:應(yīng)變量
feat_cols:變量列表
alpha:學(xué)習(xí)率,默認(rèn)0.3
num_iters:迭代次數(shù),默認(rèn)10000次
使用梯度下降法迭代權(quán)重
'''
scaler = MinMaxScaler() #最大最小值歸一化自變量
X = scaler.fit_transform(X) #歸一化
m = y.size #求出y的個(gè)數(shù)(一共多少條數(shù)據(jù))
J_history = np.zeros(num_iters) #創(chuàng)建容納每次迭代后損失值得矩陣,初始值為0
theta = np.zeros(len(feat_cols)+1) #設(shè)置默認(rèn)權(quán)重,0
for iter in np.arange(num_iters): #根據(jù)迭代次數(shù),開始迭代
t = X.dot(theta) #權(quán)重和變量點(diǎn)乘,計(jì)算出使用當(dāng)前權(quán)重時(shí)的預(yù)測值
theta = theta - alpha*(1/m)*(X.T.dot(t-y))
#對代價(jià)函數(shù)求導(dǎo),算出下降最快的方向,乘以學(xué)習(xí)率(下降的速度),再用原來的權(quán)重相減,得到新的權(quán)重
J_history[iter] = Cost_Function(X, y, theta) #求出新的權(quán)重時(shí)的損失值,存入矩陣
return(theta, J_history) #返回最終的權(quán)重和歷次迭代的損失值
這是構(gòu)造模型最為核心的部分。我們不斷迭代,尋找最優(yōu)的那條“紅線”的過程,其實(shí)是在不斷調(diào)整每個(gè)自變量的權(quán)重。而每個(gè)權(quán)重每次到底怎么調(diào)整,增大還是減?。ǚ较颍?,這就需要我們對損失函數(shù)求導(dǎo)。
如果數(shù)學(xué)不好,不理解,我們用圖來說明一下:
好比,我們站在懸崖頂端,要找到最快能達(dá)到懸崖底部的方向,那么顯而易見,你所在位置最陡峭的方向,就是正確的方向,而求導(dǎo)就是找到最陡峭的方向(切線斜率絕對值最大的點(diǎn))。
山坡是凹凸不平的,所以我們每走一步都需要重新尋找方向,這就是迭代的過程;其次,每次的步子也不能跨太大,萬一跨錯(cuò)地方了,不好糾正,所以我們又需要設(shè)置一個(gè)步子的大小——學(xué)習(xí)率。
所以梯度下降法的公式就是:
每一次更新的權(quán)重= 前一次的權(quán)重-學(xué)習(xí)率*損失函數(shù)的導(dǎo)數(shù)。
在理解了下山這個(gè)場景以后,我們就能順利的完成梯度下降法的構(gòu)建,并且通過python函數(shù)求出最后每個(gè)變量的權(quán)重和每次迭代過后的損失值。
構(gòu)建繪制損失值變化圖的函數(shù):
def plot_Cost(GD_result):
'''
繪制權(quán)重變化情況
需傳入?yún)?shù)為
GD_result:梯度下降法結(jié)果
'''
theta , Cost = GD_result #得到權(quán)重和損失值
print('theta: ',theta.ravel()) #打印權(quán)重
plt.plot(Cost) #繪制損失值變化情況
plt.title('COST change')
plt.ylabel('Cost')
plt.xlabel('Iterations')
plt.grid()
plt.show()
這個(gè)很簡單,就是通過前面梯度下降法求得的歷次迭代后的損失值,畫出變化曲線。
最后把所有函數(shù)匯總,就是我們的線性回歸模型了:
def lr_function(X,y,feat_cols):
'''
需要輸入的變量為
X:自變量
y:應(yīng)變量
feat_cols:變量列表
'''
def score(y_p,y):
'''
y_p:預(yù)測值
y:真實(shí)值
dimension:樣本數(shù)量
計(jì)算R^和調(diào)整R^
'''
aa=y_p.copy(); bb=y.copy()
if len(aa)!=len(bb):
print('not same length')
return np.nan
cc=aa-bb
wcpfh=sum(cc**2) #誤差平方和
# RR means R_Square
RR=1-sum((bb-aa)**2)/sum((bb-np.mean(bb))**2)
return RR#返回R^
X = np.c_[np.ones(X.shape[0]),X]
GD_result = GradientDescent(X, y, feat_cols)
plot_Cost(GD_result)
y_p = np.dot(X,GD_result[0])
RR = score(y_p,y)
return RR,y_p,GD_result[0] #返回R^,預(yù)測值
一般對于每個(gè)機(jī)器學(xué)習(xí)模型,都需要有一個(gè)指標(biāo)衡量其擬合程度,而線性模型我們使用的是我們所熟知的可決系數(shù)R^2。為了求出R^2,我在函數(shù)中又套用了一個(gè)簡單的求解函數(shù),具體過程不贅述了,通讀代碼就能明白。通常R^2越接近1,表示模型擬合程度越好。
模型封裝完畢,下面是見證奇跡的時(shí)刻!
model_result = lr_function(X,y,feat_cols)
print('R^2為:{}'.format(round(model_result[0],4)))
通過模型,我們求出了每個(gè)自變量的權(quán)重,圖表反應(yīng)了損失值由大變小的過程,在10000次迭代的過程中,一開始速度很快,越到后面越趨于平緩。
最后是R^2為0.70,有70%的擬合度,尚可。
————————
線性回歸模型的驗(yàn)證
————————
為了驗(yàn)證我們自己編寫的模型是否準(zhǔn)確,我們也可以使用python機(jī)器學(xué)習(xí)工具包sklearn,對同樣的數(shù)據(jù),用線性回歸模型擬合,查看最后的R^2是否一致。
先對變量標(biāo)準(zhǔn)化:
scaler = MinMaxScaler()
X = scaler.fit_transform(X)
使用LinearRegression()進(jìn)行擬合,并求出R^2:
lr = LinearRegression()
lr.fit(X,y)
R2 = lr.score(X,y)
print('R^2為:{}'.format(round(R2,4)))
R^2同樣為0.7,代表我們自己編寫的模型沒有問題。
最后,我們繪制一張真實(shí)值與預(yù)測值對比圖,可視化模型結(jié)果:
關(guān)于python中怎么推導(dǎo)線性回歸模型就分享到這里了,希望以上內(nèi)容可以對大家有一定的幫助,可以學(xué)到更多知識。如果覺得文章不錯(cuò),可以把它分享出去讓更多的人看到。
免責(zé)聲明:本站發(fā)布的內(nèi)容(圖片、視頻和文字)以原創(chuàng)、轉(zhuǎn)載和分享為主,文章觀點(diǎn)不代表本網(wǎng)站立場,如果涉及侵權(quán)請聯(lián)系站長郵箱:is@yisu.com進(jìn)行舉報(bào),并提供相關(guān)證據(jù),一經(jīng)查實(shí),將立刻刪除涉嫌侵權(quán)內(nèi)容。