Gradient Descent en Python: implementación y teoría

Introducción

Este tutorial es una introducción a una técnica de optimización simple llamada descenso de gradiente, que ha tenido una aplicación importante en modelos de Machine Learning de última generación.

Desarrollaremos una rutina de propósito general para implementar el descenso de gradientes y aplicarlo para resolver diferentes problemas, incluida la clasificación a través del aprendizaje supervisado.

En este proceso, obtendremos una idea del funcionamiento de este algoritmo y estudiaremos el efecto de varios hiperparámetros en su rendimiento. También repasaremos variantes de descenso de gradiente estocástico y por lotes como ejemplos.

¿Qué es Gradient Descent?

El descenso de gradientes es una técnica de optimización que puede encontrar el mínimo de un función objetiva. Es una técnica codiciosa que encuentra la solución óptima dando un paso en la dirección de la tasa máxima de disminución de la función.

Por el contrario, Gradient Ascent es una contraparte cercana que encuentra el máximo de una función siguiendo la dirección de la tasa máxima de aumento de la función.

Para comprender cómo funciona el descenso de gradientes, considere una función de múltiples variables (f ( textbf {w}) ), donde ( textbf w = [w_1, w_2, ldots, w_n]^ T ). Para encontrar el ( textbf {w} ) en el que esta función alcanza un mínimo, el descenso de gradiente utiliza los siguientes pasos:

• Elija un valor aleatorio inicial de ( textbf {w} )
• Elija el número máximo de iteraciones T
• Elija un valor para el tasa de aprendizaje ( eta en [a,b] )
• Repita los siguientes dos pasos hasta que (f ) no cambie o las iteraciones excedan T

  a.Cálculo: ( Delta textbf {w} = – eta nabla_ textbf {w} f ( textbf {w}) )

  segundo. actualizar ( textbf {w} ) como: ( textbf {w} leftarrow textbf {w} + Delta textbf {w} )

  El símbolo de la flecha izquierda es una forma matemáticamente correcta de escribir una declaración de asignación.

Aquí ( nabla_ textbf {w} f ) denota el gradiente de (f ) dado por:
$$
nabla_ textbf {w} f ( textbf {w}) =
begin {bmatrix}
frac { parcial f ( textbf {w})} { parcial w_1}
frac { parcial f ( textbf {w})} { parcial w_2}
vdots
frac { parcial f ( textbf {w})} { parcial w_n}
end {bmatrix}
$$

Considere una función de ejemplo de dos variables (f (w_1, w_2) = w_1 ^ 2 + w_2 ^ 2 ), luego en cada iteración ((w_1, w_2) ) se actualiza como:

$$
begin {bmatrix}
w_1 w_2
end {bmatrix} leftarrow
begin {bmatrix}
w_1 w_2
end {bmatrix} – eta
begin {bmatrix}
2w_1 2w_2
end {bmatrix}
$$

La siguiente figura muestra cómo funciona el descenso de gradientes en esta función.

Los círculos son los contornos de esta función. Si nos movemos a lo largo de un contorno, el valor de la función no cambiaría y seguiría siendo un constante.

Esto se opone a la dirección del gradiente, donde la función cambia a una velocidad máxima. Por lo tanto, la dirección del gradiente de la función en cualquier punto es normal a la tangente del contorno en ese punto.

En términos simples, el gradiente se puede tomar como una flecha que apunta en la dirección donde la función cambia más.

Seguir la dirección del gradiente negativo conduciría a puntos donde el valor de la función disminuye a una velocidad máxima. los tasa de aprendizaje, también llamado el Numero de pie, dicta qué tan rápido o lento nos movemos a lo largo de la dirección del gradiente.

Añadiendo impulso

Cuando usamos el descenso de gradientes, nos encontramos con los siguientes problemas:

• Quedar atrapado en un mínimo local, que es una consecuencia directa de que este algoritmo sea codicioso
• Sobrepasar y perder el óptimo global, esto es un resultado directo de moverse demasiado rápido a lo largo de la dirección del gradiente
• Oscilación, este es un fenómeno que ocurre cuando el valor de la función no cambia significativamente sin importar la dirección en la que avanza. Puedes pensar en ello como navegar por una meseta, estás a la misma altura sin importar a dónde vayas

