L'article traite de plusieurs façons de déterminer l'équation mathématique d'une droite de régression simple (paire).
Toutes les méthodes de résolution de l’équation discutées ici sont basées sur la méthode des moindres carrés. Notons les méthodes comme suit :
- Solution analytique
- Descente graduelle
- Descente de gradient stochastique
Pour chaque méthode de résolution de l'équation d'une droite, l'article propose diverses fonctions, qui sont principalement divisées en celles qui sont écrites sans utiliser la bibliothèque. NumPy et ceux qui utilisent pour les calculs NumPy. On pense qu'une utilisation habile NumPy réduira les coûts informatiques.
Tout le code donné dans l'article est écrit dans la langue python 2.7 en utilisant Jupyter Notebook. Le code source et le fichier contenant des exemples de données sont publiés sur
L'article s'adresse davantage aux débutants et à ceux qui ont déjà progressivement commencé à maîtriser l'étude d'une très large section de l'intelligence artificielle - l'apprentissage automatique.
Pour illustrer le matériel, nous utilisons un exemple très simple.
Exemples de conditions
Nous avons cinq valeurs qui caractérisent la dépendance Y à partir de X (Tableau n°1) :
Tableau n°1 « Exemples de conditions »
Nous supposerons que les valeurs est le mois de l'année, et — revenus ce mois-ci. En d'autres termes, les revenus dépendent du mois de l'année, et - le seul signe dont dépendent les revenus.
L'exemple est moyen, tant du point de vue de la dépendance conditionnelle des revenus au mois de l'année, que du point de vue du nombre de valeurs - il y en a très peu. Cependant, une telle simplification permettra, comme on dit, d'expliquer, pas toujours avec facilité, la matière acquise par les débutants. Et aussi la simplicité des chiffres permettra à ceux qui souhaitent résoudre l'exemple sur papier sans coûts de main d'œuvre importants.
Supposons que la dépendance donnée dans l'exemple puisse être assez bien approchée par l'équation mathématique d'une droite de régression simple (apparié) de la forme :
où est le mois au cours duquel les revenus ont été perçus, — les revenus correspondant au mois, и sont les coefficients de régression de la droite estimée.
Notez que le coefficient souvent appelée pente ou gradient de la ligne estimée ; représente le montant par lequel le lors du changement .
Évidemment, notre tâche dans l'exemple est de sélectionner de tels coefficients dans l'équation и , auquel les écarts de nos valeurs de revenus calculées par mois par rapport aux vraies réponses, c'est-à-dire les valeurs présentées dans l'échantillon seront minimes.
Méthode des moindres carrés
Selon la méthode des moindres carrés, l’écart doit être calculé en le mettant au carré. Cette technique permet d'éviter l'annulation mutuelle des écarts s'ils présentent des signes opposés. Par exemple, si dans un cas, l’écart est +5 (plus cinq), et dans l'autre -5 (moins cinq), alors la somme des écarts s’annulera et s’élèvera à 0 (zéro). Il est possible de ne pas mettre l'écart au carré, mais d'utiliser la propriété du module et alors tous les écarts seront positifs et s'accumuleront. Nous ne nous attarderons pas sur ce point en détail, mais indiquerons simplement que pour la commodité des calculs, il est d'usage de mettre l'écart au carré.
Voici à quoi ressemble la formule avec laquelle nous déterminerons la plus petite somme des écarts carrés (erreurs) :
où est une fonction d'approximation des vraies réponses (c'est-à-dire les revenus que nous avons calculés),
sont les vraies réponses (revenus fournis dans l'échantillon),
est l'indice de l'échantillon (numéro du mois au cours duquel l'écart est déterminé)
Différencions la fonction, définissons les équations aux dérivées partielles et soyons prêts à passer à la solution analytique. Mais d'abord, faisons une brève excursion sur ce qu'est la différenciation et rappelons-nous la signification géométrique de la dérivée.
Différenciation
La différenciation est l'opération consistant à trouver la dérivée d'une fonction.
A quoi sert le dérivé ? La dérivée d'une fonction caractérise le taux de changement de la fonction et nous indique sa direction. Si la dérivée en un point donné est positive, alors la fonction augmente ; sinon, la fonction diminue. Et plus la valeur de la dérivée absolue est grande, plus le taux de changement des valeurs de la fonction est élevé, ainsi que plus la pente du graphique de la fonction est forte.
Par exemple, dans les conditions d'un repère cartésien, la valeur de la dérivée au point M(0,0) est égale à +25 signifie qu'à un moment donné, lorsque la valeur est décalée à droite par une unité conventionnelle, valeur augmente de 25 unités conventionnelles. Sur le graphique, cela ressemble à une hausse assez forte des valeurs à partir d'un point donné.
Un autre exemple. La valeur dérivée est égale - 0,1 signifie que lorsqu'il est déplacé par unité conventionnelle, valeur diminue de seulement 0,1 unité conventionnelle. Dans le même temps, sur le graphique de la fonction, on peut observer une pente descendante à peine perceptible. En faisant une analogie avec une montagne, c'est comme si nous descendions très lentement une pente douce d'une montagne, contrairement à l'exemple précédent, où nous devions gravir des sommets très raides :)
Ainsi, après avoir différencié la fonction par hasard и , nous définissons des équations aux dérivées partielles du 1er ordre. Après avoir déterminé les équations, nous recevrons un système de deux équations, en résolvant lesquelles nous pourrons sélectionner de telles valeurs des coefficients и , pour lequel les valeurs des dérivées correspondantes en des points donnés changent d'une très, très petite quantité, et dans le cas d'une solution analytique ne changent pas du tout. En d'autres termes, la fonction d'erreur aux coefficients trouvés atteindra un minimum, puisque les valeurs des dérivées partielles en ces points seront égales à zéro.
Ainsi, selon les règles de différenciation, l'équation aux dérivées partielles du 1er ordre par rapport au coefficient prendra la forme :
Équation aux dérivées partielles du 1er ordre par rapport à prendra la forme :
En conséquence, nous avons reçu un système d'équations qui a une solution analytique assez simple :
begin {équation *}
commencer{cas}
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
fin{cas}
end {équation *}
Avant de résoudre l'équation, préchargeons, vérifions que le chargement est correct et formatons les données.
Chargement et formatage des données
Il est à noter que du fait que pour la solution analytique, puis pour la descente de gradient et de gradient stochastique, nous utiliserons le code selon deux variantes : en utilisant la bibliothèque NumPy et sans l'utiliser, nous aurons alors besoin d'un formatage de données approprié (voir le code).
Code de chargement et de traitement des données
# импортируем все нужные нам библиотеки
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 '********************************************'Visualisation
Maintenant, après avoir d'abord chargé les données, deuxièmement vérifié l'exactitude du chargement et enfin formaté les données, nous effectuerons la première visualisation. La méthode souvent utilisée pour cela est parcelle de couple Bibliothèque marin. Dans notre exemple, en raison du nombre limité, cela ne sert à rien d'utiliser la bibliothèque marin. Nous utiliserons la bibliothèque habituelle matplotlib et regardez simplement le nuage de points.
Code de nuage de points
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()Graphique n°1 « Dépendance des revenus au mois de l'année »
Solution analytique
Utilisons les outils les plus courants dans python et résolvez le système d’équations :
begin {équation *}
commencer{cas}
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
fin{cas}
end {équation *}
D'après la règle de Cramer on trouvera le déterminant général, ainsi que les déterminants par et , après quoi, en divisant le déterminant par au déterminant général - trouver le coefficient , de même on trouve le coefficient .
Code de solution analytique
# определим функцию для расчета коэффициентов 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)Voici ce que nous avons obtenu :
Ainsi, les valeurs des coefficients ont été trouvées, la somme des carrés des écarts a été établie. Traçons une ligne droite sur l'histogramme de diffusion conformément aux coefficients trouvés.
Code de ligne de régression
# определим функцию для формирования массива рассчетных значений выручки
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()Tableau n°2 « Réponses correctes et calculées »
Vous pouvez consulter le graphique d’écart pour chaque mois. Dans notre cas, nous n'en tirerons aucune valeur pratique significative, mais nous satisferons notre curiosité de savoir dans quelle mesure l'équation de régression linéaire simple caractérise la dépendance des revenus par rapport au mois de l'année.
Code du graphique d'écart
# определим функцию для формирования массива отклонений в процентах
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()Graphique n°3 « Écarts, % »
Pas parfait, mais nous avons accompli notre tâche.
Écrivons une fonction qui, pour déterminer les coefficients и utilise la bibliothèque NumPy, plus précisément, nous écrirons deux fonctions : l'une utilisant une matrice pseudo-inverse (non recommandée en pratique, car le processus est complexe en termes de calcul et instable), l'autre utilisant une équation matricielle.
Code de solution analytique (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_npComparons le temps passé à déterminer les coefficients и , conformément aux 3 méthodes présentées.
Code de calcul du temps de calcul
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)
Avec une petite quantité de données, une fonction « auto-écrite » apparaît en tête, qui trouve les coefficients à l'aide de la méthode de Cramer.
Vous pouvez maintenant passer à d'autres façons de trouver des coefficients и .
Descente graduelle
Tout d’abord, définissons ce qu’est un dégradé. En termes simples, le gradient est un segment qui indique la direction de croissance maximale d'une fonction. Par analogie avec l'escalade d'une montagne, là où se trouve la pente, c'est là que se trouve la montée la plus raide jusqu'au sommet de la montagne. En développant l'exemple de la montagne, nous rappelons qu'en fait nous avons besoin de la descente la plus raide pour atteindre le plus rapidement possible la plaine, c'est-à-dire le minimum - l'endroit où la fonction n'augmente ni ne diminue. À ce stade, la dérivée sera égale à zéro. Nous n’avons donc pas besoin d’un gradient, mais d’un antigradient. Pour trouver l'antigradient, il suffit de multiplier le gradient par -1 (moins un).
Faisons attention au fait qu'une fonction peut avoir plusieurs minima, et après être descendu dans l'un d'eux à l'aide de l'algorithme proposé ci-dessous, nous ne pourrons pas trouver un autre minimum, qui pourra être inférieur à celui trouvé. Détendons-nous, ce n'est pas une menace pour nous ! Dans notre cas, nous avons affaire à un seul minimum, puisque notre fonction sur le graphique se trouve une parabole régulière. Et comme nous devrions tous le savoir grâce à nos cours de mathématiques à l’école, une parabole n’a qu’un minimum.
Après avoir découvert pourquoi nous avions besoin d'un gradient, et aussi que le gradient est un segment, c'est-à-dire un vecteur avec des coordonnées données, qui sont précisément les mêmes coefficients и nous pouvons implémenter la descente de gradient.
Avant de commencer, je suggère de lire quelques phrases sur l'algorithme de descente :
- On détermine de manière pseudo-aléatoire les coordonnées des coefficients и . Dans notre exemple, nous déterminerons des coefficients proches de zéro. Il s'agit d'une pratique courante, mais chaque cas peut avoir sa propre pratique.
- De coordonnées soustraire la valeur de la dérivée partielle du 1er ordre au point . Donc, si la dérivée est positive, alors la fonction augmente. Par conséquent, en soustrayant la valeur de la dérivée, nous nous déplacerons dans le sens opposé de la croissance, c'est-à-dire dans le sens de la descente. Si la dérivée est négative, alors la fonction à ce stade diminue et en soustrayant la valeur de la dérivée, nous nous déplaçons dans le sens de la descente.
- On effectue une opération similaire avec la coordonnée : soustraire la valeur de la dérivée partielle au point .
- Afin de ne pas sauter par-dessus le minimum et voler dans l'espace lointain, il est nécessaire de régler la taille du pas dans le sens de la descente. En général, vous pourriez écrire un article entier sur la façon de régler correctement le pas et de le modifier pendant le processus de descente afin de réduire les coûts de calcul. Mais maintenant, nous avons une tâche légèrement différente qui nous attend, et nous allons établir la taille du pas en utilisant la méthode scientifique du « poke » ou, comme on dit dans le langage courant, de manière empirique.
- Une fois qu'on est à partir des coordonnées données и soustrayez les valeurs des dérivées, on obtient de nouvelles coordonnées и . Nous passons à l'étape suivante (soustraction), déjà à partir des coordonnées calculées. Et ainsi le cycle recommence encore et encore, jusqu’à ce que la convergence requise soit atteinte.
Tous! Nous sommes maintenant prêts à partir à la recherche des gorges les plus profondes de la fosse des Mariannes. Commençons.
Code pour la descente de pente
# напишем функцию градиентного спуска без использования библиотеки 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
Nous avons plongé tout en bas de la fosse des Mariannes et nous y avons trouvé toutes les mêmes valeurs de coefficient и , ce qui est exactement ce à quoi il fallait s'attendre.
Faisons une autre plongée, mais cette fois, notre véhicule hauturier sera rempli d'autres technologies, à savoir une bibliothèque NumPy.
Code pour la descente de gradient (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
Valeurs des coefficients и non modifiable.
Regardons comment l'erreur a changé pendant la descente de gradient, c'est-à-dire comment la somme des écarts au carré a changé à chaque étape.
Code pour tracer les sommes des écarts carrés
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()Graphique n°4 « Somme des écarts au carré lors de la descente de pente »
Sur le graphique on voit qu'à chaque pas l'erreur diminue, et après un certain nombre d'itérations on observe une ligne presque horizontale.
Enfin, estimons la différence de temps d’exécution du code :
Code pour déterminer le temps de calcul de la descente de pente
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)
Peut-être que nous faisons quelque chose de mal, mais encore une fois, il s’agit d’une simple fonction « écrite à la maison » qui n’utilise pas la bibliothèque. NumPy surpasse le temps de calcul d'une fonction utilisant la bibliothèque NumPy.
Mais nous ne restons pas immobiles, mais nous nous dirigeons vers l’étude d’une autre manière passionnante de résoudre l’équation de régression linéaire simple. Rencontrer!
Descente de gradient stochastique
Afin de comprendre rapidement le principe de fonctionnement de la descente de gradient stochastique, il est préférable de déterminer ses différences par rapport à la descente de gradient ordinaire. Nous, dans le cas de la descente de gradient, dans les équations des dérivées de и utilisé les sommes des valeurs de toutes les caractéristiques et les vraies réponses disponibles dans l'échantillon (c'est-à-dire les sommes de toutes и ). Dans la descente de gradient stochastique, nous n'utiliserons pas toutes les valeurs présentes dans l'échantillon, mais sélectionnerons plutôt de manière pseudo-aléatoire ce que l'on appelle l'indice d'échantillon et utiliserons ses valeurs.
Par exemple, si l'indice est déterminé comme étant le numéro 3 (trois), alors nous prenons les valeurs и , puis nous substituons les valeurs dans les équations dérivées et déterminons de nouvelles coordonnées. Ensuite, après avoir déterminé les coordonnées, nous déterminons à nouveau de manière pseudo-aléatoire l'indice de l'échantillon, substituons les valeurs correspondant à l'indice dans les équations aux dérivées partielles et déterminons les coordonnées d'une nouvelle manière и etc. jusqu'à ce que la convergence devienne verte. À première vue, cela ne semble pas fonctionner du tout, mais c’est pourtant le cas. Il est vrai qu’il convient de noter que l’erreur ne diminue pas à chaque étape, mais il existe certainement une tendance.
Quels sont les avantages de la descente de gradient stochastique par rapport à la descente de gradient conventionnelle ? Si la taille de notre échantillon est très grande et mesurée en dizaines de milliers de valeurs, il est alors beaucoup plus facile de traiter, disons, un millier de valeurs aléatoires, plutôt que l’échantillon entier. C’est là qu’intervient la descente de gradient stochastique. Dans notre cas, bien sûr, nous ne remarquerons pas beaucoup de différence.
Regardons le code.
Code pour la descente de gradient stochastique
# определим функцию стох.град.шага
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])
Nous examinons attentivement les coefficients et nous nous surprenons à nous poser la question « Comment est-ce possible ? » Nous avons d'autres valeurs de coefficient и . Peut-être que la descente de gradient stochastique a trouvé des paramètres plus optimaux pour l'équation ? Malheureusement non. Il suffit de regarder la somme des écarts carrés et de voir qu'avec de nouvelles valeurs des coefficients, l'erreur est plus grande. Nous ne sommes pas pressés de désespérer. Construisons un graphique du changement d'erreur.
Code pour tracer la somme des écarts carrés dans la descente de gradient stochastique
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()Graphique n°5 « Somme des écarts au carré lors de la descente de gradient stochastique »
En regardant le calendrier, tout se met en place et maintenant nous allons tout arranger.
Alors, qu'est-ce-qu'il s'est passé? Ce qui suit s'est produit. Lorsque nous sélectionnons un mois au hasard, c'est alors pour le mois sélectionné que notre algorithme cherche à réduire l'erreur de calcul des revenus. Ensuite, nous sélectionnons un autre mois et répétons le calcul, mais nous réduisons l'erreur pour le deuxième mois sélectionné. Rappelez-vous maintenant que les deux premiers mois s’écartent considérablement de la ligne de l’équation de régression linéaire simple. Cela signifie que lorsque l’un de ces deux mois est sélectionné, en réduisant l’erreur de chacun d’eux, notre algorithme augmente considérablement l’erreur pour l’ensemble de l’échantillon. Alors que faire? La réponse est simple : il faut réduire le pas de descente. Après tout, en réduisant le pas de descente, l'erreur cessera également de « sauter » de haut en bas. Ou plutôt, l'erreur de "saut" ne s'arrêtera pas, mais elle ne le fera pas si vite :) Vérifions.
Code pour exécuter SGD avec des incréments plus petits
# запустим функцию, уменьшив шаг в 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()
Graphique n°6 « Somme des écarts au carré lors de la descente de gradient stochastique (80 XNUMX pas) »
Les coefficients se sont améliorés, mais ne sont toujours pas idéaux. Hypothétiquement, cela peut être corrigé de cette façon. On sélectionne, par exemple, dans les 1000 dernières itérations les valeurs des coefficients avec lesquels l'erreur minimale a été commise. Certes, pour cela, nous devrons également noter les valeurs des coefficients eux-mêmes. Nous ne ferons pas cela, mais ferons plutôt attention au calendrier. Cela semble fluide et l'erreur semble diminuer uniformément. En fait, ce n'est pas vrai. Regardons les 1000 premières itérations et comparons-les avec la dernière.
Code pour le graphique SGD (1000 premières étapes)
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()Graphique n°7 « Somme des carrés des écarts SGD (1000 premiers pas) »
Graphique n°8 « Somme des carrés des écarts SGD (1000 derniers pas) »
Au tout début de la descente, on observe une diminution de l’erreur assez uniforme et forte. Dans les dernières itérations, nous voyons que l'erreur tourne autour de la valeur de 1,475 et à certains moments même égale cette valeur optimale, mais ensuite elle monte encore... Je le répète, vous pouvez noter les valeurs du coefficients и , puis sélectionnez ceux pour lesquels l'erreur est minime. Cependant, nous avons eu un problème plus grave : nous avons dû faire 80 XNUMX étapes (voir code) pour obtenir des valeurs proches de l'optimum. Et cela contredit déjà l'idée de gagner du temps de calcul avec la descente de gradient stochastique par rapport à la descente de gradient. Que peut-on corriger et améliorer ? Il n'est pas difficile de remarquer que dans les premières itérations, nous descendons avec confiance et, par conséquent, nous devrions laisser un grand pas dans les premières itérations et réduire ce pas à mesure que nous avançons. Nous ne ferons pas cela dans cet article, c'est déjà trop long. Ceux qui le souhaitent peuvent réfléchir par eux-mêmes à la manière de procéder, ce n'est pas difficile :)
Effectuons maintenant une descente de gradient stochastique à l'aide de la bibliothèque NumPy (et ne trébuchons pas sur les pierres que nous avons identifiées plus tôt)
Code pour la descente de gradient stochastique (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
Les valeurs se sont avérées presque les mêmes qu'en descendant sans utiliser NumPy. Cependant, c'est logique.
Voyons combien de temps les descentes de gradient stochastique nous ont pris.
Code pour déterminer le temps de calcul SGD (80 XNUMX étapes)
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)
Plus on s'enfonce dans la forêt, plus les nuages sont sombres : là encore, la formule « auto-écrite » donne le meilleur résultat. Tout cela suggère qu'il doit y avoir des manières encore plus subtiles d'utiliser la bibliothèque. NumPy, ce qui accélère vraiment les opérations de calcul. Dans cet article, nous n'en apprendrons pas plus à leur sujet. Il y aura quelque chose à penser pendant votre temps libre :)
Nous résumons l'
Avant de résumer, je voudrais répondre à une question qui est très probablement venue de notre cher lecteur. Pourquoi, en fait, une telle « torture » avec des descentes, pourquoi avons-nous besoin de monter et descendre la montagne (principalement en bas) pour trouver la plaine précieuse, si nous avons entre nos mains un appareil aussi puissant et simple, dans le forme d'une solution analytique, qui nous téléporte instantanément au bon endroit ?
La réponse à cette question se trouve en surface. Nous avons maintenant examiné un exemple très simple dans lequel la vraie réponse est cela dépend d'un signe . On ne voit pas cela souvent dans la vie, alors imaginons que nous ayons 2, 30, 50 signes ou plus. Ajoutons à cela des milliers, voire des dizaines de milliers de valeurs pour chaque attribut. Dans ce cas, la solution analytique peut ne pas résister au test et échouer. À son tour, la descente de gradient et ses variations nous rapprocheront lentement mais sûrement du but : le minimum de la fonction. Et ne vous inquiétez pas de la vitesse : nous chercherons probablement des moyens qui nous permettront de définir et de réguler la longueur des pas (c'est-à-dire la vitesse).
Et maintenant le bref résumé.
Premièrement, j'espère que le matériel présenté dans l'article aidera les « data scientists » débutants à comprendre comment résoudre des équations de régression linéaire simples (et pas seulement).
Deuxièmement, nous avons examiné plusieurs façons de résoudre l’équation. Désormais, en fonction de la situation, nous pouvons choisir celle qui convient le mieux pour résoudre le problème.
Troisièmement, nous avons vu la puissance des paramètres supplémentaires, à savoir la longueur du pas de descente du gradient. Ce paramètre ne peut être négligé. Comme indiqué ci-dessus, afin de réduire le coût des calculs, la longueur du pas doit être modifiée pendant la descente.
Quatrièmement, dans notre cas, les fonctions « écrites à la maison » ont montré les meilleurs résultats temporels pour les calculs. Cela est probablement dû au fait que les capacités de la bibliothèque ne sont pas utilisées de la manière la plus professionnelle possible. NumPy. Quoi qu’il en soit, la conclusion suivante s’impose. D'une part, cela vaut parfois la peine de remettre en question les opinions établies, et d'autre part, cela ne vaut pas toujours la peine de tout compliquer - au contraire, parfois une manière plus simple de résoudre un problème est plus efficace. Et puisque notre objectif était d'analyser trois approches pour résoudre une équation de régression linéaire simple, l'utilisation de fonctions « auto-écrites » nous suffisait amplement.
Littérature (ou quelque chose comme ça)
1. Régression linéaire
2. Méthode des moindres carrés
3. Dérivé
4. Dégradé
5. Descente de gradient
6. Bibliothèque NumPy
Source: habr.com