Análise de Regressão

Seleção de variáveis e validação de modelos

Prof. Dr. Juliano Bortolini

Depto. de Estatística — UFMT

Conteúdo

Métodos para construir e validar modelos de regressão (Montgomery, Cap. 10 e 11):

Seleção forward (progressiva)
Seleção backward (regressiva / eliminação)
Seleção stepwise (passo a passo)
Validação cruzada leave-one-out (LOOCV)
Validação cruzada k-fold

Exemplo-guia: os dados de cimento de Hald (Montgomery, Cap. 10).

Por que selecionar variáveis?

Ao construir um modelo com \(K\) regressores candidatos, há um compromisso:

Poucas variáveis → modelo enviesado (falta estrutura);
Muitas variáveis → variância alta, multicolinearidade, sobreajuste.

Queremos descrever a resposta com o menor número de regressores possível que ainda explique bem a variabilidade de \(y\) (princípio da parcimônia).

Os métodos stepwise avaliam poucos subconjuntos, adicionando ou removendo regressores um a um com base em critérios estatísticos.

Critérios para comparar modelos

Critério	Fórmula (modelo com \(p\) termos)	Escolha
\(R^2_p\)	\(1 - SS_{Res}(p)/SS_T\)	maior (mas sempre cresce)
\(R^2_{aj}\)	\(1 - \dfrac{SS_{Res}(p)/(n-p)}{SS_T/(n-1)}\)	maior
\(C_p\) de Mallows	\(\dfrac{SS_{Res}(p)}{\hat{\sigma}^2} - n + 2p\)	pequeno, \(C_p \approx p\)
AIC	\(n\ln\!\big(SS_{Res}/n\big) + 2p\)	menor
BIC	\(n\ln\!\big(SS_{Res}/n\big) + p\ln n\)	menor
PRESS	\(\sum_{i=1}^{n}\big(e_i/(1-h_{ii})\big)^2\)	menor (preditivo)

O step() do R usa AIC; o Minitab/olsrr usam p-valor (F-para-entrar/sair).

Os dados: cimento de Hald

Dados de cimento de Hald (Montgomery, Cap. 10 — Tabela B.21)

Calor liberado \(y\) (cal/g) na cura do cimento em função da composição:

\(x_1\): aluminato tricálcico (3CaO·Al₂O₃)
\(x_2\): silicato tricálcico (3CaO·SiO₂)
\(x_3\): aluminoferrito tetracálcico (4CaO·Al₂O₃·Fe₂O₃)
\(x_4\): silicato dicálcico (2CaO·SiO₂)

\(n = 13\) observações; \(K = 4\) regressores candidatos.

hald <- data.frame(
  x1 = c(7, 1, 11, 11, 7, 11, 3, 1, 2, 21, 1, 11, 10),
  x2 = c(26, 29, 56, 31, 52, 55, 71, 31, 54, 47, 40, 66, 68),
  x3 = c(6, 15, 8, 8, 6, 9, 17, 22, 18, 4, 23, 9, 8),
  x4 = c(60, 52, 20, 47, 33, 22, 6, 44, 22, 26, 34, 12, 12),
  y  = c(78.5, 74.3, 104.3, 87.6, 95.9, 109.2, 102.7, 72.5, 93.1, 115.9, 83.8, 113.3, 109.4)
)
head(hald)
#>   x1 x2 x3 x4     y
#> 1  7 26  6 60  78.5
#> 2  1 29 15 52  74.3
#> 3 11 56  8 20 104.3
#> 4 11 31  8 47  87.6
#> 5  7 52  6 33  95.9
#> 6 11 55  9 22 109.2

1. Seleção forward (progressiva)

Começa com o modelo nulo (só intercepto) e adiciona regressores um a um:

Entra primeiro o regressor de maior correlação simples com \(y\) (maior \(F\)-para-entrar);
A cada passo entra o de maior correlação parcial com \(y\) (dado o que já está), se seu \(F\) parcial exceder \(F_{IN}\);
Termina quando nenhum candidato exceder \(F_{IN}\) (Minitab: \(\alpha = 0{,}25\)).

Limitação: uma vez incluído, um regressor nunca é removido, mesmo que se torne redundante depois.

Forward — no R

nulo     <- lm(y ~ 1, data = hald)
completo <- lm(y ~ x1 + x2 + x3 + x4, data = hald)

# Seleção forward pelo critério AIC
forward <- step(nulo, scope = ~ x1 + x2 + x3 + x4,
                direction = "forward", trace = 0)
formula(forward)
#> y ~ x4 + x1 + x2
round(coef(forward), 4)
#> (Intercept)          x4          x1          x2 
#>     71.6483     -0.2365      1.4519      0.4161