Para combatir estos problemas, se agrega un término de impulso ( alpha ) a la expresión para ( Delta textbf {w} ) para estabilizar la tasa de aprendizaje cuando se avanza hacia el valor óptimo global.

A continuación, usamos el superíndice (i ) para denotar el número de iteración:
$$
Delta textbf {w} ^ i = – eta nabla_ textbf {w} f ( textbf {w} ^ i) + alpha textbf {w} ^ {i-1}
$$

Implementando Gradient Descent en Python

Antes de comenzar a escribir el código real para el descenso de gradientes, importemos algunas bibliotecas que utilizaremos para ayudarnos:

import numpy as np
import matplotlib
import matplotlib.pyplot as plt
import sklearn.datasets as dt
from sklearn.model_selection import train_test_split

Ahora, con eso fuera del camino, sigamos adelante y definamos un gradient_descent() función. En esta función, el ciclo termina cuando:

• El número de iteraciones supera un valor máximo
• La diferencia en los valores de la función entre dos iteraciones sucesivas cae por debajo de un cierto umbral

Los parámetros se actualizan en cada iteración de acuerdo con el gradiente de la función objetivo.

La función aceptará los siguientes parámetros:

• max_iterations: Número máximo de iteraciones para ejecutar
• threshold: Deténgase si la diferencia en los valores de la función entre dos iteraciones sucesivas cae por debajo de este umbral
• w_init: Punto inicial desde donde comenzar el descenso de gradiente
• obj_func: Referencia a la función que calcula la función objetivo
• grad_func: Referencia a la función que calcula el gradiente de la función
• extra_param: Parámetros adicionales (si es necesario) para obj_func y grad_func
• learning_rate: Tamaño de paso para descenso en pendiente. Debería estar en [0,1]
• momentum: Momentum para usar. Debería estar en [0,1]

Además, la función devolverá:

• w_history: Todos los puntos en el espacio, visitados por descenso de gradiente en los que se evaluó la función objetivo
• f_history: Valor correspondiente de la función objetivo calculada en cada punto

# Make threshold a -ve value if you want to run exactly
# max_iterations.
def gradient_descent(max_iterations,threshold,w_init,
                     obj_func,grad_func,extra_param = [],
                     learning_rate=0.05,momentum=0.8):
    
    w = w_init
    w_history = w
    f_history = obj_func(w,extra_param)
    delta_w = np.zeros(w.shape)
    i = 0
    diff = 1.0e10
    
    while  i<max_iterations and diff>threshold:
        delta_w = -learning_rate*grad_func(w,extra_param) + momentum*delta_w
        w = w+delta_w
        
        # store the history of w and f
        w_history = np.vstack((w_history,w))
        f_history = np.vstack((f_history,obj_func(w,extra_param)))
        
        # update iteration number and diff between successive values
        # of objective function
        i+=1
        diff = np.absolute(f_history[-1]-f_history[-2])
    
    return w_history,f_history

Optimización de funciones con pendiente descendente

Ahora que tenemos una implementación de propósito general del descenso de gradiente, ejecutémosla en nuestra función 2D de ejemplo (f (w_1, w_2) = w_1 ^ 2 + w_2 ^ 2 ) con contornos circulares.

La función tiene un valor mínimo de cero en el origen. Primero visualicemos la función y luego encontremos su valor mínimo.

Visualización de la función objetivo f (x)

los visualize_fw() función a continuación, genera 2500 puntos igualmente espaciados en una cuadrícula y calcula el valor de la función en cada punto.

los function_plot() La función muestra todos los puntos en diferentes colores, dependiendo del valor de (f ( textbf w) ) en ese punto. Todos los puntos en los que el valor de la función es el mismo, tienen el mismo color:

def visualize_fw():
    xcoord = np.linspace(-10.0,10.0,50)
    ycoord = np.linspace(-10.0,10.0,50)
    w1,w2 = np.meshgrid(xcoord,ycoord)
    pts = np.vstack((w1.flatten(),w2.flatten()))
    
    # All 2D points on the grid
    pts = pts.transpose()
    
    # Function value at each point
    f_vals = np.sum(pts*pts,axis=1)
    function_plot(pts,f_vals)
    plt.title('Objective Function Shown in Color')
    plt.show()
    return pts,f_vals

