Ngrampungake persamaan regresi linier prasaja

Artikel kasebut mbahas sawetara cara kanggo nemtokake persamaan matematika saka garis regresi sing prasaja (dipasangake).

Kabeh cara kanggo ngrampungake persamaan sing dibahas ing kene adhedhasar metode kuadrat paling sithik. Ayo tuduhake cara kaya ing ngisor iki:

• Solusi analitik
• Keturunan Gradien
• Penurunan gradien stokastik

Kanggo saben cara kanggo ngrampungake persamaan garis lurus, artikel kasebut nyedhiyakake macem-macem fungsi, sing utamane dipérang dadi sing ditulis tanpa nggunakake perpustakaan. NomPy lan sing digunakake kanggo petungan NomPy. Punika pitados bilih nggunakake skillful NomPy bakal ngurangi biaya komputasi.

Kabeh kode sing diwenehake ing artikel kasebut ditulis nganggo basa python-2.7 nggunakake Notebook Jupyter Kab. Kode sumber lan file kanthi data sampel dipasang ing Github

Artikel kasebut luwih ditujokake kanggo para pamula lan wong-wong sing wis mboko sithik wiwit sinau babagan bagean sing wiyar banget ing intelijen buatan - learning machine.

Kanggo nggambarake materi kasebut, kita nggunakake conto sing gampang banget.

Conto kahanan

Kita duwe limang nilai sing nggambarake katergantungan Y saka X (Tabel No. 1):

Tabel No. 1 "Contoh kahanan"

Kita bakal nganggep yen nilai kasebut iku sasi taun, lan - revenue sasi iki. Ing tembung liyane, revenue gumantung ing sasi taun, lan - mung tandha kang revenue gumantung.

Conto kaya-kaya, saka sudut pandang katergantungan kondisional saka revenue ing sasi taun, lan saka sudut pandang jumlah nilai - ana sawetara banget. Nanging, simplifikasi kasebut bakal ngidini, kaya sing dikandhakake, kanggo nerangake, ora mesthi gampang, materi sing diasimilasi para pamula. Lan uga kesederhanaan nomer kasebut bakal ngidini wong-wong sing pengin ngrampungake conto kasebut ing kertas tanpa biaya tenaga kerja sing signifikan.

Ayo kita nganggep manawa katergantungan sing diwenehake ing conto kasebut bisa dikira kanthi cukup kanthi persamaan matematika saka garis regresi prasaja (dipasangake) saka formulir:

ngendi yaiku sasi nalika entuk bathi, - revenue sing cocog karo sasi, и yaiku koefisien regresi saka garis sing dikira.

Elinga yen koefisien asring disebut slope utawa gradient saka garis kira-kira; nggantosi jumlah kang ing nalika owah-owahan .

Temenan, tugas kita ing conto yaiku milih koefisien kasebut ing persamaan и , ing ngendi panyimpangan saka nilai revenue sing diwilang miturut sasi saka jawaban sing bener, yaiku. Nilai sing ditampilake ing sampel bakal minimal.

Metode kuadrat paling sithik

Miturut cara kuadrat paling ora, panyimpangan kudu diitung kanthi kuadrat. Teknik iki ngidini sampeyan ngindhari panyimpangan bebarengan yen ana pratandha sing ngelawan. Contone, yen ing sawijining kasus, penyimpangan kasebut +5 (plus lima), lan ing liyane -5 (minus lima), banjur jumlah panyimpangan bakal mbatalake siji liyane lan jumlahe 0 (nol). Sampeyan bisa uga ora kanggo persegi panyimpangan, nanging nggunakake properti saka modulus lan banjur kabeh panyimpangan bakal positif lan bakal nglumpukake. Kita ora bakal njlentrehake babagan iki kanthi rinci, nanging mung nuduhake yen kanggo ngetung petungan, biasane kanggo kothak panyimpangan.

Iki minangka rumus sing bakal nemtokake jumlah penyimpangan kuadrat paling sithik (kesalahan):

ngendi minangka fungsi kira-kira jawaban sing bener (yaiku, revenue sing diwilang),

yaiku jawaban sing bener (asil sing diwenehake ing conto),

yaiku indeks sampel (nomer sasi ing ngendi panyimpangan ditemtokake)

Ayo mbedakake fungsi kasebut, nemtokake persamaan diferensial parsial, lan siyap nerusake menyang solusi analitis. Nanging pisanan, ayo njupuk dolan singkat babagan apa diferensiasi lan elinga makna geometris saka turunan kasebut.

