Primeira vez aqui? Seja bem vindo e cheque o FAQ!
x

Você pode me dar um exemplo de como implementar um classificador usando um modelo de resposta binária?

0 votos
186 visitas

1 Resposta

0 votos
respondida Ago 22, 2015 por danielcajueiro (5,566 pontos)  

Se você não sabe ainda o que é um modelo de resposta binária, dê uma olhada aqui.

Esses modelos já estão implementados em Python, veja por exemplo, aqui. Também já estão implementados em R ou Matlab.

A resposta que darei aqui implementa explicitamente o classificador em Python focando o entendimento e aprendizado. Outras pessoas poderiam dar uma resposta explicitando rapidez de impletamentação utilizando as ferramentas já prontas do Python, R ou Matlab.

Utilizaremos a notação para o modelo de resposta binária considerada nessa outra resposta.

Para testar nossa implementação, vamos usar dados gerados utilizando simulações Monte-Carlo. A forma mais interessante de gerar amostras de simulações monte carlo para modelos de resposta binária é:

1) Crie \(n\) valores de um vetor bidimensional de regressores com uma distribuição pre-especificada.

2) Para cada um desses \(n\) valores, calcule \(F(z_i)\), onde \(z_i=\beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2}\), usando os \(\beta\)s reais (lembre-se que estamos fazendo o papel da natureza quando usamos a simulação de Monte-Carlo) e a distribuição acumulada \(F\)
desejada (nesse caso usamos normal, mas poderiamos ter usado a distribuição logística).

3) Sorteie \(n\) valores de de uma variável aleatória uniforme \(w_i,\;\;i=1,\cdots,n\) no intervalo \([0,1]\);