# Helper function to annotate a single point
def annotate_pt(text,xy,xytext,color):
    plt.plot(xy[0],xy[1],marker="P",markersize=10,c=color)
    plt.annotate(text,xy=xy,xytext=xytext,
                 # color=color,
                 arrowprops=dict(arrowstyle="->",
                 color = color,
                 connectionstyle="arc3"))

# Plot the function
# Pts are 2D points and f_val is the corresponding function value
def function_plot(pts,f_val):
    f_plot = plt.scatter(pts[:,0],pts[:,1],
                         c=f_val,vmin=min(f_val),vmax=max(f_val),
                         cmap='RdBu_r')
    plt.colorbar(f_plot)
    # Show the optimal point
    annotate_pt('global minimum',(0,0),(-5,-7),'yellow')    

pts,f_vals = visualize_fw()

Ejecución de Gradient Descent con diferentes hiperparámetros

Ahora es el momento de ejecutar el descenso de gradiente para minimizar nuestra función objetivo. Llamar gradient_descent(), definimos dos funciones:

• f(): Calcula la función objetivo en cualquier punto w
• grad(): Calcula el gradiente en cualquier punto w

Para comprender el efecto de varios hiperparámetros en el descenso de gradientes, la función solve_fw() llamadas gradient_descent() con 5 iteraciones para diferentes valores de velocidad e impulso de aprendizaje.

La función visualize_learning(), traza los valores de ((w_1, w_2) ), con los valores de la función mostrados en diferentes colores. Las flechas en el gráfico facilitan el seguimiento de qué punto se actualizó desde el último:

# Objective function
def f(w,extra=[]):
    return np.sum(w*w)

# Function to compute the gradient
def grad(w,extra=[]):
    return 2*w

# Function to plot the objective function
# and learning history annotated by arrows
# to show how learning proceeded
def visualize_learning(w_history):  
    
    # Make the function plot
    function_plot(pts,f_vals)
    
    # Plot the history
    plt.plot(w_history[:,0],w_history[:,1],marker="o",c="magenta") 
    
    # Annotate the point found at last iteration
    annotate_pt('minimum found',
                (w_history[-1,0],w_history[-1,1]),
                (-1,7),'green')
    iter = w_history.shape[0]
    for w,i in zip(w_history,range(iter-1)):
        # Annotate with arrows to show history
        plt.annotate("",
                    xy=w, xycoords="data",
                    xytext=w_history[i+1,:], textcoords="data",
                    arrowprops=dict(arrowstyle="<-",
                            connectionstyle="angle3"))     
    
