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

Como implementar o algoritmo de seleção de variáveis conhecido como Forwards (Stepwise) Selection em python?

+1 voto
167 visitas
perguntada Jul 1, 2016 em Aprendizagem de Máquinas por Matheus Costa (26 pontos)  

Considere a questão 5.18 do livro "Modern multivariate statistical techniques" - Alan Julian Izenman que versa sobre o uso da base de dados "Boston housing data"e variable selection algorithms. Interprete os resultados e documente cuidadosamente o seu material e implementação. Vamos considerar essa questão dividida em partes:

(a) Como implementar o algoritmo de seleção de variáveis em Python conhecido como Forwards (Stepwise) Selection em python? Compare esses resultados com as estimativas de OLS usando os dados descritos acima.

Compartilhe

1 Resposta

0 votos
respondida Jul 1, 2016 por Matheus Costa (26 pontos)  

O método Forward Selection de seleção de variáveis para regressão se baseia em incluir uma a uma as variáveis independentes disponíveis na ordem de magnitude do teste F entre a regressão sem tal variável e a que à inclui. O processo é repetido até que se atinja um critério mínimo predeterminado para a inclusão.

O teste F é realizado através da comparação entre a soma dos quadrados dos resíduos de duas regressões, uma incluindo o conjunto total de variáveis independentes e outra com um conjunto restrito. Para aplicação de Forward Selection, estaremos comparando sempre regressões com apenas uma variável de diferença. O valor do teste F é dado por:

\[ F = \frac{RSS_0 - RSS_1}{RSS_1}\frac{df_1}{df_0-df_1} \]

Onde \(RSS_0\) é a soma dos quadrados dos resíduos do modelo original, com \(df_0\) graus de liberdade, e \(RSS_1\) é a do modelo estendido, ou seja, incluindo a nova variável, com \(df_1\) graus de liberdade. Note que \(df_0-df_1 = 1\) para todos os testes que realizaremos (dado que incluímos uma variável de cada vez).

Nós partimos do estágio inicial que inclui apenas o intercepto e gradativamente incluímos as novas variáveis. Como critério mínimo para inclusão utilizaremos o valor de F de 6,635, valor da estatística F para 1% de significância.

O conjunto de dados utilizados é o "Boston House Prices", que elenca as medianas dos valores de casas ocupadas por região de Boston e diversas variáveis com potencial explicativo para estes preços. Obtivemos os dados através da biblioteca scikit-learn. As variáveis explicativas disponíveis são:

    1- CRIM     per capita crime rate by town
    2- ZN       proportion of residential land zoned for 
                   lots over 25,000 sq.ft.
    3- INDUS    proportion of non-retail business acres per town
    4- CHAS     Charles River dummy variable 
                    (= 1 if tract bounds river; 0 otherwise)
    5- NOX      nitric oxides concentration (parts per 10 million)
    6- RM       average number of rooms per dwelling
    7- AGE      proportion of owner-occupied units built 
                     prior to 1940
    8- DIS      weighted distances to five Boston 
                     employment centres
    9- RAD      index of accessibility to radial highways
  10- TAX      full-value property-tax rate per $10,000
  11- PTRATIO  pupil-teacher ratio by town
  12- B        1000(Bk - 0.63)^2 where Bk is the proportion 
                     of blacks by town
  13- LSTAT    % lower status of the population
  14- MEDV     Median value of owner-occupied 
                     homes in $1000's

Sendo a última, MEDV, nossa variável dependente. Antes de implementar o método vamos definir funções para realizar o teste F e as regressões lineares (MQO) que serão necessárias.

import numpy as np
from numpy import *
from sklearn.datasets import load_boston
boston = load_boston()

y= boston.target
X= np.array(boston.data)
i = np.array([[1] for j in range(len(X))])
X = np.concatenate((X,i), axis = 1)  
print boston.DESCR
print X.shape

SST = 0
ym = mean(y)
for i in range(len(X)):
        SST = SST + y[i]**2 - ym**2
print SST

def Ftest(a,b,c):
    return ((a-b)/b)*(c)