Diferensiasi

Diferensiasi yaiku operasi nemokake turunan saka sawijining fungsi.

Apa derivatif digunakake kanggo? Turunan saka sawijining fungsi menehi ciri tingkat owah-owahan fungsi lan ngandhani arahe. Yen turunan ing titik tartamtu positif, banjur fungsi mundhak; Lan luwih gedhe nilai turunan absolut, luwih dhuwur tingkat owah-owahan nilai fungsi, uga kemiringan grafik fungsi sing luwih dhuwur.

Contone, miturut kondisi sistem koordinat Cartesian, nilai turunan ing titik M(0,0) padha karo + 25 tegese ing titik tartamtu, nalika Nilai wis pindah sisih tengen dening unit conventional, Nilai mundhak dening 25 Unit conventional. Ing grafik katon kaya mundhak cukup tajem ing nilai saka titik tartamtu.

Conto liyane. Nilai turunan padha -0,1 tegese nalika dipindhah saben unit konvensional, nilai nyuda mung 0,1 unit konvensional. Ing wektu sing padha, ing grafik fungsi kasebut, kita bisa mirsani slope mudhun sing meh ora katon. Nggambar analogi karo gunung, kaya-kaya kita mudhun alon-alon saka gunung, ora kaya conto sadurunge, ing ngendi kita kudu mendaki puncak sing curam banget :)

Mangkono, sawise mbedakake fungsi dening rintangan и , kita nemtokake persamaan diferensial parsial urutan 1. Sawise nemtokake persamaan, kita bakal nampa sistem rong persamaan, kanthi ngrampungake sing bakal bisa milih nilai koefisien kasebut. и , sing nilai turunan sing cocog ing titik tartamtu diganti kanthi jumlah sing cilik banget, lan ing kasus solusi analitis ora owah. Ing tembung liya, fungsi kesalahan ing koefisien sing ditemokake bakal tekan minimal, amarga nilai turunan parsial ing titik kasebut bakal padha karo nol.

Dadi, miturut aturan diferensiasi, persamaan derivatif parsial saka urutan 1 babagan koefisien. bakal njupuk formulir:

Persamaan turunan parsial urutan 1 babagan bakal njupuk formulir:

Akibaté, kita nampa sistem persamaan sing nduweni solusi analitis sing cukup prasaja:

miwiti {persamaan*}
miwiti {kasus}
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
pungkasan {kasus}
pungkasan {persamaan*}

Sadurunge ngrampungake persamaan, ayo preload, priksa manawa loading wis bener, lan format data.

Loading lan format data

Perlu dicathet menawa amarga kasunyatan manawa kanggo solusi analitik, lan sabanjure kanggo gradien lan gradient stokastik, kita bakal nggunakake kode kasebut ing rong variasi: nggunakake perpustakaan NomPy lan tanpa nggunakake, banjur kita kudu format data cocok (ndeleng kode).

Data loading lan kode pangolahan

# импортируем все нужные нам библиотеки
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 '********************************************'

Visualisasi

Saiki, sawise kita wis, pisanan, dimuat data, sareh, mriksa bener loading lan pungkasanipun format data, kita bakal nindakake visualisasi pisanan. Cara sing asring digunakake kanggo iki yaiku pasangan perpustakaan segara. Ing conto kita, amarga jumlah winates, ora ana gunane nggunakake perpustakaan segara. Kita bakal nggunakake perpustakaan biasa matplotlib lan mung katon ing scatterplot.

Kode Scatterplot

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()

Bagan No. 1 "Ketergantungan revenue ing sasi ing taun"

Solusi analitik

Ayo nggunakake alat sing paling umum ing python lan ngrampungake sistem persamaan:

miwiti {persamaan*}
miwiti {kasus}
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
pungkasan {kasus}
pungkasan {persamaan*}

Miturut aturan Cramer kita bakal nemokake determinant umum, uga determinants dening lan dening , sawise kang, dibagi determinant dening menyang determinan umum - golek koefisien , uga kita nemokake koefisien .

Kode solusi analitik

# определим функцию для расчета коэффициентов 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)

Punika ingkang kita pikantuk:

Dadi, nilai-nilai koefisien wis ditemokake, jumlah penyimpangan kuadrat wis ditemtokake. Ayo tarik garis lurus ing histogram hamburan sesuai karo koefisien sing ditemokake.