def solve_fw():
    # Setting up
    rand = np.random.RandomState(19)
    w_init = rand.uniform(-10,10,2)
    fig, ax = plt.subplots(nrows=4, ncols=4, figsize=(18, 12))
    learning_rates = [0.05,0.2,0.5,0.8]
    momentum = [0,0.5,0.9]
    ind = 1
    
    # Iteration through all possible parameter combinations
    for alpha in momentum:
        for eta,col in zip(learning_rates,[0,1,2,3]):
            plt.subplot(3,4,ind)        
            w_history,f_history = gradient_descent(5,-1,w_init, f,grad,[],eta,alpha)
            
            visualize_learning(w_history)
            ind = ind+1
            plt.text(-9, 12,'Learning Rate="+str(eta),fontsize=13)
            if col==1:
                plt.text(10,15,"momentum = ' + str(alpha),fontsize=20)

    fig.subplots_adjust(hspace=0.5, wspace=.3)
    plt.show()

Corramos solve_fw() y vea cómo la tasa de aprendizaje y el impulso afectan el descenso del gradiente:

solve_fw()

Este ejemplo aclara el papel tanto del impulso como de la tasa de aprendizaje.

En la primera gráfica, con impulso cero y tasa de aprendizaje establecida en 0.05, el aprendizaje es lento y el algoritmo no alcanza el mínimo global. Aumentar el impulso acelera el aprendizaje, como podemos ver en los gráficos de la primera columna. El otro extremo es la última columna, donde la tasa de aprendizaje se mantiene alta. Esto provoca oscilaciones, que pueden controlarse hasta cierto punto añadiendo impulso.

La pauta general para el descenso de gradientes es utilizar valores pequeños de tasa de aprendizaje y valores más altos de impulso.

Descenso de gradiente para minimizar el error cuadrático medio

El descenso de gradientes es una técnica sencilla y agradable para minimizar el error cuadrático medio en un problema de regresión o clasificación supervisada.

Supongamos que se nos dan (m ) ejemplos de entrenamiento ([x_{ij}]) con (i = 1 ldots m ), donde cada ejemplo tiene (n ) características, es decir, (j = 1 ldots n ). Si los valores objetivo y de salida correspondientes para cada ejemplo son (t_i ) y (o_i ) respectivamente, entonces la función de error cuadrático medio (E ) (en este caso nuestra función de objeto) se define como:

$$
E = frac {1} {m} Sigma_ {i = 1} ^ m (t_i – o_i) ^ 2
$$

Donde la salida (o_i ) está determinada por una combinación lineal ponderada de entradas, dada por:

$$
o_i = w_0 + w_1 x_ {i1} + w_2 x_ {i2} + ldots + w_n x_ {in}
$$

El parámetro desconocido en la ecuación anterior es el vector de peso ( textbf w = [w_0,w_1,ldots,w_n]^ T ).

La función objetivo en este caso es el error cuadrático medio con un gradiente dado por:

$$
nabla _ { textbf w} E ( textbf w) = – Sigma_ {i = 1} ^ {m} (t_i – o_i) textbf {x} _i
$$

Donde (x_ {i} ) es el i-ésimo ejemplo. o una serie de características de tamaño n.

Todo lo que necesitamos ahora es una función para calcular el gradiente y una función para calcular el error cuadrático medio.

los gradient_descent() La función se puede utilizar tal cual. Tenga en cuenta que todos los ejemplos de entrenamiento se procesan juntos al calcular el gradiente. Por lo tanto, esta versión de descenso de gradiente para actualizar pesos se conoce como actualización por lotes o aprendizaje por lotes:

# Input argument is weight and a tuple (train_data, target)
def grad_mse(w,xy):
    (x,y) = xy
    (rows,cols) = x.shape
    
    # Compute the output
    o = np.sum(x*w,axis=1)
    diff = y-o
    diff = diff.reshape((rows,1))    
    diff = np.tile(diff, (1, cols))
    grad = diff*x
    grad = -np.sum(grad,axis=0)
    return grad

# Input argument is weight and a tuple (train_data, target)
def mse(w,xy):
    (x,y) = xy
    
    # Compute output
    # keep in mind that wer're using mse and not mse/m
    # because it would be relevant to the end result
    o = np.sum(x*w,axis=1)
    mse = np.sum((y-o)*(y-o))
    mse = mse/2
    return mse    

Ejecución de Gradient Descent en OCR

Para ilustrar el descenso de gradiente en un problema de clasificación, hemos elegido los conjuntos de datos de dígitos incluidos en sklearn.datasets.

Para simplificar las cosas, hagamos una prueba de descenso de gradiente en un problema de dos clases (dígito 0 frente a dígito 1). El siguiente código carga los dígitos y muestra los primeros 10 dígitos. Esto nos da una idea de la naturaleza de los puntos de entrenamiento:

# Load the digits dataset with two classes
digits,target = dt.load_digits(n_class=2,return_X_y=True)
fig,ax = plt.subplots(nrows=1, ncols=10,figsize=(12,4),subplot_kw=dict(xticks=[], yticks=[]))

# Plot some images of digits
for i in np.arange(10):
    ax[i].imshow(digits[i,:].reshape(8,8),cmap=plt.cm.gray)   
plt.show()

También necesitamos el método train_test_split desde sklearn.model_selection para dividir los datos de entrenamiento en un tren y un conjunto de prueba. El siguiente código ejecuta el descenso de gradiente en el conjunto de entrenamiento, aprende los pesos y traza el error cuadrático medio en diferentes iteraciones.

Al ejecutar el descenso de gradiente, mantendremos la velocidad de aprendizaje y el impulso muy pequeños, ya que las entradas no están normalizadas ni estandarizadas. Además, la versión por lotes del descenso de gradientes requiere una tasa de aprendizaje menor:

# Split into train and test set
x_train, x_test, y_train, y_test = train_test_split(
                        digits, target, test_size=0.2, random_state=10)

# Add a column of ones to account for bias in train and test
x_train = np.hstack((np.ones((y_train.size,1)),x_train))
x_test  = np.hstack((np.ones((y_test.size,1)),x_test))

# Initialize the weights and call gradient descent
rand = np.random.RandomState(19)
w_init = rand.uniform(-1,1,x_train.shape[1])*.000001
w_history,mse_history = gradient_descent(100,0.1,w_init,
                              mse,grad_mse,(x_train,y_train),
                             learning_rate=1e-6,momentum=0.7)

# Plot the MSE
plt.plot(np.arange(mse_history.size),mse_history)
plt.xlabel('Iteration No.')
plt.ylabel('Mean Square Error')
plt.title('Gradient Descent on Digits Data (Batch Version)')
plt.show()

¡Esto luce genial! Comprobemos la tasa de error de nuestro OCR en los datos de entrenamiento y prueba. A continuación se muestra una pequeña función para calcular la tasa de error de clasificación, que se llama en el conjunto de entrenamiento y prueba:

# Returns error rate of classifier
# total miclassifications/total*100
def error(w,xy):
    (x,y) = xy
    o = np.sum(x*w,axis=1)
    
    #map the output values to 0/1 class labels
    ind_1 = np.where(o>0.5)
    ind_0 = np.where(o<=0.5)
    o[ind_1] = 1
    o[ind_0] = 0
    return np.sum((o-y)*(o-y))/y.size*100
    
train_error = error(w_history[-1],(x_train,y_train))
test_error = error(w_history[-1],(x_test,y_test))

print("Train Error Rate: " + "{:.2f}".format(train_error))
print("Test Error Rate: " + "{:.2f}".format(test_error))

Train Error Rate: 0.69
Test Error Rate: 1.39

Descenso de gradiente estocástico en Python

En la sección anterior, usamos el esquema de actualización por lotes para el descenso de gradientes.

Otra versión del descenso de gradientes es el en línea o estocástico esquema de actualización, donde cada ejemplo de entrenamiento se toma uno a la vez para actualizar los pesos.

Una vez que se completan todos los ejemplos de entrenamiento, decimos que se completa una época. Los ejemplos de entrenamiento se barajan antes de cada época para obtener mejores resultados.

El siguiente fragmento de código es una pequeña modificación del gradient_descent() función para incorporar su contraparte estocástica. Esta función toma el (conjunto de entrenamiento, objetivo) como parámetro en lugar del parámetro adicional. El término ‘iteraciones’ ha sido renombrado a ‘épocas’:

# (xy) is the (training_set,target) pair
def stochastic_gradient_descent(max_epochs,threshold,w_init,
                                obj_func,grad_func,xy,
                                learning_rate=0.05,momentum=0.8):
    (x_train,y_train) = xy
    w = w_init
    w_history = w
    f_history = obj_func(w,xy)
    delta_w = np.zeros(w.shape)
    i = 0
    diff = 1.0e10
    rows = x_train.shape[0]
    
    # Run epochs
    while  i<max_epochs and diff>threshold:
        # Shuffle rows using a fixed seed to reproduce the results
        np.random.seed(i)
        p = np.random.permutation(rows)
        
        # Run for each instance/example in training set
        for x,y in zip(x_train[p,:],y_train[p]):
            delta_w = -learning_rate*grad_func(w,(np.array([x]),y)) + momentum*delta_w
            w = w+delta_w
            
        i+=1
        w_history = np.vstack((w_history,w))
        f_history = np.vstack((f_history,obj_func(w,xy)))
        diff = np.absolute(f_history[-1]-f_history[-2])
        
    return w_history,f_history

Ejecutemos el código para ver cómo son los resultados para la versión estocástica del descenso de gradiente:

rand = np.random.RandomState(19)
w_init = rand.uniform(-1,1,x_train.shape[1])*.000001
w_history_stoch,mse_history_stoch = stochastic_gradient_descent(
                                100,0.1,w_init,
                              mse,grad_mse,(x_train,y_train),
                             learning_rate=1e-6,momentum=0.7)

# Plot the MSE
plt.plot(np.arange(mse_history_stoch.size),mse_history_stoch)
plt.xlabel('Iteration No.')
plt.ylabel('Mean Square Error')
plt.title('Gradient Descent on Digits Data (Stochastic Version)')

plt.show()

Revisemos también la tasa de error:

train_error_stochastic = error(w_history_stoch[-1],(x_train,y_train))
test_error_stochastic = error(w_history_stoch[-1],(x_test,y_test))

print("Train Error rate with Stochastic Gradient Descent: " + 
      "{:.2f}".format(train_error_stochastic))
print("Test Error rate with Stochastic Gradient Descent: "  
      + "{:.2f}".format(test_error_stochastic))

Train Error rate with Stochastic Gradient Descent: 0.35
Test Error rate with Stochastic Gradient Descent: 1.39

Comparación de versiones por lotes y estocásticas

Comparemos ahora las versiones por lotes y estocásticas del descenso de gradiente.

Fijaremos la tasa de aprendizaje para ambas versiones en el mismo valor y variaremos el impulso para ver qué tan rápido convergen ambas. Los pesos iniciales y los criterios de detención para ambos algoritmos siguen siendo los mismos:

fig, ax = plt.subplots(nrows=3, ncols=1, figsize=(10,3))

rand = np.random.RandomState(11)
w_init = rand.uniform(-1,1,x_train.shape[1])*.000001
eta = 1e-6
for alpha,ind in zip([0,0.5,0.9],[1,2,3]):
	w_history,mse_history = gradient_descent(
                                100,0.01,w_init,
                              mse,grad_mse,(x_train,y_train),
                             learning_rate=eta,momentum=alpha)

    w_history_stoch,mse_history_stoch = stochastic_gradient_descent(
                                100,0.01,w_init,
                              mse,grad_mse,(x_train,y_train),
                             learning_rate=eta,momentum=alpha)
    
    # Plot the MSE
    plt.subplot(130+ind)
    plt.plot(np.arange(mse_history.size),mse_history,color="green")
    plt.plot(np.arange(mse_history_stoch.size),mse_history_stoch,color="blue")
    plt.legend(['batch','stochastic'])
    
    # Display total iterations
    plt.text(3,-30,'Batch: Iterations="+
             str(mse_history.size) )
    plt.text(3,-45,"Stochastic: Iterations="+
             str(mse_history_stoch.size))
    plt.title("Momentum = ' + str(alpha))   
    
    # Display the error rates
    train_error = error(w_history[-1],(x_train,y_train))
    test_error = error(w_history[-1],(x_test,y_test))
    
    train_error_stochastic = error(w_history_stoch[-1],(x_train,y_train))
    test_error_stochastic = error(w_history_stoch[-1],(x_test,y_test))
    
    print ('Momentum = '+str(alpha))
    
    print ('tBatch:')
    print ('ttTrain error: ' + "{:.2f}".format(train_error) )
    print ('ttTest error: ' + "{:.2f}".format(test_error) )
    
    print ('tStochastic:')
    print ('ttTrain error: ' + "{:.2f}".format(train_error_stochastic) )
    print ('ttTest error: ' + "{:.2f}".format(test_error_stochastic) )
    
        
plt.show()

Momentum = 0
	Batch:
		Train error: 0.35
		Test error: 1.39
	Stochastic:
		Train error: 0.35
		Test error: 1.39
Momentum = 0.5
	Batch:
		Train error: 0.00
		Test error: 1.39
	Stochastic:
		Train error: 0.35
		Test error: 1.39
Momentum = 0.9
	Batch:
		Train error: 0.00
		Test error: 1.39
	Stochastic:
		Train error: 0.00
		Test error: 1.39

Si bien no hay una diferencia significativa en la precisión entre las dos versiones del clasificador, la versión estocástica es un claro ganador cuando se trata de la velocidad de convergencia. Se necesitan menos iteraciones para lograr el mismo resultado que su contraparte por lotes.

Conclusiones

El descenso de gradientes es una técnica simple y fácil de implementar.

En este tutorial, ilustramos el descenso de gradientes en función de dos variables con contornos circulares. Luego ampliamos nuestro ejemplo para minimizar el error cuadrático medio en un problema de clasificación y construimos un sistema OCR simple. También discutimos la versión estocástica del descenso de gradiente.

En este tutorial se desarrolló una función de propósito general para implementar el descenso de gradientes. Animamos a los lectores a utilizar esta función en diferentes problemas de regresión y clasificación, con distintos hiperparámetros, para una mejor comprensión de su funcionamiento.