一元線性回歸的Python實現

来源:https://www.cnblogs.com/hznudmh/archive/2022/06/17/16382786.html
-Advertisement-
Play Games

1 問題的提出 對於給定的數據集 \(D = \{(x_1,y_1),(x_2,y_2),\cdots,(x_m,y_m)\}\),線性回歸 (linear regression) 試圖學得一個線性模型以儘可能準確地預測是指輸出標記. 2 原理 設給定的數據集 \(D = \{(x_i,y_i)\} ...


目錄

1 問題的提出

對於給定的數據集 \(D = \{(x_1,y_1),(x_2,y_2),\cdots,(x_m,y_m)\}\),線性回歸 (linear regression) 試圖學得一個線性模型以儘可能準確地預測是指輸出標記.

2 原理

設給定的數據集 \(D = \{(x_i,y_i)\}_{i=1}^m, \ x_i,y_i \in \mathcal{R}\). 對於離散屬性,如果屬性值間存在“序”(order)的關係,可通過連續化將其轉化為連續值,例如二值屬性“身高”的取值“高”“矮”可轉化為 \(\{1.0,0.0\}\),三值屬性“高度”的取值“高”“中”“低”可轉化為 \(\{1.0,0.5,0.0\}\);若屬性之間不存在有序關係,假定有 \(k\) 個屬性值,則通常轉化為 \(k\) 維向量,例如屬性“瓜類”的取值“西瓜”“南瓜”“黃瓜”可轉化為 \((0,0,1),(0,1,0),(1,0,0)\).

線性回歸試圖學得

\[f(x_i) = wx_i + b, \ s.t. f(x_i) \simeq y_i \tag{1} \]

2.1 代價函數

我們需要確定 \(w\)\(b\) 的值使得 \(f(x)\)\(y\) 之間的差別儘可能小. 因此我們引入均方誤差,均方誤差是回歸問題中最常用的性能度量工具,我們可以試圖讓均方誤差最小化,即

\[\begin{aligned} (w^*,b^*) &= \arg\min\limits_{(w,b)} \sum_{i=1}^m(f(x_i)-y_i)^2\\ &=\arg\min\limits_{(w,b)} \sum_{i=1}^m(y_i-wx_i-b)^2 \end{aligned} \tag{2} \]

\(\arg\min\limits_{\theta}f(x;\theta)\) 意思就是找出一個 \(\theta\) 使得 \(f(x;\theta)\) 的值最小,即他的返回值是 \(f(x;\theta)\) 的最小值所對應的 \(\theta\) 的值.

均方誤差有很好的幾何意義,它對應了常用的“歐式距離”(Euclidean distance),基於均方誤差最小化的進行模型求解的方法稱為“最小二乘法”(least square method). 線上性回歸中,最小二乘法就是試圖尋找一條直線,使得所有樣本點到直線的歐氏距離之和最小.

在均方誤差的基礎上進一步構造代價函數

\[J(w,b) = \frac{1}{2m}\sum_{i=1}^m(f(x_i)-y_i)^2 = \frac{1}{2m}\sum_{i=1}^m(wx_i + b - y_i)^2 \]

這裡分母的 \(2\) 是為了後續求導的方便

求解 \(w\)\(b\) 使 \(J(w,b)\) 最小化的過程,稱為線性回歸模型的最小二乘“參數估計”(parameter estimation). 我們可將 \(J(w,b)\) 分別對 \(w\)\(b\) 求導,得到

\[\begin{aligned} &\frac{\partial J(w,b)}{\partial w} = \frac{1}{m}\sum_{i=1}^m(wx_i+b-y_i)x_i\\ &\frac{\partial J(w,b)}{\partial b} = \frac{1}{m}\sum_{i=1}^m(wx_i+b-y_i) \end{aligned} \]

令上兩式等於0,解得

\[w = \frac{\sum_{i=1}^my_i(x_i-\bar{x})}{\sum_{i=1}^mx_i^2-\frac{1}{m}(\sum_{i=1}^mx_i)^2}\\ b = \frac{1}{m}\sum_{i=1}^m(y_i-wx_i) \]

\(J(w,b)\) 是關於 \(w,b\) 的凸函數,根據凸函數的性質,其偏導為 \(0\) 時就是 \(w\)\(b\) 的最優解.

其中 \(\bar{x} = \frac{1}{m}\sum_\limits{i=1}^mx_i\)\(x\) 的均值.

2.2 模型的評價

2.2.1 皮爾遜相關係數

使用相關係數衡量線性相關性的強弱,皮爾遜相關係數的公式如下:

\[r_{xy} = \frac{COV(X,Y)}{\sqrt{Var(X)Var(Y)}} \]

相關度越高,皮爾遜相關係數的值就趨於 1 或 -1 (趨於 1 表示它們呈正相關, 趨於 -1 表示它們呈負相關);如果相關係數等於0,表明它們之間不存線上性相關關係.

2.2.2 決定繫數

決定繫數 \(R^2\) 也稱擬合優度,反應了 \(y\) 的波動有多少百分比能被 \(x\) 的波動所描述. 決定繫數越接近 1 ,說明擬合程度越好.

總平方和

\[SST = \sum_{i=1}^n(y_i-\bar{y})^2 \]

回歸平方和

\[SSR = \sum_{i=1}^n(\hat{y}_i-\bar{y})^2 \]

殘差平方和

\[SSE = \sum_{i=1}^n(y_i-\hat{y}_i)^2 \]

其中

\[SST = SSR + SSE\\ R^2 = \frac{SSR}{SST} = 1- \frac{SSE}{SST} \]

3 Python 實現

3.1 不調sklearn庫

Step1. 調庫

import numpy as np
from numpy import genfromtxt
import matplotlib.pyplot as plt

Step2. 數據導入並繪製散點圖

data = genfromtxt("data.csv", delimiter = ",")
x = data[:, 0, np.newaxis]
y = data[:, 1, np.newaxis]
plt.scatter(x, y)
plt.show()

Step3. 求回歸繫數
根據先前的推導,已經知道

\[w = \frac{\sum_{i=1}^my_i(x_i-\bar{x})}{\sum_{i=1}^mx_i^2-\frac{1}{m}(\sum_{i=1}^mx_i)^2}\\ b = \frac{1}{m}\sum_{i=1}^m(y_i-wx_i) \]

m = len(x)
x_bar = np.mean(x)
w = (np.sum((x - x_bar)*y))/(np.sum(x**2)-(1/m)*(np.sum(x))**2)
b = np.mean(y-w*x)
print(w,b)

1.3224310227553517 7.991020982270779

Step4. 擬合圖像

plt.plot(x, y, 'b.')
plt.plot(x, w*x+b, 'r')
plt.show()

Step5. 計算相關係數和決定繫數

COVxy = np.cov(x.T,y.T)
r = COVxy[0,1]/(x.std()*y.std())
print(r)

0.78154393928063

y_hat = w*x+b
SSR = np.sum((y_hat-np.mean(y))**2)
SST = np.sum((y-np.mean(y))**2)
R2 = SSR/SST
print(R2)

0.5986557915386548

3.2 調 sklearn 庫

建模:

from sklearn.linear_model import LinearRegression
model = LinearRegression()
model.fit(x, y)

LinearRegression()

擬合圖像的得出

plt.plot(x, y, 'b.')
plt.plot(x, model.predict(x), 'r')
plt.show()

回歸繫數

w = model.intercept_
b = model.coef_
print("截距為 {0}, 回歸繫數為 {1}".format(w[0], b[0][0]))

截距為 7.991020982270385, 回歸繫數為 1.32243102275536

決定繫數

model.score(x, y)

0.598655791538662

4 梯度下降法

4.1 原理

由於代價函數是凸函數,因此只有全局最小值,梯度下降法的原理是先定一個初始值,然後利用導數就像階梯一樣慢慢逼近全局最小值

圖片出處:https://zhuanlan.zhihu.com/p/36564434

已知代價函數 \(J(w,b)\),我們需要找一組 \(w,b\) 使得 \(J(w,b)\) 最小,給定一個演算法:

\[\begin{aligned} &給定 \ w,b \ 的初始值 \ w_0,b_0\\ \\ &repeat until convergence \ \{\\ &\quad\quad temp0 = w - \alpha \frac{\partial}{\partial w}J(w,b)\\ &\quad\quad temp1 = b - \alpha \frac{\partial}{\partial b}J(w,b)\\ &\quad\quad w = temp0\\ &\} \end{aligned} \]

其中 \(\alpha\) 為學習率,學習率不能太大也不能太小,可以多次嘗試 \(0.1,0.03,0.01,0.003,0.001,0.0003,\cdots\).

根據已知條件,有

\[J(w,b) = \frac{1}{2m}\sum_{i=1}^m(wx_i + b - y_i)^2\\ \frac{\partial J(w,b)}{\partial w} = \frac{1}{m}\sum_{i=1}^m(wx_i+b-y_i)x_i\\ \frac{\partial J(w,b)}{\partial b} = \frac{1}{m}\sum_{i=1}^m(wx_i+b-y_i) \]

於是

\[w = w - \alpha \frac{1}{m}\sum_{i=1}^m(wx_i+b-y_i)\\ b = b - \alpha \frac{1}{m}\sum_{i=1}^m(wx_i+b-y_i)x_i \]

4.2 Python實現

# 學習率learning rate
lr = 0.0001
# 截距初值
b = 0
# 斜率初值
w = 0
# 最大迭代次數
epochs = 50
# 最小二乘法
def compute_error(b, w, x, y):
  totalError = 0
  for i in range(0, len(x)):
    totalError += (y[i] - (w * x[i] + b)) ** 2
  return totalError / float(len(x)) / 2.0

def gradient_descent_runner(x, y, b, w, lr, epochs):
  # 計算總數據量
  m = float(len(x))
  # 迴圈epochs次
  for i in range(epochs):
    b_grad = 0
    w_grad = 0
    # 計算梯度的總和再求平均
    for j in range(0, len(x)):
      b_grad += -(1/m) * (y[j] - ((w * x[j]) + b))
      w_grad += -(1/m) * x[j] * (y[j] - ((w * x[j]) + b))
     
    # 更新 b 和 w
    b = b - (lr * b_grad)
    w = w - (lr * w_grad)
    
    # 每迭代5次,輸出一次圖像
    if i % 5 == 0:
      print('epochs:', i)
      plt.plot(x, y, 'b.')
      plt.plot(x, w*x + b, 'r')
      plt.show()
  return b,w
print('Starting b = {0}, w = {1}, error = {2}'.format(b, w, compute_error(b, w, x, y)))
print('Running')
b, w = gradient_descent_runner(x, y, b, w, lr, epochs)
print('After {0} iterations b = {1}, w = {2}, error = {3}'.format(epochs, b, w, compute_error(b,w,x,y)))
# 畫圖
plt.plot(x, y, 'b.')
plt.plot(x, w * x + b, 'r')
plt.show()

迭代50次後得到 \(b = 0.03207192, w = 1.47886174\),最小二乘誤差為 56.3244305.

參考

[1]. 周志華.《機器學習》.清華大學出版社
[2]. https://www.bilibili.com/video/BV1Rt411q7WJ?p=4&vd_source=08d2535b05740d396ec0fc720c52f36a


您的分享是我們最大的動力!

-Advertisement-
Play Games
更多相關文章
  • 這裡給大家分享我在OpenLayers 地圖開發工作中總結出的一下代碼和註意點,希望對大家有所幫助 效果如下: 核心代碼展示:附帶講解註釋 var map = new ol.Map({ // 初始化地圖 target: 'map',// 選擇地圖對象 layers: [ new ol.layer.T ...
  • springboot整合ueditor實現圖片上傳和文件上傳功能 寫在前面: 在閱讀本篇之前,請先按照我的這篇隨筆完成對ueditor的前期配置工作: springboot+layui 整合百度富文本編輯器ueditor入門使用教程(踩過的坑)https://www.cnblogs.com/rain ...
  • 1、前言 作為嵌入式軟體開發,可能經常會使用命令行或者顯示屏等設備實現人機交互的功能,功能中通常情況都包含 UI 菜單設計;很多開發人員都會有自己的菜單框架模塊,防止重覆造輪子,網上有很多這種菜單框架的代碼,但是大多耦合性太強,無法獨立出來適配不同的菜單設計。 本文介紹一個降低了耦合性,完全獨立的菜 ...
  • 目錄 一.簡介 二.效果演示 三.源碼下載 四.猜你喜歡 零基礎 OpenGL (ES) 學習路線推薦 : OpenGL (ES) 學習目錄 >> OpenGL ES 基礎 零基礎 OpenGL (ES) 學習路線推薦 : OpenGL (ES) 學習目錄 >> OpenGL ES 轉場 零基礎 O ...
  • 分類方式 按參數分: 有參構造(預設構造) & 無參構造 按類型分: 普通構造 & 拷貝構造 調用方式 括弧法 顯示法 隱式轉換法 PS:下方所有文本均以此代碼為基礎 1 class Person { 2 public: 3 //無參構造函數 4 Person() { 5 std::cout << ...
  • 引言:今天工作遇到了一個需要按行讀取txt文件數據的需求,查詢了一下getline()函數,發現這竟然是一個C++的標準庫函數,而且設計的很好,特地做一下記錄。getline本質是一個定界流輸入截取函數,預設是換行符‘/n’ 個人技術博客(文章整理+源碼): https://zobolblog.gi ...
  • 大佬的理解->《IO流和File》 1、File類 File類是IO包中唯一代表磁碟文件本身的對象,File類定義了一些與平臺無關的方法來操作文件。通過調用File類提供的各種方法,能夠完成創建、刪除文件、重命名文件、判斷文件的讀寫許可權許可權是否存在、設置和查詢文件的最近修改時間等操作。 ​ File ...
  • 前言 大部分情況下,地理繪圖可使用 Arcgis 等工具實現。但正版的 Arcgis 並非所有人可以承受。本文基於 Python 的 cartopy 和 matplotlib 等庫,為地理空間繪圖的代碼實現提供參考。 所有所需庫如下: gma、cartopy、matplotlib、numpy 更多內 ...
一周排行
    -Advertisement-
    Play Games
  • Timer是什麼 Timer 是一種用於創建定期粒度行為的機制。 與標準的 .NET System.Threading.Timer 類相似,Orleans 的 Timer 允許在一段時間後執行特定的操作,或者在特定的時間間隔內重覆執行操作。 它在分散式系統中具有重要作用,特別是在處理需要周期性執行的 ...
  • 前言 相信很多做WPF開發的小伙伴都遇到過表格類的需求,雖然現有的Grid控制項也能實現,但是使用起來的體驗感並不好,比如要實現一個Excel中的表格效果,估計你能想到的第一個方法就是套Border控制項,用這種方法你需要控制每個Border的邊框,並且在一堆Bordr中找到Grid.Row,Grid. ...
  • .NET C#程式啟動閃退,目錄導致的問題 這是第2次踩這個坑了,很小的編程細節,容易忽略,所以寫個博客,分享給大家。 1.第一次坑:是windows 系統把程式運行成服務,找不到配置文件,原因是以服務運行它的工作目錄是在C:\Windows\System32 2.本次坑:WPF桌面程式通過註冊表設 ...
  • 在分散式系統中,數據的持久化是至關重要的一環。 Orleans 7 引入了強大的持久化功能,使得在分散式環境下管理數據變得更加輕鬆和可靠。 本文將介紹什麼是 Orleans 7 的持久化,如何設置它以及相應的代碼示例。 什麼是 Orleans 7 的持久化? Orleans 7 的持久化是指將 Or ...
  • 前言 .NET Feature Management 是一個用於管理應用程式功能的庫,它可以幫助開發人員在應用程式中輕鬆地添加、移除和管理功能。使用 Feature Management,開發人員可以根據不同用戶、環境或其他條件來動態地控制應用程式中的功能。這使得開發人員可以更靈活地管理應用程式的功 ...
  • 在 WPF 應用程式中,拖放操作是實現用戶交互的重要組成部分。通過拖放操作,用戶可以輕鬆地將數據從一個位置移動到另一個位置,或者將控制項從一個容器移動到另一個容器。然而,WPF 中預設的拖放操作可能並不是那麼好用。為瞭解決這個問題,我們可以自定義一個 Panel 來實現更簡單的拖拽操作。 自定義 Pa ...
  • 在實際使用中,由於涉及到不同編程語言之間互相調用,導致C++ 中的OpenCV與C#中的OpenCvSharp 圖像數據在不同編程語言之間難以有效傳遞。在本文中我們將結合OpenCvSharp源碼實現原理,探究兩種數據之間的通信方式。 ...
  • 一、前言 這是一篇搭建許可權管理系統的系列文章。 隨著網路的發展,信息安全對應任何企業來說都越發的重要,而本系列文章將和大家一起一步一步搭建一個全新的許可權管理系統。 說明:由於搭建一個全新的項目過於繁瑣,所有作者將挑選核心代碼和核心思路進行分享。 二、技術選擇 三、開始設計 1、自主搭建vue前端和. ...
  • Csharper中的表達式樹 這節課來瞭解一下表示式樹是什麼? 在C#中,表達式樹是一種數據結構,它可以表示一些代碼塊,如Lambda表達式或查詢表達式。表達式樹使你能夠查看和操作數據,就像你可以查看和操作代碼一樣。它們通常用於創建動態查詢和解析表達式。 一、認識表達式樹 為什麼要這樣說?它和委托有 ...
  • 在使用Django等框架來操作MySQL時,實際上底層還是通過Python來操作的,首先需要安裝一個驅動程式,在Python3中,驅動程式有多種選擇,比如有pymysql以及mysqlclient等。使用pip命令安裝mysqlclient失敗應如何解決? 安裝的python版本說明 機器同時安裝了 ...