4) Gere os \(n\) valores de \(y_i\) fazendo o seguinte teste: Se \(w_i\in [0,F(x_i'\beta)]\) então \(y_i=1\). Em caso contrário, \(y_i=0,\;\;i=1,\cdots,n\).

Entenda os passos (3) e (4) com cuidado. Na vida real, muitas vezes temos características que nos sugerem que vamos tomar um determinado tipo de atitude e tomamos uma atitude diferente daquela prevista (erros de modelagem). Então esses passos nos dizem que quanto maior a probabilidade medida por \(F\), maior a chance de eu me comportar como esperado.

Apenas para comparação vamos também gerar valores de \(y_i\) que refletem perfeitamente as características individuais. Substitua os passos (3) e (4) pelo passo (3\(^\prime\)):

3\(^\prime\)) Gerar os $n$ valores de \(y_i\) fazendo o seguinte teste: Se \(F(x_i'\beta)>1/2\) então \(y_i=1\). Em caso contrário, \(y_i=0,\;\;i=1,\cdots,n\).

A figura abaixo apresenta essas duas simulações acima. Se usarmos o passo (3\(^\prime\)) geramos dados que são perfeitamente separáveis e em caso contrário geramos dados que refletem mais situações reais e não são perfeitamente separáveis.

A imagem será apresentada aqui.

A reta usada para classificação é exatamente

\[z=\beta_0 + \beta_1x_1 +\beta_2 x_2=0\]

Note que a reta contínua utiliza os parâmetros reais e a reta tracejada usa os parâmetros estimados.

Note também que podemos provar geometricamente que a distância de um ponto \((x_1,x_2)\) a reta usada para classificação nos dá uma medida de quanto ela atende bem as características daquela região do classificador.

O código utilizado para essa implementação em Python está aqui (se você ainda não sabe Python e gostaria de aprender, siga por aqui):

import matplotlib.pyplot as plt
from scipy import stats
from scipy.optimize import fmin
import numpy as np
import math as math
np.random.seed(50)

plt.figure(num=None, figsize=(8, 4), dpi=80, facecolor='w', edgecolor='k')

def geraClassificador(beta0,beta1,beta2,x1,parameter):
    minx1=np.min(x1)
    maxx1=np.max(x1)    

    minx2=(beta0+beta1*minx1)/(-beta2)
    maxx2=(beta0+beta1*maxx1)/(-beta2)
    if(parameter):    
        plt.plot([minx1,maxx1],[minx2,maxx2],'k-')  
    else:
        plt.plot([minx1,maxx1],[minx2,maxx2],'y--')  

def geracaoDados(n,beta0,beta1,beta2):


    x1=np.random.normal(1,0.64,n) # Primeiro regressor
    x2=np.random.normal(1,0.36,n) # Segundo regressor
    z=beta0*np.ones(n)+beta1*x1+beta2*x2 # Variavel latente

    w=np.random.uniform(0,1,n)
    y=np.empty([n])
    yP=np.empty([n])    


    for i in range(n):
        if(stats.norm.cdf(z[i])>0.5):
            yP[i]=1
        else:
            yP[i]=0    


        if(w[i]< stats.norm.cdf(z[i])):
            y[i]=1
        else:
            y[i]=0    

    return yP,y,x1,x2


def probitEst(betas,*data):
    (beta0,beta1,beta2)=betas
    y,x1,x2=data
    n=np.size(y)
    mLike=0;
    for i in range(n):
        z=beta0+beta1*x1[i]+beta2*x2[i]
        F=stats.norm(0, 1).cdf(z)
        if((F != 0) and (F != 1)):
            mLike=mLike+ y[i]*math.log(F) + (1-y[i])*math.log(1-F);
    return -mLike

def estimacaoRegBin(y,x1,x2):
    data=(y,x1,x2)
    beta00=0.5
    beta10=0.5
    beta20=0.5
    betas = fmin(probitEst, [beta00,beta10,beta20], data,xtol=1e-4)
    return betas


    beta, alpha, r_value, p_value, std_err = stats.linregress(x,y)
    return alpha,beta

if __name__ == '__main__':
    n=100 # numero de observacoes
    beta0=2.5
    beta1=-1.5
    beta2=-1

    [yP,y,x1,x2]=geracaoDados(n,beta0,beta1,beta2)        

    axes=plt.subplot(121)

    plt.hold(True)
    for i in range(n):
        if(yP[i]==1):
            plt.plot(x1[i],x2[i],u'ro',markersize=5, markeredgewidth=0)
        else:
            plt.plot(x1[i],x2[i],u'bo',markersize=5, markeredgewidth=0)   

    axes.set_xlabel('$x_1$')
    axes.set_ylabel('$x_2$')          
    axes.set_title('Perfeitamente Separavel')
    geraClassificador(beta0,beta1,beta2,x1,True)
    [hatBeta0,hatBeta1,hatBeta2]=estimacaoRegBin(yP,x1,x2)
    print hatBeta0,hatBeta1,hatBeta2
    geraClassificador(hatBeta0,hatBeta1,hatBeta2,x1,False)    

    axes=plt.subplot(122)

    plt.hold(True)
    for i in range(n):
        if(y[i]==1):
            plt.plot(x1[i],x2[i],u'ro',markersize=5, markeredgewidth=0)
        else:
            plt.plot(x1[i],x2[i],u'bo',markersize=5, markeredgewidth=0)   
    axes.set_xlabel('$x_1$')
    axes.set_ylabel('$x_2$')  
    axes.set_title('Nao Perfeitamente Separavel')    
    geraClassificador(beta0,beta1,beta2,x1,True)
    [hatBeta0,hatBeta1,hatBeta2]=estimacaoRegBin(y,x1,x2)
    print hatBeta0,hatBeta1,hatBeta2
    geraClassificador(hatBeta0,hatBeta1,hatBeta2,x1,False)    

Talvez alguém possa contribuir com o mesmo código em R ou Matlab em uma outra resposta. Mais detalhes sobre modelos de resposta binária e essa implementação você pode encontrar nessa referência.

comentou Fev 7 por Everton (1 ponto)  
Boa tarde professor Daniel, testando o exemplo não estava plotando, ai eu acrescentei
plt.show()
ao final de todo o código para poder exibir.
comentou Fev 7 por danielcajueiro (5,566 pontos)  
Ótimo! Acho que isso tem a ver com versão e plataforma. Eu escrevi esses códigos no Python 2.7 e deveria ter atualizado pois os uso como material de um curso, mas não tive mais tempo. Se vc estiver usando o Python 3.5, vc deve ter tido que colocar parenteses no "print" também. Obrigado.
...