Aula_Regularizacao.knit

class: inverse, middle

## Seleção de variáveis e Regularização

.large[#### Agatha Rodrigues
]

---
class: middle, inverse

# Contexto

---
class: middle

# Contexto

- Em um cenário de regressão, o modelo linear padrão é:
`\begin{align*}
  Y = \beta_0 + \beta_1X_1 + \ldots + \beta_pX_p + \epsilon  
\end{align*}`

- Alternativas para selecionar as variáveis e reduzir dimensão:

-- **Seleção de subconjuntos**: Nessa abordagem selecionamos um subconjunto das `$p$` covariáveis.

-- **Encurtamento** (Shrinkage): Essa abordagem envolve ajustar um modelo com as `$p$` covariáveis. As covariáveis vistas como não importantes são encurtadas para zero. Esse encurtamento (regularização) tem o efeito de reduzir a variância. Dependendo do tipo de encurtamento realizado, alguns coeficientes podem ter estimativas exatamente iguais a zero. Métodos de regularização podem também realizar seleção de variáveis.

-- **Redução de dimensão**: Essa abordagem envolve em projetar os `$p$` preditores em um subespaço `$M$`-dimensional, com `$M<p$`. Isso é obtido por calcular `$M$` combinações lineares, ou projeções, das variáveis. Essas `$M$` projeções são usadas como preditores para ajustar o modelo de regressão.

---
class: middle

# Por que não utilizar todas as covariáveis?

1) Várias variáveis importam pouco;

2) Muitos coeficientes para estimar

Lembrando que: EQM = Viés^2 + Variância + Erro irredutível

Variância: erro por ter que estimar muitos coeficientes.

---
class:middle

---
class: middle, inverse

# Seleção de subconjuntos

---
class: middle

# Seleção de subconjuntos

- Para selecionar o melhor subconjunto de covariáveis, ajustamos modelo de regressão para cada combinação das `$p$` covariáveis.

- Ajustamos todos os `$p$` modelos com só um preditor, `$p(p-1)/2$` modelos que contem exatamente dois preditores e assim por diante.

- O problema de selecionar o melhor modelo dentre as `$2^p$` possibilidades consideradas pelo melhor subconjunto não é trivial. Por isso, a sugestão a seguir de realizar em dois estágios:

1) Seja `$\mathcal{M}_0$` o modelo sem preditores.
  
  2) Para `$k=1,2,\ldots,p$`:
  
   a) Ajuste de todos os `$\binom{p}{k}$` modelos que contém exatamente `$k$` covariáveis;
  
   b) Selecione o melhor entre os `$\binom{p}{k}$` e o chame de  `$\mathcal{M}_k$`. Esse melhor modelo pode ser aquele com menor SQR (ou equivalente a ter maior `$R^2$`). 
  
  3) Selecione o melhor dentre os `$\mathcal{M}_0, \mathcal{M}_1,\ldots,\mathcal{M}_p$` modelos, usando um critério de escolha que penalize o número de preditores.

---
class: middle

# Critérios `$C_p$`, AIC, BIC e `$R^2a$`

- O Erro Quadrático Médio (EQM) é `$SQR/n$`. Essas medidas e o `$R^2$` não podem ser utilizados para escolher modelos com diferentes números de covariáveis.

- Outros critérios podem ser considerados para selecionar o melhor modelo dentro de um conjunto com modelos com diferentes números de preditores.

## `$C_p$`

Para um modelo contendo `$d$` preditores, o `$C_p$` estimado de uma amostra de validação com `$m$` observações é dado por:

`\begin{align*}
  C_p = \frac{1}{m}(SQR + 2d\widehat{\sigma}^2),
\end{align*}`
em que `$\widehat{\sigma}^2$` é uma estimativa da variância do `$\epsilon$`, usualmente usando o modelo com os `$p$` preditores.

---
class: middle

## AIC

O critério de AIC de um modelo com `$d$` preditores é dado por:
`\begin{align*}
  AIC=2d-2\ln({\hat {L}}),
\end{align*}`
em que `$\hat {L}$` é o valor da função de verossimilhança avaliada nas estimativas de máxima verossimilhança.

No caso do modelo de regressão normal, máxima verossimilhança e mínimos quadrados são a mesma coisa.  Por esse motivo, `$C_p$` e `$AIC$` são proporcionais, isto é:
`\begin{align*}
  AIC \approx \frac{1}{m}(SQR + 2d\widehat{\sigma}^2).
\end{align*}`