Para reproduzir a seleção por p-valor do livro (como no Minitab), use o pacote olsrr: ols_step_forward_p(completo, penter = 0.25).

Forward — modelo final

Resultado (Montgomery, Exemplo 10.2) — entram \(x_4\), depois \(x_1\), depois \(x_2\):

\[ \hat{y} = 71{,}65 + 1{,}452\,x_1 + 0{,}416\,x_2 - 0{,}237\,x_4 \]

Regressores selecionados: \(x_1, x_2, x_4\).

Comentário: \(x_4\) entra primeiro (maior correlação com \(y\)), mas depois de entrar \(x_2\) ele se torna quase redundante — como o forward não remove, \(x_4\) permanece no modelo final.

2. Seleção backward (eliminação)

Começa com o modelo completo (todos os \(K\) regressores) e remove um a um:

Calcula o \(F\) parcial (ou \(t\)) de cada regressor, como se fosse o último a entrar;
Remove o de menor \(F\) parcial, se estiver abaixo de \(F_{OUT}\);
Reajusta e repete; termina quando o menor \(F\) parcial já não for menor que \(F_{OUT}\) (Minitab: \(\alpha = 0{,}10\)).

Muito usada por analistas que gostam de ver o efeito de todos os regressores antes de eliminar.

Backward — no R

backward <- step(completo, direction = "backward", trace = 0)
formula(backward)
#> y ~ x1 + x2 + x4
round(coef(backward), 4)
#> (Intercept)          x1          x2          x4 
#>     71.6483      1.4519      0.4161     -0.2365

Resultado (Montgomery, Exemplo 10.3) — sai \(x_3\) (t = 0,14), depois \(x_4\):

\[ \hat{y} = 52{,}58 + 1{,}468\,x_1 + 0{,}662\,x_2 \]

Regressores selecionados: \(x_1, x_2\).

Obs.: por p-valor (livro) o backward para em \(\{x_1,x_2\}\); por AIC o step() mantém \(x_4\) — critérios diferentes podem levar a modelos diferentes.

3. Seleção stepwise (passo a passo)

Combina forward e backward (algoritmo de Efroymson, 1960):

É um forward em que, a cada passo, todos os regressores já incluídos são reavaliados pelo seu \(F\) (ou \(t\)) parcial;
Um regressor incluído antes pode se tornar redundante e ser removido;
Usa dois cortes: \(F_{IN}\) (entrar) e \(F_{OUT}\) (sair), em geral \(F_{IN} \ge F_{OUT}\).

Corrige a principal limitação do forward (não conseguir remover).

Stepwise — no R

stepwise <- step(nulo, scope = ~ x1 + x2 + x3 + x4,
                 direction = "both", trace = 0)
formula(stepwise)
#> y ~ x4 + x1 + x2
round(coef(stepwise), 4)
#> (Intercept)          x4          x1          x2 
#>     71.6483     -0.2365      1.4519      0.4161

Resultado (Montgomery, Exemplo 10.4) — entra \(x_4\), \(x_1\), \(x_2\) e depois remove \(x_4\):

\[ \hat{y} = 52{,}58 + 1{,}468\,x_1 + 0{,}662\,x_2 \]

Regressores selecionados: \(x_1, x_2\) (igual ao backward e ao all possible regressions).

Obs.: como no backward, o step() por AIC mantém \(x_4\); a seleção por p-valor (livro/olsrr) resulta em \(\{x_1, x_2\}\).

Comentário: os três métodos concordam?

Para os dados de Hald (Montgomery, p. 366):

Método	Selecionados
Forward	\(x_1,\ x_2,\ x_4\)
Backward	\(x_1,\ x_2\)
Stepwise	\(x_1,\ x_2\)

Os métodos nem sempre concordam — a ordem de entrada/saída depende das intercorrelações entre regressores;
\(x_4\) é fortemente correlacionado com \(x_2\) (multicolinearidade: \(VIF_4 \approx 18{,}9\));
A ordem de entrada não indica ordem de importância;
Use os métodos com cautela: em geral há vários modelos igualmente bons.

4–5. Validação cruzada (Cap. 11)

Como o modelo se comporta ao prever novas observações? A validação cruzada (uma forma de data splitting) separa os dados em treino e teste.

Ideia — Montgomery (Cap. 11) mostra que a estatística PRESS já é uma forma de validação cruzada:

\[ PRESS = \sum_{i=1}^{n}\left(\frac{e_i}{1-h_{ii}}\right)^2, \qquad R^2_{\text{pred}} = 1 - \frac{PRESS}{SS_T}. \]

