Regressão Linear com uma variável

Regress˜ao Linear com uma vari´avel Manuel Leiria July 31, 2016

Abstract A regress˜ ao linear prevˆe um valor real baseado num valor de input. Vamos usar como exemplo o pre¸co de casas para fazer previs˜oes de pre¸cos de casas que n˜ao estejam nos dados. Ser´a introduzida a no¸c˜ ao de fun¸c˜ ao de custo (cost function) e o m´etodo do gradiente descendente (gradient descent) para a aprendizagem. Este algoritmo de regress˜ao linear faz parte do conjunto de algoritmos que se chamam aprendizagem supervisionada (supervised learning). Aprendizagem supervisionada significa que os dados que tenho para treinar o algoritmo apresentam as respostas correctas.

1

Introdu¸c˜ ao

Vamos come¸car com um exemplo pr´atico e a partir da´ı mostrar como a regress˜ao linear pode ajudar-nos a resolver o problema. Suponhamos que eu tenho um apartamento para vender em determinada cidade. Conhe¸co os pre¸c de alguns apartamentos parecidos na zona onde habito e quero ter uma ideia de quanto ´e que o meu vale. Este ´e um caso de aprendizagem supervisionada porque eu tenho as ”respostas correctas”, ou seja, tenho o tamanho do apartamento e tenho o respectivo pre¸co.

2

Dados Dispon´ıveis

Os dados tenho dispon´ıveis que tenho s˜ao:

1

Table 1: Pre¸co dos apartamentos em Lisboa. Tamanho (m2 ) Pre¸co (euros) 80

250.000

83

260.000

70

240.000

100

300.000

50

150.000

35

100.000

85

270.000

82

250.000

70

230.000

Visualiza¸c˜ ao dos dados: NOTA: em algumas tarefas,que se querem r´apidas usarei c´odigo Python (como por exemplo mostrar um gr´ afico) mas este artigo e o desenvolvimento do algoritmo para regress˜ ao linear ser´a em Java. C´odigo Python import matplotlib.pyplot as plt import pylab

x = [80, 83, 70, 100, 50, 35, 85, 82, 70] y = [250000, 260000, 240000, 300000, 150000, 100000, 270000, 250000, 230000]

plt.scatter(x,y) plt.xlabel(’Tamanho em m^2’) plt.ylabel(’Preco em euros’) plt.show() A ideia da regress˜ ao linear ´e tra¸car um linha recta entre os pontos de forma a que consigamos calcular o pre¸co de um apartamento que n˜ao perten¸ca aos pontos (por exemplo, um

2

Figure 1: Pre¸co do apartamento em fun¸ca˜o do tamanho apartamento com 65 m2 ). Nota¸ c˜ ao: m: N´ umero de exemplos de treino xs: Vari´avel de input (features) ys: vari´avel de output (x, y): Um exemplo de treino (uma linha da tabela) (x(i) , y (i) ): Exemplo de treino da linha i da tabela Por exemplo, em rela¸c˜ ao aos dados, (x(3) , y (3) ) = (70, 240000)

3

H´ıp´ otese

N´os queremos encontrar uma fun¸c˜ ao h que receba como input o tamanho da casa (x) e retorne o seu pre¸co (y). Esquematicamente :

3

Dados de Treino

Algoritmo de Aprendizagem

Tamanho do apartamento (x)

h

Preco Estimado (y)

Portanto, n´ os queremos uma fun¸c˜ ao h que mapei valores de x em y. Representamos h como:

hθ (x) = θ0 + θ1 x

(1)

Esta ´e uma fun¸c˜ ao linear. Graficamente e para alguns valores de θ0 e θ1 , temos: C´odigo Python import matplotlib.pyplot as plt import pylab import numpy as np

x = [80, 83, 70, 100, 50, 35, 85, 82, 70] y = [250000, 260000, 240000, 300000, 150000, 100000, 270000, 250000, 230000]

x1 = np.arange(0, 150, 1)

plt.scatter(x,y) plt.plot(3000*x1) plt.plot(3300*x1) plt.plot(3500*x1) plt.plot(3700*x1) plt.xlabel(’Tamanho em m^2’) plt.ylabel(’Preco em euros’) plt.show()

4

A quest˜ao que se coloca ´e saber quais os valores de θ que fazem a linha mais aproximada dos valores reais y. A isto chama-se regress˜ ao linear a uma vari´avel. Podemos aplicar a regress˜ao linear a uma vari´avel quando os dados para os quais queremos fazer previs˜oes podem ser aproximados por uma linha recta. No entanto, e veremos mais ´a frente que facilmente podemos generalizar esta teoria para polin´ omios de ordens mais elevadas.

