El artículo analiza varias formas de determinar la ecuación matemática de una recta de regresión simple (emparejada).
Todos los métodos para resolver la ecuación discutida aquí se basan en el método de mínimos cuadrados. Denotemos los métodos de la siguiente manera:
- Solucion analitica
- Descenso de gradiente
- Descenso de gradiente estocástico
Para cada método para resolver la ecuación de una línea recta, el artículo proporciona varias funciones, que se dividen principalmente en aquellas que se escriben sin usar la biblioteca. NumPy y los que se utilizan para los cálculos. NumPy. Se cree que el uso hábil. NumPy reducirá los costos de computación.
Todo el código proporcionado en el artículo está escrito en el idioma. python 2.7 utilizando Cuaderno Jupyter. El código fuente y el archivo con datos de muestra se publican en
El artículo está más dirigido tanto a principiantes como a aquellos que ya han comenzado gradualmente a dominar el estudio de una sección muy amplia de la inteligencia artificial: el aprendizaje automático.
Para ilustrar el material, utilizamos un ejemplo muy sencillo.
Condiciones de ejemplo
Tenemos cinco valores que caracterizan la dependencia. Y de X (Tabla No. 1):
Cuadro No. 1 “Condiciones de ejemplo”
Supondremos que los valores es el mes del año y - ingresos este mes. En otras palabras, los ingresos dependen del mes del año, y - el único signo del que dependen los ingresos.
El ejemplo es regular, tanto desde el punto de vista de la dependencia condicional de los ingresos del mes del año como desde el punto de vista del número de valores: hay muy pocos. Sin embargo, tal simplificación permitirá, como suele decirse, explicar, no siempre con facilidad, el material que asimilan los principiantes. Y además la sencillez de los números permitirá a quienes quieran resolver el ejemplo en papel sin costes laborales importantes.
Supongamos que la dependencia dada en el ejemplo puede aproximarse bastante bien mediante la ecuación matemática de una recta de regresión simple (emparejada) de la forma:
donde es el mes en el que se recibieron los ingresos, — ingresos correspondientes al mes, и son los coeficientes de regresión de la recta estimada.
Tenga en cuenta que el coeficiente a menudo se le llama pendiente o gradiente de la línea estimada; representa la cantidad por la cual al cambiar .
Obviamente, nuestra tarea en el ejemplo es seleccionar dichos coeficientes en la ecuación и , en el que las desviaciones de nuestros valores de ingresos calculados por mes de las respuestas verdaderas, es decir Los valores presentados en la muestra serán mínimos.
Método de mínimos cuadrados
Según el método de mínimos cuadrados, la desviación debe calcularse elevándola al cuadrado. Esta técnica permite evitar la cancelación mutua de desviaciones si tienen signos opuestos. Por ejemplo, si en un caso, la desviación es +5 (más cinco), y en el otro -5 (menos cinco), entonces la suma de las desviaciones se cancelará entre sí y ascenderá a 0 (cero). Es posible no elevar al cuadrado la desviación, sino utilizar la propiedad del módulo y entonces todas las desviaciones serán positivas y se acumularán. No nos detendremos en este punto en detalle, simplemente indicaremos que para facilitar los cálculos, se acostumbra elevar al cuadrado la desviación.
Así luce la fórmula con la que determinaremos la mínima suma de desviaciones al cuadrado (errores):
donde es una función de aproximación de las respuestas verdaderas (es decir, los ingresos que calculamos),
son las respuestas verdaderas (ingresos proporcionados en la muestra),
es el índice de la muestra (número del mes en el que se determina la desviación)
Diferenciamos la función, definamos las ecuaciones diferenciales parciales y estemos listos para pasar a la solución analítica. Pero primero, hagamos un breve recorrido sobre qué es la diferenciación y recordemos el significado geométrico de la derivada.
Diferenciación
La derivación es la operación de encontrar la derivada de una función.
¿Para qué se utiliza la derivada? La derivada de una función caracteriza la tasa de cambio de la función y nos dice su dirección. Si la derivada en un punto dado es positiva, entonces la función aumenta; en caso contrario, la función disminuye. Y cuanto mayor sea el valor de la derivada absoluta, mayor será la tasa de cambio de los valores de la función, así como más pronunciada la pendiente de la gráfica de la función.
Por ejemplo, bajo las condiciones de un sistema de coordenadas cartesiano, el valor de la derivada en el punto M(0,0) es igual a +25 significa que en un punto dado, cuando el valor se desplaza a la derecha por una unidad convencional, valor aumenta en 25 unidades convencionales. En el gráfico parece un aumento bastante pronunciado de los valores. desde un punto dado.
Otro ejemplo. El valor de la derivada es igual. -0,1 significa que cuando se desplaza por unidad convencional, valor disminuye en sólo 0,1 unidad convencional. Al mismo tiempo, en la gráfica de la función podemos observar una pendiente descendente apenas perceptible. Haciendo una analogía con una montaña, es como si estuviéramos descendiendo muy lentamente por una suave pendiente de una montaña, a diferencia del ejemplo anterior, donde tuvimos que subir picos muy empinados :)
Así, después de derivar la función por probabilidades и , definimos ecuaciones diferenciales parciales de primer orden. Después de determinar las ecuaciones, obtendremos un sistema de dos ecuaciones, al resolver el cual podremos seleccionar dichos valores de los coeficientes. и , para el cual los valores de las derivadas correspondientes en puntos dados cambian en una cantidad muy, muy pequeña, y en el caso de una solución analítica no cambian en absoluto. En otras palabras, la función de error en los coeficientes encontrados alcanzará un mínimo, ya que los valores de las derivadas parciales en estos puntos serán iguales a cero.
Entonces, de acuerdo con las reglas de diferenciación, la ecuación derivada parcial de primer orden con respecto al coeficiente tomará la forma:
Ecuación derivada parcial de primer orden con respecto a tomará la forma:
Como resultado, obtuvimos un sistema de ecuaciones que tiene una solución analítica bastante simple:
comenzar {ecuación*}
comenzar {casos}
na + bsumlimits_{i=1}^nx_i — sumalimits_{i=1}^ny_i = 0
sumalimits_{i=1}^nx_i(a +bsumlimits_{i=1}^nx_i - sumalimits_{i=1}^ny_i) = 0
fin {casos}
fin {ecuación*}
Antes de resolver la ecuación, realicemos una precarga, verifiquemos que la carga sea correcta y formateemos los datos.
Cargando y formateando datos
Cabe señalar que debido al hecho de que para la solución analítica, y posteriormente para el gradiente y el descenso de gradiente estocástico, usaremos el código en dos variaciones: usando la biblioteca NumPy y sin usarlo, necesitaremos un formato de datos adecuado (ver código).
Código de carga y procesamiento de datos.
# импортируем все нужные нам библиотеки
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 '********************************************'Visualización
Ahora, después de haber cargado, en primer lugar, los datos, en segundo lugar, haber comprobado la exactitud de la carga y finalmente haber formateado los datos, realizaremos la primera visualización. El método más utilizado para esto es parcela Biblioteca nacido en el mar. En nuestro ejemplo, debido al número limitado, no tiene sentido utilizar la biblioteca. nacido en el mar. Usaremos la biblioteca normal. matplotlib y solo mira el diagrama de dispersión.
código de diagrama de dispersión
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()Cuadro No. 1 “Dependencia de los ingresos del mes del año”
Solucion analitica
Utilicemos las herramientas más comunes en pitón y resuelve el sistema de ecuaciones:
comenzar {ecuación*}
comenzar {casos}
na + bsumlimits_{i=1}^nx_i — sumalimits_{i=1}^ny_i = 0
sumalimits_{i=1}^nx_i(a +bsumlimits_{i=1}^nx_i - sumalimits_{i=1}^ny_i) = 0
fin {casos}
fin {ecuación*}
Según la regla de Cramer encontraremos el determinante general, así como los determinantes por y , después de lo cual, dividiendo el determinante por al determinante general - encuentre el coeficiente , de manera similar encontramos el coeficiente .
Código de solución analítica
# определим функцию для расчета коэффициентов 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)Esto es lo que tenemos:
Entonces, se encontraron los valores de los coeficientes, se estableció la suma de las desviaciones al cuadrado. Dibujemos una línea recta en el histograma de dispersión de acuerdo con los coeficientes encontrados.
Código de línea de regresión
# определим функцию для формирования массива рассчетных значений выручки
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()Cuadro No. 2 “Respuestas correctas y calculadas”
Puedes consultar el gráfico de desviación de cada mes. En nuestro caso, no obtendremos ningún valor práctico significativo de ello, pero satisfaremos nuestra curiosidad sobre qué tan bien la ecuación de regresión lineal simple caracteriza la dependencia de los ingresos del mes del año.
Código del gráfico de desviación
# определим функцию для формирования массива отклонений в процентах
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()Cuadro No. 3 “Desviaciones, %”
No es perfecto, pero completamos nuestra tarea.
Escribamos una función que, para determinar los coeficientes. и usa la biblioteca NumPy, más precisamente, escribiremos dos funciones: una usando una matriz pseudoinversa (no recomendada en la práctica, ya que el proceso es computacionalmente complejo e inestable), la otra usando una ecuación matricial.
Código de solución analítica (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_npComparemos el tiempo dedicado a determinar los coeficientes. и , de acuerdo con los 3 métodos presentados.
Código para calcular el tiempo de cálculo.
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)
Con una pequeña cantidad de datos, sale adelante una función “autoescrita” que encuentra los coeficientes utilizando el método de Cramer.
Ahora puedes pasar a otras formas de encontrar coeficientes. и .
Descenso de gradiente
Primero, definamos qué es un gradiente. En pocas palabras, el gradiente es un segmento que indica la dirección de máximo crecimiento de una función. Por analogía con escalar una montaña, donde se enfrenta la pendiente es donde se encuentra la subida más empinada hasta la cima de la montaña. Desarrollando el ejemplo con la montaña, recordamos que en realidad necesitamos el descenso más pronunciado para llegar a la tierra baja lo más rápido posible, es decir, el mínimo, el lugar donde la función no aumenta ni disminuye. En este punto la derivada será igual a cero. Por tanto, no necesitamos un gradiente, sino un antigradiente. Para encontrar el antigradiente solo necesitas multiplicar el gradiente por -1 (menos uno).
Prestemos atención al hecho de que una función puede tener varios mínimos, y habiendo descendido a uno de ellos mediante el algoritmo propuesto a continuación, no podremos encontrar otro mínimo, que puede ser menor que el encontrado. ¡Relajémonos, esto no es una amenaza para nosotros! En nuestro caso estamos ante un mínimo único, ya que nuestra función en la gráfica hay una parábola regular. Y como todos deberíamos saber muy bien en nuestro curso de matemáticas escolar, una parábola tiene solo un mínimo.
Después descubrimos por qué necesitábamos un gradiente, y también que el gradiente es un segmento, es decir, un vector con coordenadas dadas, que son exactamente los mismos coeficientes. и Podemos implementar el descenso de gradiente.
Antes de comenzar, sugiero leer sólo unas frases sobre el algoritmo de descenso:
- Determinamos de forma pseudoaleatoria las coordenadas de los coeficientes. и . En nuestro ejemplo, determinaremos coeficientes cercanos a cero. Esta es una práctica común, pero cada caso puede tener su propia práctica.
- Desde coordenada restar el valor de la derivada parcial de primer orden en el punto . Entonces, si la derivada es positiva, entonces la función aumenta. Por tanto, al restar el valor de la derivada, nos moveremos en el sentido contrario al crecimiento, es decir, en el sentido de descenso. Si la derivada es negativa, entonces la función en este punto disminuye y restando el valor de la derivada nos movemos en la dirección de descenso.
- Realizamos una operación similar con la coordenada. : restar el valor de la derivada parcial en el punto .
- Para no saltar por encima del mínimo y volar al espacio profundo, es necesario establecer el tamaño del paso en la dirección de descenso. En general, se podría escribir un artículo completo sobre cómo configurar el paso correctamente y cómo cambiarlo durante el proceso de descenso para reducir los costos computacionales. Pero ahora tenemos una tarea ligeramente diferente y estableceremos el tamaño del paso utilizando el método científico de "empujar" o, como dicen en el lenguaje común, empíricamente.
- Una vez que estemos en las coordenadas dadas и restamos los valores de las derivadas, obtenemos nuevas coordenadas и . Damos el siguiente paso (resta), ya a partir de las coordenadas calculadas. Y así el ciclo comienza una y otra vez, hasta que se logra la convergencia requerida.
¡Todo! Ahora estamos listos para salir en busca del desfiladero más profundo de la Fosa de las Marianas. Empecemos.
Código para descenso de gradiente
# напишем функцию градиентного спуска без использования библиотеки 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
Nos sumergimos hasta el fondo de la Fosa de las Marianas y allí encontramos los mismos valores de coeficientes. и , que es exactamente lo que se esperaba.
Hagamos otra inmersión, solo que esta vez nuestro vehículo de aguas profundas estará lleno de otras tecnologías, en particular una biblioteca. NumPy.
Código para descenso de gradiente (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
Valores de coeficiente и inmutable.
Veamos cómo cambió el error durante el descenso del gradiente, es decir, cómo cambió la suma de las desviaciones al cuadrado con cada paso.
Código para trazar sumas de desviaciones al cuadrado.
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()Gráfico No. 4 “Suma de desviaciones al cuadrado durante el descenso del gradiente”
En el gráfico vemos que con cada paso el error disminuye, y después de un cierto número de iteraciones observamos una línea casi horizontal.
Finalmente, estimemos la diferencia en el tiempo de ejecución del código:
Código para determinar el tiempo de cálculo del descenso de gradiente
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)
Quizás estemos haciendo algo mal, pero nuevamente es una función simple "escrita en casa" que no utiliza la biblioteca. NumPy supera el tiempo de cálculo de una función usando la biblioteca NumPy.
Pero no nos quedamos quietos, sino que avanzamos hacia el estudio de otra forma interesante de resolver la ecuación de regresión lineal simple. ¡Encontrarse!
Descenso de gradiente estocástico
Para comprender rápidamente el principio de funcionamiento del descenso de gradiente estocástico, es mejor determinar sus diferencias con el descenso de gradiente ordinario. Nosotros, en el caso del descenso de gradiente, en las ecuaciones de derivadas de и utilizó las sumas de los valores de todas las características y respuestas verdaderas disponibles en la muestra (es decir, las sumas de todas и ). En el descenso de gradiente estocástico, no usaremos todos los valores presentes en la muestra, sino que seleccionaremos pseudoaleatoriamente el llamado índice de muestra y usaremos sus valores.
Por ejemplo, si se determina que el índice es el número 3 (tres), entonces tomamos los valores и , luego sustituimos los valores en las ecuaciones derivadas y determinamos nuevas coordenadas. Luego, una vez determinadas las coordenadas, nuevamente determinamos pseudoaleatoriamente el índice de muestra, sustituimos los valores correspondientes al índice en las ecuaciones diferenciales parciales y determinamos las coordenadas de una nueva manera. и etc. hasta que la convergencia se vuelva verde. A primera vista, puede parecer que esto no podría funcionar en absoluto, pero funciona. Es cierto que cabe destacar que el error no disminuye con cada paso, pero ciertamente hay una tendencia.
¿Cuáles son las ventajas del descenso de gradiente estocástico sobre el convencional? Si el tamaño de nuestra muestra es muy grande y se mide en decenas de miles de valores, entonces es mucho más fácil procesar, digamos, miles de ellos al azar, en lugar de toda la muestra. Aquí es donde entra en juego el descenso de gradiente estocástico. En nuestro caso, por supuesto, no notaremos mucha diferencia.
Miremos el código.
Código para descenso de gradiente estocástico
# определим функцию стох.град.шага
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])
Observamos detenidamente los coeficientes y nos sorprendemos preguntándonos "¿Cómo puede ser esto?" Tenemos otros valores de coeficientes. и . ¿Quizás el descenso de gradiente estocástico haya encontrado parámetros más óptimos para la ecuación? Lamentablemente no. Basta mirar la suma de las desviaciones al cuadrado y ver que con nuevos valores de los coeficientes, el error es mayor. No tenemos prisa por desesperarnos. Construyamos un gráfico del cambio de error.
Código para trazar la suma de las desviaciones al cuadrado en el descenso de gradiente estocástico
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()Gráfico No. 5 “Suma de desviaciones al cuadrado durante el descenso del gradiente estocástico”
Mirando el cronograma, todo encaja y ahora lo arreglaremos todo.
¿Entonces qué pasó? Sucedió lo siguiente. Cuando seleccionamos un mes al azar, es para el mes seleccionado que nuestro algoritmo busca reducir el error en el cálculo de los ingresos. Luego seleccionamos otro mes y repetimos el cálculo, pero reducimos el error para el segundo mes seleccionado. Ahora recordemos que los primeros dos meses se desvían significativamente de la recta de la ecuación de regresión lineal simple. Esto significa que cuando se selecciona cualquiera de estos dos meses, al reducir el error de cada uno de ellos, nuestro algoritmo aumenta seriamente el error de toda la muestra. ¿Entonces lo que hay que hacer? La respuesta es simple: es necesario reducir el paso de descenso. Después de todo, al reducir el paso de descenso, el error también dejará de “saltar” hacia arriba y hacia abajo. O mejor dicho, el error de “salto” no se detendrá, pero no lo hará tan rápido :) Comprobemos.
Código para ejecutar SGD con incrementos más pequeños
# запустим функцию, уменьшив шаг в 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()
Gráfico No. 6 “Suma de desviaciones al cuadrado durante el descenso del gradiente estocástico (80 mil pasos)”
Los coeficientes han mejorado, pero todavía no son ideales. Hipotéticamente, esto se puede corregir de esta manera. Seleccionamos, por ejemplo, en las últimas 1000 iteraciones los valores de los coeficientes con los que se cometió el error mínimo. Es cierto que para ello también tendremos que anotar los valores de los propios coeficientes. No haremos esto, sino que prestaremos atención al cronograma. Parece suave y el error parece disminuir uniformemente. Actualmente, esto no es verdad. Veamos las primeras 1000 iteraciones y compárelas con la última.
Código para gráfico SGD (primeros 1000 pasos)
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()Gráfico No. 7 “Suma de desviaciones al cuadrado SGD (primeros 1000 pasos)”
Gráfico No. 8 “Suma de desviaciones al cuadrado SGD (últimos 1000 pasos)”
Al comienzo del descenso, observamos una disminución bastante uniforme y pronunciada del error. En las últimas iteraciones vemos que el error da vueltas y vueltas alrededor del valor de 1,475 y en algunos momentos incluso iguala este valor óptimo, pero luego sigue subiendo... Repito, puedes anotar los valores de coeficientes и y luego seleccione aquellos para los cuales el error es mínimo. Sin embargo, tuvimos un problema más grave: tuvimos que dar 80 mil pasos (ver código) para conseguir valores cercanos al óptimo. Y esto ya contradice la idea de ahorrar tiempo de cálculo con el descenso de gradiente estocástico en relación con el descenso de gradiente. ¿Qué se puede corregir y mejorar? No es difícil notar que en las primeras iteraciones vamos bajando con confianza y, por lo tanto, debemos dejar un gran paso en las primeras iteraciones y reducir el paso a medida que avanzamos. No haremos esto en este artículo; ya es demasiado largo. Quien lo desee puede pensar por sí mismo cómo hacerlo, no es difícil :)
Ahora realicemos un descenso de gradiente estocástico usando la biblioteca. NumPy (y no tropecemos con las piedras que identificamos anteriormente)
Código para el descenso del gradiente estocástico (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
Los valores resultaron ser casi los mismos que al descender sin usar. NumPy. Sin embargo, esto es lógico.
Averigüemos cuánto tiempo nos llevaron los descensos de gradiente estocástico.
Código para determinar el tiempo de cálculo de SGD (80 mil pasos)
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)
Cuanto más se adentra en el bosque, más oscuras son las nubes: una vez más, la fórmula "autoescrita" muestra el mejor resultado. Todo esto sugiere que debe haber formas aún más sutiles de utilizar la biblioteca. NumPy, que realmente aceleran las operaciones de cálculo. En este artículo no aprenderemos sobre ellos. Habrá algo en qué pensar en tu tiempo libre :)
Vamos a resumir
Antes de resumir, me gustaría responder a una pregunta que muy probablemente le surgió a nuestro querido lector. ¿Por qué, de hecho, tanta “tortura” con los descensos, por qué tenemos que subir y bajar la montaña (principalmente hacia abajo) para encontrar las preciadas tierras bajas, si tenemos en nuestras manos un dispositivo tan poderoso y simple, en el forma de solución analítica, que instantáneamente nos teletransporta al lugar correcto?
La respuesta a esta pregunta está en la superficie. Ahora hemos visto un ejemplo muy simple, en el que la verdadera respuesta es depende de un signo . Esto no se ve a menudo en la vida, así que imaginemos que tenemos 2, 30, 50 o más signos. Agreguemos a esto miles, o incluso decenas de miles de valores para cada atributo. En este caso, es posible que la solución analítica no resista la prueba y falle. A su vez, el descenso del gradiente y sus variaciones nos acercarán lenta pero seguramente al objetivo: el mínimo de la función. Y no se preocupe por la velocidad: probablemente buscaremos formas que nos permitan establecer y regular la longitud del paso (es decir, la velocidad).
Y ahora el breve resumen real.
En primer lugar, espero que el material presentado en el artículo ayude a los "científicos de datos" principiantes a comprender cómo resolver ecuaciones de regresión lineal simples (y no solo).
En segundo lugar, analizamos varias formas de resolver la ecuación. Ahora bien, dependiendo de la situación, podremos elegir el que mejor se adapte a solucionar el problema.
En tercer lugar, vimos el poder de las configuraciones adicionales, es decir, la longitud del paso de descenso del gradiente. Este parámetro no se puede descuidar. Como se señaló anteriormente, para reducir el costo de los cálculos, se debe cambiar la longitud del paso durante el descenso.
En cuarto lugar, en nuestro caso, las funciones "escritas en casa" mostraron los mejores resultados de tiempo para los cálculos. Probablemente esto se deba a que el uso de las capacidades de la biblioteca no es el más profesional. NumPy. Pero sea como fuere, se sugiere la siguiente conclusión. Por un lado, a veces vale la pena cuestionar las opiniones establecidas y, por otro, no siempre vale la pena complicarlo todo; por el contrario, a veces una forma más sencilla de resolver un problema es más eficaz. Y como nuestro objetivo era analizar tres enfoques para resolver una ecuación de regresión lineal simple, el uso de funciones "autoescritas" fue suficiente para nosotros.
Literatura (o algo así)
1. Regresión lineal
2. Método de mínimos cuadrados
3. Derivado
4. Gradiente
5. Descenso de gradiente
6. Biblioteca NumPy
Fuente: habr.com