Kode garis regresi

# определим функцию для формирования массива рассчетных значений выручки
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()

Bagan No. 2 “Jawaban sing bener lan diwilang”

Sampeyan bisa ndeleng grafik panyimpangan saben wulan. Ing kasus kita, kita ora bakal entuk nilai praktis sing signifikan saka iku, nanging kita bakal nyukupi rasa penasaran babagan carane persamaan regresi linier sing prasaja nggambarake katergantungan revenue ing sasi taun.

Kode grafik panyimpangan

# определим функцию для формирования массива отклонений в процентах
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()

Bagan No. 3 "Penyimpangan, %"

Ora sampurna, nanging kita ngrampungake tugas kita.

Ayo dadi nulis fungsi sing, kanggo nemtokake koefisien и nggunakake perpustakaan NomPy, luwih tepat, kita bakal nulis rong fungsi: siji nggunakake matriks pseudoinverse (ora dianjurake ing laku, amarga proses komputasi rumit lan ora stabil), liyane nggunakake persamaan matriks.

Kode Solusi Analitik (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

Ayo mbandhingake wektu sing ditindakake kanggo nemtokake koefisien и , sesuai karo 3 cara sing diwenehake.

Kode kanggo ngitung wektu pitungan

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)

Kanthi jumlah data sing cilik, fungsi "ditulis dhewe" metu ing ngarep, sing nemokake koefisien nggunakake metode Cramer.

Saiki sampeyan bisa pindhah menyang cara liya kanggo nemokake koefisien и .

Keturunan Gradien

Pisanan, ayo nemtokake apa gradien. Cukup, gradien minangka bagean sing nuduhake arah pertumbuhan maksimal sawijining fungsi. Kanthi analogi karo pendakian gunung, ing endi gradien ngadhepi ing ngendi pendakian sing paling curam menyang puncak gunung. Ngembangake conto karo gunung, kita elinga yen sejatine kita butuh turunan sing paling curam supaya bisa tekan dataran rendah kanthi cepet, yaiku minimal - papan sing fungsine ora mundhak utawa nyuda. Ing titik iki turunan bakal padha karo nol. Mulane, kita ora butuh gradien, nanging antigradient. Kanggo nemokake antigradient sampeyan mung kudu Multiply gradien dening -1 (minus siji).

Ayo kita mbayar manungsa waé kanggo kasunyatan sing fungsi bisa duwe sawetara minima, lan sawise mudhun menyang salah siji saka wong-wong mau nggunakake algoritma ngajokaken ing ngisor iki, kita ora bakal bisa kanggo golek minimal liyane, kang bisa uga luwih murah tinimbang ditemokaké. Ayo santai, iki dudu ancaman kanggo kita! Ing kasus kita dealing karo minimal siji, wiwit fungsi kita ing grafik punika parabola biasa. Lan kita kabeh kudu ngerti banget babagan kursus matematika sekolah, parabola mung duwe minimal siji.

Sawise kita ngerteni kenapa kita butuh kecerunan, lan uga manawa kecerunan kasebut minangka bagean, yaiku, vektor kanthi koordinat sing diwenehake, sing padha karo koefisien. и kita bisa ngleksanakake keturunan gradien.

Sadurunge miwiti, aku saranake maca sawetara ukara babagan algoritma keturunan:

• Kita nemtokake kanthi cara pseudo-acak koordinat saka koefisien и . Ing conto kita, kita bakal nemtokake koefisien cedhak nol. Iki minangka praktik umum, nanging saben kasus bisa uga duwe praktik dhewe.
• Saka koordinat nyuda nilai turunan parsial urutan 1 ing titik kasebut . Dadi, yen turunan positif, fungsi kasebut mundhak. Mulane, kanthi nyuda nilai turunan, kita bakal pindhah menyang arah ngelawan saka wutah, yaiku, ing arah keturunan. Yen turunan kasebut negatif, mula fungsi ing titik iki mudhun lan kanthi nyuda nilai turunan kasebut, kita pindhah menyang arah keturunan.
• Kita nindakake operasi sing padha karo koordinat : nyuda nilai turunan parsial ing titik kasebut .
• Supaya ora mlumpat kanthi minimal lan mabur menyang papan sing jero, perlu nyetel ukuran langkah ing arah mudhun. Umumé, sampeyan bisa nulis artikel lengkap babagan carane nyetel langkah kanthi bener lan carane ngganti nalika proses mudhun kanggo ngurangi biaya komputasi. Nanging saiki kita duwe tugas rada beda ahead saka kita, lan kita bakal netepake ukuran langkah nggunakake metode ilmiah "poke" utawa, lagi ngomong ing parlance umum, empiris.
• Sawise kita saka koordinat diwenehi и nyuda nilai turunan, kita entuk koordinat anyar и . Kita njupuk langkah sabanjure (subtraction), wis saka koordinat sing diwilang. Dadi siklus kasebut diwiwiti maneh lan maneh, nganti konvergensi sing dibutuhake bisa ditindakake.