4

Fun¸c˜ ao de custo

Vamos definir a fun¸c˜ ao de custo, que nos vai ajudar a construir a linha que mais se aproxima dos nossos dados reais. A nossa hip´ otese ´e hθ (x) = θ0 + θ1 x que nos vai ajudar a fazer previs˜oes. Chamamos aos θ, parˆametros e o grande problema passa a ser como escolher os parˆametros. Tal como referido anteriormente, queremos escolher os parˆametros θ0 e θ1 de forma que hθ (x) seja o mais pr´oximo poss´ıvel de y para os nossos dados de treino (x, y). Formalizando, n´ os queremos minimizar a m´edia da diferen¸ca do quadrado do erro em θ0 e θ1 para todos os dados de teste: m

1 X (hθ (x(i) − y (i) )2 minimizar 2m

(2)

i=1

. Mais formalmente, quero minimizar a fun¸c˜ao de custo (ou fun¸c˜ao do quadrado do erro) J(θ0 , θ1 ):

m

1 X J(θ0 , θ1 ) = (hθ (x(i) − y (i) )2 2m i=1

5

(3)

. Resumindo: Hip´otese: hθ (x) = θ0 + θ1 x Parˆametros: θ0 , θ1 Fun¸c˜ao de custo: J(θ0 , θ1 ) =

1 2m

Pm

i=1 (hθ (x

(i)

− y (i) )2

Objectivo: minimizar em θ0 , θ1 a fun¸c˜ao J(θ0 , θ1 ) Note-se que para J(θ0 , θ1 ) fixos, hθ (x) ´e uma fun¸c˜ao de x.

5

Gradiente Descendente

Resumindo mas tamb´em generalizando um pouco: Temos uma fun¸c˜ ao J(θ0 , θ1 , θ2 , ..., θn ) e queremos minimizar em fun¸c˜ao dos parˆametros θ. Come¸camos com uns valores iniciais θ0 , θ1 , θ2 , ..., θn . Vamos mudando esses valores para reduzir J(θ0 , θ1 , θ2 , ..., θn ) at´e encontrarmos o m´ınimo global. H´ a algumas considera¸c˜ oes matem´aticas que se devem ter em conta, como por exemplo o problema de alcan¸car um m´ınimo local em vez do m´ınimo global mas tais quest˜ oes est˜ao fora do ˆ ambito deste artigo. Em termos matem´ aticos, encontrar o m´ınimo de uma fun¸c˜ao corresponde a igualar a zero a derivada dessa fun¸c˜ ao.

5.0.1

Algoritmo do Gradiente Descendente

O que temos de fazer ´e: repetir 0 ,θ1 ) para j = 0 e j = 1 θj := θj − α ∂J(θ ∂θj

at´e convergir. Aten¸c˜ao que os updates em j tˆem de ser simultˆaneos. Na equa¸c˜a supra, o α ´e um n´ umero que se chama taxa de aprendizagem (learning rate) e basicamente o que faz ´e controlar o tamanho do passo na descida do gradiente. Se o valor de α for muito pequeno, a descida do gradiente pode ser bastante lenta. Se α for muito grande, a descida do gradiente ´e r´apida mas pode n˜ ao convergir ou at´e mesmo divergir (overshoot the minimum.

6

5.0.2

Algoritmo do Gradiente Descendente para Regress˜ ao Linear

Resumindo, temos: Algoritmo do Gradiente Descendente Repetir at´e convergir: 0 ,θ1 ) θj := θj − α ∂J(θ ∂θj

para j = 0 e j = 1 Modelo de Regress˜ ao Linear hθ (x) = θ0 + θ1 x (hip´ otese linear) P m 1 (i) − y (i) )2 (fun¸ J(θ0 , θ1 ) = 2m c˜ao de custo do quadrado do erro. i=1 (hθ (x O que vamos fazer ´e aplicar o algoritmo do gradiente descendente ´a fun¸c˜ao de custo. Saltando o c´alculo das derivadas parciais (quem souber fazer, tudo bem, quem n˜ao souber, tudo bem na mesma porque n˜ ao ´e relevante para o caso), obtemos: m

m

i=1

i=1

∂J(θ0 , θ1 ) ∂ 1 X ∂ 1 X = (hθ (x(i) ) − y (i) )2 = (θ0 + θ1 x(i) − y (i) )2 ∂θj ∂θj 2m ∂θj 2m

(4)

Calculando as derivadas parciais, obtemos: 1 Pm 0 ,θ1 ) (i) (i) =m j = 0 : ∂J(θ i=1 (hθ (x ) − y ) ∂θ0 1 Pm 0 ,θ1 ) (i) (i) (i) j = 1 : ∂J(θ =m i=1 (hθ (x ) − y )x ∂θ1 Portanto, voltando ao algoritmo do gradiente descendente, repetir 1 Pm (i) (i) θ0 := θ0 − α m i=1 (hθ (x ) − y ) 1 Pm (i) (i) (i) θ1 := θ1 − α m i=1 (hθ (x ) − y )x at´e convergir. A fu¸c˜ao de custo para a regress˜ ao linear ´e uma fu¸c˜ao convexa o que significa que n˜ao tem m´ınimos locais. H´ a s´ o um m´ınimo global. De notar que este algoritmo em cada passo (step) utiliza todos os dados de treino. Por u ´ltimo e antes de avan¸car para a implementa¸c˜ao em Java, resta referir que, para resolver este tipo de problemas em que a regress˜ao linear se aplica, h´a um outro m´etodo chamado Equa¸c˜ao normalmas n˜ ao escala bem para conjuntos de dados muito grandes. Outro ponto que vale a pena mencionar ´e que esta implementa¸c˜ao n˜ao est´a restringida a uma vari´avel (feature) e pode ser generalizada (por exemplo se eu tiver dados de pre¸cos dos apartamentos em fun¸c˜ ao n˜ ao s´ o do tamanho em metros quadrados mas tamb´em do n´ umero de assoalhadas. Por exemplo:

7

Table 2: Pre¸co dos apartamentos em Lisboa. Tamanho (m2 ) Num assoalhadas

Pre¸co (euros)

80

3

250.000

83

4

260.000

70

3

240.000

100

5

300.000

50

2

150.000

35

1

100.000

85

3

270.000

82

3

250.000

70

3

230.000

Neste caso a nossa fun¸c˜ ao de hip´otese fica: hθ (x) = θ0 + θ1 x1 + θ2 x2 , onde x1 corresponde ao tamanho, x2 ´e o n´ umero de assoalhadas e y, tal como anteriormente, ´e o pre¸co.Nesta nota¸c˜ ao, por exemplo x32 = 3. Corresponde ao valor da segunda coluna, terceira linha. Na implementa¸c˜ ao, e isto aplica-se para qualquer linguagem de programa¸c{ao, podemos reescrever a hip´ otese como: hθ (x) = θ0 x0 + θ1 x1 + θ2 x2 , onde x0 = 1, sempre. Isto permite-nos trabalhar em formato matricial. Primeiro come¸co por apresentar a class Matrix.java, que representa, tal como o nome indica, uma matriz. Penso que o c´ odigo est´ a simples e ´e auto explicativo. A classe ´e um pouco grande porque tˆem v´ arios m´etodos gerais para trabalhar com matrizes:

package pt.mleiria.machinelearning.matrixAlgebra;

/** * * @author manuel

8

*/ public class Matrix {

protected double[][] components;

/** * * @param a */ public Matrix(final double[][] a) { components = a; }

/** * Creates a zero matrix of given dimensions * * @param n number of rows * @param m number of columns */ public Matrix(final int n, final int m) { if (n columns()) { throw new IllegalArgumentException("Column Index provided is greater than matrix number of columns"); } double[] colArray = new double[n]; for (int i = 0; i < n; i++) { colArray[i] = component(i, colIndex); } return colArray; } } Suponhamos que carrego os dados numa matriz nxm, ou seja com n colunas e m linhas. A u ´ltima  coluna representa  o y Tomando os dados da tabela 2 como exemplo, seria:  80 3 250000     83 4 260000       70 3 240000      100 5 300.000      a=  50 2 150.000    35 1 100.000      85 3 270.000        82 3 250.000   70 3 230.000 Uma vez carregada esta estrutura no objecto Matrix, vamos separar em duas matrizes: uma com os valores de x e outra (que ser´ a uma matriz s´o com uma coluna, ou vector) com os valores de y: /*

23

alpha value. you can try several values (e.g., 0.05, 0.1,...) and see how it goes */ double alpha = 0.01; /* Number of iterations */ int numIter = 1500; /* or until some precision is reached */ double precision = 0.00001; final Matrix[] splitedM = a.split(1); final Matrix featuresX = splitedM[0]; final Matrix outputY = splitedM[1]; Com estas matrizes vamos alimentar uma nova classe chamada GradientDescent.java. Esta classe extende uma classe abstracta chamada IteratorProcessor.java que usada para definir os m´etodos que todas as classes que usam processos iterativos devem implementar:

package pt.mleiria.machinelearning.regression;

import pt.mleiria.machinelearning.iterations.IteratorProcessor; import pt.mleiria.machinelearning.matrixAlgebra.Matrix;

/** * * @author manuel */ public class GradientDescent extends IteratorProcessor {

private Matrix featuresX;

24

private Matrix outputY; private Matrix theta; private double alpha; private final int dataSize; private double[] costHistory; /** * * @param featuresX * @param outputY * @param alpha */ public GradientDescent(Matrix featuresX, Matrix outputY, double alpha) { this.featuresX = featuresX; this.outputY = outputY; this.alpha = alpha; this.dataSize = outputY.rows(); } /** * */ @Override public void initializeIterations() { costHistory = new double[getMaximumIterations() + 1]; setFeaturesX(featuresX.addOnes()); setTheta(new Matrix(featuresX.columns(), 1)); computeCost(); } /** * */ private void computeCost() { double coeff = 1.0 / (2.0 * dataSize);

25

final Matrix a = featuresX.multiply(theta); final Matrix b = a.subtract(outputY); final Matrix c = b.transpose(); final Matrix d = c.multiply(b); final Matrix result = d.multiply(coeff); costHistory[iterations] = result.component(0, 0); }

/** * theta = theta - (alpha / m) * X’ * (X * theta - y); * * b = X * theta * c = b - y * d = X’ * c * e = (alpha/m) * d * * @return */ @Override public double evaluateIteration() { double coeff = (alpha / (double) dataSize); Matrix b = featuresX.multiply(theta); Matrix c = b.subtract(outputY); Matrix d = featuresX.transpose(); Matrix e = d.multiply(c); Matrix f = e.multiply(coeff); theta = theta.subtract(f); computeCost(); double precision = costHistory[iterations] - costHistory[iterations 1]; return Math.abs(precision); }

26

/** * @return the featuresX */ public Matrix getFeaturesX() { return featuresX; }

/** * @param featuresX the featuresX to set */ public void setFeaturesX(Matrix featuresX) { this.featuresX = featuresX; }

/** * @return the outputY */ public Matrix getOutputY() { return outputY; }

/** * @param outputY the outputY to set */ public void setOutputY(Matrix outputY) { this.outputY = outputY; }

/** * @return the theta */

27

public Matrix getTheta() { return theta; }

/** * @param theta the theta to set */ public void setTheta(Matrix theta) { this.theta = theta; }

/** * @return the alpha */ public double getAlpha() { return alpha; }

/** * @param alpha the alpha to set */ public void setAlpha(double alpha) { this.alpha = alpha; }

public double[] getCostHistory() { return costHistory; }

}

/*

28

* A general structure for managing iterations */ package pt.mleiria.machinelearning.iterations;

import org.apache.log4j.Logger; import pt.mleiria.LogTypes; /** * * @author manuel */ public abstract class IteratorProcessor { protected static final Logger log = Logger.getLogger(LogTypes.MLEARNING_LOG);

/** * Number of iterations performed. */ protected int iterations; /** * Maximum allowed number of iterations. */ private int maximumIterations = 5000; /** * Desired precision. */ private double desiredPrecision; /** * Achieved precision. */ private double precision;

/**

29

* Generic constructor. */ public IteratorProcessor() { }

/** * Performs the iterative process */ public void evaluate() { iterations = 0; initializeIterations(); while (iterations++ < maximumIterations) { precision = evaluateIteration(); //LOG.log(Level.INFO, "Precision:{0}", precision); if (hasConverged()) { log.info("Converged at iteration: [" + iterations + "]"); log.info("Converged with precision:[" + precision + "]"); break; } } finalizeIterations(); }

/** * Evaluate the result of the current iteration. * * @return the estimated precision of the result. */ abstract public double evaluateIteration();

/** * Clean-up operations

30

*/ public void finalizeIterations() { log.info("Iterations:[" + iterations +"]" ); }

/** * @return desired precision. */ public double getDesiredPrecision() { return desiredPrecision; }

/** * @return the number of iterations performed. */ public int getIterations() { return iterations; }

/** * @return the maximum allowed number of iterations. */ public int getMaximumIterations() { return maximumIterations; }

/** * @return the attained precision. */ public double getPrecision() { return precision; }

31

/** * Check to see if the result has been attained. * * @return boolean */ public boolean hasConverged() { return precision < desiredPrecision; }

/** * Initializes internal parameters to start the iterative process. */ public void initializeIterations() { }

/** * Defines the desired precision. * * @param prec */ public void setDesiredPrecision(double prec) throws IllegalArgumentException { if (prec