def linreg(X,y):
                # Implementacao da regressao linear
    n=len(X)
    Xt = X.T
    XXI=linalg.inv(matmul(Xt,X))
    Xy = matmul(Xt,y)
    b = matmul(XXI,Xy)
    res = [0 for i in range(n)]
    SSR = 0
    for i in range (n):
        res[i] = y[i] - matmul(X,b)[i]
        SSR = SSR + res[i]**2
    R2 = 1 - SSR/SST
    return b, SSR, R2

Passando agora ao método de fato. X é o vetor de dados, já com o vetor unitário para cálculo do intercepto. y é a variável independente. Cr será o critério de parada, ou seja, o mínimo valor do teste F para que a variável seja incluída. O algoritmo se utiliza de um procedimento recursivo, incluindo a cada iteração a variável que expressa o maior valor no teste F com relação à regressão anterior. O vetor "a" incorpora as variáveis até então incluídas e utilizaremos um vetor "z" para elencar as variáveis que permanecem disponíveis.

z=[]
def StepF(X,y,Cr, a=[]):
    global z
    n=len(X)   # Numero de observacoes
    b=[z,[]]   # 'b' elenca na primeira linha as variaveis
               # nao incluidas e na segunda o valor
               # do teste F obtido com a
               # inclusao desta variavel.
    k=len(a)+1
    global R0  # R0 sera o a soma dos quadrados 
               # dos residuos da regressao restrita.
    if len(a)==0:
        R0 = SST
        z=[i for i in range(len(X[0]))]
    for i in z:
               # Usamos a funcao 'for' para obter o valor
               # do teste F com cada uma das variaveis ainda 
               # nao incluidas em cada iteracao.
        ai = a[:]
        ai.append(i)
        Z = np.array(X[:,ai])
        Ri = linreg(Z, y)[1] 
               # 'linreg(Z, y)[1]' eh o valor do SSR da regressao
        b[1].append(Ftest(R0, Ri, n-k-1))
    if len(z)==0:
        Xi=np.array(X[:,a])
        return linreg(Xi,y), a
    else:
        j = np.argmax(b[1])
        if b[1][j] >= Cr:
            a.append(z[j])
            Xi = np.array(X[:,a])
            R0 = linreg(Xi,y)[1]
            z.remove(z[j])
            return StepF(X,y,Cr,a)
                # Caso a variavel com o maior valor de F
                # satisfaca o criterio de inclusao, a 
                # adicionamos ao vetor 'a' e chamamos 
                # de novo a a funcao StepF com o novo vetor 'a'
        else:
            Xi=np.array(X[:,a])
            return linreg(Xi,y), 'a =', a, 'z = ', z
                # Caso nenhuma das variaveis satisfaca o criterio 
                # de inclusao, encerramos a funcao e retornamos os
                # os coeficientes da regressao obtida, seu SSR, R2
                # o vetor 'a' de variaveis incluidas e
                # o vetor z de variaveis excluidas.

print StepF(X,y,6.635)
print linreg(X,y)

Como resultado, obtemos que as variáveis incluídas em ordem de relevância foram:

    5- NOX      nitric oxides concentration (parts per 10 million)
   12- B        1000(Bk - 0.63)^2 where Bk is the proportion 
                     of blacks by town
   10- TAX      full-value property-tax rate per $10,000
   11- PTRATIO  pupil-teacher ratio by town
    7- AGE      proportion of owner-occupied units built 
                     prior to 1940
   13- LSTAT    % lower status of the population
    4- CHAS     Charles River dummy variable 
                    (= 1 if tract bounds river; 0 otherwise)
    3- INDUS    proportion of non-retail business acres per town
    1- CRIM     per capita crime rate by town

Notamos entre elas que a mais importante é a concentração óxidos de nitrogênio. Entre as outras, no entanto, podemos perceber algumas variáveis cuja endogeneidade não pode ser facilmente descartada, como o valor da taxa sobre propriedade (10) e a porcentagem de negros entre a população (11), esta última diretamente relacionada com o nível de riqueza da população. Mas não faz parte do escopo da pergunta levantar este tipo de discussão.

Comparando com a regressão linear irrestrita, com todas variáveis, percebemos que apesar da exclusão de 5 das variáveis disponíveis, o método de Forward Selection atingiu um nível de previsão muito próximo ao do modelo original, provando sua eficácia.

\[R^2_{FS} = 0.7266\]
\[R^2_{Ols} = 0.7406\]

...