Kabeh! Saiki kita wis siyap kanggo nggoleki jurang paling jero saka Palung Mariana. Ayo dadi miwiti.

Kode kanggo keturunan gradien

# напишем функцию градиентного спуска без использования библиотеки 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

Kita nyilem menyang paling ngisor Mariana Trench lan ing kono kita nemokake kabeh nilai koefisien sing padha и , sing persis kaya sing dikarepake.

Ayo nyilem maneh, mung wektu iki, kendaraan laut jero bakal diisi karo teknologi liyane, yaiku perpustakaan. NomPy.

Kode keturunan gradien (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

Nilai koefisien и ora bisa diganti.

Ayo goleki carane kesalahan diganti nalika keturunan gradien, yaiku, carane jumlah penyimpangan kuadrat diganti karo saben langkah.

Kode kanggo ngrancang jumlah penyimpangan kuadrat

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()

Grafik No. 4 "Jumlah penyimpangan kuadrat sajrone turunan gradien"

Ing grafik, kita weruh yen saben langkah kesalahane suda, lan sawise sawetara iterasi kita mirsani garis sing meh horisontal.

Pungkasan, ayo ngira prabédan ing wektu eksekusi kode:

Kode kanggo nemtokake wektu pitungan gradien keturunan

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)

Bisa uga ana sing salah, nanging maneh minangka fungsi "ditulis ing omah" sing ora nggunakake perpustakaan. NomPy outperforms wektu pitungan fungsi nggunakake perpustakaan NomPy.

Nanging kita ora mandheg, nanging terus sinau babagan cara sing nyenengake kanggo ngatasi persamaan regresi linier sing prasaja. Ketemu!

Penurunan gradien stokastik

Supaya cepet ngerti prinsip operasi turunan gradient stokastik, luwih becik nemtokake bedane saka turunan gradien biasa. Kita, ing kasus keturunan gradien, ing persamaan turunan saka и nggunakake jumlah nilai kabeh fitur lan jawaban sing bener sing kasedhiya ing sampel (yaiku, jumlah kabeh и ). Ing turunan kecerunan stokastik, kita ora bakal nggunakake kabeh nilai sing ana ing sampel, nanging pseudo-acak milih indeks sampel sing diarani lan nggunakake nilai kasebut.

Contone, yen indeks ditemtokake dadi nomer 3 (telu), banjur njupuk nilai и , banjur kita ngganti nilai kasebut menyang persamaan turunan lan nemtokake koordinat anyar. Banjur, sawise nemtokake koordinat, kita maneh pseudo-acak nemtokake indeks sampel, ngganti nilai sing cocog karo indeks menyang persamaan diferensial parsial, lan nemtokake koordinat kanthi cara anyar. и lsp. nganti konvergensi dadi ijo. Sepisanan, bisa uga ora kaya ngono, nanging bisa ditindakake. Sampeyan kudu nyatet sing kesalahan ora suda saben langkah, nanging mesthi ana kecenderungan.

Apa kaluwihan saka turunan gradient stokastik tinimbang sing konvensional? Yen ukuran sampel kita gedhe banget lan diukur ing puluhan ewu nilai, mula luwih gampang diproses, ucapake, sewu acak, tinimbang kabeh sampel. Iki ngendi keturunan stochastic gradient teka menyang muter. Ing kasus kita, mesthi, kita ora bakal sok dong mirsani akeh prabédan.

Ayo katon ing kode.

Kode kanggo turunan gradient stokastik

# определим функцию стох.град.шага
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])