Vamos comparar dois modelos (Montgomery, Tabela 10.4): (1) \(x_1, x_2, x_4\) e (2) \(x_1, x_2\).

Leave-one-out (LOOCV)

Cada observação é usada uma vez como teste, treinando nas outras \(n-1\); o erro é a média dos \((y_i - \hat{y}_{(i)})^2\). Para MQO há forma fechada:

\[ y_i - \hat{y}_{(i)} = \frac{e_i}{1-h_{ii}} \;\;\Rightarrow\;\; \text{MSE}_{LOOCV} = \frac{1}{n}\sum_{i=1}^{n}\left(\frac{e_i}{1-h_{ii}}\right)^2 = \frac{PRESS}{n}. \]

loocv_mse <- function(mod) mean( (residuals(mod) / (1 - hatvalues(mod)))^2 )

m1 <- lm(y ~ x1 + x2 + x4, data = hald)   # 3 regressoras
m2 <- lm(y ~ x1 + x2,      data = hald)   # 2 regressoras
c(m1_x1x2x4 = loocv_mse(m1), m2_x1x2 = loocv_mse(m2))
#> m1_x1x2x4   m2_x1x2 
#>    6.5655    7.2217

LOOCV — com o pacote `boot`

cv.glm() com K = n faz o LOOCV automaticamente:

library(boot)
g1 <- glm(y ~ x1 + x2 + x4, data = hald)
g2 <- glm(y ~ x1 + x2,      data = hald)

c(m1_x1x2x4 = cv.glm(hald, g1)$delta[1],
  m2_x1x2   = cv.glm(hald, g2)$delta[1])
#> m1_x1x2x4   m2_x1x2 
#>    6.5655    7.2217

Conclusão (LOOCV) — o modelo (1) \(x_1,x_2,x_4\) tem erro menor (\(\approx 6{,}57\)) que o (2) \(x_1,x_2\) (\(\approx 7{,}22\)). Equivale a comparar \(PRESS\): \(85{,}35 < 93{,}88\).

k-fold cross-validation

Divide os dados em \(k\) grupos (folds): treina em \(k-1\) e testa no restante, repetindo \(k\) vezes; o erro é a média dos \(k\) erros. LOOCV é o caso \(k=n\).

set.seed(1)               # o resultado depende da partição aleatória
c(m1_x1x2x4 = cv.glm(hald, g1, K = 5)$delta[1],
  m2_x1x2   = cv.glm(hald, g2, K = 5)$delta[1])
#> m1_x1x2x4   m2_x1x2 
#>    6.1104    5.8000

Com \(k = 5\), o modelo (1) \(x_1,x_2,x_4\) também apresenta erro de previsão ligeiramente menor — a conclusão é a mesma do LOOCV (use set.seed() para reprodutibilidade, pois os folds são sorteados).

Conclusão dos itens 4–5

Tanto o LOOCV quanto o k-fold indicam erro de previsão um pouco menor para o modelo com três regressoras \(\{x_1,x_2,x_4\}\) do que com duas \(\{x_1,x_2\}\) (coerente com \(PRESS\): 85,35 vs 93,88 e \(R^2_{\text{pred}}\): 0,968 vs 0,965).

Mas há um compromisso: o ganho preditivo é pequeno e \(x_4\) introduz forte multicolinearidade (\(VIF \approx 18{,}9\)), o que torna os coeficientes instáveis. Por parcimônia e estabilidade, o modelo \(\{x_1,x_2\}\) costuma ser preferido — foi o escolhido por backward, stepwise e all possible regressions.

Síntese

Forward / backward / stepwise: adicionam/removem regressores um a um por critério estatístico (p-valor ou AIC); podem levar a modelos diferentes;
Avalie os candidatos por \(R^2_{aj}\), \(C_p\), AIC, BIC e PRESS;
Validação cruzada estima o erro de previsão: LOOCV (\(=PRESS/n\) no MQO) e k-fold;
Nos dados de Hald, \(\{x_1,x_2,x_4\}\) prevê um pouco melhor, mas \(\{x_1,x_2\}\) é mais parcimonioso e sem multicolinearidade.

Seleção automática é um ponto de partida, não um substituto para a análise crítica do modelo (resíduos, influência, multicolinearidade e contexto).

Referências

Montgomery, D. C.; Peck, E. A.; Vining, G. G. Introduction to Linear Regression Analysis, 6ª ed. Wiley, 2021 — Cap. 10 (Variable Selection and Model Building) e Cap. 11 (Validation of Regression Models).

Prof. Dr. Juliano Bortolini · www.julianobortolini.com.br

@profjulianobortolini · linkedin · github