Régression linéaire en Python : 3 méthodes

 La régression linéaire est un des piliers du machine learning.

Je vais présenter 3 méthodes pour trouver l'équation de la droite qui résume au plus près un nuage de point..

Le chargement des données. 

Les données sont issues d'un jeu de donnée des valeurs foncières en opendata (disponible sur github) . 

Les données seront chargées dans un dataframe de Pandas.

import pandas as pd
data = pd.read_csv('parcelles_ext.csv', sep=';')
import numpy as np

data.replace('None',np.nan , inplace = True )
data2 = data[['valeur_fonciere','surface_reelle_bati','adresse_code_voie', 'nombre_pieces_principales']]
data2 = data2.dropna()
y= data2['valeur_fonciere'].to_numpy()
X=data2['surface_reelle_bati'].to_numpy()
X = np.reshape(X,(-1,1))
y = np.reshape(y,(-1,1))


Seules deux colonnes seront utilisées.

Méthode par un batch de descente en gradient.

Deux fonctions sont nécessaires: une qui calcule le 'coût' (l'erreur)  et une autre qui va répéter autant de fois que necessaire le calcul de cout puis le gradient et va se rapprocher petit à petit du minimum de la fonction de coût.

def  cost_function(theta,X,y):   
    m = len(y)
    predictions = X.dot(theta)
    cost = (1/2*m) * np.sum(np.square(predictions-y))
    return cost

Le batch

def gradient_descent(X,y,theta, alpha=0.01,iterations=100):
    m = len(y)
    cost_history = np.zeros(iterations)
    theta_history = np.zeros((iterations,2))
    for it in range(iterations):   
        prediction = np.dot(X,theta)
        theta = theta -(1/m)*alpha*( X.T.dot(prediction - y))
        theta_history[it,:] =theta.T
        cost_history[it]  = cost_function(theta,X,y)     
    return theta, cost_history, theta_history

Pour maintenant activer le batch , on a besoin d'ajouter une colonne  'biais' aux données en entrée.

Et de choisir un taux d'apprentissage et un nombre d'itération.

lr =0.0001
n_iter = 100000*8
theta = np.random.randn(2,1)
initial_theta = theta
X_b = np.c_[np.ones((len(X),1)),X]

Le 'alpha' (learn rate) est ici très bas et le nombre d'itération très élevés.

La raison principale est l'ordre de grandeur des données: 

 [  1. 113.]
 [  1.  91.]
 [  1.  61.]
 [  1. 150.]]

113 , 150 etc , les données ne sont pas à l'echelle (voir la suite) mais on arrive à un résultat.

Méthode par une fonction de recherche de minimum.

au lieu de devoir faire un batch à la main, il est possible d'utiliser une fonction intégrée qui va se charger de trouver le minimum d'une fonction en optimisant les itérations par l'emploi de transformation. 

On va individualiser la fonction gradient:

def gradient(theta,X,y):
    #print(X.shape)
    #print(theta.shape)
    prediction = np.dot(X,theta)
    t_gradient = 1/len(y)*(X.T.dot(prediction-y))
    return t_gradient

Et appeler la fonction de recherche de minimum (fmin_tnc de scipy) en lui passant le nom des fonctions de cout et de gradient.

Les deux méthodes arrivent à des résultats comparables.

Méthode par l'équation  normale.

Il est possible de trouver directement l'équation de la droite normale mais avec un niveau de complexité qui explose avec le nombre de donnée.

 

A retenir : toujours normaliser ou mettre à l'echelle les données. 

On peut utiliser des méthodes intégrées:

ou le faire à la main:

L'avantage est d'avoir en retour les valeurs de mu et sigma pour l'appliquer aux prévisions.

Le learn rate et le  nombre d'itération  seront moins extrêmes.

Le code est sous forme de notebook pour jupyter à retrouver sur mon github.