Kita ndeleng kanthi teliti ing koefisien lan nyekel dhéwé takon pitakonan "Kepiye iki bisa?" Kita entuk nilai koefisien liyane и . Mungkin keturunan kecerunan stokastik wis nemokake paramèter sing luwih optimal kanggo persamaan kasebut? Sayange ora. Cukup kanggo ndeleng jumlah penyimpangan kuadrat lan ndeleng manawa kanthi nilai koefisien anyar, kesalahane luwih gedhe. Kita ora kesusu putus asa. Ayo nggawe grafik saka owah-owahan kesalahan.

Kode kanggo ngrancang jumlah penyimpangan kuadrat ing turunan gradien stokastik

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()

Grafik No. 5 "Jumlah penyimpangan kuadrat sajrone turunan gradien stokastik"

Deleng jadwal, kabeh wis ana lan saiki bakal ndandani kabeh.

Dadi apa sing kedadeyan? Ing ngisor iki kedadeyan. Nalika kita milih sasi kanthi acak, banjur kanggo sasi sing dipilih algoritma kita ngupaya nyuda kesalahan ing ngitung revenue. Banjur kita pilih sasi liyane lan baleni pitungan, nanging kita nyuda kesalahan kanggo sasi kapindho sing dipilih. Saiki elinga yen rong sasi kapisan nyimpang sacara signifikan saka garis persamaan regresi linier sing prasaja. Iki tegese yen salah siji saka rong sasi iki dipilih, kanthi ngurangi kesalahan saben wong, algoritma kita akeh nambah kesalahan kanggo kabeh sampel. Dadi apa sing kudu ditindakake? Jawaban iki prasaja: sampeyan kudu ngurangi langkah keturunan. Sawise kabeh, kanthi ngurangi langkah mudhun, kesalahan uga bakal mandheg "mlumpat" munggah lan mudhun. Utawa, kesalahan "mlumpat" ora bakal mandheg, nanging ora bakal ditindakake kanthi cepet :) Ayo priksa.

Kode kanggo mbukak SGD kanthi tambahan sing luwih cilik

# запустим функцию, уменьшив шаг в 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()

Grafik No. 6 "Jumlah penyimpangan kuadrat sajrone turunan gradien stokastik (80 ewu langkah)"

Koefisien wis apik, nanging isih ora becik. Kanthi hipotesis, iki bisa didandani kanthi cara iki. Kita milih, contone, ing 1000 iterasi pungkasan nilai koefisien sing kesalahan minimal digawe. Bener, kanggo iki kita uga kudu nulis nilai koefisien dhewe. Kita ora bakal nindakake iki, nanging menehi perhatian marang jadwal. Katon lancar lan kesalahan katon mudhun kanthi rata. Bener iki ora bener. Ayo katon ing 1000 iterasi pisanan lan mbandhingake karo sing pungkasan.

Kode kanggo grafik SGD (1000 langkah pisanan)

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()

Grafik No. 7 "Jumlah sisihan kuadrat SGD (1000 langkah pisanan)"

Grafik No. 8 "Jumlah penyimpangan kuadrat SGD (1000 langkah pungkasan)"

Ing awal banget mudhun, kita mirsani nyuda cukup seragam lan tajem ing kesalahan. Ing pengulangan pungkasan, kita bisa ndeleng manawa kesalahan kasebut ngubengi lan ngubengi nilai 1,475 lan ing sawetara wektu malah padha karo nilai optimal iki, nanging banjur isih munggah ... Aku mbaleni, sampeyan bisa nulis mudhun nilai saka koefisien и , banjur pilih sing kesalahane minimal. Nanging, kita duwe masalah sing luwih serius: kita kudu njupuk 80 ewu langkah (deleng kode) kanggo entuk nilai sing paling optimal. Lan iki wis mbantah ide kanggo ngirit wektu komputasi kanthi turunan gradient stokastik relatif marang keturunan gradien. Apa sing bisa didandani lan didandani? Iku ora angel kanggo sok dong mirsani sing ing iterasi pisanan kita manteb ing ati mudhun lan, mulane, kita kudu ninggalake langkah gedhe ing iterasi pisanan lan ngurangi langkah nalika kita maju. Kita ora bakal nindakake iki ing artikel iki - wis suwe banget. Sing pengin bisa mikir dhewe carane nindakake iki, ora angel :)

Saiki ayo nindakake turunan gradient stokastik nggunakake perpustakaan NomPy (lan aja nganti kesandhung ing watu sing wis dingerteni sadurunge)

