O artigo analiza varias formas de determinar a ecuación matemática dunha recta de regresión simple (parellada).
Todos os métodos para resolver a ecuación aquí discutida baséanse no método dos mínimos cadrados. Denotamos os métodos do seguinte xeito:
- Solución analítica
- Descenso Gradiente
- Descenso de gradiente estocástico
Para cada método de resolución da ecuación dunha liña recta, o artigo proporciona varias funcións, que se dividen principalmente en aquelas que se escriben sen utilizar a biblioteca. numpy e os que empregan para os cálculos numpy. Crese que o uso hábil numpy reducirá os custos informáticos.
Todo o código indicado no artigo está escrito no idioma pitón 2.7 usando Caderno Jupyter. O código fonte e o ficheiro con datos de mostra están publicados en
O artigo está máis dirixido tanto a principiantes como a aqueles que xa comezaron a dominar pouco a pouco o estudo dunha sección moi ampla en intelixencia artificial: a aprendizaxe automática.
Para ilustrar o material, utilizamos un exemplo moi sinxelo.
Condicións de exemplo
Temos cinco valores que caracterizan a dependencia Y de X (Táboa no 1):
Táboa no 1 “Condicións exemplo”
Asumiremos que os valores é o mes do ano, e - ingresos deste mes. Noutras palabras, os ingresos dependen do mes do ano, e - o único sinal do que dependen os ingresos.
O exemplo é así, tanto desde o punto de vista da dependencia condicional dos ingresos do mes do ano, como desde o punto de vista do número de valores, hai moi poucos. Non obstante, tal simplificación permitirá, como din, explicar, non sempre con facilidade, o material que asimilan os principiantes. E tamén a sinxeleza dos números permitirá a quen desexen resolver o exemplo en papel sen custos laborais significativos.
Supoñamos que a dependencia dada no exemplo pódese aproximar bastante ben mediante a ecuación matemática dunha recta de regresión simple (parellada) da forma:
onde é o mes no que se percibiron os ingresos, - Ingresos correspondentes ao mes, и son os coeficientes de regresión da recta estimada.
Teña en conta que o coeficiente moitas veces chamada pendente ou pendente da liña estimada; representa a cantidade en que o cando cambia .
Obviamente, a nosa tarefa no exemplo é seleccionar tales coeficientes na ecuación и , no que as desviacións dos nosos ingresos calculados valoran por mes a partir das verdadeiras respostas, é dicir. os valores presentados na mostra serán mínimos.
Método de mínimos cadrados
Segundo o método dos mínimos cadrados, a desviación debe calcularse cadrando. Esta técnica permítelle evitar a cancelación mutua das desviacións se teñen signos opostos. Por exemplo, se nun caso, a desviación é +5 (máis cinco), e no outro -5 (menos cinco), entón a suma das desviacións anularase entre si e ascenderá a 0 (cero). Non tes que cadrar a desviación, pero usa a propiedade do módulo e entón todas as desviacións serán positivas e acumularanse. Non nos deteremos neste punto en detalle, senón que simplemente indicaremos que, para a comodidade dos cálculos, é costume cadrar a desviación.
Este é o aspecto da fórmula coa que determinaremos a mínima suma de desviacións ao cadrado (erros):
onde é unha función da aproximación de respostas verdadeiras (é dicir, os ingresos que calculamos),
son as verdadeiras respostas (ingresos proporcionados na mostra),
é o índice mostral (número do mes no que se determina a desviación)
Diferenciamos a función, definamos as ecuacións en derivadas parciais e esteamos preparados para pasar á solución analítica. Pero primeiro, fagamos unha pequena excursión sobre o que é a diferenciación e lembremos o significado xeométrico da derivada.
Diferenciación
A diferenciación é a operación de atopar a derivada dunha función.
Para que serve a derivada? A derivada dunha función caracteriza a taxa de cambio da función e indícanos a súa dirección. Se a derivada nun punto dado é positiva, entón a función aumenta; se non, a función diminúe. E canto maior sexa o valor da derivada absoluta, maior será a taxa de cambio dos valores da función, así como máis inclinada será a gráfica da función.
Por exemplo, nas condicións dun sistema de coordenadas cartesianas, o valor da derivada no punto M(0,0) é igual a + 25 significa que nun punto dado, cando o valor se despraza á dereita por unha unidade convencional, valor aumenta en 25 unidades convencionais. Na gráfica parece un aumento bastante pronunciado dos valores dende un punto dado.
Outro exemplo. O valor da derivada é igual -0,1 significa que cando se despraza valor por unidade convencional diminúe só en 0,1 unidade convencional. Ao mesmo tempo, na gráfica da función, podemos observar unha pendente descendente apenas perceptible. Facendo unha analoxía cunha montaña, é coma se baixamos moi lentamente unha suave pendente dende unha montaña, a diferenza do exemplo anterior, onde tivemos que subir picos moi empinados :)
Así, despois de diferenciar a función por probabilidades и , definimos ecuacións en derivadas parciais de 1a orde. Despois de determinar as ecuacións, recibiremos un sistema de dúas ecuacións, resolvendo as cales poderemos seleccionar tales valores dos coeficientes и , no que os valores das derivadas correspondentes en determinados puntos cambian nunha cantidade moi, moi pequena, e no caso dunha solución analítica non cambia en absoluto. Noutras palabras, a función de erro nos coeficientes atopados alcanzará un mínimo, xa que os valores das derivadas parciais nestes puntos serán iguais a cero.
Entón, segundo as regras de diferenciación, a ecuación en derivada parcial de 1a orde con respecto ao coeficiente adoptará a forma:
Ecuación en derivada parcial de 1a orde con respecto a adoptará a forma:
Como resultado, recibimos un sistema de ecuacións que ten unha solución analítica bastante sinxela:
comezar{ecuación*}
comezar{casos}
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
rematar{casos}
fin {ecuación*}
Antes de resolver a ecuación, imos cargar previamente, comprobar que a carga é correcta e formatear os datos.
Cargando e formateando datos
Cómpre sinalar que, debido ao feito de que para a solución analítica, e posteriormente para o descenso de gradientes e gradientes estocásticos, utilizaremos o código en dúas variantes: usando a biblioteca. numpy e sen usalo, entón necesitaremos un formato de datos axeitado (ver código).
Código de carga e procesamento 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
Agora, despois de que, en primeiro lugar, cargamos os datos, en segundo lugar, comprobamos a corrección da carga e, finalmente, formateamos os datos, realizaremos a primeira visualización. O método usado a miúdo para iso é trama de parellas bibliotecas Nado no mar. No noso exemplo, debido ao número limitado, non ten sentido usar a biblioteca Nado no mar. Utilizaremos a biblioteca habitual matplotlib e só mira o 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()Gráfico no 1 “Dependencia dos ingresos do mes do ano”
Solución analítica
Usemos as ferramentas máis comúns en python e resolve o sistema de ecuacións:
comezar{ecuación*}
comezar{casos}
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
rematar{casos}
fin {ecuación*}
Segundo a regra de Cramer atoparemos o determinante xeral, así como os determinantes por e por , despois de que, dividindo o determinante entre ao determinante xeral - atopar o coeficiente , do mesmo xeito atopamos o 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)Aquí tes o que temos:
Así, atopáronse os valores dos coeficientes, estableceuse a suma das desviacións ao cadrado. Debuxemos unha liña recta no histograma de dispersión segundo os coeficientes atopados.
Código de liña 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()Cadro n.o 2 “Respostas correctas e calculadas”
Podes ver o gráfico de desviación de cada mes. No noso caso, non derivaremos del ningún valor práctico significativo, pero satisfaceremos a nosa curiosidade sobre o ben que a ecuación de regresión lineal simple caracteriza a dependencia dos ingresos do mes do ano.
Código do gráfico de desviacións
# определим функцию для формирования массива отклонений в процентах
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()Gráfico no 3 "Desviacións, %"
Non é perfecto, pero rematamos a nosa tarefa.
Escribamos unha función que, para determinar os coeficientes и usa a biblioteca numpy, máis precisamente, escribiremos dúas funcións: unha utilizando unha matriz pseudoinversa (non recomendada na práctica, xa que o proceso é computacionalmente complexo e inestable), a outra mediante unha 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 o tempo empregado na determinación dos coeficientes и , de acordo cos 3 métodos presentados.
Código para calcular o tempo 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)
Cunha pequena cantidade de datos, aparece unha función "autoescrita", que atopa os coeficientes mediante o método de Cramer.
Agora podes pasar a outras formas de atopar coeficientes и .
Descenso Gradiente
En primeiro lugar, imos definir o que é un gradiente. Simplemente, o gradiente é un segmento que indica a dirección de crecemento máximo dunha función. Por analoxía coa escalada dunha montaña, onde as caras do desnivel é onde se atopa a subida máis pronunciada ata o cumio da montaña. Desenvolvendo o exemplo coa montaña, lembramos que de feito necesitamos o descenso máis pronunciado para chegar á terra baixa o máis rápido posible, é dicir, o mínimo -o lugar onde a función non aumenta nin diminúe-. Neste punto a derivada será igual a cero. Polo tanto, non necesitamos un gradiente, senón un antigradiente. Para atopar o antigradiente só tes que multiplicar o gradiente por -1 (menos un).
Prestemos atención a que unha función pode ter varios mínimos, e despois de descender a un deles mediante o algoritmo proposto a continuación, non poderemos atopar outro mínimo, que pode ser inferior ao atopado. Relaxémonos, isto non é unha ameaza para nós! No noso caso estamos ante un único mínimo, xa que a nosa función na gráfica hai unha parábola regular. E como todos debemos saber moi ben polo curso de matemáticas do noso colexio, unha parábola só ten un mínimo.
Despois de descubrir por que necesitabamos un gradiente, e tamén que o gradiente é un segmento, é dicir, un vector con coordenadas dadas, que son precisamente os mesmos coeficientes и podemos implementar o descenso de gradientes.
Antes de comezar, suxiro ler só algunhas frases sobre o algoritmo de descenso:
- Determinamos de forma pseudoaleatoria as coordenadas dos coeficientes и . No noso exemplo, definiremos coeficientes próximos a cero. Esta é unha práctica común, pero cada caso pode ter a súa propia práctica.
- Desde coordenadas resta o valor da derivada parcial de 1a orde no punto . Entón, se a derivada é positiva, entón a función aumenta. Polo tanto, restando o valor da derivada, moverémonos no sentido contrario do crecemento, é dicir, no sentido do descenso. Se a derivada é negativa, entón a función neste punto diminúe e restando o valor da derivada movémonos na dirección do descenso.
- Realizamos unha operación similar coa coordenada : resta o valor da derivada parcial no punto .
- Para non saltar o mínimo e voar ao espazo profundo, é necesario establecer o tamaño do paso na dirección do descenso. En xeral, podes escribir un artigo completo sobre como configurar o paso correctamente e como cambialo durante o proceso de descenso para reducir os custos computacionais. Pero agora temos unha tarefa lixeiramente diferente por diante, e estableceremos o tamaño do paso mediante o método científico do "poke" ou, como se di na linguaxe común, empíricamente.
- Unha vez que somos das coordenadas dadas и restar os valores das derivadas, obtemos novas coordenadas и . Damos o seguinte paso (resta), xa a partir das coordenadas calculadas. E así o ciclo comeza unha e outra vez, ata conseguir a converxencia requirida.
Todo! Agora estamos preparados para ir á procura do desfiladeiro máis profundo da Fosa das Marianas. Imos comezar.
Código para o descenso de gradientes
# напишем функцию градиентного спуска без использования библиотеки 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
Mergullámonos ata o fondo da fosa das Marianas e alí atopamos todos os mesmos valores de coeficiente и , que é exactamente o que era de esperar.
Imos facer outro mergullo, só que esta vez, o noso vehículo de profundidade estará cheo doutras tecnoloxías, é dicir, unha biblioteca numpy.
Código para o descenso de gradientes (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 do coeficiente и inmutable.
Vexamos como cambiou o erro durante o descenso do gradiente, é dicir, como cambiou a suma das desviacións ao cadrado con cada paso.
Código para representar sumas de desviacións ao cadrado
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 das desviacións ao cadrado durante o descenso do gradiente”
Na gráfica vemos que a cada paso o erro diminúe, e despois dun certo número de iteracións observamos unha liña case horizontal.
Finalmente, imos estimar a diferenza no tempo de execución do código:
Código para determinar o tempo de cálculo do descenso do 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)
Quizais esteamos facendo algo mal, pero de novo é unha simple función "escrita na casa" que non usa a biblioteca numpy supera o tempo de cálculo dunha función usando a biblioteca numpy.
Pero non estamos parados, senón que avanzamos cara a estudar outro xeito interesante de resolver a ecuación de regresión lineal simple. Coñece!
Descenso de gradiente estocástico
Para comprender rapidamente o principio de funcionamento do descenso do gradiente estocástico, é mellor determinar as súas diferenzas co descenso do gradiente ordinario. Nós, no caso de descenso de gradientes, nas ecuacións de derivadas de и utilizou as sumas dos valores de todas as características e as respostas verdadeiras dispoñibles na mostra (é dicir, as sumas de todas as и ). No descenso de gradientes estocásticos, non usaremos todos os valores presentes na mostra, senón que seleccionaremos pseudoaleatoriamente o chamado índice de mostra e utilizaremos os seus valores.
Por exemplo, se se determina que o índice é o número 3 (tres), tomamos os valores и , entón substituímos os valores nas ecuacións derivadas e determinamos novas coordenadas. Despois, unha vez determinadas as coordenadas, determinamos de novo pseudoaleatoriamente o índice de mostra, substituímos os valores correspondentes ao índice nas ecuacións en derivadas parciais e determinamos as coordenadas dun xeito novo. и etc. ata que a converxencia se volva verde. A primeira vista, pode parecer que isto non pode funcionar en absoluto, pero é así. É certo que convén sinalar que o erro non diminúe a cada paso, pero hai certamente unha tendencia.
Cales son as vantaxes do descenso do gradiente estocástico fronte ao convencional? Se o tamaño da nosa mostra é moi grande e se mide en decenas de miles de valores, entón é moito máis doado procesar, por exemplo, un milleiro deles ao azar, en lugar da mostra completa. Aquí é onde entra en xogo o descenso de gradientes estocásticos. No noso caso, por suposto, non notaremos moita diferenza.
Vexamos o código.
Código para o descenso de gradientes estocásticos
# определим функцию стох.град.шага
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 coidadosamente os coeficientes e captámonos facendo a pregunta "Como pode ser isto?" Obtivemos outros valores de coeficiente и . Quizais o descenso do gradiente estocástico atopou parámetros máis óptimos para a ecuación? Desgraciadamente non. Basta con mirar a suma das desviacións ao cadrado e ver que con novos valores dos coeficientes, o erro é maior. Non temos présa por desesperarnos. Imos construír unha gráfica do cambio de erro.
Código para representar a suma das desviacións ao cadrado no descenso do 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 das desviacións ao cadrado durante o descenso do gradiente estocástico"
Mirando o horario, todo cae no seu sitio e agora arranxaremos todo.
Entón, que pasou? O seguinte pasou. Cando seleccionamos un mes ao azar, entón o noso algoritmo busca reducir o erro ao calcular os ingresos para o mes seleccionado. Despois seleccionamos outro mes e repetimos o cálculo, pero reducimos o erro para o segundo mes seleccionado. Agora lembra que os dous primeiros meses desvían significativamente da liña da ecuación de regresión lineal simple. Isto significa que cando se selecciona algún destes dous meses, ao reducir o erro de cada un deles, o noso algoritmo aumenta seriamente o erro para toda a mostra. Entón, que facer? A resposta é sinxela: cómpre reducir o paso de baixada. Despois de todo, ao reducir o paso de descenso, o erro tamén deixará de "saltar" cara arriba e abaixo. Ou mellor dito, o erro de "salto" non parará, pero non o fará tan rápido :) Comprobamos.
Código para executar SGD con incrementos máis pequenos
# запустим функцию, уменьшив шаг в 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 desviacións ao cadrado durante o descenso do gradiente estocástico (80 mil pasos)"
Os coeficientes melloraron, pero aínda non son ideais. Hipotéticamente, isto pódese corrixir deste xeito. Seleccionamos, por exemplo, nas últimas 1000 iteracións os valores dos coeficientes cos que se cometeu o erro mínimo. É certo, para iso tamén teremos que anotar os propios valores dos coeficientes. Non o faremos, senón que prestaremos atención ao horario. Parece suave e o erro parece diminuír uniformemente. En realidade, isto non é certo. Vexamos as 1000 primeiras iteracións e comparámolas coas últimas.
Código para o gráfico SGD (primeiros 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 desviacións ao cadrado SGD (primeiros 1000 pasos)"
Gráfico no 8 "Suma de desviacións ao cadrado SGD (últimos 1000 pasos)"
Ao comezo mesmo do descenso observamos unha diminución do erro bastante uniforme e pronunciada. Nas últimas iteracións, vemos que o erro dá a volta ao valor de 1,475 e nalgúns momentos mesmo iguala este valor óptimo, pero despois aínda aumenta... Repito, podes anotar os valores do coeficientes и e, a continuación, seleccione aqueles para os que o erro é mínimo. Non obstante, tivemos un problema máis grave: tivemos que dar 80 mil pasos (ver código) para conseguir valores próximos ao óptimo. E isto xa contradí a idea de aforrar tempo de cálculo co descenso do gradiente estocástico en relación ao descenso do gradiente. Que se pode corrixir e mellorar? Non é difícil notar que nas primeiras iteracións imos baixando con confianza e, polo tanto, debemos deixar un gran paso nas primeiras iteracións e ir reducindo o paso a medida que avanzamos. Non o faremos neste artigo: xa é demasiado longo. Os que o desexen poden pensar por si mesmos como facelo, non é difícil :)
Agora imos realizar un descenso de gradientes estocásticos usando a biblioteca numpy (e non tropecemos coas pedras que identificamos antes)
Código para Descenso de 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
Os valores resultaron ser case os mesmos que ao descender sen usar numpy. Non obstante, isto é lóxico.
Descubramos canto tempo nos levaron os descensos en gradiente estocástico.
Código para determinar o tempo 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)
Canto máis dentro do bosque, máis escuras son as nubes: de novo, a fórmula "escrita por si mesmo" mostra o mellor resultado. Todo isto suxire que debe haber formas aínda máis sutís de usar a biblioteca numpy, que realmente aceleran as operacións de cálculo. Neste artigo non aprenderemos sobre eles. Haberá algo no que pensar no teu tempo libre :)
Resumamos
Antes de resumir, gustaríame responder a unha pregunta que moi probablemente xurdiu do noso querido lector. Por que, de feito, tal “tortura” con baixadas, por que necesitamos andar arriba e abaixo da montaña (principalmente abaixo) para atopar a preciada terra baixa, se temos nas nosas mans un dispositivo tan poderoso e sinxelo, no forma dunha solución analítica, que nos teletransporta instantáneamente ao lugar correcto?
A resposta a esta pregunta está na superficie. Agora miramos un exemplo moi sinxelo, no que a verdadeira resposta é depende dun signo . Non ves isto a miúdo na vida, así que imaxinemos que temos 2, 30, 50 ou máis signos. Imos engadir a isto miles, ou incluso decenas de miles de valores para cada atributo. Neste caso, a solución analítica pode non soportar a proba e fallar. Á súa vez, o descenso en gradiente e as súas variacións achegaranos lenta pero seguramente ao obxectivo: o mínimo da función. E non te preocupes pola velocidade: probablemente busquemos formas que nos permitan establecer e regular a lonxitude do paso (é dicir, a velocidade).
E agora o breve resumo real.
En primeiro lugar, espero que o material presentado no artigo axude aos "científicos de datos" iniciados a comprender como resolver ecuacións de regresión lineal sinxelas (e non só).
En segundo lugar, analizamos varias formas de resolver a ecuación. Agora, segundo a situación, podemos elixir a que máis se adapte para resolver o problema.
En terceiro lugar, vimos o poder das configuracións adicionais, é dicir, a lonxitude do paso de baixada do gradiente. Este parámetro non se pode descoidar. Como se indicou anteriormente, para reducir o custo dos cálculos, a lonxitude do paso debe cambiarse durante o descenso.
En cuarto lugar, no noso caso, as funcións "escritas na casa" mostraron os mellores resultados de tempo para os cálculos. Isto probablemente débese ao uso máis profesional das capacidades da biblioteca numpy. Pero sexa como for, suxire a seguinte conclusión. Por unha banda, ás veces paga a pena cuestionar opinións establecidas e, por outra banda, non sempre paga a pena complicalo todo, pola contra, ás veces é máis eficaz unha forma máis sinxela de resolver un problema. E dado que o noso obxectivo era analizar tres enfoques para resolver unha ecuación de regresión lineal sinxela, o uso de funcións "escritas por conta propia" abondaba para nós.
Literatura (ou algo así)
1. Regresión lineal
2. Método dos mínimos cadrados
3. Derivada
4 Gradiente
5. Descenso en pendiente
6. Biblioteca NumPy
Fonte: www.habr.com