本文討論了確定簡單（配對）迴歸線的數學方程式的幾種方法。
這裡討論的所有求解方程式的方法都基於最小平方法。 讓我們將這些方法表示如下：
- 解析解
- 梯度下降
- 隨機梯度下降
對於每種求解直線方程式的方法，文章都提供了各種函數，主要分為不使用函式庫編寫的函數 數字貨幣 以及那些用於計算的 數字貨幣。 相信熟練運用 數字貨幣 將降低計算成本。
文章中給出的所有程式碼都是用該語言編寫的 python 2.7 同 Jupyter筆記本。 原始碼和包含範例資料的檔案發佈在
本文更針對初學者和那些已經逐漸開始掌握人工智慧非常廣泛的部分（機器學習）的研究的人。
為了說明本資料，我們使用一個非常簡單的例子。
條件範例
我們有五個表徵依賴的價值觀 Y 從 X （表 1）：
表 1“範例條件”
我們假設這些值 是一年中的月份，並且 ——本月收入。 換句話說，收入取決於一年中的月份，並且 - 收入所依賴的唯一標誌。
這個例子馬馬虎虎，無論是從收入對一年中月份的條件依賴的角度來看，還是從值的數量來看——它們都很少。 然而，正如他們所說，這種簡化將有可能解釋初學者所吸收的材料，但並不總是那麼容易。 而且數字的簡單性將使那些希望在紙上解決範例的人無需花費大量的勞動力。
讓我們假設範例中給出的相關性可以透過以下形式的簡單（成對）迴歸線的數學方程式很好地近似：
哪裡 是收到收入的月份， — 對應月份的收入， и 是估計線的迴歸係數。
注意係數 通常稱為估計線的斜率或坡度； 代表的金額 當它改變時 .
顯然，我們在範例中的任務是在方程中選擇這樣的係數 и ，其中我們按月計算的收入值與真實答案的偏差，即樣本中呈現的值將是最小的。
最小平方法
根據最小二乘法，應透過對其進行平方來計算偏差。 如果偏差具有相反的符號，則此技術可讓您避免相互抵消。 例如，如果在一種情況下，偏差是 +5 （加五），以及在其他 -5 （負五），那麼偏差總和將相互抵消並等於 0（零）。 可以不對偏差進行平方，而是利用模數的性質，然後所有偏差都將是正值並且會累積。 關於這一點，我們不做詳細闡述，只是簡單說明，為了計算方便，習慣上對偏差平方。
這就是我們用來確定最小偏差平方和（誤差）的公式：
哪裡 是真實答案的近似函數（即我們計算的收入），
是真實答案（樣本中提供的收入），
是樣本索引（確定偏差的月份數）
讓我們對函數進行微分，定義偏微分方程，並準備好繼續求解解析解。 但首先，讓我們先簡單了解微分是什麼，並記住導數的幾何意義。
差異化
微分是求函數導數的運算。
衍生性商品有什麼用？ 函數的導數表徵了函數的變化率並告訴我們它的方向。 如果給定點的導數為正，則函數會增大；否則，函數會減少。 且絕對導數的值越大，函數值的變化率就越高，函數圖形的斜率就越陡。
例如，在直角座標系條件下，M(0,0)點的導數值等於 + 25 意味著在給定點，當值發生變化時 向右按常規單位，值 增加25個常規單位。 在圖表上，數值看起來相當陡峭的上升 從給定點。
另一個例子。 導數值相等 -0,1 意味著當移位時 每一個常規單位，價值 僅減少0,1個常規單位。 同時，在函數圖上，我們可以觀察到幾乎不明顯的向下斜率。 與一座山進行類比，就好像我們從一座山上慢慢地走下緩坡，不像前面的例子，我們必須爬上非常陡峭的山峰:)
因此，對函數求導後 按賠率 и ，我們定義一階偏微分方程。 確定方程式後，我們將得到一個由兩個方程式組成的系統，透過求解該方程組，我們將能夠選擇係數的值 и ，給定點處對應導數的值變化非常非常小，而在解析解的情況下根本不變化。 換句話說，找到的係數處的誤差函數將達到最小值，因為這些點的偏導數值將等於零。
因此，根據微分規則，關於係數的一階偏導方程 將採用以下形式：
一階偏導方程式關於 將採用以下形式：
結果，我們得到了一個具有相當簡單解析解的方程組：
開始{方程式*}
開始{案例}
na + bsumlimits_{i=1}^nx_i — sumlimits_{i=1}^ny_i = 0
sumlimits_{i=1}^nx_i(a +bsumlimits_{i=1}^nx_i — sumlimits_{i=1}^ny_i) = 0
結束{案例}
結束{方程式*}
在解方程式之前，我們先預先加載，檢查加載是否正確，並格式化資料。
載入和格式化數據
應該注意的是，由於對於解析解以及隨後的梯度和隨機梯度下降，我們將使用兩種變體的程式碼：使用庫 數字貨幣 如果不使用它，那麼我們將需要適當的資料格式（請參閱程式碼）。
資料載入和處理程式碼
# импортируем все нужные нам библиотеки
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import math
import pylab as pl
import random
# графики отобразим в Jupyter
%matplotlib inline
# укажем размер графиков
from pylab import rcParams
rcParams['figure.figsize'] = 12, 6
# отключим предупреждения Anaconda
import warnings
warnings.simplefilter('ignore')
# загрузим значения
table_zero = pd.read_csv('data_example.txt', header=0, sep='t')
# посмотрим информацию о таблице и на саму таблицу
print table_zero.info()
print '********************************************'
print table_zero
print '********************************************'
# подготовим данные без использования NumPy
x_us = []
[x_us.append(float(i)) for i in table_zero['x']]
print x_us
print type(x_us)
print '********************************************'
y_us = []
[y_us.append(float(i)) for i in table_zero['y']]
print y_us
print type(y_us)
print '********************************************'
# подготовим данные с использованием NumPy
x_np = table_zero[['x']].values
print x_np
print type(x_np)
print x_np.shape
print '********************************************'
y_np = table_zero[['y']].values
print y_np
print type(y_np)
print y_np.shape
print '********************************************'可視化
現在，在我們首先載入數據，其次檢查載入的正確性並最後格式化資料之後，我們將進行第一次視覺化。 經常使用的方法是 配對圖 圖書館 海生。 在我們的範例中，由於數量有限，使用該庫沒有意義 海生。 我們將使用常規庫 Matplotlib 只需查看散點圖即可。
散點圖程式碼
print 'График №1 "Зависимость выручки от месяца года"'
plt.plot(x_us,y_us,'o',color='green',markersize=16)
plt.xlabel('$Months$', size=16)
plt.ylabel('$Sales$', size=16)
plt.show()圖1“收入對月份的依賴性”
解析解
讓我們使用最常用的工具 蟒蛇 並求解方程組：
開始{方程式*}
開始{案例}
na + bsumlimits_{i=1}^nx_i — sumlimits_{i=1}^ny_i = 0
sumlimits_{i=1}^nx_i(a +bsumlimits_{i=1}^nx_i — sumlimits_{i=1}^ny_i) = 0
結束{案例}
結束{方程式*}
根據克萊默法則 我們將找到一般行列式以及以下行列式 並由 ，然後將行列式除以 一般行列式 - 求係數 ，類似地我們求出係數 .
解析解程式碼
# определим функцию для расчета коэффициентов a и b по правилу Крамера
def Kramer_method (x,y):
# сумма значений (все месяца)
sx = sum(x)
# сумма истинных ответов (выручка за весь период)
sy = sum(y)
# сумма произведения значений на истинные ответы
list_xy = []
[list_xy.append(x[i]*y[i]) for i in range(len(x))]
sxy = sum(list_xy)
# сумма квадратов значений
list_x_sq = []
[list_x_sq.append(x[i]**2) for i in range(len(x))]
sx_sq = sum(list_x_sq)
# количество значений
n = len(x)
# общий определитель
det = sx_sq*n - sx*sx
# определитель по a
det_a = sx_sq*sy - sx*sxy
# искомый параметр a
a = (det_a / det)
# определитель по b
det_b = sxy*n - sy*sx
# искомый параметр b
b = (det_b / det)
# контрольные значения (прооверка)
check1 = (n*b + a*sx - sy)
check2 = (b*sx + a*sx_sq - sxy)
return [round(a,4), round(b,4)]
# запустим функцию и запишем правильные ответы
ab_us = Kramer_method(x_us,y_us)
a_us = ab_us[0]
b_us = ab_us[1]
print ' 33[1m' + ' 33[4m' + "Оптимальные значения коэффициентов a и b:" + ' 33[0m'
print 'a =', a_us
print 'b =', b_us
print
# определим функцию для подсчета суммы квадратов ошибок
def errors_sq_Kramer_method(answers,x,y):
list_errors_sq = []
for i in range(len(x)):
err = (answers[0] + answers[1]*x[i] - y[i])**2
list_errors_sq.append(err)
return sum(list_errors_sq)
# запустим функцию и запишем значение ошибки
error_sq = errors_sq_Kramer_method(ab_us,x_us,y_us)
print ' 33[1m' + ' 33[4m' + "Сумма квадратов отклонений" + ' 33[0m'
print error_sq
print
# замерим время расчета
# print ' 33[1m' + ' 33[4m' + "Время выполнения расчета суммы квадратов отклонений:" + ' 33[0m'
# % timeit error_sq = errors_sq_Kramer_method(ab,x_us,y_us)這是我們得到的：
因此，係數的值已經找到，偏差平方和也已經確定。 讓我們根據找到的係數在散射直方圖上畫一條直線。
迴歸線程式碼
# определим функцию для формирования массива рассчетных значений выручки
def sales_count(ab,x,y):
line_answers = []
[line_answers.append(ab[0]+ab[1]*x[i]) for i in range(len(x))]
return line_answers
# построим графики
print 'Грфик№2 "Правильные и расчетные ответы"'
plt.plot(x_us,y_us,'o',color='green',markersize=16, label = '$True$ $answers$')
plt.plot(x_us, sales_count(ab_us,x_us,y_us), color='red',lw=4,
label='$Function: a + bx,$ $where$ $a='+str(round(ab_us[0],2))+',$ $b='+str(round(ab_us[1],2))+'$')
plt.xlabel('$Months$', size=16)
plt.ylabel('$Sales$', size=16)
plt.legend(loc=1, prop={'size': 16})
plt.show()圖 2“正確和計算出的答案”
您可以查看每個月的偏差圖。 在我們的例子中，我們不會從中獲得任何重要的實用價值，但我們將滿足我們的好奇心，即簡單的線性迴歸方程式如何很好地描述收入對一年中月份的依賴性。
偏差圖程式碼
# определим функцию для формирования массива отклонений в процентах
def error_per_month(ab,x,y):
sales_c = sales_count(ab,x,y)
errors_percent = []
for i in range(len(x)):
errors_percent.append(100*(sales_c[i]-y[i])/y[i])
return errors_percent
# построим график
print 'График№3 "Отклонения по-месячно, %"'
plt.gca().bar(x_us, error_per_month(ab_us,x_us,y_us), color='brown')
plt.xlabel('Months', size=16)
plt.ylabel('Calculation error, %', size=16)
plt.show()第 3 號圖表“偏差，%”
並不完美，但我們完成了任務。
讓我們寫一個函數來確定係數 и 使用圖書館 數字貨幣，更準確地說，我們將編寫兩個函數：一個使用偽逆矩陣（在實踐中不推薦，因為該過程計算複雜且不穩定），另一個使用矩陣方程式。
分析解決方案代碼 (NumPy)
# для начала добавим столбец с не изменяющимся значением в 1.
# Данный столбец нужен для того, чтобы не обрабатывать отдельно коэффицент a
vector_1 = np.ones((x_np.shape[0],1))
x_np = table_zero[['x']].values # на всякий случай приведем в первичный формат вектор x_np
x_np = np.hstack((vector_1,x_np))
# проверим то, что все сделали правильно
print vector_1[0:3]
print x_np[0:3]
print '***************************************'
print
# напишем функцию, которая определяет значения коэффициентов a и b с использованием псевдообратной матрицы
def pseudoinverse_matrix(X, y):
# задаем явный формат матрицы признаков
X = np.matrix(X)
# определяем транспонированную матрицу
XT = X.T
# определяем квадратную матрицу
XTX = XT*X
# определяем псевдообратную матрицу
inv = np.linalg.pinv(XTX)
# задаем явный формат матрицы ответов
y = np.matrix(y)
# находим вектор весов
return (inv*XT)*y
# запустим функцию
ab_np = pseudoinverse_matrix(x_np, y_np)
print ab_np
print '***************************************'
print
# напишем функцию, которая использует для решения матричное уравнение
def matrix_equation(X,y):
a = np.dot(X.T, X)
b = np.dot(X.T, y)
return np.linalg.solve(a, b)
# запустим функцию
ab_np = matrix_equation(x_np,y_np)
print ab_np讓我們比較一下確定係數所花費的時間 и ，依照所提出的 3 種方法。
計算計算時間的程式碼
print ' 33[1m' + ' 33[4m' + "Время выполнения расчета коэффициентов без использования библиотеки NumPy:" + ' 33[0m'
% timeit ab_us = Kramer_method(x_us,y_us)
print '***************************************'
print
print ' 33[1m' + ' 33[4m' + "Время выполнения расчета коэффициентов с использованием псевдообратной матрицы:" + ' 33[0m'
%timeit ab_np = pseudoinverse_matrix(x_np, y_np)
print '***************************************'
print
print ' 33[1m' + ' 33[4m' + "Время выполнения расчета коэффициентов с использованием матричного уравнения:" + ' 33[0m'
%timeit ab_np = matrix_equation(x_np, y_np)
對於少量數據，「自寫」函數出現在前面，它使用 Cramer 方法找到係數。
現在您可以繼續使用其他方法來尋找係數 и .
梯度下降
首先，我們來定義什麼是漸層。 簡單來說，梯度就是表示函數最大成長方向的一段。 類比爬山，坡度面向的地方就是登頂山頂最陡的地方。 以山為例，我們記得實際上我們需要最陡的下降才能盡快到達低地，即最小值 - 函數不增加或減少的地方。 此時導數將等於零。 因此，我們不需要梯度，而是反梯度。 要找到反梯度，您只需將梯度乘以 -1 （減一）。
讓我們注意這樣一個事實：一個函數可以有多個最小值，並且使用下面提出的演算法下降到其中一個最小值時，我們將無法找到另一個最小值，可能低於找到的最小值。 放心吧，這對我們沒有威脅！ 在我們的例子中，我們正在處理一個最小值，因為我們的函數 圖上是一條正規拋物線。 正如我們從學校數學課程中應該知道的那樣，拋物線只有一個最小值。
當我們發現為什麼需要梯度，並且梯度是一個段，即給定座標的向量，它們是完全相同的係數之後 и 我們可以實現梯度下降。
在開始之前，我建議先閱讀幾句話有關下降演算法的內容：
- 我們以偽隨機方式決定係數的座標 и 。 在我們的範例中，我們將確定接近零的係數。 這是一種常見的做法，但每種情況可能都有自己的做法。
- 從座標 減去該點的一階偏導數的值 。 因此，如果導數為正，則函數會增大。 因此，透過減去導數的值，我們將向成長的相反方向移動，即向下降的方向移動。 如果導數為負，則此時的函數會減少，透過減去導數的值，我們將向下降方向移動。
- 我們用座標進行類似的操作 ：減去該點的偏導數值 .
- 為了不跳過最小值而飛入深空，需要設定下降方向的步長。 一般來說，您可以寫整篇文章來介紹如何正確設定步長以及如何在下降過程中更改步長以減少計算成本。 但現在我們面前的任務略有不同，我們將使用「戳」的科學方法，或者正如他們通常所說的那樣，憑經驗確定步長。
- 一旦我們從給定的座標開始 и 減去導數的值，我們得到新的座標 и 。 我們已經根據計算出的座標進行下一步（減法）。 如此循環一次又一次地開始，直到達到所需的收斂。
全部！ 現在我們準備出發尋找馬裡亞納海溝最深的峽谷。 讓我們開始吧。
梯度下降的程式碼
# напишем функцию градиентного спуска без использования библиотеки NumPy.
# Функция на вход принимает диапазоны значений x,y, длину шага (по умолчанию=0,1), допустимую погрешность(tolerance)
def gradient_descent_usual(x_us,y_us,l=0.1,tolerance=0.000000000001):
# сумма значений (все месяца)
sx = sum(x_us)
# сумма истинных ответов (выручка за весь период)
sy = sum(y_us)
# сумма произведения значений на истинные ответы
list_xy = []
[list_xy.append(x_us[i]*y_us[i]) for i in range(len(x_us))]
sxy = sum(list_xy)
# сумма квадратов значений
list_x_sq = []
[list_x_sq.append(x_us[i]**2) for i in range(len(x_us))]
sx_sq = sum(list_x_sq)
# количество значений
num = len(x_us)
# начальные значения коэффициентов, определенные псевдослучайным образом
a = float(random.uniform(-0.5, 0.5))
b = float(random.uniform(-0.5, 0.5))
# создаем массив с ошибками, для старта используем значения 1 и 0
# после завершения спуска стартовые значения удалим
errors = [1,0]
# запускаем цикл спуска
# цикл работает до тех пор, пока отклонение последней ошибки суммы квадратов от предыдущей, не будет меньше tolerance
while abs(errors[-1]-errors[-2]) > tolerance:
a_step = a - l*(num*a + b*sx - sy)/num
b_step = b - l*(a*sx + b*sx_sq - sxy)/num
a = a_step
b = b_step
ab = [a,b]
errors.append(errors_sq_Kramer_method(ab,x_us,y_us))
return (ab),(errors[2:])
# запишем массив значений
list_parametres_gradient_descence = gradient_descent_usual(x_us,y_us,l=0.1,tolerance=0.000000000001)
print ' 33[1m' + ' 33[4m' + "Значения коэффициентов a и b:" + ' 33[0m'
print 'a =', round(list_parametres_gradient_descence[0][0],3)
print 'b =', round(list_parametres_gradient_descence[0][1],3)
print
print ' 33[1m' + ' 33[4m' + "Сумма квадратов отклонений:" + ' 33[0m'
print round(list_parametres_gradient_descence[1][-1],3)
print
print ' 33[1m' + ' 33[4m' + "Количество итераций в градиентном спуске:" + ' 33[0m'
print len(list_parametres_gradient_descence[1])
print
我們潛入馬裡亞納海溝的最底部，在那裡我們發現了所有相同的係數值 и ，這正是我們所期望的。
讓我們再潛水一次，只是這次，我們的深海飛行器將充滿其他技術，即圖書館 數字貨幣.
梯度下降代碼 (NumPy)
# перед тем определить функцию для градиентного спуска с использованием библиотеки NumPy,
# напишем функцию определения суммы квадратов отклонений также с использованием NumPy
def error_square_numpy(ab,x_np,y_np):
y_pred = np.dot(x_np,ab)
error = y_pred - y_np
return sum((error)**2)
# напишем функцию градиентного спуска с использованием библиотеки NumPy.
# Функция на вход принимает диапазоны значений x,y, длину шага (по умолчанию=0,1), допустимую погрешность(tolerance)
def gradient_descent_numpy(x_np,y_np,l=0.1,tolerance=0.000000000001):
# сумма значений (все месяца)
sx = float(sum(x_np[:,1]))
# сумма истинных ответов (выручка за весь период)
sy = float(sum(y_np))
# сумма произведения значений на истинные ответы
sxy = x_np*y_np
sxy = float(sum(sxy[:,1]))
# сумма квадратов значений
sx_sq = float(sum(x_np[:,1]**2))
# количество значений
num = float(x_np.shape[0])
# начальные значения коэффициентов, определенные псевдослучайным образом
a = float(random.uniform(-0.5, 0.5))
b = float(random.uniform(-0.5, 0.5))
# создаем массив с ошибками, для старта используем значения 1 и 0
# после завершения спуска стартовые значения удалим
errors = [1,0]
# запускаем цикл спуска
# цикл работает до тех пор, пока отклонение последней ошибки суммы квадратов от предыдущей, не будет меньше tolerance
while abs(errors[-1]-errors[-2]) > tolerance:
a_step = a - l*(num*a + b*sx - sy)/num
b_step = b - l*(a*sx + b*sx_sq - sxy)/num
a = a_step
b = b_step
ab = np.array([[a],[b]])
errors.append(error_square_numpy(ab,x_np,y_np))
return (ab),(errors[2:])
# запишем массив значений
list_parametres_gradient_descence = gradient_descent_numpy(x_np,y_np,l=0.1,tolerance=0.000000000001)
print ' 33[1m' + ' 33[4m' + "Значения коэффициентов a и b:" + ' 33[0m'
print 'a =', round(list_parametres_gradient_descence[0][0],3)
print 'b =', round(list_parametres_gradient_descence[0][1],3)
print
print ' 33[1m' + ' 33[4m' + "Сумма квадратов отклонений:" + ' 33[0m'
print round(list_parametres_gradient_descence[1][-1],3)
print
print ' 33[1m' + ' 33[4m' + "Количество итераций в градиентном спуске:" + ' 33[0m'
print len(list_parametres_gradient_descence[1])
print
係數值 и 不可改變的。
讓我們來看看梯度下降過程中誤差如何變化，即偏差平方和如何隨著每一步而變化。
繪製偏差平方和的程式碼
print 'График№4 "Сумма квадратов отклонений по-шагово"'
plt.plot(range(len(list_parametres_gradient_descence[1])), list_parametres_gradient_descence[1], color='red', lw=3)
plt.xlabel('Steps (Iteration)', size=16)
plt.ylabel('Sum of squared deviations', size=16)
plt.show()圖 4“梯度下降期間偏差平方和”
在圖表上我們看到，每一步誤差都會減小，經過一定次數的迭代後，我們觀察到一條幾乎水平的線。
最後，我們來估計一下程式碼執行時間的差異：
確定梯度下降計算時間的程式碼
print ' 33[1m' + ' 33[4m' + "Время выполнения градиентного спуска без использования библиотеки NumPy:" + ' 33[0m'
%timeit list_parametres_gradient_descence = gradient_descent_usual(x_us,y_us,l=0.1,tolerance=0.000000000001)
print '***************************************'
print
print ' 33[1m' + ' 33[4m' + "Время выполнения градиентного спуска с использованием библиотеки NumPy:" + ' 33[0m'
%timeit list_parametres_gradient_descence = gradient_descent_numpy(x_np,y_np,l=0.1,tolerance=0.000000000001)
也許我們做錯了什麼，但這又是一個簡單的「自製」函數，不使用該函式庫 數字貨幣 優於使用函式庫的函數的計算時間 數字貨幣.
但我們並沒有停滯不前，而是正在研究另一種令人興奮的方法來求解簡單的線性迴歸方程式。 見面！
隨機梯度下降
為了快速理解隨機梯度下降的運作原理，最好確定它與普通梯度下降的差異。 我們，在梯度下降的情況下，在導數方程中 и 使用樣本中可用的所有特徵和真實答案的值的總和（即所有特徵的總和 и ）。 在隨機梯度下降中，我們不會使用樣本中存在的所有值，而是偽隨機地選擇所謂的樣本索引並使用其值。
例如，如果索引確定為數字 3（三），則我們取值 и ，然後我們將這些值代入導數方程式並確定新的座標。 然後，確定了座標後，我們再次偽隨機決定樣本索引，將索引對應的值代入偏微分方程中，以新的方式決定座標 и ETC。 直到收斂變成綠色。 乍一看，這似乎根本行不通，但確實行得通。 確實，值得注意的是，誤差並不會隨著每一步而減少，但肯定有一種趨勢。
隨機梯度下降與傳統梯度下降相比有哪些優點？ 如果我們的樣本量非常大並且以數萬個值進行測量，那麼處理隨機的數千個值而不是整個樣本要容易得多。 這就是隨機梯度下降發揮作用的地方。 當然，就我們而言，我們不會注意到太大的差異。
讓我們看一下程式碼。
隨機梯度下降的程式碼
# определим функцию стох.град.шага
def stoch_grad_step_usual(vector_init, x_us, ind, y_us, l):
# выбираем значение икс, которое соответствует случайному значению параметра ind
# (см.ф-цию stoch_grad_descent_usual)
x = x_us[ind]
# рассчитывыаем значение y (выручку), которая соответствует выбранному значению x
y_pred = vector_init[0] + vector_init[1]*x_us[ind]
# вычисляем ошибку расчетной выручки относительно представленной в выборке
error = y_pred - y_us[ind]
# определяем первую координату градиента ab
grad_a = error
# определяем вторую координату ab
grad_b = x_us[ind]*error
# вычисляем новый вектор коэффициентов
vector_new = [vector_init[0]-l*grad_a, vector_init[1]-l*grad_b]
return vector_new
# определим функцию стох.град.спуска
def stoch_grad_descent_usual(x_us, y_us, l=0.1, steps = 800):
# для самого начала работы функции зададим начальные значения коэффициентов
vector_init = [float(random.uniform(-0.5, 0.5)), float(random.uniform(-0.5, 0.5))]
errors = []
# запустим цикл спуска
# цикл расчитан на определенное количество шагов (steps)
for i in range(steps):
ind = random.choice(range(len(x_us)))
new_vector = stoch_grad_step_usual(vector_init, x_us, ind, y_us, l)
vector_init = new_vector
errors.append(errors_sq_Kramer_method(vector_init,x_us,y_us))
return (vector_init),(errors)
# запишем массив значений
list_parametres_stoch_gradient_descence = stoch_grad_descent_usual(x_us, y_us, l=0.1, steps = 800)
print ' 33[1m' + ' 33[4m' + "Значения коэффициентов a и b:" + ' 33[0m'
print 'a =', round(list_parametres_stoch_gradient_descence[0][0],3)
print 'b =', round(list_parametres_stoch_gradient_descence[0][1],3)
print
print ' 33[1m' + ' 33[4m' + "Сумма квадратов отклонений:" + ' 33[0m'
print round(list_parametres_stoch_gradient_descence[1][-1],3)
print
print ' 33[1m' + ' 33[4m' + "Количество итераций в стохастическом градиентном спуске:" + ' 33[0m'
print len(list_parametres_stoch_gradient_descence[1])
我們仔細觀察這些係數，不禁問自己“這怎麼可能？” 我們得到了其他係數值 и 。 也許隨機梯度下降已經為方程式找到了更優化的參數？ 很不幸的是，不行。 只需查看偏差平方和即可發現，使用新的係數值，誤差會更大。 我們並不急於絕望。 讓我們建立一個誤差變化圖。
繪製隨機梯度下降中偏差平方和的程式碼
print 'График №5 "Сумма квадратов отклонений по-шагово"'
plt.plot(range(len(list_parametres_stoch_gradient_descence[1])), list_parametres_stoch_gradient_descence[1], color='red', lw=2)
plt.xlabel('Steps (Iteration)', size=16)
plt.ylabel('Sum of squared deviations', size=16)
plt.show()圖 5“隨機梯度下降過程中的偏差平方和”
看看日程安排，一切都已就位，現在我們將解決所有問題。
所以發生了什麼事？ 發生了以下情況。 當我們隨機選擇一個月時，我們的演算法會針對所選月份尋求減少計算收入的誤差。 然後我們選擇另一個月份並重複計算，但我們減少了第二個選定月份的誤差。 現在請記住，前兩個月明顯偏離簡單線性迴歸方程式的直線。 這意味著，當選擇這兩個月中的任何一個月時，透過減少其中每個月的誤差，我們的演算法會嚴重增加整個樣本的誤差。 那該怎麼辦？ 答案很簡單：你需要減少下降步長。 畢竟，透過減少下降步長，誤差也會停止上下「跳躍」。 或者更確切地說，「跳躍」錯誤不會停止，但不會那麼快發生:)讓我們檢查一下。
以較小增量執行 SGD 的程式碼
# запустим функцию, уменьшив шаг в 100 раз и увеличив количество шагов соответсвующе
list_parametres_stoch_gradient_descence = stoch_grad_descent_usual(x_us, y_us, l=0.001, steps = 80000)
print ' 33[1m' + ' 33[4m' + "Значения коэффициентов a и b:" + ' 33[0m'
print 'a =', round(list_parametres_stoch_gradient_descence[0][0],3)
print 'b =', round(list_parametres_stoch_gradient_descence[0][1],3)
print
print ' 33[1m' + ' 33[4m' + "Сумма квадратов отклонений:" + ' 33[0m'
print round(list_parametres_stoch_gradient_descence[1][-1],3)
print
print ' 33[1m' + ' 33[4m' + "Количество итераций в стохастическом градиентном спуске:" + ' 33[0m'
print len(list_parametres_stoch_gradient_descence[1])
print 'График №6 "Сумма квадратов отклонений по-шагово"'
plt.plot(range(len(list_parametres_stoch_gradient_descence[1])), list_parametres_stoch_gradient_descence[1], color='red', lw=2)
plt.xlabel('Steps (Iteration)', size=16)
plt.ylabel('Sum of squared deviations', size=16)
plt.show()
圖 6“隨機梯度下降期間的偏差平方和（80 萬步）”
係數有所改善，但仍不理想。 假設，這可以透過這種方式來糾正。 例如，我們在最後 1000 次迭代中選擇誤差最小的係數值。 確實，為此我們還必須寫下係數本身的值。 我們不會這樣做，而是要注意日程。 看起來很平滑，誤差似乎均勻減少。 其實這不是真的。 讓我們來看看前 1000 次迭代，並將它們與最後一次進行比較。
SGD 圖表代碼（前 1000 步）
print 'График №7 "Сумма квадратов отклонений по-шагово. Первые 1000 итераций"'
plt.plot(range(len(list_parametres_stoch_gradient_descence[1][:1000])),
list_parametres_stoch_gradient_descence[1][:1000], color='red', lw=2)
plt.xlabel('Steps (Iteration)', size=16)
plt.ylabel('Sum of squared deviations', size=16)
plt.show()
print 'График №7 "Сумма квадратов отклонений по-шагово. Последние 1000 итераций"'
plt.plot(range(len(list_parametres_stoch_gradient_descence[1][-1000:])),
list_parametres_stoch_gradient_descence[1][-1000:], color='red', lw=2)
plt.xlabel('Steps (Iteration)', size=16)
plt.ylabel('Sum of squared deviations', size=16)
plt.show()圖 7“偏差平方和 SGD（前 1000 步）”
圖 8“偏差平方和 SGD（最後 1000 步）”
在下降的一開始，我們觀察到誤差相當均勻且急劇下降。 在最後一次迭代中，我們看到誤差在 1,475 的值附近徘徊，在某些時刻甚至等於這個最佳值，但隨後它仍然上升......我再說一遍，你可以寫下 的值係數 и ，然後選擇誤差最小的。 然而，我們遇到了一個更嚴重的問題：我們必須執行 80 萬步驟（參見程式碼）才能獲得接近最優的值。 而這已經與相對於梯度下降用隨機梯度下降節省計算時間的想法相矛盾了。 哪些方面可以修正和改進？ 不難注意到，在第一次迭代中，我們自信地向下走，因此，我們應該在第一次迭代中留下一個大的步驟，並在前進時減少步驟。 我們不會在本文中這樣做 - 它已經太長了。 有興趣的可以自己想想怎麼做，並不難：）
現在讓我們使用該庫執行隨機梯度下降 數字貨幣 （我們不要被我們之前確定的石頭絆倒）
隨機梯度下降代碼 (NumPy)
# для начала напишем функцию градиентного шага
def stoch_grad_step_numpy(vector_init, X, ind, y, l):
x = X[ind]
y_pred = np.dot(x,vector_init)
err = y_pred - y[ind]
grad_a = err
grad_b = x[1]*err
return vector_init - l*np.array([grad_a, grad_b])
# определим функцию стохастического градиентного спуска
def stoch_grad_descent_numpy(X, y, l=0.1, steps = 800):
vector_init = np.array([[np.random.randint(X.shape[0])], [np.random.randint(X.shape[0])]])
errors = []
for i in range(steps):
ind = np.random.randint(X.shape[0])
new_vector = stoch_grad_step_numpy(vector_init, X, ind, y, l)
vector_init = new_vector
errors.append(error_square_numpy(vector_init,X,y))
return (vector_init), (errors)
# запишем массив значений
list_parametres_stoch_gradient_descence = stoch_grad_descent_numpy(x_np, y_np, l=0.001, steps = 80000)
print ' 33[1m' + ' 33[4m' + "Значения коэффициентов a и b:" + ' 33[0m'
print 'a =', round(list_parametres_stoch_gradient_descence[0][0],3)
print 'b =', round(list_parametres_stoch_gradient_descence[0][1],3)
print
print ' 33[1m' + ' 33[4m' + "Сумма квадратов отклонений:" + ' 33[0m'
print round(list_parametres_stoch_gradient_descence[1][-1],3)
print
print ' 33[1m' + ' 33[4m' + "Количество итераций в стохастическом градиентном спуске:" + ' 33[0m'
print len(list_parametres_stoch_gradient_descence[1])
print
結果與不使用下降時的數值幾乎相同 數字貨幣。 然而，這是合乎邏輯的。
讓我們看看隨機梯度下降花了多長時間。
決定SGD計算時間的代碼（80萬步）
print ' 33[1m' + ' 33[4m' +
"Время выполнения стохастического градиентного спуска без использования библиотеки NumPy:"
+ ' 33[0m'
%timeit list_parametres_stoch_gradient_descence = stoch_grad_descent_usual(x_us, y_us, l=0.001, steps = 80000)
print '***************************************'
print
print ' 33[1m' + ' 33[4m' +
"Время выполнения стохастического градиентного спуска с использованием библиотеки NumPy:"
+ ' 33[0m'
%timeit list_parametres_stoch_gradient_descence = stoch_grad_descent_numpy(x_np, y_np, l=0.001, steps = 80000)
越深入森林，雲層就越暗：「自寫」公式再次顯示了最好的結果。 所有這些都表明必須有更微妙的方式來使用該庫 數字貨幣，這確實加快了計算操作。 在本文中我們不會了解它們。 閒暇時會有一些事情可以思考:)
我們總結一下
在進行總結之前，我想回答一個很可能是我們親愛的讀者提出的問題。 事實上，如果我們手中有如此強大而簡單的設備，為什麼我們需要在山上走上走下（大部分是下山）才能找到珍貴的低地？分析解決方案的形式，它立即將我們傳送到正確的地方？
這個問題的答案就在表面上。 現在我們來看一個非常簡單的例子，其中真正的答案是 取決於一個標誌 。 這種情況在生活中並不常見，所以讓我們想像一下我們有 2 個、30 個、50 個或更多的標誌。 讓我們為每個屬性添加數千甚至數萬個值。 在這種情況下，解析解可能無法經得起考驗而失敗。 反過來，梯度下降及其變化將緩慢但肯定地使我們更接近目標——函數的最小值。 不要擔心速度 - 我們可能會尋找允許我們設定和調節步長（即速度）的方法。
現在是實際的簡要總結。
首先，我希望本文中提供的材料能夠幫助剛入門的「資料科學家」了解如何求解簡單（而不僅僅是）線性迴歸方程式。
其次，我們研究了求解方程式的幾種方法。 現在，我們可以根據情況選擇最適合解決問題的方案。
第三，我們看到了附加設定的威力，即梯度下降步長。 這個參數不能忽略。 如上所述，為了減少計算成本，應該在下降過程中改變步長。
第四，在我們的例子中，「自製」函數顯示了最佳的計算時間結果。 這可能是由於沒有最專業地使用圖書館的功能 數字貨幣。 但無論如何，以下結論不言自明。 一方面，有時值得質疑既定的觀點，另一方面，並不總是值得讓一切變得複雜——相反，有時解決問題的更簡單的方法更有效。 由於我們的目標是分析解決簡單線性迴歸方程式的三種方法，因此使用「自編寫」函數對我們來說就足夠了。
文學（或類似的東西）
1. 線性迴歸
2.最小平方法
3. 衍生性商品
4。 梯度
5.梯度下降
6.NumPy庫
來源： www.habr.com