Kode kanggo Stochastic Gradient Descent (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

Nilai-nilai kasebut meh padha karo nalika mudhun tanpa nggunakake NomPy. Nanging, iki logis.

Ayo ngerteni sepira suwene penurunan gradien stokastik.

Kode kanggo nemtokake wektu pitungan SGD (80 ewu langkah)

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)

Luwih menyang alas, awan sing luwih peteng: maneh, rumus "ditulis dhewe" nuduhake asil sing paling apik. Kabeh iki nuduhake yen kudu ana cara sing luwih halus kanggo nggunakake perpustakaan NomPy, sing tenan nyepetake operasi komputasi. Ing artikel iki kita ora bakal sinau babagan. Bakal ana sing kudu dipikirake ing wektu luang :)

Ringkesan

Sadurunge ngringkes, aku pengin mangsuli pitakon sing paling mungkin muncul saka pembaca sing dikasihi. Yagene, nyatane, "nyiksa" kasebut kanthi turun-tumurun, kenapa kita kudu mlaku munggah-mudhun gunung (biasane mudhun) kanggo nemokake dataran rendah sing ditresnani, yen ana ing tangan kita piranti sing kuat lan prasaja, ing wangun solusi analitis, kang enggal teleports kita Panggonan tengen?

Jawaban kanggo pitakonan iki dumunung ing permukaan. Saiki kita wis katon ing conto banget prasaja, kang jawaban bener gumantung ing siji tandha . Sampeyan ora weruh iki asring ing urip, supaya ayo kang mbayangno sing kita duwe 2, 30, 50 utawa liyane pratandha. Ayo ditambahake menyang ewu, utawa malah puluhan ewu nilai kanggo saben atribut. Ing kasus iki, solusi analitis bisa uga ora tahan uji lan gagal. Sabanjure, turunan gradien lan variasi bakal alon-alon nanging mesthi nggawa kita nyedhaki tujuan - minimal fungsi kasebut. Lan aja kuwatir babagan kacepetan - kita bakal bisa ndeleng cara sing bakal ngidini kita nyetel lan ngatur dawa langkah (yaiku kacepetan).

Lan saiki ringkesan ringkes nyata.

Kaping pisanan, muga-muga materi sing diwenehake ing artikel kasebut bakal mbantu miwiti "ilmuwan data" kanggo mangerteni carane ngatasi persamaan regresi linier sing prasaja (lan ora mung).

Kapindho, kita ndeleng sawetara cara kanggo ngatasi persamaan kasebut. Saiki, gumantung saka kahanan, kita bisa milih sing paling cocog kanggo ngatasi masalah kasebut.

Katelu, kita weruh kekuwatan setelan tambahan, yaiku dawa langkah keturunan gradien. Parameter iki ora bisa diabaikan. Kaya kasebut ing ndhuwur, kanggo ngurangi biaya petungan, dawa langkah kudu diganti nalika mudhun.

Kaping papat, ing kasus kita, fungsi "ditulis ing omah" nuduhake asil wektu paling apik kanggo petungan. Iki mbokmenawa amarga ora nggunakake paling profesional kapabilitas perpustakaan kang NomPy. Nanging apa wae, kesimpulan ing ngisor iki nyaranake dhewe. Ing tangan siji, kadhangkala iku worth pitakonan panemu mantep, lan ing tangan liyane, iku ora tansah worth complicating kabeh - ing nalisir, kadhangkala cara prasaja kanggo ngatasi masalah luwih efektif. Lan amarga tujuane kanggo nganalisa telung pendekatan kanggo ngrampungake persamaan regresi linier sing sederhana, panggunaan fungsi "ditulis dhewe" cukup kanggo kita.

Sastra (utawa sapanunggalane)

1. Regresi linier

http://statistica.ru/theory/osnovy-lineynoy-regressii/

2. Metode kuadrat paling sithik

mathprofi.ru/metod_naimenshih_kvadratov.html

3. Turunan

www.mathprofi.ru/chastnye_proizvodnye_primery.html

4. Gradien

mathprofi.ru/proizvodnaja_po_napravleniju_i_gradient.html

5. Keturunan gradien

habr.com/en/post/471458

habr.com/en/post/307312

artemarakcheev.com//2017-12-31/linear_regression

6. Pustaka NumPy

docs.scipy.org/doc/numpy-1.10.1/reference/generated/numpy.linalg.solve.html

docs.scipy.org/doc/numpy-1.10.0/reference/generated/numpy.linalg.pinv.html

pythonworld.ru/numpy/2.html

Source: www.habr.com