discesa del gradiente usando pitone e NumPy

Question

def gradient(X_norm,y,theta,alpha,m,n,num_it): 
    temp=np.array(np.zeros_like(theta,float)) 
    for i in range(0,num_it): 
     h=np.dot(X_norm,theta) 
     #temp[j]=theta[j]-(alpha/m)*( np.sum((h-y)*X_norm[:,j][np.newaxis,:]) ) 
     temp[0]=theta[0]-(alpha/m)*(np.sum(h-y)) 
     temp[1]=theta[1]-(alpha/m)*(np.sum((h-y)*X_norm[:,1])) 
     theta=temp 
    return theta 



X_norm,mean,std=featureScale(X) 
#length of X (number of rows) 
m=len(X) 
X_norm=np.array([np.ones(m),X_norm]) 
n,m=np.shape(X_norm) 
num_it=1500 
alpha=0.01 
theta=np.zeros(n,float)[:,np.newaxis] 
X_norm=X_norm.transpose() 
theta=gradient(X_norm,y,theta,alpha,m,n,num_it) 
print theta 

mio theta dal codice precedente è 100.2 100.2, ma dovrebbe essere 100.2 61.09 in MATLAB che è corretto.

Answer 1

Penso che il tuo codice sia un po 'troppo complicato e ha bisogno di più struttura, perché altrimenti ti perderai in tutte le equazioni e operazioni. Alla fine questa regressione riduce a quattro operazioni:

1. calcolare l'ipotesi h = X * theta
2. Calcolare la perdita = h - y e forse il costo quadrato (perdita^2)/2m
3. calcolare il gradiente = X' * perdita/m
4. aggiornare i parametri di teta = theta - alfa * gradient

Nel tuo caso, credo che hai confuso con mn. Qui m indica il numero di esempi nel set di allenamento, non il numero di funzioni.

Diamo uno sguardo alla mia variante del codice:

import numpy as np 
import random 

# m denotes the number of examples here, not the number of features 
def gradientDescent(x, y, theta, alpha, m, numIterations): 
    xTrans = x.transpose() 
    for i in range(0, numIterations): 
     hypothesis = np.dot(x, theta) 
     loss = hypothesis - y 
     # avg cost per example (the 2 in 2*m doesn't really matter here. 
     # But to be consistent with the gradient, I include it) 
     cost = np.sum(loss ** 2)/(2 * m) 
     print("Iteration %d | Cost: %f" % (i, cost)) 
     # avg gradient per example 
     gradient = np.dot(xTrans, loss)/m 
     # update 
     theta = theta - alpha * gradient 
    return theta 


def genData(numPoints, bias, variance): 
    x = np.zeros(shape=(numPoints, 2)) 
    y = np.zeros(shape=numPoints) 
    # basically a straight line 
    for i in range(0, numPoints): 
     # bias feature 
     x[i][0] = 1 
     x[i][1] = i 
     # our target variable 
     y[i] = (i + bias) + random.uniform(0, 1) * variance 
    return x, y 

# gen 100 points with a bias of 25 and 10 variance as a bit of noise 
x, y = genData(100, 25, 10) 
m, n = np.shape(x) 
numIterations= 100000 
alpha = 0.0005 
theta = np.ones(n) 
theta = gradientDescent(x, y, theta, alpha, m, numIterations) 
print(theta) 

In un primo momento ho creato un piccolo insieme di dati casuali che dovrebbe assomigliare a questo:

Come potete vedere ho ha anche aggiunto la linea di regressione generata e la formula che è stata calcolata da Excel.

È necessario prestare attenzione all'intuizione della regressione utilizzando la discesa del gradiente. Mentre esegui un passaggio completo sui tuoi dati X, devi ridurre le perdite m di ogni esempio in un singolo aggiornamento del peso. In questo caso, questa è la media della somma sopra i gradienti, quindi la divisione di m.

La prossima cosa che devi fare attenzione è tracciare la convergenza e regolare il tasso di apprendimento. Del resto dovresti sempre tenere traccia del costo di ogni iterazione, magari anche tracciarlo.

Se si esegue il mio esempio, il theta restituito sarà simile a questa:

Iteration 99997 | Cost: 47883.706462 
Iteration 99998 | Cost: 47883.706462 
Iteration 99999 | Cost: 47883.706462 
[ 29.25567368 1.01108458] 

che in realtà è abbastanza vicino alla equazione che è stato calcolato da Excel (y = x + 30). Si noti che quando abbiamo passato il bias nella prima colonna, il primo valore theta indica il peso di bias.

Answer 2

Di seguito è possibile trovare la mia implementazione della discesa del gradiente per il problema di regressione lineare.

In un primo momento, si calcola il gradiente come X.T * (X * w - y)/N e si aggiorna contemporaneamente il theta corrente con questo gradiente.

• X: Funzione matrice
• y: bersaglio Valori
• w: pesi/valori
• N: dimensione della formazione impostata

Ecco il codice pitone:

import pandas as pd 
import numpy as np 
from matplotlib import pyplot as plt 
import random 

def generateSample(N, variance=100): 
    X = np.matrix(range(N)).T + 1 
    Y = np.matrix([random.random() * variance + i * 10 + 900 for i in range(len(X))]).T 
    return X, Y 

def fitModel_gradient(x, y): 
    N = len(x) 
    w = np.zeros((x.shape[1], 1)) 
    eta = 0.0001 

    maxIteration = 100000 
    for i in range(maxIteration): 
     error = x * w - y 
     gradient = x.T * error/N 
     w = w - eta * gradient 
    return w 

def plotModel(x, y, w): 
    plt.plot(x[:,1], y, "x") 
    plt.plot(x[:,1], x * w, "r-") 
    plt.show() 

def test(N, variance, modelFunction): 
    X, Y = generateSample(N, variance) 
    X = np.hstack([np.matrix(np.ones(len(X))).T, X]) 
    w = modelFunction(X, Y) 
    plotModel(X, Y, w) 


test(50, 600, fitModel_gradient) 
test(50, 1000, fitModel_gradient) 
test(100, 200, fitModel_gradient) 

Answer 3

So che questa domanda già stato risposta, ma ho fatto qualche aggiornamento alla funzione GD:

### COST FUNCTION 

def cost(theta,X,y): 
    ### Evaluate half MSE (Mean square error) 
    error = np.dot(X,theta) - y 
    J = np.sum(error ** 2)/(2*m) 
    return J 

cost(theta,X,y) 



def GD(X,y,theta,alpha): 

    cost_histo = [0] 
    theta_histo = [0] 

    # an arbitrary gradient, to pass the initial while() check 
    delta = [np.repeat(1,len(X))] 
    # Initial theta 
    old_cost = cost(theta,X,y) 

    while (np.max(np.abs(delta)) > 1e-6): 
     error = np.dot(X,theta) - y 
     delta = np.dot(np.transpose(X),error)/len(y) 
     trial_theta = theta - alpha * delta 
     trial_cost = cost(trial_theta,X,y) 
     while (trial_cost >= old_cost): 
      trial_theta = (theta +trial_theta)/2 
      trial_cost = cost(trial_theta,X,y) 
      cost_histo = cost_histo + trial_cost 
      theta_histo = theta_histo + trial_theta 
     old_cost = trial_cost 
     theta = trial_theta 
    Intercept = theta[0] 
    Slope = theta[1] 
    return [Intercept,Slope] 

res = GD(X,y,theta,alpha) 

Questa funzione riduce l'alfa sopra l'iterazione far funzionare troppo convergere velocemente vedere Estimating linear regression with Gradient Descent (Steepest Descent) per un esempio in R. Applico la stessa logica ma in Python.