---
class: middle

## BIC

O critério BIC é dado por (ignorando constantes irrelevantes):
`\begin{align*}
  BIC = \frac{1}{m}(SQR + \log(m)d\widehat{\sigma}^2).
\end{align*}`

- Note que o critério BIC substitui o `$2d\widehat{\sigma}^2$` do `$C_p$` por `$\log(m)d\widehat{\sigma}^2$`.

- Como `$\log(m) >2$` para algum `$m> 7$`, o critério BIC impõe uma penalidade mais pesada aos modelos com muitas variáveis, resultando na seleção de modelos menores que o `$C_p$`.

---
class: middle

## `$R^2$` ajustado

- O `$R^2$` é definido como `$1- SQR/SQT$`, em que `$SQT= \sum(y_i-\bar{y})^2$` é a soma dos quadrados total.

- Como SQR sempre decresce conforme mais preditores são adicionados no modelo, o `$R^2$` sempre cresce quando mais covariáveis são adicionadas.

- O `$R^2$` ajustado ($R^2a$) de um modelo com `$p$` covariáveis é dado por:
`\begin{align*}
  R^2a = 1 - \frac{SQR/(m-d-1)}{SQT/(m-1)}.
\end{align*}`

- Maximizar o `$R^2$` ajustado é equivalente a minimizar `$SQR/(m-d-1)$`.

---
class: middle

## Seleção Stepwise

- Por questões computacionais, a abordagem de selecionar o melhor subconjunto não pode ser considerado em cenários com `$p$` muito grande. 
Quanto maior o espaço de buscas, maior a chance de encontrar modelos que parecem bons na amostra treinamento e não terem bom poder preditivo. Assim, uma busca muito, muito grande pode levar ao sobreajuste.

- Enquanto a abordagem de selecionar o melhor subconjunto ajusta `$2^p$` modelos, a seleção stepwise considera um conjunto menor de modelos.

- O método **forward** começa com um modelo sem nenhum preditor e vai adicionando preditores, um de cada vez, até chegar no modelo com `$p$` preditores.

- O método **backward** começa com o modelo completo (com os `$p$` preditores) e vai removendo preditores, um de cada vez, até chegar no modelo sem preditores.

---
class: middle

## Seleção Forward stepwise

**Algoritmo:**

1. Seja `$\mathcal{M}_0$` o modelo sem nenhum preditor;

2. Para `$k=0,\ldots,p-1$`:

a) Considere todos os `$p-k$` modelos que aumentam os preditores em `$\mathcal{M}_k$` com um preditor adicional.
   
   b) Escolha o melhor entre os `$p-k$` modelos e o chame de `$\mathcal{M}_{k+1}$`. Aqui o melhor modelo é definido como tendo a menor SQR (ou maior `$R^2$`). 
   
3. Selecione o melhor modelo dentre `$\mathcal{M}_0, \mathcal{M}_1,\ldots,\mathcal{M}_p$` usando `$C_p$` (AIC), BIC ou `$R^2a$`.

**Observação 1:** A abordagem do melhor subconjunto e a seleção forward podem não selecionar o mesmo modelo.

**Observação 2:** Na seleção forward, são ajustados no total `$1+ p + (p-1) + (p-2) + \ldots + 2 + 1 = 1 + p(p+1)/2$`, menor que `$2^p$` modelos ajustados necessários no método da seleção do melhor subconjunto.

---
class: middle

## Seleção Backward stepwise

**Algoritmo:**

1. Seja `$\mathcal{M}_p$` o modelo com todos os preditores;

2. Para `$k=p,p-1,\ldots,1$`:

a) Considere todos os `$k$` modelos com todos preditores, com exceção do preditor `$k$` em `$\mathcal{M}_k$`, para um total de `$k-1$` preditores. 
   
   b) Escolha o melhor entre os `$k$` modelos e o chame de `$\mathcal{M}_{k-1}$`. Aqui o melhor modelo é definido como tendo a menor SQR (ou maior `$R^2$`). 
   
3. Selecione o melhor modelo dentre `$\mathcal{M}_0, \mathcal{M}_1,\ldots,\mathcal{M}_p$` usando `$C_p$` (AIC), BIC ou `$R^2a$`.

