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

Como fazer classificação usando "Stepwise logistic regression" usando o python?

+1 voto
67 visitas
perguntada Jul 5 em Estatística por IgorNascimento (51 pontos)  
republicada Jul 5 por IgorNascimento
Compartilhe

1 Resposta

+1 voto
respondida Jul 5 por IgorNascimento (51 pontos)  
editado Jul 6 por IgorNascimento

O modelo logit é utilizado para medir probabilidade da existência de uma determinada características em um unidade \(i\), representado por \(Y_i=1\), utilizando informações \(X\):

\(P(Y_i = 1 | X) = \frac{e^{\beta_0 + \sum_i^{n} \beta_i X_i} }{1 + e^{\beta_0 + \sum_i^{n} \beta_i X_i}}\)

Dessa forma, esse modelo tem como objetivo determinar a probabilidade da existência ou não de terminada características a partir de variáveis explicativas.

A definição de quais informações fazem parte do \(X\) pode seguir um processo sistemático, útil, principalmente, há um grande número de variáveis candidatas a explicativas.

O método Stepwise permite a escolha do conjunto regressor de forma iterativa, admitindo a variação foward, partindo o modelo naive são incluída variáveis para o completo (todas as variáveis), ou backward, partindo do modelo completo são excluídas variáveis para um modelo mais parcimonioso.

Considerando o modelo Stepwise forward, segundo Izenman(2008), o processo de decisão entre inclusão/exclusão de variáveis se dá por meio da comparação da razão F com o valor 2 ou 4:
\(F = \frac{\frac{SS0 - SS1}{df0 - df1}}{\frac{SS1}{df1}}\)
sendo SS0 a soma de resíduos e df0 graus de liberdade do modelo de comparação, SS1 a soma de resíduos e df0 graus de liberdade do modelo mais completo.

Como exemplo, será utilizada a base de dados que identificam se um e-mail é spam com base conteúdo lexical do corpo do texto. Mais informações da base, tais como descrição e download podem ser acessadas no link.

Dessa forma, a variável \( Y\) identifica se a mensagem é spam ou não será modelada utilizando o modelo logit com seleção stepwise forward de \(X\), formado por 57 variáveis que determinam a intensidade de algumas palavras no corpo da mensagem.

Para isso, as bibliotecas utilizadas foram:

from scipy import stats
from scipy.optimize import fmin
import numpy as np
import math as math
import copy

A primeira função definida é estabelece o cálculo da log-verossimilhança do modelo logit:

def logitEst(betas):
    mLike=0;
    z = np.matrix(np.zeros((len(y),1),dtype = float))
    for i in range(n):
        for j in range(np.shape(xs)[1]):
            z[i,0] = z[i,0] + xs[i,j] * betas[j]
    for i in range(n):
        try:
            prob = (math.exp(z[i,0])) / (math.exp(z[i,0]) + 1)
        except OverflowError:
            prob = float(.9999999999)
        logprob = math.log(stats.binom.pmf(y[i],1,prob))
        mLike = mLike +  logprob
    return(-mLike)

Essa função é importante, pois será utilizada para o processo de otimização utilizando o pacote scipy.optimize. O chute inicial para o vetor de \(\beta\) do modelo é de 0.5 para todos.

def estimacaoRegBin(y,xs):

    betas = np.transpose(np.matrix(np.ones(( 1 ,np.shape(xs)[1])) / 2 ))    
    betas = np.matrix(fmin(logitEst, betas ,xtol=1e-4))
    return(betas)

Em seguida a base de dados baixa é importada por meio da seguinte programação:

    martiz_resultado = np.matrix(np.empty((58,3),float))
        ## Importando a base
    results = []
    with open('C:\\Users\\Igor\\Documents\\Doutorado\\Métodos Computacionais\\Estrela4\\base.txt') as inputfile:
        for line in inputfile:
            results.append(line.strip().split(','))

    base = np.array(results)
    b = np.empty((4601,58),float)
    for i in range(0,4601): 
        for j in range(0,58):
            b[i,j] = float(base[i,j])

    y  = b[:,57]
    n = len(y)

