Na aprendizagem de máquina supervisionada, quando a resposta que precisa ser prevista é qualitativa (discreta), é usado o processo chamado de Classificação. São problemas de classificação, por exemplo, decidir se um e-mail recebido é spam ou não, ou decidir se uma foto contém a imagem de um gato, cachorro ou um pássaro.
Para funcionar, um classificador precisa ser treinado usando dados de treinamento. Cada exemplo nos dados de treinamento possui um vetor de características $(x_1, x_2, ..., x_n)$ e uma resposta associada $y$ que é a classificação do exemplo.
Apesar do nome, Regressão Logística é um algoritmo de classificação que pode ser usado para resolver problemas com duas ou mais possíveis classes. Entretanto, para facilitar o entendimento, consideraremos que $y$ pode ter apenas dois valores, 0 ou 1, sendo 0 uma classe e 1 a outra.
Com as respostas variando apenas entre dois valores, 0 e 1, usar regressão linear não é adequado, pois ela pode prever qualquer valor contínuo que pode ir muito além dos limites de 0 e 1. Assim, mudaremos a função $f(x)$, que representa a hipótese para explicar os dados, para uma função que varie apenas entre 0 e 1.
Na regressão logística usaremos a função logística da seguinte forma:
$$f(x) = g(w^Tx)$$$$g(z) = \frac{1}{1+e^{-z}}$$Perceba que continua-se usando $w^Tx$ como na regressão linear, mas agora aplica-se o resultado na função logística $g(z)$.
A seguir temos a implementação das duas funções acima, $f(x)$ e $g(z)$. $f(x)$ recebe o vetor de coeficientes $w$ e as características $x$, já $g(z)$ recebe apenas um valor que é aplicado a função logística.
import math
def f(x, w):
z = 0
for i in range(len(x)):
z += x[i]*w[i]
return g(z)
def g(z):
return 1 / (1 + math.exp(-z))
A função logística tem uma forma de "S" como pode ser vista no gráfico seguir:
import numpy as np
from bokeh.io import show, output_notebook
from bokeh.plotting import figure
output_notebook()
x_points = np.linspace(-10, 10, 100)
y_points = [g(i) for i in x_points]
p = figure(plot_width=700, plot_height=250)
p.xaxis.axis_label = "z"
p.yaxis.axis_label = "g(z)"
p.line(x_points, y_points, line_width=2)
show(p, notebook_handle=True)
Note que quanto menor o valor de $z$, ou seja, o valor de $w^Tx$, mais próximo de zero é o resultado da função. Por outro lado, quanto maior o valor de $z$, mais próximo de 1 é o resultado. Dessa forma, o resultado nunca ultrapassará o intervalo entre 0 e 1.
Na regressão logística, interpretamos o resultado de $f(x)$ como a probabilidade da resposta ser 1 para o vetor de características $x$. Assim, se $f(x) = 0,8$, então $x$ tem probabilidade $0,8$ de ser da classe 1, e $0,2$ de ser da classe 0.
Para obter uma classificação discreta, é necessário definir um limiar em relação à probabilidade. Por exemplo, se $f(x) \geq 0,5 \Rightarrow y = 1$. Ou seja, se a probabilidade de ser 1 for maior ou igual 50%, então classificamos como 1. Adicionalmente, se $f(x) < 0,5 \Rightarrow y = 0$.
Na função logística $g(z)$, se $z \geq 0$, então $g(z) \geq 0,5$. Como a hipótese na regressão logística é $f(x) = g(w^Tx)$, então, se $w^Tx \geq 0$, isto resultará em $g(z) \geq 0,5$.
Com isso, $w^Tx$ especifica uma fronteira de decisão que separa os pontos que são classificados como 0 ou 1. Definimos matematicamente a fronteira de decisão da seguinte forma:
$$w^Tx = 0 $$Assim,
$$w^Tx \geq 0 \Rightarrow y = 1$$$$w^Tx < 0 \Rightarrow y = 0$$Dependendo do tipo de função definida por $w^Tx$ a fronteira de decisão pode ter diferentes formas.
Por exemplo, usando a função $f(x) = w_0 + w_1x_1 + w_2x_2$, definimos a fronteira como:
$$w_0 + w_1x_1 + w_2x_2 = 0$$Essa equação define uma reta.
Abaixo temos o gráfico mostrando a fronteira usando a equação da reta anterior para $w^T = \begin{bmatrix} -3 & 0,8 & 1 \end{bmatrix}$. Ou seja, a equação é: $0,8x_1 + x_2 - 3 = 0$.
Note que no gráfico, os pontos são separados em dois grupos pela reta, sendo os de cor azul correspondentes a $y=0$ e os de cor vermelha a $y=1$.
def boundary(x):
return 3-0.8*x
x = [0.5, 1, 1.2, 1.1, 1.5, 2.6, 2.3, 3, 2.8, 2.4]
y = [0.7, 1.2, 0.5, 1.8, 1, 2.5, 2.8, 2.5, 3, 2.2]
colors = ["blue", "blue", "blue", "blue", "blue", "red", "red", "red", "red", "red"]
p = figure(plot_width=700, plot_height=250)
p.circle(x, y, radius=0.025, color=colors)
p.line([0, 3.5], [boundary(0), boundary(3.5)], line_width=2)
show(p, notebook_handle=True)
<Bokeh Notebook handle for In[3]>
Usando funções diferentes, podemos ter curvas mais complexas. No exemplo a seguir usamos a função $f(x) = -1,2 + 3x_1^2 + 3x_2^2$, criando a fronteira $3x_1^2 + 3x_2^2 - 1,2 = 0$. Esta equação define uma elipse.
import math
import itertools
def boundary(x1):
return math.sqrt((1.2-3*x1**2)/0.3)
x = [0.5, -0.1, 0.2, -0.3, -0.3, 0, 0.7, 1, 0.8, 0, -0.5, -0.8, -0.5]
y = [0.2, -0.8, 0.5, -0.1, 0.7, 3.5, 2.0, 0.5, -1.5, -3.5, -2, 0, 2]
colors = ["blue", "blue", "blue", "blue", "blue", "red", "red", "red", "red", "red", "red", "red", "red"]
curve_x_points = np.linspace(-0.8, 0.8, 100000)
curve_x_points = [i for i in curve_x_points if ((1.2-3*i**2)/0.3) >= 0]
curve_y_points = [boundary(i) for i in curve_x_points]
curve_neg_y_points = [-boundary(i) for i in curve_x_points]
p = figure(plot_width=700, plot_height=250)
p.circle(x, y, radius=0.015, color=colors)
p.line(curve_x_points, curve_y_points, line_width=2)
p.line(curve_x_points, curve_neg_y_points, line_width=2)
show(p, notebook_handle=True)
<Bokeh Notebook handle for In[4]>
Dois grupos foram criados pela fronteira. Os pontos que estão fora da elipse possuem classificação $y=1$ e os que estão dentro possuem classificação $y=0$.
Para descobrir os valores dos coeficientes em $w$ vamos, assim como nas outras regressões, minimizar uma função de erro através da descida de gradiente.
A função de erro é definida assim:
$$J(w) = \frac{1}{m} \sum_{i=1}^{m} Cost(f(x^{(i)}), y^{(i)})$$Onde $m$ é a quantidade de exemplos de treinamento.
$Cost(f(x), y)$ varia de acordo com o valor de $y$. Para $y=1$, temos:
$$Cost(f(x), y) = -\log(f(x))$$Assim, quanto mais próximo de 1 for $f(x)$, mais próximo de 0 é o valor de $Cost$. Caso contrário, quanto mais próximo de 0 for $f(x)$, mais próximo de infinito é o valor de $Cost$.
Para $y=0$:
$$Cost(f(x), y) = -\log(1 - f(x))$$Nesse caso, temos o oposto, quanto mais próximo de 0 for $f(x)$, mais próximo de 0 é o valor de $Cost$ e quanto mais próximo de 1 for $f(x)$, mais próximo de infinito é o valor de $Cost$.
Com isso, se o modelo atribui probabilidades altas de forma correta paras as classificações dos exemplos, então o valor de $J(w)$ é baixo. Caso contrário, se ele atribui de forma errada probabilidades altas, o valor de $J(w)$ tende a infinito.
Os gráficos a seguir mostram como os valores de $Cost$ variam para $y=1$ (esquerda) e $y=0$ (direita).
from bokeh.layouts import row
x_points = np.linspace(0.01, 0.99, 10000)
y1_points = [-math.log(i) for i in x_points]
y0_points = [-math.log(1-i) for i in x_points]
p1 = figure(plot_width=350, plot_height=250)
p1.xaxis.axis_label = "f(x)"
p1.yaxis.axis_label = "Cost(f(x), 1)"
p1.line(x_points, y1_points, line_width=2)
p0 = figure(plot_width=350, plot_height=250)
p0.xaxis.axis_label = "f(x)"
p0.yaxis.axis_label = "Cost(f(x), 0)"
p0.line(x_points, y0_points, line_width=2)
show(row(p1, p0))
A função $Cost$ pode ser condensada em apenas uma definição:
$$Cost(f(x), y) = -y\log(f(x)) - (1-y)\log(1-f(x))$$Note que quando $y=1$, a função se resume a $Cost(f(x), y) = -\log(f(x))$ e quando $y=0$ ela tem a forma $Cost(f(x), y) = -\log(1-f(x))$. Ou seja, as mesmas definições mostradas anteriormente.
Com isso, usaremos essa definição diretamente na função de erro $J(w)$:
$$J(w) = \frac{1}{m} \sum_{i=1}^{m} [-y^{(i)}\log(f(x^{(i)})) - (1-y^{(i)})\log(1-f(x^{(i)}))]$$Para simplificar, multiplicamos todo o somatório por $-1$, assim trocando os sinais internos do somatório.
$$J(w) = -\frac{1}{m} \sum_{i=1}^{m} [y^{(i)}\log(f(x^{(i)})) + (1-y^{(i)})\log(1-f(x^{(i)}))]$$A seguir é mostrada a implementação da função de erro $J(w)$. Ela recebe as características dos dados de treinamento (parâmetro x) e as suas respostas (parâmetro y), além dos coeficientes (parâmetro w).
O parâmetro x é uma matriz e cada linha dessa matriz corresponde ao vetor de características de um exemplo.
Outro ponto importante é lembrar que o primeiro elemento do vetor de características é sempre 1, por isso é feita a concatenação com [1]
na primeira linha dentro do for
.
def cost(x, y, w):
costs = 0
for i in range(len(x)):
xi = np.concatenate([[1], x[i]])
yi = y[i]
f_value = f(xi, w)
if yi == 1:
if (f_value == 0):
return math.inf
else:
costs -= math.log(f_value)
else:
if (f_value == 1):
return math.inf
else:
costs -= math.log(1 - f_value)
return (1/(2*len(x)))*costs
Vamos aplicar a regressão logística no conhecido conjunto de dados sobre as flores iris. No código abaixo, carregamos e visualizamos dois tipos de flores (azul e vermelho) de acordo com duas características ($x_1$ e $x_2$).
from sklearn import datasets
iris_data = datasets.load_iris()
x_chart = [i[0] for i in iris_data.data[:100]]
y_chart = [i[1] for i in iris_data.data[:100]]
colors = (["blue"] * 50) + (["red"] * 50)
p = figure(plot_width=700, plot_height=250)
p.xaxis.axis_label = "x1"
p.yaxis.axis_label = "x2"
p.circle(x_chart, y_chart, radius=0.015, color=colors)
show(p, notebook_handle=True)
<Bokeh Notebook handle for In[7]>
Para minimizar $J(w)$ usando a descida de gradiente é necessário calcular a derivada parcial da função:
\begin{equation*} \frac{\partial}{\partial w_j} J(w) = \frac{1}{m} \sum_{i=1}^{m} (f(x^{(i)}) - y^{(i)})x_j^{(i)} \end{equation*}As derivadas são calculadas pelo código abaixo. A função recebe as características dos dados de treinamento (parâmetro x), as suas respostas (parâmetro y), os coeficientes (parâmetro w) e o índice $j$ (parâmetro j) da equação da derivada parcial que especifica qual derivada parcial está sendo calculada.
def d_error(x, y, w, j):
total = 0
for i in range(len(x)):
xi = np.concatenate([[1], x[i]])
yi = y[i]
fxi = f(xi, w)
error = (fxi - yi)
total += error*xi[j]
return total/(len(x))
Com isso, atualizamos os coeficientes iterativamente assim:
\begin{equation*} w_j := w_j - \alpha \frac{\partial}{\partial w_j} J(w) \end{equation*}ou seja
\begin{equation*} w_j := w_j - \alpha \frac{1}{m} \sum_{i=1}^{m} (f(x^{(i)}) - y^{(i)})x_j^{(i)} \end{equation*}onde $j$ varia de 0 até $N$, sendo $N+1$ o total de coeficientes.
Vamos ajustar uma reta como fronteira de decisão, ela tem a forma: $w_0 + w_1x_1 + w_2x_2 = 0$.
x = [[i[0], i[1]] for i in iris_data.data[:100]]
y = iris_data.target[:100]
w = [0, 0, 0]
alpha = 0.1
for n in range(10000):
w_temp = [0, 0, 0]
for i in range(len(w)):
w_temp[i] = w[i] - alpha*d_error(x, y, w, i)
w = w_temp
Os coeficientes encontrados foram:
print(w)
[-3.3872040867695046, 6.0550889793412503, -9.5045216464262339]
Dessa forma, foi encontrada a reta: $-3,3872 + 6,055x_1 -9,5045x_2 = 0$. Traçando ela no gráfico, percebemos claramente a divisão entre os dois tipos de flores.
p = figure(plot_width=700, plot_height=250)
p.xaxis.axis_label = "x1"
p.yaxis.axis_label = "x2"
p.circle(x_chart, y_chart, radius=0.015, color=colors)
p.line([4, 8], [(-3.3872+6.055*4)/9.5045, (-3.3872+6.055*8)/9.5045], line_width=2)
show(p, notebook_handle=True)
<Bokeh Notebook handle for In[11]>
Vamos testar o classificador para os pontos (6, 2) e (4,3).
print("Probabilidade do ponto (6,2) ter a classificação 1 (vermelho): {}".format(f([1, 6, 2], w)))
print("Probabilidade do ponto (4,3) ter a classificação 1 (vermelho): {}".format(f([1, 4, 3], w)))
Probabilidade do ponto (6,2) ter a classificação 1 (vermelho): 0.999999111994042 Probabilidade do ponto (4,3) ter a classificação 1 (vermelho): 0.00046157094623368804
Podemos também usar uma fronteira definida por um polinômio. Por exemplo: $w_0 + w_1x_1 + w_2x_2^2 = 0$.
Para isso, escreveremos a função f
e a função para calcular a sua derivada de forma específica para esse exemplo, pois a implementação anterior não tratava expoentes nas variáveis.
def f(x, w):
z = w[0]*x[0] + w[1]*x[1] + w[2]*(x[2]**2)
return g(z)
def d_error(x, y, w, j):
total = 0
for i in range(len(x)):
xi = np.concatenate([[1], x[i]])
yi = y[i]
fxi = f(xi, w)
error = (fxi - yi)
if (j == 2):
total += error*xi[j]**2
else:
total += error*xi[j]
return total/(len(x))
Executando a descida de gradiente:
x = [[i[0], i[1]] for i in iris_data.data[:100]]
y = iris_data.target[:100]
w = [0, 0, 0]
alpha = 0.01
for n in range(1000000):
w_temp = [0, 0, 0]
for i in range(len(w)):
w_temp[i] = w[i] - alpha*d_error(x, y, w, i)
w = w_temp
Os coeficientes encontrados foram:
print(w)
[-27.651233755257426, 8.0858324530145236, -1.6503621951007472]
Traçando a curva no gráfico, temos:
x_points = np.linspace(4, 8, 100)
y_points = [math.sqrt(((w[0] + w[1]*i))/(-w[2])) for i in x_points]
p = figure(plot_width=700, plot_height=250)
p.xaxis.axis_label = "x1"
p.yaxis.axis_label = "x2"
p.circle(x_chart, y_chart, radius=0.015, color=colors)
p.line(x_points, y_points, line_width=2)
show(p, notebook_handle=True)
<Bokeh Notebook handle for In[16]>
Vamos testar o classificador para os pontos (6, 2) e (4,3).
print("Probabilidade do ponto (6,2) ter a classificação 1 (vermelho): {}".format(f([1, 6, 2], w)))
print("Probabilidade do ponto (4,3) ter a classificação 1 (vermelho): {}".format(f([1, 4, 3], w)))
Probabilidade do ponto (6,2) ter a classificação 1 (vermelho): 0.9999993603296184 Probabilidade do ponto (4,3) ter a classificação 1 (vermelho): 3.8640780844487264e-05