**Observação:**  A seleção backward requer que o número de observações na amostra `$n$` seja maior que o número de preditores `$p$` (assim o modelo completo pode ser ajustado). Em contraste, a seleção forward pode ser usada mesmo quando `$n<p$`, e é o único método de subconjunto viável quando p é muito grande.

---
class: middle

## Abordagem híbrica (forward e backward)

Nessa abordagem, as variáveis são adicionadas sequencialmente no modelo, uma analogia à seleção forward. No entanto, após adicionar cada nova variável, o método também pode remover algumas variáveis que não adicionam nenhuma melhora ao ajuste.

Essa abordagem tenta imitar a seleção do melhor subconjunto, enquanto conserva as vantagens computacionais das seleções forward e backward.

---
class: middle, inverse

# Regularização

---
class: middle

## Regularização

O termo regularização refere-se a um conjunto de técnicas utilizadas para especificar modelos que ajustem a um conjunto de dados evitando o sobreajuste. Em termos gerais, essas técnicas servem para ajustar modelos de regressão baseados em uma função de perda que contém um termo de penalização. Esse termo tem a finalidade de reduzir a influência de coeficientes responsáveis por flutuações excessivas.

Consideramos aqui três técnicas:

- Ridge (regularização `$L_2$`);

- Lasso (regularização `$L_1$`);

- Elastic Net (mistura das duas anteriores).

---
class: middle

Antes de ir para os métodos de regularização:

- A busca pelo melhor subconjunto como explicado anteriormente é exatamente encontrar os coeficientes `$\beta$`'s que:
`\begin{align*}
 argmin_{\beta}\left[ \frac{1}{m}SQR + K\sum_{j=1}^p {I}(\beta_j \neq 0)\right]
\end{align*}`

Lembre que no modelo linear, EQM = SQR/m e se o critério de seleção for o AIC, `$K=2\widehat{\sigma}^2/m$`.

- O problema dessa equação é que não tem um jeito rápido de resolvê-la, uma vez que tem a somatória de variáveis indicadoras. Na verdade, a única maneira de resolver esse problema é investigando todos `$2^p$` modelos e escolher entre todos eles aquele que me dá o menor risco estimado.

- Nessa forma, a magnitude dos dados não importa! Ele apenas analisa se é ou não igual a 0.

---
class: middle

## Regularização `$L_2$` (Ridge)

- A técnica de regressão Ridge foi introduzida inicialmente para tratar problemas de multicolinearidade, mas também pode ser utilizada para corrigir problemas ligados ao sobreajuste.

- O modelo linear:

`\begin{align*}
  y_i = \beta_0 + \beta_1x_{i1} + \ldots + \beta_px_{ip} + \epsilon_i,
\end{align*}`
para `$i=1,\ldots,n$`, e com as `$p$` covariáveis reunidas no vetor $ = (1, x_{i1}, \ldots, x_{ip})^{\top}$ e sendo `$\beta = (\beta_{0}, \beta_{1},\ldots, \beta_{p})^{\top}$` o vetor de parâmetros a serem estimados.

- Os **estimadores de mínimos quadrados penalizados** correspondem à solução de 
`\begin{align*}
  \widehat{\beta}_{Ridge}(\lambda) = argmin_{\beta}\left[\sum_{i=1}^n(y_i-\beta^{\top}{\bf x}_i)^2 + \lambda\sum_{j=1}^p\beta_j^2\right],
\end{align*}`
em que `$\lambda\geq 0$` é o **coeficiente de regularização**, que controla a importância relativa entre a minimização da soma dos quadrados dos erros e o termo de penalização.

-- Se `$\lambda = \infty$`, não há variáveis a serem incluídas no modelo.

-- Se `$\lambda = 0$`, obtemos os estimadores de mínimos quadrados.

---
class: middle

## Regularização `$L_2$` (Ridge)

Obter o mínimo da expressão do slide anterior é equivalente a minimizar a soma de quadrados não regularizada sujeita à restrição
`\begin{align*}
  \sum_{j=1}^p\beta_j^2 \leq s,
\end{align*}`
para algum valor apropriado de `$s$`, ou seja, é um problema de otimização com **multiplicadores de Lagrange**. Em outras palavras, para cada valor de `$\lambda$`, existe algum `$s$` tal que as duas soluções retornarão as mesmas estimativas dos coeficientes Ridge.

O **estimador Ridge** pode ser expresso como
`\begin{align*}
   \widehat{\beta}_{Ridge}(\lambda) =({\bf X}^{\top}{\bf X} + \lambda {\bf I})^{-1}{\bf X}^{\top}{\bf y},
\end{align*}`
com `${\bf X}$` sendo a natriz de especificação do modelo e `${\bf I}$` é a matriz identidade.