Em seguida, são definidas as variáveis iniciais para o processo stepwise, sendo vecvars o vetor que armazena as variáveis utilizadas pelo modelo.

    vecvars = []
    medida_ajuste = float(0)
    RSS0 = 10**6
    df0 = n

Por fim, o código a seguir realiza as chamadas das funções para o processo sequencial de inclusão de variáveis e armazena os resultados.

    for j in range(-1,np.shape(b)[1]-1): # scanning variables
        if j == -1:
             xs = np.matrix(np.ones((len(y),1),dtype = float))
             vecvars_temp =  copy.deepcopy(vecvars)
        else:
            vecvars_temp =  copy.deepcopy(vecvars)
            vecvars_temp.append(j)
            nvars = len(vecvars_temp)
            xs = np.matrix(np.ones((len(y),nvars+1),dtype = float))
            for i in vecvars_temp :
                xs[:,i +1] = np.transpose(np.matrix(b[:,i]))


        resultado = estimacaoRegBin(y,xs)
        z = xs * np.transpose(resultado)
        martiz_probz = np.matrix(np.empty((4601,3),float))
        for i in range(0,n):
            try:
                prob = (math.exp(z[i,0])) / (math.exp(z[i,0]) + 1)
                martiz_probz[i,0] = prob
            except OverflowError:
                prob = float(.9999999999)
                martiz_probz[i,0] = prob

            if (prob >= 0.5) :
                z[i,0] = 1
            else:
                z[i,0] = 0

        diag = 0 
        sum_erro = 0 
        for i in range(n):
            if z[i,0] == float(y[i]):
                diag = diag + 1
            sum_erro = sum_erro + (z[i,0] - y[i])**2

        df1      = n - np.shape(resultado)[1]   
        ## Utilizando aproximação normal
        Fvalue = ((RSS0 - sum_erro)/(df0 - df1))/(sum_erro/df1)

        print("rss0= ",RSS0," e rss1=",sum_erro," com valor F=",Fvalue)
        print("vf= ",vecvars)
        print("vt= ",vecvars_temp)
        print("df0=",df0," e df1=",df1)





        ## testando o quantidade de valores corretos
        medida_ajuste_temp = float(diag / n )
        if medida_ajuste_temp > medida_ajuste:
            medida_ajuste = copy.deepcopy(medida_ajuste_temp)

        ## testando o valor F
        if Fvalue > 2:
            RSS0 = copy.deepcopy(sum_erro)
            df0 = copy.deepcopy(df1)
            vecvars=  copy.deepcopy(vecvars_temp)

        #armazenando SS
        martiz_resultado[j+1,0] = sum_erro
        martiz_resultado[j+1,1] = Fvalue
        martiz_resultado[j+1,2] = medida_ajuste_temp

As variáveis que estão foram selecionadas para compor o modelo estão armazenadas no objeto vecvars.

A imagem a seguir representa o processo, sendo cada barra uma variável e seu tamanho a proporções da diagonal (acurácia) da matriz de confusão. As palavras representadas pelas barras cinzas não entram no modelo. Ao todo são 40 barras vermelhas, isto é, são´incluídas 40 variáveis das 57, sendo a acurácia final do modelo de 92.37%.

A imagem será apresentada aqui.

Referências

Izenman, A. (2008). Modern Multivariate Statistical Techniques. Assessment (Vol. 65, p. 733). Springer New York. https://doi.org/10.1007/978-0-387-78189-1

Lichman, M. (2013). UCI Machine Learning Repository [http://archive.ics.uci.edu/ml]. Irvine, CA: University of California, School of Information and Computer Science.

comentou Jul 10 por Carlos Eduardo Véras (11 pontos)  
Explicação simples e direta! Parabéns! O código funcionou corretamente. Foi bacana o fato de você ter utilizado o tratamento de exceções no código.  Só sugiro colocar o arquivo de dados no mesmo diretório do código de modo que não seja necessário digitar todo o caminho quando da sua leitura.
...