---
class: middle

## Regularização `$L_2$` (Ridge)

Se as covariáveis são ortonormais `$(X^{\top}X=I)$`, temos:
`\begin{align*}
   \widehat{\beta}_{Ridge}  = \widehat{\beta}_{MQ}/(1+\lambda) 
\end{align*}`

Assim, 
`\begin{align*}
   Var[\widehat{\beta_j}_{Ridge}]  = \frac{Var[\widehat{\beta_j}_{MQ}] }{(1+\lambda) ^2} 
\end{align*}`

-- Encolhimento dos coeficientes

---
class: middle

## Regularização `$L_1$` (Lasso)

- O **estimador Lasso** é obtido por 
`\begin{align*}
  \widehat{\beta}_{Lasso}(\lambda) = argmin_{\beta}\left[\sum_{i=1}^n(y_i-\beta^{\top}{\bf x}_i)^2 + \lambda\sum_{j=1}^p|\beta_j|\right]. 
\end{align*}`

Nesse caso, a restrição é 
`\begin{align*}
  \sum_{j=1}^p|\beta_j| \leq s.
\end{align*}`

Algumas propriedades do estimador Lasso:

i) Parâmetros que correspondem a preditores redundantes são encolhidos para zero.

ii) Sob certas condições, o estimador Lasso descarta as variáveis irrelevantes do modelo atribuindo pesos nulos aos respectivos coeficientes -- esparsidade dos coeficientes.

iii) Por isso, a técnica Lasso pode ser considerada para seleção de variáveis.

---
class: middle

## Regularizações Ridge e Lasso

<img src="img/representacao_lasso_ridge.jpg" width="75%" style="display: block; margin: auto;" />
Representação da restrição Ridge à esquerda e da restrição Lasso à direita

---
class: middle

## Regularização Elastic Net

- O **estimador Elastic Net** (EN) é obtido por 
`\begin{align*}
  \widehat{\beta}_{EN}(\lambda) = argmin_{\beta}\left[\sum_{i=1}^n(y_i-\beta^{\top}{\bf x}_i)^2 + \lambda_1\sum_{j=1}^p\beta_j^2  + \lambda_2\sum_{j=1}^p|\beta_j|\right]. 
\end{align*}`

Nesse caso, a restrição é 
`\begin{align*}
  \lambda_1\sum_{j=1}^p\beta_j^2  + \lambda_2\sum_{j=1}^p|\beta_j|\leq s.
\end{align*}`

---
class: middle

## Regularização Elastic Net
 
 Na figura abaixo apresentamos esquematicamente uma região delimitada pela restrição 
`\begin{align*}
  J(\beta) = (1-\alpha)\sum_{j=1}^p\beta_j^2 + \alpha\sum_{j=1}^p|\beta_j| \leq s
\end{align*}`

<img src="img/representacao_elastic-net.jpg" width="55%" style="display: block; margin: auto;" />
Geometria das restrições Elastic Net (curva contínua), Ridge (curva tracejada) e Lasso (curva pontilhada).

---
class: middle

## Comentários sobre a Regressão Logística Lasso e Ridge

Seja `$l(\beta)$` o logaritmo da função de verossimilhança para `$\beta$`. O valor de `$\beta$` que maximiza `$l(beta)$` é o mesmo que minimiza o `$l^*(\beta)= -l(\beta)$`.

O estimador Lasso de `$\beta$` é o que minimiza
`\begin{align*}
  l^*(\beta) + \lambda\sum_{j=1}^p |\beta_j|
\end{align*}`

O estimador Ridge de `$\beta$` é o que minimiza
`\begin{align*}
  l^*(\beta) + \lambda\sum_{j=1}^p \beta_j^2
\end{align*}`

O estimador Elastic Net de `$\beta$` é o que minimiza
`\begin{align*}
  l^*(\beta) + \lambda\left[\alpha \sum_{j=1}^p |\beta_j| + (1-\alpha)\sum_{j=1}^p \beta_j^2\right]
\end{align*}`

---
class: middle

## Observações finais

**Como escolher `$\lambda$`?**

Em um grid de possíveis valores de `$\lambda$`, selecionar aquele com o menor risco estimado na amostra de validação.

**Pacote no R glmnet**