4.2. Regressão Linear Múltipla
Regressão múltipla é uma coleção de técnicas estatísticas para construir modelos que descrevem de maneira razoável relações entre várias variáveis explicativas de um determinado processo. A diferença entre a regressão linear simples e a múltipla é que na múltipla são tratadas duas ou mais variáveis explicativas.
Apresentamos a teoria para analisar e validar Modelos de Regressão Linear Múltipla (MRLM), focando nos seguinte tópicos:
-
Modelo estatístico;
-
Estimação dos parâmetros;
-
Propriedades dos estimadores;
-
Testes de hipóteses e intervalos de confiança para os parâmetros do modelo;
-
Análise de diagnóstico do modelo;
-
Aplicações.
Motivação 2
O ganho de um transistor consiste na diferença entre o emissor e o coletor. A variável Ganho (em hFE) pode ser controlada no processo de deposição de íons por meio das variáveis Tempo de emissão (em minutos) e Dose de íons ($\times 10^{14}$). Os dados encontram-se na Tabela 4.2.1.
| Observação | Tempo (min) | Ganho (hFe) | Dose de íons ($\times10^{14}$) |
|---|---|---|---|
| 1 | 195 | 1004 | 4 |
| 2 | 255 | 1636 | 4 |
| 3 | 195 | 852 | 4,6 |
| 4 | 255 | 1506 | 4,6 |
| 5 | 225 | 1272 | 4,2 |
| 6 | 225 | 1270 | 4,1 |
| 7 | 225 | 1269 | 4,6 |
| 8 | 195 | 903 | 4,3 |
| 9 | 255 | 1555 | 4,3 |
| 10 | 225 | 1260 | 4 |
| 11 | 225 | 1146 | 4,7 |
| 12 | 225 | 1276 | 4,3 |
| 13 | 225 | 1225 | 4,72 |
| 14 | 230 | 1321 | 4,3 |
Tabela 4.2.1: Dados do ganho em um conjunto de pistões à diferentes níveis de emissão e coleta.
Ilustração do Modelo
Figura 4.2.1: Ganho vs Tempo vs Dose.
A análise gráfica (Figura 4.2.1) nos mostra que existe relação entre o ganho dos transmissores e as duas variáveis de entrada (Tempo e Dose). Assim, é razoável supor uma relação linear entre a variável Ganho e as variáveis Tempo e Dose.
2.1 Modelo Estatístico
Como visto na “Motivação 2”, supor a construção de um modelo para relacionar a variável ganho (em íons) com as variáveis explicativas emissor de tempo e emissor da dose é razoável. Assim, definimos o modelo de regressão linear múltipla dado por
$$Y=\beta_0+\beta_1x_1+\beta_2x_2+\varepsilon \tag{2.1}$$
em que $Y$ representa a variável resposta (o ganho em íons), $x_1$ e $x_2$ representam as variáveis explicativas (o emissor de tempo e o emissor de dose, respectivamente) e $\varepsilon$ representa o erro experimental. Esse é um modelo de regressão linear múltipla com duas variáveis independentes ou explicativas ($x_1$ e $x_2$). O termo linear é usado pois a equação (2.1) é uma função linear de parâmetros desconhecidos $\beta_0,\beta_1$ e $\beta_2,$ denominados coeficientes da regressão.
Interpretação dos parâmetros do modelo
- O parâmetro $\beta_0$ corresponde ao intercepto do plano com o eixo z. Se $x=(x1,x2)=(0,0)$ o parâmetro $\beta_0$ fornece a resposta média nesse ponto. Caso contrário, não é possível interpretar o parâmetro $\beta_0$.
- O parâmetro $\beta_1$ indica uma mudança na resposta média a cada unidade de mudança em $x1$, quando as demais variáveis são mantidas fixas.
- De forma semelhante é a interpretação para o parâmetro $\beta_2$, que indica uma mudança na resposta média a cada unidade de mudança em $x2$, quando $x1$ é mantido constante.
Figura 4.2.2: Ilustração do modelo
Supondo $E(\varepsilon)=0$, temos $E(Y|x)=\beta_0+\beta_1x_1+\beta_2x_2, $ que descreve um plano bidimensional, denominado superfície de resposta.
De maneira geral, a variável resposta $Y$ pode ser relacionada a um número $p$ de variáveis de entrada. O modelo de regressão linear múltipla (MRLM) com $p$ variáveis explicativas é dado por
$$Y_i=\beta_{0}+\beta_{1}x_{i1}+\beta_{2}x_{i2}+…+\beta_{p}x_{ip}+\epsilon_i, \qquad i=1,…,n \tag{2.2}$$
em que
- $x_{i1},x_{i2},…,x_{ip}$ são valores das variáveis explicativas, constantes conhecidas;
- $\beta_{0},\beta_{1},\beta_{2},…,\beta_{p}$ são parâmetros ou coeficientes da regressão;
- $\epsilon_i$ são erros aleatórios independentes.
Este modelo descreve um hiperplano p-dimensional referente às variáveis explicativas.
Efeito das interações
Modelos mais complexos do que o “Modelo 2.2” também são analisados usando técnicas de regressão linear múltipla. Consideremos o modelo de regressão linear múltipla com duas variáveis regressoras, $x_1$ e $x_2$, dado por
$$ Y=\beta_0+\beta_1 x_{1}+\beta_2 x_{2}+\beta_{3}\underbrace{x_{1} x_{2}}_{\hbox{interação}} + \varepsilon$$
Neste caso, $x_1x_2$ representa a interação existente entre as variáveis $x_1$ e $x_2$. Se a interação está presente e é significativa, o efeito de $x_{1}$ na resposta média depende do nível de $x_{2}$ e analogamente o efeito de $x_{2}$ na resposta média depende do nível de $x_{1}.$
Sabendo que $E(\varepsilon)=0$, tem-se que
$$E(Y|x)=\beta_0+\beta_1x_1+\beta_2x_2+\beta_3x_1x_2$$
A interpretação para os parâmetros $\beta_1$ e $\beta_2,$ no modelo com interação, não é o mesmo visto anteriormente.
Suposições para o modelo
As suposições necessárias para o Modelo de Regressão Linear Múltipla são:
i) O erro tem média zero e variância $\sigma^2$, desconhecida;
ii) Os erros são não correlacionados;
iii) Os erros têm distribuição normal;
iv) As variáveis regressoras $x_{1}, x_{2}, \ldots, x_{p}$ assumem valores fixos.
As suposições (i)-(iii), simbolicamente, podem ser representadas por
$$\varepsilon_{i} \overset{\hbox{iid}}{\sim} N(0, \sigma^2)$$
Se as suposições do MRLM se verificam, então a variável $Y$ tem distribuição normal com variância $\sigma^2$ e média
$$E(Y\mid x)=\beta_0+\beta_1 x_{1}+\beta_2 x_{2}+\ldots+\beta_p x_{p}$$
Neste caso, os parâmetros $\beta_j,$ $j=1,\dots,p$ representam a variação (média) esperada na variável resposta ($Y$) quando a variável $x_j$ sofre um acréscimo unitário, enquanto todas as outras variáveis $x_i \quad (i\neq j)$ são mantidas constantes. Por esse motivo os $\beta_j$ são chamados de coeficientes parciais.
Se os valores de $x_j$ incluem os valores $x_j=0,j=1,\dots,p$ então $\beta_0$ é a média de $Y$ quando $x_j=0.$ Em caso contrário não existe interpretação prática.
2.2 Estimação dos Parâmetros do Modelo
Suponha que temos $n$ observações $(n> p)$ da variável resposta e das $p$ variáveis explicativas. Assim, $y_i$ é o valor da variável resposta na i-ésima observação enquanto que $x_{ij}$ é o valor da variável $x_j$ na i-ésima observação, $j=1,\dots,p.$ Os dados de um MRLM podem ser representados da seguinte forma:
| $y$ | $x_1$ | $x_2$ | $\ldots$ | $x_p$ |
|---|---|---|---|---|
| $y_1$ | $x_{11}$ | $x_{12}$ | $\ldots$ | $x_{1p}$ |
| $y_2$ | $x_{21}$ | $x_{22}$ | $\ldots$ | $x_{2p}$ |
| $\vdots$ | $\vdots$ | $\vdots$ | $\vdots$ | $\vdots$ |
| $y_n$ | $x_{n1}$ | $x_{n2}$ | $\ldots$ | $x_{np}$ |
Tabela 4.2.2: Representação dos dados.
em que cada observação satisfaz $$\begin{equation*}Y_i=\beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \ldots + \beta_p x_{ip} + \varepsilon_i, \quad i=1,\ldots, n.\end{equation*}$$
2.2.1 Método dos Mínimos Quadrados
O objetivo é minimizar a função $$L=\sum_{i=1}^n\varepsilon_i^2=\sum_{i=1}^n\left(Y_i-\beta_0-\beta_1 x_{i1}-\beta_2 x_{i2}-\ldots -\beta_px_{ip}\right)^2.$$
Derivando L em função dos $\beta$’s obtemos $$\dfrac{\partial L}{\partial \beta_0}=-2\sum_{i=1}^n[Y_i-\beta_0-\beta_1x_{i1}-\beta_2x_{i2}-\dots-\beta_px_{ip}],$$
$$\dfrac{\partial L}{\partial \beta_j}=-2\sum_{i=1}^n[Y_i-\beta_0-\beta_1x_{i1}-\beta_2x_{i2}-\dots-\beta_px_{ip}]x_{ji},\quad j=1,2,\dots,p.$$
Igualando as derivadas parciais a zero e substituindo $\beta_0,\beta_1,\dots,\beta_p$ por $\widehat{\beta} _0,\widehat{\beta} _1,\dots,\widehat{\beta} _p$, temos o sistema de equações
$$ \begin{cases} n \widehat{\beta} _{0} + \widehat{\beta} _{1} \displaystyle\sum_{i=1}^n x_{i1} + \widehat{\beta} _{2} \displaystyle\sum_{i=1}^n x_{i2} + \ldots + \widehat{\beta} _{p} \displaystyle\sum_{i=1}^n x_{ip} = \displaystyle\sum_{i=1}^n Y_i \cr \cr \widehat{\beta} _{0} \displaystyle\sum_{i=1}^n x_{i1} + \widehat{\beta} _{1} \displaystyle\sum_{i=1}^n x_{i1}^2 + \widehat{\beta} _{2} \displaystyle\sum_{i=1}^n x_{i1}x_{i2} + \ldots + \widehat{\beta} _{p} \displaystyle\sum_{i=1}^n x_{i1}x_{ip} = \displaystyle\sum_{i=1}^n x_{i1}Y_i \cr \cr \qquad \vdots \qquad \qquad \qquad \vdots \cr \cr \widehat{\beta} _{0} \displaystyle\sum_{i=1}^n x_{ip} + \widehat{\beta}_{1} \displaystyle\sum_{i=1}^n x_{ip}x_{i1} + \widehat{\beta} _{2} \displaystyle\sum_{i=1}^n x_{ip}x_{i2} + \ldots + \widehat{\beta} _{p} \displaystyle\sum_{i=1}^n x_{ip}^2 = \displaystyle\sum_{i=1}^n x_{ip} Y_i \end{cases} $$
Resolvendo este sistema, obtemos os estimadores de mínimos quadrados $\widehat{\beta} _0,\dots,\widehat{\beta}_p$ dos parâmetros do modelo em questão.
2.2.2 Representação matricial do MRLM
Notemos que os estimadores de mínimos quadrados dos parâmetros do “Modelo 2.2” podem ser facilmente encontrados considerando a notação matricial dos dados, que é de fácil manipulação. Desta forma, considerando a entrada de dados apresentada na Tabela 2.2.1, o modelo de Regressão Linear Múltipla pode ser escrito como $$\begin{equation*}Y=X\beta+\varepsilon,\end{equation*}$$
com
$$Y=\begin{bmatrix} Y_1 \cr Y_2 \cr \vdots \cr Y_n \end{bmatrix} \quad , \quad X=\begin{bmatrix} \ x_{11} \ x_{12}~\ldots~x_{1p} \cr 1~x_{21}~x_{22}~\ldots~x_{2p} \cr \quad ~\vdots \quad ~\vdots \quad \ddots \ \vdots \cr 1~x_{n1}~x_{n2}~\ldots~x_{np} \end{bmatrix} \quad , \quad \beta =\begin{bmatrix} \beta_0 \cr \beta_1 \cr \vdots \cr \beta_p \end{bmatrix} \hbox{e} \quad \varepsilon= \begin{bmatrix}\varepsilon_1 \cr \varepsilon_2 \cr \vdots \cr \varepsilon_n \end{bmatrix} ,$$
em que
-
$Y$ é um vetor $n\times 1$ cujos componentes corresponde às n respostas;
-
$X$ é uma matriz de dimensão $n\times (p+1)$ denominada matriz do modelo;
-
$\varepsilon$ é um vetor de dimensão $n\times 1$ cujos componentes são os erros e
-
$\beta$ é um vetor $(p+1)\times 1$ cujos elementos são os coeficientes de regressão.
O método de mínimos quadrados tem como objetivo encontrar o vetor $\widehat{\beta}$ que minimiza
$$L=\sum_{i=1}^n\varepsilon_i^2=\varepsilon^\prime \varepsilon= (Y-X\beta)^\prime (Y-X\beta)=$$
$$=Y^\prime Y-Y^\prime X\beta-\beta^\prime X^\prime Y+\beta^\prime X^\prime X\beta=Y^\prime Y-2\beta^\prime X^\prime Y+\beta^\prime X^\prime X\beta,$$
sendo que $Y^\prime X\beta=\beta^\prime X^\prime Y$ pois o produto resulta em um escalar. A notação $X^\prime $ representa o transposto da matriz $X$ enquanto que $Y^\prime$ e $\beta^\prime$ representam os transpostos dos vetores $Y$ e $\beta$, respectivamente. Usando a técnica de derivação (em termos matriciais) obtemos
$$\dfrac{\partial L}{\partial\beta}=-2X^\prime Y+2X^\prime X\beta.$$
Igualando a zero e substituindo o vetor $\beta$ pelo vetor $\widehat{\beta}$, temos
$$(X^\prime X)\widehat{\beta}=X^\prime Y.$$
Em geral, a matriz $(X^\prime X)$ é não singular, ou seja, tem determinante diferente de zero, e portanto é invertível. Desta forma, conclui-se que os estimadores para os parâmetros $\beta_j,\quad j=0,\dots,p$ são dados pelo vetor
$$\begin{equation*}\widehat{\beta}=(X^\prime X)^{-1} X^\prime Y.\end{equation*}$$
Portanto, o modelo de regressão linear ajustado e o vetor de resíduos são respectivamente
$$\begin{equation*}\widehat{Y}=X\widehat{\beta}\quad {e}\quad e=Y-\widehat{Y}=Y - \widehat{Y}.\end{equation*}$$
Ao substituirmos os estimadores de mínimos quadrados, obtemos que $\widehat{Y} = HY$ no qual $H=X(X^\prime X)^{-1} X^\prime$ é a matriz chapéu, ou matriz de projeção do vetor de respostas $Y$ no vetor de respostas ajustadas $\widehat{Y}$.
Exemplo 2.2.1
Com os dados do exemplo na “Motivação 2”, obter as estimativas dos parâmetros do “Modelo 2.2”.
Solução:
Sejam a variável resposta Ganho (Y) e as variáveis explicativas Tempo ($x_1$) e Dose ($x_2$).
Temos que $n=14,$ $\quad \sum\limits_{i=1}^{14}Y_i= 17.495$; $\quad\sum\limits_{i=1}^{14}x_{i1}= 3.155$ e $\quad \sum\limits_{i=1}^{14} x_{i2}= 60,72.$ Além disso,
$\sum\limits_{i=1}^{14} x_{i1}^2= 716.425$; $\quad\sum\limits_{i=1}^{14}x_{i2}^2= 264,2584$; $\quad \sum\limits_{i=1}^{14} x_{i1} x_{i2}= 13.683,50$; $\quad \sum\limits_{i=1}^{14} x_{i1} Y_i= 4.001.120$ e $\quad \sum\limits_{i=1}^{14} x_{i2}Y_i= 75.738,30.$
As equações normais serão
$$\left\lbrace \begin{array}{c} n \widehat{\beta} _0+\widehat{\beta} _1 \displaystyle\sum_{i=1}^n x_{i1}+\widehat{\beta} _2\displaystyle\sum_{i=1}^n x_{i2}=\displaystyle\sum_{i=1}^n Y_i \cr \cr \widehat{\beta} _0 \displaystyle\sum_{i=1}^n x_{i1}+\widehat{\beta} _1 \displaystyle\sum_{i=1}^n x_{i1}^2 +\widehat{\beta} _2\displaystyle\sum_{i=1}^n x_{i1} x_{i2}=\displaystyle\sum_{i=1}^n x_{i1} Y_i \cr \cr \widehat{\beta} _0 \displaystyle\sum_{i=1}^n x_{i2}+\widehat{\beta} _1 \displaystyle\sum_{i=1}^n x_{i1} x_{i2}+\widehat{\beta} _2\displaystyle\sum_{i=1}^n x_{i2}^2=\sum_{i=1}^n x_{i2} Y_i\end{array}\right.$$
Substituindo os valores para esse exemplo, temos
$$\left\lbrace \begin{array}{c}14 \widehat{\beta} _0 + 3.155 \widehat{\beta} _1 + 60,72 \widehat{\beta}_2 = 17.495 \cr 3.155 \widehat{\beta} _0 + 716.425 \widehat{\beta} _1 + 13.683,50 \widehat{\beta}_2 = 4.001.120 \cr 60,72 \widehat{\beta} _0 + 13.683,50 \widehat{\beta} _1 + 264,2584 \widehat{\beta} _2 = 75.738,30\end{array}\right.$$
Resolvendo o sistema, obtemos
$\widehat{\beta} _0= -520,08,$ $\widehat{\beta} _1 = 10,78 $ e $\widehat{\beta} _2= -152,15.$
Na representação matricial temos
$$ (X’X)= \begin{bmatrix} 14,00 & 3.155,00 & 60,72 \cr 3.155,00 & 716.425,00 & 13.683,50\cr 60,72 & 13.683,50 & 264,26 \end{bmatrix} \quad\Rightarrow\quad (X’X)^{-1}= \begin{bmatrix} 30,24760 & -0,04172 & -4,78994 \cr -0,04172 & 0,00018 & 0,00004 \cr -4,78994 & 0,00004 & 1,10244 \end{bmatrix} $$
$$ X’y = \begin{bmatrix} 17.495,0 \cr 4.001.120,0\cr 75.738,3 \end{bmatrix} $$
Logo, as estimativas $\widehat{\beta}$ são dadas por
$$ \widehat{\beta} = (X´X)^{-1} X’y = \begin{bmatrix} -520,08 \cr 10,78 \cr -152,15 \end{bmatrix} $$
e portanto, $\widehat{y} = -520,08 + 10,78 x_1 -152,15 x_2.$
Não podemos afirmar que a variável $x_1$ aumenta o ganho, pois as variáveis $x_1$ e $x_2$ estão em unidades diferentes. No exemplo 2.2.3 abaixo, as variáveis foram transformadas para que fiquem na mesma unidade.
Exemplo 2.2.2
Interpretar os resultados do “Exemplo 2.2.1”.
Solução:
Com os dados do exemplo na “Motivação 2”, obtemos o modelo
$$E(Y\mid x)= -520,08 - 152,15 \hbox{Dose} + 10,78 \ \hbox{Tempo}.$$
Se o $\hbox{Tempo} = 225$ (constante), então $$E(Y\mid x)= 1.905,42 - 152,15 \ \hbox{Dose}.$$
Assim,
$\beta_1=-152,15$ indica que a cada acréscimo de uma unidade na $Dose$ a resposta média decrescerá $152,15$ unidades;
$\beta_2=10,78$ indica um acréscimo na resposta média de $10,78$ unidades para cada acréscimo de uma unidade na variável $\hbox{Dose}$.
2.2.3 Transformação de dados
Em determinadas situações é usual transformarmos as variáveis explicativas dos dados originais para que fiquem na mesma unidade (escala), facilitando a interpretação dos resultados. Suponha que $\xi_{1},\xi_{2},\dots,\xi_{p}$ são as variáveis explicativas originais do modelo. A expressão para transformar as variáveis explicativas $\xi$ é dada por
$$x_{ij}=\dfrac{\xi_{ij}-\dfrac{[max(\xi_{ij})+min(\xi_{ij})]}{2}}{\dfrac{[max(\xi_{ij})-min(\xi_{ij})]}{2}},$$
em que j indica a variável que está sendo transformada, $j=1,\dots,p$. Desta forma, temos que os valores das variáveis explicativas transformadas estão todos entre -1 e 1.
Notemos que os valores das estimativas dos parâmetros do modelo de regressão linear múltipla não são os mesmos considerando as variáveis originais e as variáveis transformadas.
Exemplo 2.2.3
Considerando os dados do Exemplo 2.2.1, transformar as variáveis explicativas para que fiquem na mesma unidade (escala). Então, ajustar o modelo de regressão linear múltipla aos dados transformados.
Considerando que $\xi_{1}$ e $\xi_{2}$ são as variáveis explicativas no conjunto de dados original, temos que as variáveis explicativas transformadas, denotadas por $x_1$ e $x_2$ são dadas por
- Para $x_{i1}$:
$$x_{i1}=\dfrac{\xi_{i1} - 225}{30}.$$
- Para $x_{i2}$:
$$x_{i2}=\frac{\xi_{i2} -4,36}{0,36}.$$
| Observação | Tempo (min) $(\xi_{1})$ | Dose de íons $10^{14}$ $(\xi_{2})$ | $x_1$ | $x_2$ | Ganho (y) |
|---|---|---|---|---|---|
| 1 | 195 | 4 | -1 | -1 | 1.004 |
| 2 | 255 | 4 | 1 | -1 | 1.636 |
| 3 | 195 | 4,6 | -1 | 0,6667 | 852 |
| 4 | 255 | 4,6 | 1 | 0,6667 | 1.506 |
| 5 | 225 | 4,2 | 0 | -0,4444 | 1.272 |
| 6 | 225 | 4,1 | 0 | -0,7222 | 1.270 |
| 7 | 225 | 4,6 | 0 | 0,6667 | 1.269 |
| 8 | 195 | 4,3 | -1 | -0,1667 | 903 |
| 9 | 255 | 4,3 | 1 | -0,1667 | 1.555 |
| 10 | 225 | 4 | 0 | -1 | 1.260 |
| 11 | 225 | 4,7 | 0 | 0,9444 | 1.146 |
| 12 | 225 | 4,3 | 0 | -0,1667 | 1.276 |
| 13 | 225 | 4,72 | 0 | 1 | 1.225 |
| 14 | 230 | 4,3 | 0,1667 | -0,1667 | 1.321 |
Tabela 4.2.3: Dados Transformados.
Considerando os dados transformados, temos que a matriz $X$ e o vetor $y$ são dados respectivamente por
$$ X = \begin{bmatrix} 1 & -1 & -1 \cr 1 & 1 & -1 \cr 1 & -1 & 0,6667 \cr 1 & 1 & 0,6667 \cr 1 & 0 & -0,4444 \cr 1 & 0 & -0,7222 \cr 1 & 0 & 0,6667 \cr 1 & -1 & -0,1667 \cr 1 & 1 & -0,1667 \cr 1 & 0 & -1 \cr 1 & 0 & 0,9444 \cr 1 & 0 & -0,1667 \cr 1 & 0 & 1 \cr 1 & 0,1667 & -0,1667 \end{bmatrix}, \quad y = \begin{bmatrix} 1004 \cr 1636 \cr 852 \cr 1506 \cr 1272 \cr 1270 \cr 1269 \cr 903 \cr 1555 \cr 1260 \cr 1146 \cr 1276 \cr 1225 \cr 1321 \end{bmatrix}. $$
Assim, a matriz $X^\prime X$ é
$$(X’X) = \begin{bmatrix} 14 & 0,1667 & -0,8889 \cr 0,1667 & 6,0278 & -0,0278 \cr -0,8889 & -0,0278 & 7,0556 \end{bmatrix} \Rightarrow (X’X)^{-1} = \begin{bmatrix} 0,0720 & -0,0019 & 0,0091 \cr -0,0019 & 0,1660 & 0,0004 \cr 0,0091 & 0,0004 & 0,1429 \end{bmatrix}$$
$$\hbox{e} \qquad X’y = \begin{bmatrix} 17495 \cr 2158,167 \cr -1499,722 \end{bmatrix}. $$
Portanto, as estimativas $\widehat{\beta}$ são
$$ \widehat{\beta} = (X´X)^{-1} X’y = \begin{bmatrix} 1242,31 \cr 323,43 \cr 54,77 \end{bmatrix} $$
Assim, a equação da regressão é dada por
$$\widehat{y}= 1.242,31 + 323,43~x_1 - 54,77~x_2.$$
Usando o software Action temos os seguintes resultados:
- Considerando a escala original dos dados
| Estimativa | Desvio Padrão | Estat.t | P-valor | |
|---|---|---|---|---|
| Intercepto | -520.0766769 | 192.1070916 | -2.7072227 | 0.0203918 |
| Tempo | 10.7811579 | 0.4743196 | 22.7297315 | 0 |
| Dose de ions | -152.1488747 | 36.675439 | -4.1485223 | 0.0016205 |
Tabela 4.2.4: Coeficientes com a escala original
- Considerando a transformação
| Estimativa | Desvio Padrão | Estat.t | P-valor | |
|---|---|---|---|---|
| Intercepto | 1242.3139816 | 9.3744517 | 132.5212414 | 0 |
| x1 | 323.4344987 | 14.2295393 | 22.7297942 | 0 |
| x2 | -54.7730076 | 13.2031035 | -4.1484949 | 0.0016206 |
Tabela 4.2.5: Coeficientes com a transformação
2.3 Propriedades dos Estimadores
2.3.1 Propriedades dos estimadores de mínimos quadrados e do estimador para $\sigma^2$
Consideremos o “Modelo 2.2” na forma matricial. Pelo Teorema de Gauss-Markov temos que o estimador de mínimos quadrados $\widehat{\beta}$ é não viciado e tem variância mínima entre todos os estimadores não viciados que são combinações lineares dos $Y_i$. Assim,
1. Valor esperado (média) de $\widehat{\beta}$:
$$E(\widehat{\beta})= E[(X^\prime X)^{-1}X^\prime Y]=E[(X^\prime X)^{-1}X^\prime (X\beta+\varepsilon)]=E[(X^\prime X)^{-1}X^\prime X\beta+(X^\prime X)^{-1}X^\prime \varepsilon]$$ $$=E[I\beta]+E[(X^\prime X)^{-1}X^\prime \varepsilon]=\beta+(X^\prime X)^{-1}X^\prime E[\varepsilon]=\beta,$$
em que $E[\varepsilon]=0$ e $(X^\prime X)^{-1}X^\prime X=I$ (matriz identidade).
2. Matriz de covariâncias de $\widehat{\beta}$:
Para calcular a variância de $\widehat{\beta},$ vamos primeiramente destacar a definição de variância no caso matricial, ou seja, se $W$ é um vetor de variáveis aleatórias, então a matriz de covariâncias de W é dado por $$Cov(W)= E \left[ W W^\prime\right] - E[W] E(W)^\prime,$$
que na forma matricial é escrita como
$$Cov[W] = \begin{bmatrix} Cov[W_1,W_1] \quad Cov[W_1,W_2] \quad Cov[W_1,W_3] \quad \ldots \quad Cov[W_1,W_n] \cr Cov[W_2,W_1] \quad Cov[W_2,W_2] \quad Cov[W_2,W_3] \quad \ldots \quad Cov[W_2,W_n] \cr Cov[W_3,W_1] \quad Cov[W_3,W_2] \quad Cov[W_3,W_3] \quad \ldots \quad Cov[W_3,W_n] \cr \qquad \vdots \qquad \qquad \qquad \vdots \quad \quad \quad \quad \quad \quad \vdots \quad \quad \ddots \quad \quad \quad \quad \vdots \cr Cov[W_n,W_1] \quad Cov[W_n,W_2] \quad Cov[W_n,W_3] \quad \ldots \quad Cov[W_n,W_n] \end{bmatrix}$$
Com isso, a matriz de covariâncias $\widehat{\beta}$ é
$$Cov(\widehat{\beta})= E \left[\widehat{\beta}\widehat{\beta}^\prime \right] - E[\widehat{\beta}] E(\widehat{\beta})^\prime = E \left\lbrace \left[(X^\prime X)^{-1} X^\prime Y\right] \ \left[(X^\prime X)^{-1} X^\prime Y\right]^\prime \right\rbrace - \beta \beta^\prime$$
$$ = (X’X)^{-1} X' E(YY') X (X’X)^{-1} - \beta \beta' = (X’X)^{-1} X' \left[ \operatorname{Cov}(Y) + E(Y)E(Y)' \right] X (X’X)^{-1} - \beta \beta'$$
$$=(X^\prime X)^{-1} X^\prime ~ Cov(Y) ~ X(X^\prime X)^{-1} + (X^\prime X)^{-1}X^\prime ~ E(Y)E(Y)^\prime X(X^\prime X)^{-1} - \beta \beta^\prime.$$
Fazendo $Cov(Y)=\sigma^2~I$ e também que $E(Y)=X \beta$
$$Cov(\widehat{\beta})=\sigma^2 (X^\prime X)^{-1}\overbrace{ X^\prime I X(X^\prime X)^{-1}}^{=~I} + (X^\prime X)^{-1}X^\prime(X\beta)(X\beta)^\prime X(X^\prime X)^{-1} - \beta \beta^\prime$$
$$ = \sigma^2 (X^\prime X)^{-1} + \overbrace{(X^\prime X)^{-1}X^\prime X}^{=~I}\beta ~ \beta^\prime \overbrace{X^\prime X(X^\prime X)^{-1}}^{=~I} ~ - ~ \beta \beta^\prime$$
$$ = \sigma^2 (X^\prime X)^{-1} + \beta \beta^\prime - \beta \beta^\prime =\sigma^2(X^\prime X)^{-1}.$$
3. Estimador não viciado para $\sigma^2$
Consideremos a soma de quadrados dos resíduos dada por Anchor$$SQE= \sum_{i=1}^n e^2_i = e^\prime e = (Y-\widehat{Y})^\prime (Y-\widehat{Y})= (Y-X\widehat{\beta})^\prime (Y-X\widehat{\beta})$$
$$= Y^\prime Y - Y^\prime X\widehat{\beta} - \widehat{\beta}^\prime X^\prime Y + \widehat{\beta}^\prime X^\prime X\widehat{\beta}= Y^\prime Y - 2\widehat{\beta}^\prime X^\prime Y + \widehat{\beta}^\prime X^\prime X\widehat{\beta}.$$
Desde que $X^\prime X\widehat{\beta}=X^\prime Y,$ segue que
$$SQE = Y^\prime Y - 2\widehat{\beta}^\prime X^\prime Y + \widehat{\beta}^\prime X^\prime Y = Y^\prime Y - \widehat{\beta}^\prime X^\prime Y \quad \hbox{ou ainda,}$$
$$= Y^\prime Y - Y^\prime X (X^\prime X)^{-1}X^\prime Y = Y^\prime (I - X(X^\prime X)^{-1}X^\prime)Y.$$
Portanto, $$SQE = Y^\prime(I - X(X^\prime X)^{-1}X^\prime)Y.$$
Veremos a SQE com mais detalhes em “Análise de Variância”.
Teorema - Distribuição de forma quadrática: Se $Y \sim N_p(\mu;\Sigma)$, então, $Y^\prime AY \sim \chi^2_{r(A),\delta}$ (Qui-quadrado não central) se, e somente se, $A\Sigma$ é idempotente, em que
-
r(A): representa o rank da matriz A, ou seja, o número de colunas linearmente independentes da matriz A.
-
$\delta=\dfrac{1}{2}\mu^\prime A \mu$: representa o parâmetro de não centralidade.
-
Idempotente: $A\Sigma \times A\Sigma = A\Sigma$
Como assumimos que o vetor de erro $\varepsilon~\sim~N_p(0;\sigma^2I_p)$, segue que $Y\sim N_p(X\beta;\sigma^2I_p)$ e então, $$\dfrac{Y}{\sigma}\sim N_p\left(\dfrac{X\beta}{\sigma};I_p\right).$$
Desta forma, utilizando o teorema obtemos que
$$\dfrac{SQE}{\sigma^2}=\dfrac{Y^\prime }{\sigma} [I - X (X^\prime X)^{-1}X^\prime]\dfrac{Y}{\sigma}\sim\chi^{2}_{r[I -X(X^\prime X)^{-1}X^\prime ];\delta}$$
já que a matriz $(I - X (X^\prime X)^{-1}X^\prime )$ é idempotente. Como $$\delta=\dfrac{1}{2}\dfrac{\beta^\prime X^\prime (I-X(X^\prime X)^{-1}X^\prime )X\beta}{\sigma^2}=0 \quad {e}$$ $$r(I-X(X^\prime X)^{-1}X^\prime )=n-(p+1),$$
então $$\dfrac{SQE}{\sigma^2}\sim\chi^{2}_{n-(p+1)}.$$
Portanto, um estimador não viciado para $\sigma^2$ é dado por
$$\widehat{\sigma}^2=QME=\dfrac{SQE}{n-p-1}.$$
4. Matriz de covariâncias estimada de $\widehat{\beta}$:
$$ \widehat{\hbox{Cov}}(\widehat{\beta}) = \hat{\sigma}^2 (X^\prime X)^{-1} = \begin{bmatrix} \hat{\sigma}^2(\hat{\beta}_0) & \widehat{\hbox{Cov}}(\hat{\beta}_0,\hat{\beta}_1) & \widehat{\hbox{Cov}}(\hat{\beta}_0,\hat{\beta}_2) & \ldots & \widehat{\hbox{Cov}}(\hat{\beta}_0,\hat{\beta}_p) \cr \widehat{\hbox{Cov}}(\hat{\beta}_1,\hat{\beta}_0) & \hat{\sigma}^2(\hat{\beta}_1) & \widehat{\hbox{Cov}}(\hat{\beta}_1,\hat{\beta}_2) & \ldots & \widehat{\hbox{Cov}}(\hat{\beta}_1,\hat{\beta}_p) \cr \vdots & \vdots & \vdots & \ddots & \vdots \cr \widehat{\hbox{Cov}}(\hat{\beta}_p,\hat{\beta}_0) & \widehat{\hbox{Cov}}(\hat{\beta}_p,\hat{\beta}_1) & \widehat{\hbox{Cov}}(\hat{\beta}_p,\hat{\beta}_2) & \ldots & \hat{\sigma}^2(\hat{\beta}_p) \end{bmatrix} $$
em que $\widehat{Cov}(\bf{\widehat{\beta}})$ é uma matriz $(p+1)\times (p+1)$, sendo $p$ o número de variáveis explicativas do modelo.
Sendo $C$ a diagonal da matriz ${(X^\prime X)}^{-1}$, isto é
$$C=\begin{bmatrix}C_{00} \quad \quad \quad \quad \quad \quad \quad ~ \cr \quad \quad ~C_{11} \quad \quad \quad \quad \quad \cr \quad \quad \quad \quad \quad \ddots \quad \quad ~ \cr \quad \quad \quad \quad \quad \quad \quad ~C_{pp}\end{bmatrix},$$
podemos escrever a variância estimada dos $\widehat{\beta}_{j}$ como
$$\widehat{\sigma}^{2}(\widehat{\beta_j})=\widehat{\sigma}^2C_{jj},\quad j=0,1,\dots,p.$$
Exemplo 2.3.1
Calcular a matriz de covariâncias estimada considerando os dados transformados do “Exemplo 2.2.3”.
A matriz ${(X^\prime X)}$ neste caso é dada por
$$ (X^\prime X) = \begin{bmatrix} 14,00 & 0,17 & -0,89 \cr 0,17 & 6,03 & -0,03 \cr -0,89 & -0,03 & 7,06 \end{bmatrix} \Rightarrow (X^\prime X)^{-1} = \begin{bmatrix} 0,0720 & -0,0019 & 0,0091 \cr -0,0019 & 0,1660 & 0,0004 \cr 0,0091 & 0,0004 & 0,1429 \end{bmatrix} $$
Temos também que
$$ \hat{\sigma}^2 = QME = \frac{SQE}{n-p-1} = \frac{y^\prime y - \hat{\beta}^\prime X^\prime y}{14 - 2 - 1} = \frac{22.527.889 - 22.514.467,9}{11} = \frac{13.421,2}{11} = 1.220,1. $$
Logo,
$$ \widehat{\text{Cov}}(\widehat{\beta}) = \hat{\sigma}^2 (X^\prime X)^{-1} = \begin{bmatrix} 87,881 & -2,379 & 11,062 \cr -2,379 & 202,481 & 0,497 \cr 11,062 & 0,497 & 174,323 \end{bmatrix}. $$
Podemos também utilizar a matriz $C$, dada por
$$C=\begin{bmatrix} 0,0720 \quad \quad \quad \quad \quad \quad \quad ~ \cr \quad \quad 0,1660 \quad \quad \quad \quad \quad \cr \quad \quad \quad \quad \quad \ddots \quad \quad ~ \cr \quad \quad \quad \quad \quad \quad \quad 0,1429\end{bmatrix},$$
Neste caso, a variância estimada dos estimadores $\widehat{\beta}_{j}\quad j=0,1,2$ é
$$\widehat{\sigma}^{2}(\widehat{\beta_{0}})=\widehat{\sigma}^{2}C_{00}=1.220,1(0,0720)=87,881.$$
$$\widehat{\sigma}^{2}(\widehat{\beta_{1}})=\widehat{\sigma}^{2}C_{11}=1.220,1(0,1660)=202,481.$$
$$\widehat{\sigma}^{2}(\widehat{\beta_{2}})=\widehat{\sigma}^{2}C_{22}=1.220,1(0,1429)=174,323.$$
O desvio padrão dos estimadores é
$$\widehat{\sigma}(\widehat{\beta_{0}})=\sqrt{87,881}=9,374.$$
$$\widehat{\sigma}(\widehat{\beta_{1}})=\sqrt{202,481}=14,230.$$
$$\widehat{\sigma}(\widehat{\beta_{2}})=\sqrt{174,323}=13,203.$$
Usando o software Action temos os seguintes resultados:
| Estimativa | Desvio Padrão | Estat.t | P-valor | |
|---|---|---|---|---|
| Intercepto | 1242.3139816 | 9.3744517 | 132.5212414 | 0 |
| x1 | 323.4344987 | 14.2295393 | 22.7297942 | 0 |
| x2 | -54.7730076 | 13.2031035 | -4.1484949 | 0.0016206 |
Tabela 4.2.6: Estimativas dos coeficientes transformados
2.4 Análise de Variância (Teste F) - Medidas de Associação
2.4.1 Notação Matricial das Somas de Quadrados
Como visto na “Análise de Variância” no Modelo de Regressão Linear Simples, podemos decompor a variabilidade total na variabilidade do modelo mais a variabilidade dos erros.
A soma de quadrados total (SQT), considerando a notação matricial do modelo 2.2, é dada por
$$SQT=\displaystyle\sum\limits_{i=1}^n(Y_i-\overline{Y})^2=Y^\prime Y-\dfrac{Y^\prime J Y}{n}=Y^\prime\left(I-\dfrac{J}{n}\right)Y,$$
em que
$$J =\begin{bmatrix}1 \quad 1 \quad \ldots \quad 1 \cr 1 \quad 1 \quad \ldots \quad 1 \cr \vdots \quad \vdots \quad \ddots \quad \vdots \cr 1 \quad 1 \quad \ldots \quad 1 \cr \end{bmatrix}_{n \times n}$$
Além disso, de “Propriedades dos Estimadores” temos que a soma de quadrados dos erros (dos resíduos) é dada por
$$SQE=Y^\prime Y-\widehat{\beta}^\prime X^\prime Y=Y^\prime (I-X (X^\prime X )^{-1}X^\prime )Y=Y^\prime(I-H)Y.$$
A matriz $I$ é a matriz identidade com dimensão n x n e a matriz $H=X [X^\prime X]^{-1}X^\prime$ é chamada matriz chapéu que transforma o vetor de respostas Y no vetor de valores ajustados $\widehat{Y}$, pois
$$\widehat{Y}=X\widehat{\beta}=X[X^\prime X]^{-1}X^\prime Y=HY.$$
Desta forma, obtemos que a soma de quadrados da regressão é dada por
$$SQR=SQT-SQE=\left(Y^\prime Y-\dfrac{Y^\prime JY}{n}\right)-(Y^\prime Y-\widehat{\beta}^\prime X^\prime Y)=\widehat{\beta}^\prime X^\prime Y-\dfrac{Y^\prime JY}{n}$$ $$=Y^\prime X(X^\prime X)^{-1}X^\prime Y-\dfrac{Y^\prime JY}{n}=Y^\prime \left(X(X^\prime X)^{-1}X^\prime -\dfrac{J}{n}\right)Y=Y^\prime \left(H-\dfrac{J}{n}\right)Y.$$
Notemos que as somas de quadrados da Análise de Variância no caso do MRLM são representadas na forma quadrática $Y’AY$, em que $A$ é uma matriz simétrica.
Vale ressaltar que
-
$H$ é quadrada de dimensão $n$ x $n$ e envolve apenas $X$ (as variáveis independentes assumidas como constantes).
-
$H$ é idempotente, isto é, $HH=H$.
-
As matrizes $\left(I-\dfrac{J}{n}\right)$, $\left(H-\dfrac{J}{n}\right)$ e $(I-H)$ são idempotentes.
Com os valores das somas de quadrados, podemos obter a Tabela Anova.
| Fonte | Soma de quadrados | GL | Quadrado Médio |
|---|---|---|---|
| Regressão | $SQR$ | $p$ | $\dfrac{SQR}{p}$ |
| Erro (Resíduo) | $SQE$ | $n-p-1$ | $\dfrac{SQE}{n-p-1}$ |
| Total | $SQT$ | $n-1$ |
Tabela 4.2.7: Tabela da ANOVA
2.4.2 Teste F
Em problemas de regressão linear múltipla, certos testes de hipóteses sobre os parâmetros do modelo são úteis para verificar a “adequabilidade” do modelo.
2.4.2.1 Teste para significância da regressão
O teste para significância da regressão é um teste para determinar se há uma relação linear entre a variável resposta $Y$ e algumas das variáveis regressora $x_1,x_2,\dots,x_p$. Consideremos as hipóteses
$$\begin{cases}H_0 : \quad \beta_1=\beta_2=\ldots=\beta_p=0 \cr H_1: \quad \beta_j\neq 0 \quad \hbox{ para qualquer}~j~=~1, \cdots, p \cr \end{cases}$$
Se rejeitamos $H_0$, temos que ao menos uma variável explicativa $x_1,x_2,\dots,x_p$ contribui significativamente para o modelo.
Sob $H_0,$ temos pelo “Teorema - Distribuição de forma quadrática” que
$$\dfrac{SQR}{\sigma^2} \sim \chi^2_{(p)} \quad \quad {e~que} \quad \quad \dfrac{SQE}{\sigma^2} \sim \chi^2_{(n-p-1)}.$$
Além disso, temos que $SQR$ e $SQE$ são independentes. Logo, concluímos sob $H_0$ que
$$F_0=\dfrac{\dfrac{SQR}{p}}{\dfrac{SQE}{n-p-1}}= \dfrac{QMR}{QME}~\sim ~F_{(p ; n-p-1)}.$$
Portanto, rejeitamos $H_0$ se $F_0 > F_{(1-\alpha ; p ; n-p-1)}$ e se $p-valor=P[F_{p;n-p-1} > F_0]<\alpha,$ em que $\alpha$ é o nível de significância considerado. Geralmente adotamos $\alpha=5$ %.
A Tabela Anova com a estatística F é dada por
| Fonte | Soma de quadrados | GL | Quadrado Médio | $F_0$ |
|---|---|---|---|---|
| Regressão | $SQR$ | $p$ | $QMR=\dfrac{SQR}{p}$ | $F_0=\dfrac{QMR}{QME}$ |
| Erro (Resíduo) | $SQE$ | $n-p-1$ | $QME=\dfrac{SQE}{n-p-1}$ | |
| Total | $SQT$ | $n-1$ |
Tabela 2.4.8: Tabela da ANOVA e estatística de $F_0$
Exemplo 2.4.1
Construir a tabela da ANOVA considerando os dados transformados no “Exemplo 2.2.3”.
Solução:
Temos as hipóteses:
$$\begin{cases}H_0: \quad \beta_1=\beta_2=0 \cr H_1: \quad \beta_j\neq 0 \quad\hbox{para qualquer}~j~=~1, \cdots, p\end{cases}$$
As somas de quadrados são
$$SQT=y^\prime y-\dfrac{\left( \displaystyle\sum\limits_{i=1}^{14} y_i\right)^2}{14}=22.527.889 - \dfrac{(17.495)^2}{14}$$
$$= 665.387,2.$$
$$SQR= \widehat{\beta}^\prime X^\prime y - \dfrac{\left(\displaystyle\sum\limits_{i=1}^{14}y_i\right)^2}{14} = 22.514.467,9 - 21.862.501,8$$
$$= 651.966,1.$$
$$SQE= SQT-SQR = 665.387,2 - 651.966,1$$
$$= 13.421,1.$$
Resumidamente temos,
| Fonte | Soma de quadrados | GL | Quadrado Médio | $F_0$ | $P-Valor$ |
|---|---|---|---|---|---|
| Regressão | $651.966,1$ | $~2~$ | $\dfrac{651.966,1}{2}=325.983,0$ | $267,18$ | $0,00$ |
| Erro (Resíduo) | $13.421,2$ | $~11~$ | $\dfrac{13.421,2}{11}=1.220,1$ | ||
| Total | $655.387,2$ | $~13~$ |
Tabela 4.2.9: Tabela da ANOVA do Exemplo 2.2.3
Para $\alpha=0,05$, temos que $F_0 = 267,18 > F_{0,95;2;11}=3,98$. Analisando o p-valor, temos que p_valor $= P[F_{2;11} > F_0] = 0,000$. Assim, rejeitamos $H_0$ com um nível de confiança de 95%.
2.4.3 Medidas de Associação
2.4.3.1 Coeficiente de Determinação Múltiplo - $R^2$
O coeficiente de determinação múltiplo é dado por $$R^2 = \frac{SQR}{SQT} = 1-\frac{SQE}{SQT}.$$
Ele representa a proporção da variabilidade de Y explicada pelas variáveis regressoras. Assim, quanto mais próximo $R^2$ estiver de 1, maior é a explicação da variável resposta pelo modelo ajustado.
2.4.3.2 Coeficiente de Determinação Ajustado - $R^2_a$
O coeficiente de determinação ajustado é definido como $$R^2_a=1-\left(\frac{n-1}{n-p}\right)(1-R^2).$$
Este coeficiente ajustado pode ser menor quando outra variável X entra no modelo, pois a diminuição na SQE pode ser compensada pela perda de 1 grau de liberdade no denominador n-p.
Exemplo 2.4.2
Considerando o exemplo na “Motivação 2”, calcular o coeficiente de determinação $(R^2)$ e o coeficiente de determinação ajustado $(R^2_{a})$.
Solução:
$$R^2 = \frac{SQR}{SQT} = \frac{651.966,1}{665.387,2}=0,9798$$
$$R^2_{a}=1-\left(\frac{13}{11}\right)(1-0,9798)= 0,9762.$$
Usando o software Action temos os seguintes resultados:
| G.L. | Soma de Quadrados | Quadrado Médio | Estat.t | P-valor | |
|---|---|---|---|---|---|
| x1 | 1 | 630968.32383546 | 630968.32383546 | 517.146904974 | 0 |
| x2 | 1 | 20997.846071104 | 20997.846071104 | 17.210009911 | 0.001620574 |
| Resíduos | 11 | 13421.04437915 | 1220.094943559 |
Tabela 4.2.10: Tabela da ANOVA
| Desvio Padrão dos Resíduos | Graus de Liberdade | R^2 | R^2 Ajustado |
|---|---|---|---|
| 34.92985748 | 11 | 0.979829723 | 0.9761624 |
Tabela 4.2.11: Medida Descritiva da Qualidade do Ajuste
O Action constrói a tabela da ANOVA de forma diferente da calculada no exemplo 2.4.1. Para obtermos os mesmos resultados do exemplo 2.4.1, basta proceder da seguinte forma:
-
Soma de quadrados da regressão: $SQR=SQ_{x_1}+SQ_{x_2}=630976,866+20998,23=651966,17$
-
Graus de liberdade: $GL=1+1=2$
-
Quadrado médio da regressão: $QMR=\frac{SQR}{2}=\frac{651966,1}{2}=325983$
-
Estatística F: $F_0=\frac{QMR}{QME}=\frac{325983}{1220,1}=267,18$
2.5 Testes Individuais e Intervalos de Confiança para os Parâmetros
2.5.1 Testes individuais para os coeficientes da regressão
Testes de hipóteses individuais para os coeficientes da regressão são fundamentais para se determinar se cada variável explicativa é importante para o modelo de regressão. Por exemplo, o modelo pode ser mais eficaz com a inclusão ou com a exclusão de novas variáveis.
Adicionar uma variável ao modelo de regressão sempre causa um aumento na soma dos quadrados da regressão e um decréscimo na soma dos quadrados do erro. Entretanto, a adição de variáveis regressoras também aumenta a variância do valor ajustado $\widehat{Y}$. Por isso, devemos ter cuidado para incluir somente variáveis regressoras que realmente explicam a variável resposta.
As hipóteses para testar a significância de qualquer coeficiente de regressão individualmente são dadas por,
$$\begin{cases}H_0:\beta_j =0 \cr H_1:\beta_j \neq0 \cr \end{cases} \quad j = 0, 1, \ldots, p.$$
Se $H_0$ ($\beta_j =0$) não é rejeitada, então podemos retirar $x_j$ do modelo já que esta variável não influencia a resposta de forma significativa.
Sabemos que $Y\sim N_p(X\beta;\sigma^2I_p)$ e que $\widehat{\beta}=(X^\prime X)^{-1}X^\prime Y$. Como $\widehat{\beta}$ é combinação linear de distribuições normais, segue que $\widehat{\beta}$ também é normal, isto é
$$\widehat{\beta}\sim N_p(\beta;\sigma^2C),$$
em que $C=(X^\prime X)^{-1}$. Logo, temos que $\widehat{\beta} _j \sim N(\beta_j;\sigma^2C_{jj})$ com $C_{jj}$ sendo o $j$-ésimo elemento da diagonal de $(X^\prime X)^{-1}$, $j=0,1,\dots,p$. Portanto, obtemos
$$N_0=\dfrac{\widehat{\beta} _ j-\beta_j}{\sqrt{\sigma^2 C_{jj}}} \sim N(0,1).$$
Temos também que $$\dfrac{(n-p-1)\widehat{\sigma}^2}{\sigma^2} \sim \chi ^2_{(n-p-1)},$$
independente de $N_0$. Logo, sob $H_0$ temos que a estatística de teste é dada por
$$t_0= \dfrac{\dfrac{\widehat{\beta_j}}{\sqrt{\sigma^2 C_{jj}}}}{\sqrt{\dfrac{(n-p-1)\widehat{\sigma}^2}{\sigma^2}}}=\dfrac{\widehat{\beta}_j}{\sqrt{\widehat{\sigma}^2C_{jj}}} \sim t_{(n-p-1)}, \tag{2.5.1.1}$$
A hipótese nula $H_0:\beta_j =0$ é rejeitada se
$$\mid t_0\mid > t_{(1-\frac{\alpha}{2}, n-p-1)}.$$
Considerando o p-valor, dado por meio da expressão
$$2*P\left( t_{n-p-1}> \mid t_0 \mid \right),$$
rejeitamos $H_0$ se p_valor $<~\alpha$.
O denominador é frequentemente chamado de erro padrão de $\widehat{\beta}_j$ e denotado por
$$se_{(\widehat{\beta}_j)}=\sqrt{\widehat{\sigma}^2C_{jj}}.$$
2.5.2 Intervalo de confiança para os coeficientes da regressão
Considerando a estatística dada em (2.5.1.1), um intervalo com $100(1-\alpha)\char37$ de confiança para os coeficientes da regressão $\beta_j, \quad ~ j=0,1,2,\ldots,p,$ é dado por
$$\left[\widehat{\beta}_j-t_{\left(\frac{\alpha}{2},n-p-1\right)}{\sqrt{\widehat{\sigma}^2 C_{jj}}} ; \widehat{\beta}_j+t_{\left(\frac{\alpha}{2},n-p-1\right)}{\sqrt{\widehat{\sigma}^2 C_{jj}}} \right].$$
2.5.3 Testes parciais para os coeficientes da regressão (Teste F Parcial)
Em alguma situações estamos interessados em saber se um subconjunto de variaveis explicativas são importantes para o modelo de regressão. Considere o modelo escrito na forma matricial com k covariaveis
$$y=X\beta + \varepsilon$$
em que, dado $p=k+1$
- $y$ é um vetor representando a variável resposta com dimensão $n\times 1$;
- $X$ é uma matriz com dimensão $n \times p$
- $\beta$ é o vetor de parâmetros com dimensão $n \times 1$;
- $\varepsilon$ é um vetor com os erros do modelo de dimensão $n \times 1$
Queremos determinar se um subconjunto de covariaveis são significantes ao modelo de regressão. Para isso, considere a seguinte partição do vetor de parâmetros $\beta$:
$$\beta = \begin{bmatrix} \beta_1 \cr - \cr \beta_2 \cr \end{bmatrix} $$
em que $\beta_1$ é um vetor de parâmetros com dimensão $(p-r)\times 1$ e $\beta_2$ é o complemento de $\beta_1$, ou seja, um vetor de parâmetros com dimensão $r \times 1$. Portanto, estamos interessados em testar
$$\begin{cases} H_0: \beta_2 = 0 \cr H_1: \beta_2 \neq 0 \end{cases} $$
Utilizando a partição de $\beta$ podemos reescrever o modelo
$$Y= X\beta + \varepsilon = X_1\beta_1 + X_2\beta_2 + \varepsilon, \tag{2.5.3.1}$$
em que $X_1$ é uma matriz com dimensão $n\times (p-r)$ que representa as colunas de $X$ associadas aos parametros do vetor $\beta_1$ e $X_2$ é uma matriz com dimensão $n \times r$ que representa as colunas da matrix $X$ associadas ao vetor $\beta_2$. Denominamos o modelo 2.5.3.1 por modelo completo.
Para o modelo completo, temos que
$$SQR(\beta)=\widehat{\beta}^{\prime} X^{\prime} y-\dfrac{Y^{\prime}JY}{n}$$
e
$$QME=\dfrac{y^{\prime} y - \widehat{\beta}^{\prime} X^{\prime} y}{n-p}$$
O modelo sob a hipótese nula $H_0:\beta_2 = 0$ verdadeira é dado por
$$Y = X_1\beta_1 + \varepsilon, \tag{2.5.3.2} $$
Denominamos o modelo 2.5.3.2 como modelo reduzido. Portanto temos
$$SQR(\beta_1)=\widehat{\beta} _1^{\prime} X_1^{\prime} y-\dfrac{Y^{\prime}JY}{n}$$
Desta forma, a soma dos quadrados da regressão referente a $\beta_2$ dado que $\beta_1$ ja esta no modelo pode ser determinada por
$$SQR(\beta_2|\beta_1) = SQR(\beta) - SQR(\beta_1)$$
em que $SQR(\beta_2|\beta_1)$ possui $p-(p-r)=r$ graus de liberdade. Essa soma de quadrados da regressão parcial representa a quantidade adicional que teríamos na soma de quadrados da regressão ao adicionar $r$ covariaveis no modelo reduzido. Então, podemos testar a hipótese nula $H_0:\beta_2=0$ utilizando a estatística
$$F_0 = \dfrac{SQR(\beta_2|\beta_1)/r}{QME} \sim F_{r, (n-p)} \tag{2.5.3.3}$$
Portanto se $F_0 > F_{\alpha, r, (n-p)}$ rejeitamos a hipótese nula ao nível de significância $\alpha$ e concluímos que pelo menos um dos parâmetros contidos no vetor $\beta_2$ é diferente de zero. Na literatura o teste 2.5.3.3 é conhecido como Teste F Parcial.
Exemplo 2.5.1
Para ilustrar o uso da estatística $t$, utilizamos novamente os dados transformados no “Exemplo 2.2.3”. Vamos agora construir a estatística $t$ para as hipóteses: $H_0:\beta_0=0$, $H_0:\beta_1=0$ e $H_0:\beta_2=0$.
Solução:
Do “Exemplo 2.2.3” temos que
$$ (X’X)= \begin{bmatrix} 14,00 & 0,17 & -0,89 \cr 0,17 & 6,03 & -0,03 \cr -0,89 & -0,03 & 7,06 \end{bmatrix} \quad\Rightarrow\quad (X’X)^{-1}= \begin{bmatrix} 0,0720 & -0,0019 & 0,0091 \cr -0,0019 & 0,1660 & 0,00004 \cr 0,0091 & 0,0004 & 0,1429 \end{bmatrix} $$
Assim, temos que $C_{00}=0,0720$, $C_{11}=0,1660$ e $C_{22}=0,1429$. Além disso, pelo “Exemplo 2.3.1” segue que $\widehat{\sigma}^2=1.220,10$. Logo, as estatísticas $t$ são dadas por
Para $H_0:\beta_0=0,$ $$t_0=\dfrac{\widehat{\beta} _0}{\sqrt{\widehat{\sigma}^2~C_{00}}}$$
$$=\dfrac{1.242,31}{\sqrt{(1.220,1)(0,0720)}}=\dfrac{1.242,31}{9,3745}=132,521.$$
Para $H_0:\beta_1=0,$ $$t_0=\dfrac{\widehat{\beta} _1}{\sqrt{\widehat{\sigma}^2~C_{11}}}$$
$$=\dfrac{323,43}{\sqrt{(1.220,1)(0,1660)}}=\dfrac{323,43}{14,2296}=22,730.$$
Para $H_0:\beta_2=0,$ $$t_0=\dfrac{\widehat{\beta_2}}{\sqrt{\widehat{\sigma}^2~C_{22}}}$$
$$=\dfrac{-54,77}{\sqrt{(1.220,1)(0,1429)}}=\dfrac{-54,77}{13,203}=-4,149.$$
Como todos os valores absolutos destas estatísticas são maiores do que o valor crítico $t_{(0,975;11)}=2,201,$ as hipóteses $H_0$ são rejeitadas nos três casos. Desta forma, as variáveis tempo e dose de íons contribuem significativamente para o modelo.
Analisando o P-valor, temos que
Para $H_0:\beta_0=0$, $$2*P\left( t_{(11)}>\mid 132,521 \mid \right)=0,000.$$
Para $H_0:\beta_1=0$, $$2*P\left( t_{(11)}>\mid 22,730 \mid \right)=0,000.$$
Para $H_0:\beta_2=0$, $$2*P\left( t_{(11)}>\mid -4,149 \mid \right)=0,00162.$$
Com isso, rejeitamos $H_0$ para $\beta_0$, $\beta_1$ e $\beta_2$ pois os respectivos p-valores são menores do que $\alpha$.
Exemplo 2.5.2
Construir um intervalo de confiança com $95\char37$ para o parâmetro $\beta_1$, considerando os dados do “Exemplo 2.2.3”.
Solução:
Lembramos que $\widehat{\beta} _1=323,43$, $\widehat{\sigma}^2=1.220,1$ e que $C_{11}=0,1660$. Assim, o intervalo com 95\char37 de confiança para $\beta_1$ é dado por
$$\left[\widehat{\beta} _1-t_{(0,025;11)}\sqrt{\widehat{\sigma}^2 C_{11}} ; \widehat{\beta} _1+t_{(0,025;11)}\sqrt{\widehat{\sigma}^2C_{11}}\right]$$
$$\left[323,43-2,201\sqrt{(1.220,1)(0,1660)} ; 323,43 + 2,201\sqrt{(1.220,1)(0,1660)}\right]$$
$$\left[323,43-2,201(14,23) ; 323,43+2,201(14,23)\right]$$ $$\left[292,10 ; 354,75\right]$$
Portanto, $292,10 \leq \beta_1 \leq 354,75.$
Usando o software Action Stat temos os seguintes resultados:
| G.L. | Soma de Quadrados | Quadrado Médio | Estat.t | P-valor | |
|---|---|---|---|---|---|
| x1 | 1 | 630968.32383546 | 630968.32383546 | 517.146904974 | 0 |
| x2 | 1 | 20997.846071104 | 20997.846071104 | 17.210009911 | 0.001620574 |
| Resíduos | 11 | 13421.04437915 | 1220.094943559 |
Tabela 4.2.12: Tabela da ANOVA
| Mínimo | 1Q | Mediana | Média | 3Q | Máximo |
|---|---|---|---|---|---|
| -44.58635322 | -26.348187101 | -3.263076957 | 0 | 26.003395891 | 63.203182565 |
Tabela 4.2.13: Análise Exploratória (resíduos)
| Estimativa | Desvio Padrão | Estat.t | P-valor | |
|---|---|---|---|---|
| Intercepto | 1242.313981613 | 9.374451738 | 132.521241389 | 0 |
| x1 | 323.434498746 | 14.229539255 | 22.729794195 | 0 |
| x2 | -54.773007617 | 13.203103529 | -4.148494897 | 0.001620574 |
Tabela 4.2.14: Coeficientes
| Desvio Padrão dos Resíduos | Graus de Liberdade | $R^2$ | $R^2$ Ajustado |
|---|---|---|---|
| 34.92985748 | 11 | 0.979829723 | 0.9761624 |
Tabela 4.2.15: Medida Descritiva da Qualidade do Ajuste
| Intervalo de confiança para os parâmetros | 2.5 % | 97.5 % |
|---|---|---|
| Intercepto | 1221.680952455 | 1262.947010772 |
| x1 | 292.115494011 | 354.75350348 |
| x2 | -83.832842553 | -25.713172682 |
Tabela 4.2.16: Intervalo de confiança para os parâmetros
2.6 Intervalo de Confiança para Resposta Média e Predição
2.6.1 Intervalo de confiança para a Predição
Podemos obter o intervalo de confiança para a resposta média em um ponto particular da amostra. Seja o vetor $x_0^\prime =[1, x_{01}, x_{02}, \ldots, x_{0p}]^\prime $ um ponto da amostra.
A resposta média para este ponto é
$$E(Y\mid x_0)=\mu_{Y\mid x_0}=\beta_0+\beta_1x_{01}+\beta_2x_{02}+\ldots+\beta_px_{0p}=x^\prime _0\beta.$$
O estimador da média neste ponto é
$$\widehat{\mu}_{Y\mid x_0}=\widehat{Y}_0=x^\prime _0\widehat{\beta}.$$
Notemos que o estimador da resposta média é não viciado. De fato,
$$E(\widehat{\mu}_{Y \mid x_0})= E(x^\prime_0 \widehat{\beta}) = x ^\prime _0 \beta=\mu _{Y\mid x_0} \quad \hbox{(não viciado)}.$$
Além disso, temos que a variância de $\widehat{\mu}_{Y\mid x_0}$ é dada por
$$Var(\widehat{\mu}_{Y\mid x_0})=\sigma^2 x^\prime _0(X^\prime X)^{-1}x_0.$$
Assim, temos que
$$\dfrac{\widehat{\mu} _{Y\mid x_0}-\mu_{Y\mid x_0}}{\sqrt{\sigma^2 x^\prime _0 (X^\prime X)^{-1} x_0}}\sim N(0,1) \quad \hbox{e que}$$
$$\dfrac{(n-p-1)\widehat{\sigma}^2}{\sigma^2}\sim\chi^2_{(n-p-1)} \quad \hbox{independentes entre si}.$$
Logo, segue que
$$\dfrac{\widehat{\mu}_{Y\mid x_0}-\mu_{Y\mid x_0}}{\sqrt{\widehat{\sigma}^2x^\prime _0(X^\prime X)^{-1}x_0}}\sim t_{(n-p-1)}.$$
Portanto, o intervalo com confiança de $100(1-\alpha)\char37$ para a resposta média ($\mu_{Y\mid x_0}$) no ponto $x_{0}$ é dado por
$$\left[\widehat{\mu}_{Y\mid x_0}-t_{(\frac{\alpha}{2},n-p-1)}\sqrt{\widehat{\sigma}^2 x^\prime _0(X^\prime X)^{-1}x_0} \quad ; \quad \widehat{\mu}_{Y\mid x_0}+t_{(\frac{\alpha}{2},n-p-1)}\sqrt{\widehat{\sigma}^2 x^\prime _0(X^\prime X)^{-1}x_0}\right].$$
Exemplo 2.6.1
Determinar um intervalo de confiança para a resposta média no ponto $x_{01}=225$ min e $x_{02}=4,36 \times10^{14}$ íons.
Solução:
Em termos das variáveis codificadas, temos que $x_{01}=0$ e $x_{02}=0$, em que $x_{01}$ é a variável tempo e $x_{02}$ é a variável dose. Desta forma, seja o ponto $x^\prime_0=[1,0,0].$ A estimativa da resposta média neste ponto é
$$\widehat{\mu}_{y\mid x_0}=x_0^\prime\widehat{\beta}=[1~0~0]\left[\begin{array}{c}1.242,31 \cr 323,43 \cr -54,77 \cr \end{array}\right]=1.242,31.$$
Temos também que
$\hbox{Var}(\widehat{\mu}_{y\mid x_0}) = \sigma^2 x^\prime_0 (X^\prime X)^{-1} x_0 $
$ \qquad \qquad = \sigma^2 [1~0~0] \begin{bmatrix} 0,0720 & -0,0019 & 0,0091 \cr -0,0019 & 0,1660 & 0,0004 \cr 0,0091 & 0,0004 & 0,1429 \cr \end{bmatrix} \begin{bmatrix} 1 \cr 0 \cr 0 \end{bmatrix}$
$ \qquad \qquad = \sigma^2 (0,072)$
Sabemos também que o estimador de $\sigma^2$ é $\widehat{\sigma}^2=QME =1.220,1.$
Portanto, temos que o intervalo de confiança é dado por
$$\left[\widehat{\mu}_{y\mid x_0}-t_{(\frac{\alpha}{2},n-p-1)}\sqrt{\sigma^2 x^\prime _0(X^\prime X)^{-1}x_0} ; \widehat{\mu}_{y\mid x_0}+t_{(\frac{\alpha}{2},n-p-1)}\sqrt{\sigma^2 x^\prime _0(X^\prime X)^{-1}x_0}\right]$$
$$\left[1.242,31-2,201\sqrt{1.220,1 \times 0,072} ; 1.242,31+2,201 \sqrt{1.220,1\times0,072}\right]$$
$$\left[1.242,31-20,63 ; 1.242,31+20,63\right]$$ $$\left[1.221,68 ; 1.262,94\right]$$
Portanto, $1.221,68\leq\mu_{y\mid x_0}\leq 1.262,94.$
2.6.2 Previsões para novas observações
O modelo de regressão pode ser usado para predizer futuras observações na variável resposta $Y$ correspondente a valores particulares das variáveis explicativas $x_{01},x_{02},\dots,x_{0p}$. Aqui assumimos que os valores de $x_{01},x_{02,\dots,x_{0p}}$ são não observados. Se $x^\prime _0=[1,x_{01},x_{02},\ldots,x_{0p}]$, o valor predito para a nova observação é dada por
$$\widehat{y}(x_0)=x^\prime _0\widehat{\beta}$$
Um intervalo de previsão $100(1-\alpha)\char37$ para futuras observações ($\mu_{y\mid x_0}$) é
$$\left[\widehat{y}(x_0)-t_{(\frac{\alpha}{2},n-p-1)}\sqrt{\widehat{\sigma}^2(1+x^\prime _0(X^\prime X)^{-1}x_0)} \quad ; \quad \widehat{y}(x_0)+t_{(\frac{\alpha}{2},n-p-1)}\sqrt{\widehat{\sigma}^2(1+x^\prime _0(X^\prime X)^{-1}x_0)}\right].$$
É preciso cuidado para estimar a variável resposta para pontos que estão além da região da amostra. É bem provável que um modelo ajustado dentro da região dos dados possa não estar mais ajustado quando usarmos observações fora dessa região. Em regressão múltipla é fácil extrapolar novos dados mas para tal devemos ter certeza de que os níveis das variáveis definam uma região conjunta que contenha os dados.
Como exemplo considere a Figura 4.2.3 que ilustra a região citada para um modelo de duas variáveis.
Figura 4.2.3: Um exemplo de extrapolação na regressão múltipla.
Exemplo 2.6.2
Encontrar um intervalo de predição com $95\char37$ de confiança para o ganho do transistor, no ponto de $x_{01}=225$ min e $x_{02}= 4,36 \times 10^{14}$ íons.
Solução:
Com as variáveis codificadas, estes pontos passam a ser $x_{01}=0$ e $x_{02}=0$, sendo $x_{01}$ a variável tempo e $x_{02}$ a variável dose. Assim, o valor predito para o ganho no ponto $x^\prime _0=[1,0,0]$ é dado por
$$\widehat{y}(x_0)=x_0^\prime \widehat{\beta}=[1, 0, 0] \left[\begin{array}{c}1.242,31 \cr 323,43 \cr -54,77 \cr \end{array}\right]=1.242,31.$$
Do “Exemplo 2.6.1” temos que $$x_0^\prime (X^\prime X)^{-1}x_0 = 0,072.$$
Então o intervalo de predição será
$$\left[\widehat{y}(x_0)-t_{(\frac{\alpha}{2},n-p-1)}\sqrt{\widehat{\sigma}^2(1+x^\prime _0(X^\prime X)^{-1}x_0)} ; \widehat{y}(x_0)+t_{(\frac{\alpha}{2},n-p-1)}\sqrt{\widehat{\sigma}^2(1+x^\prime _0(X^\prime X)^{-1}x_0)}\right]$$
$$\left[1.242,31-2,201 \sqrt{1.220,1(1 + 0,072)} \quad ; \quad 1.242,31+2,201 \sqrt{1.220,1 (1 + 0,072)}\right]$$
$$\left[1.242,31-79,60 \quad ; \quad 1.242,31+79,60\right]$$
$$\left[1.162,71 ; 1.321,91\right]$$
Portanto, $1.162,71 \leq \mu_{y\mid x_0} \leq 1.321,91.$
Observemos que se compararmos o intervalo de predição com o intervalo de confiança para a resposta média no mesmo ponto para os dois exemplos, poderemos constatar que o intervalo de predição é bem mais amplo. Isto reflete o fato de que existe dificuldade em predizer um valor individual futuro.
Usando o software Action temos os seguintes resultados:
| Ganho | x1 | x2 | Valor Ajustado | Limite Inferior | Limite Superior | Desvio Padrão |
|---|---|---|---|---|---|---|
| 1 | 1004 | -1 | -1 | 973.652490485 | 927.050213313 | 1020.254767657 |
| 2 | 1636 | 1 | -1 | 1620.521487976 | 1574.521246208 | 1666.521729744 |
| 3 | 852 | -1 | 0.6667 | 882.362318689 | 839.081201282 | 925.643436097 |
| 4 | 1506 | 1 | 0.6667 | 1529.231316181 | 1486.410999704 | 1572.051632657 |
| 5 | 1272 | 0 | -0.4444 | 1266.655106199 | 1243.312698635 | 1289.997513762 |
| 6 | 1270 | 0 | -0.7222 | 1281.871047715 | 1253.785960402 | 1309.956135027 |
| 7 | 1269 | 0 | 0.6667 | 1205.796817435 | 1176.258075823 | 1235.335559047 |
| 8 | 903 | -1 | -0.1667 | 928.010143238 | 890.114775471 | 965.905511004 |
| 9 | 1555 | 1 | -0.1667 | 1574.879140729 | 1537.618630694 | 1612.139650764 |
| 10 | 1260 | 0 | -1 | 1297.086989231 | 1262.983951497 | 1331.190026964 |
| 11 | 1146 | 0 | 0.9444 | 1190.58635322 | 1154.807580996 | 1226.365125444 |
| 12 | 1276 | 0 | -0.1667 | 1251.444641983 | 1230.676344134 | 1272.212939833 |
| 13 | 1225 | 0 | 1 | 1187.540973996 | 1150.427976073 | 1224.653971919 |
| 14 | 1321 | 0.1667 | -0.1667 | 1305.361172924 | 1284.039755034 | 1326.682590814 |
Tabela 4.2.17: Intervalo de predição para as observações com 95% de confiança
2.7 Seleção de Variáveis
Em modelos de regressão múltipla é necessário determinar um subconjunto de variáveis independentes que melhor explique a variável resposta, isto é, dentre todas as variáveis explicativas disponíveis, devemos encontrar um subconjunto de variáveis importantes para o modelo.
Construir um modelo que inclui apenas um subconjunto de variáveis explicativas envolve dois objetivos conflitantes:
-
Obter o máximo de informação por meio de um modelo com tantas variáveis independentes possíveis;
-
Diminuir a variância da estimativa e o custo da coleta por meio de um modelo com menor número possível de variáveis.
Desta forma, obter um equilíbrio entre esses dois compromissos é de interesse. Para isto, utilizamos uma técnica, denominada de seleção de variáveis.
Existem duas principais estratégias no processo de seleção de variáveis:
- Todos os modelos possíveis: considera todos os subconjuntos possíveis de variáveis explicativas, e considerando critérios de avaliação, seleciona o melhor deles.
- Seleção Automática: faz uma busca do melhor subconjunto de variáveis explicativas sem considerar todos os possíveis subconjuntos.
Na prática, assumimos que a correta especificação funcional das variáveis explicativas é conhecida (por exemplo, $1/x_1$, $ln~x_2$) e que não há outliers ou pontos influentes e então, aplicamos a técnica de seleção de variáveis. Entretanto, o ideal seria inicialmente,
- Identificar outliers e pontos influentes,
- Identificar eventuais colinearidade e heteroscedasticidade,
- Realizar quaisquer transformações que sejam necessárias,
e então, aplicar seleção de variáveis.
A seguir apresentamos detalhadamente as estratégias Todos os modelos possíveis e Seleção Automática.
2.7.1 Seleção Todos os Modelos Possíveis
Considere o modelo de regressão linear múltipla
$$Y=\beta_0+\beta_1 x_1+\beta_2 x_2+…+\beta_p x_p+\epsilon,$$
e suas suposições. O método de todos os modelos possíveis possibilita a análise do ajuste de todos os submodelos compostos pelos possíveis subconjuntos das p variáveis e identifica os melhores desses subconjuntos, segundo critérios de avaliação.
Exemplo 2.7.1
Considere o exemplo na “Motivação 2”.
Temos duas variáveis explicativas e por isso os modelos considerando todas as combinações possíveis são:
Modelo 1: $\hbox{Ganho}=\beta_0+\beta_2\hbox{Tempo}.$
Modelo 2: $\hbox{Ganho}=\beta_0+\beta_1\hbox{Dose}.$
Modelo 3: $\hbox{Ganho}=\beta_0+\beta_1\hbox{Dose}+\beta_2\hbox{Tempo}.$
Alguns critérios para avaliar os modelos são: $R_p^2$, $R_a^2$, $QME$ (Quadrado Médio do Erro), $C_p$ de Mallows, $AIC$, $BIC$ e $PRESS_p$. A seguir uma abordagem de cada um deles.
2.7.1.1 Coeficiente de Determinação Múltipla
Seja $R_p^2$ a notação do coeficiente de determinação múltipla de um modelo com p variáveis explicativas, isto é, p coeficientes e o intercepto $\beta_0$.
$$R^2_p= \frac{SQR}{SQT}=1-\frac{SQE}{SQT},$$
em que SQR, SQE e SQT são a soma dos quadrados do modelo, soma dos quadrados dos resíduos (erros) e soma dos quadrados total, respectivamente.
O critério utilizado nesse método é que se adicionarmos uma variável insignificante teremos um aumento mínimo (pequeno) de $R_p^2$. Assim, ele é mais usado para julgar quando parar de adicionar variáveis do que para encontrar o melhor modelo já que $R_p^2$ nunca diminui quando $p$ aumenta.
Exemplo 2.7.1.1
Para o exemplo na “Motivação 2”, temos os valores do coeficiente de determinação múltipla para cada modelo possível:
- Modelo 1:
O modelo 1, como descrevemos no “Exemplo 2.7.1”, tem apenas a variável explicativa Tempo.
Com o ajuste do modelo, as somas de quadrados são:
$SQR=630.968$, $SQE=34.419$, $SQT=665.387$.
$$R^2_p= \frac{630.968}{665.387}=1-\frac{34.419}{665.387}=0,948.$$
- Modelo 2:
O modelo 2 tem apenas como variável explicativa a Dose de íons.
As somas dos quadrados do ajuste do modelo são:
$SQR=21.612$, $SQE=643.775$, $SQT=665.387$.
$$R^2_p= \frac{21.612}{665.387}=1-\frac{643.775}{665.387}=0,032.$$
- Modelo 3:
O modelo 3 tem como variáveis explicativas o Tempo e Dose de íons.
As somas de quadrados do modelo são:
$SQR=651.966$, $SQE=13.421$, $SQT=665.387$.
$$R^2_p= \frac{651.966}{665.387}=1-\frac{13.421}{665.387}=0,9798.$$
Usando o software Action temos os seguintes resultados:
| Variáveis | R² | CP | R² Ajustado | AIC | BIC | PRESS | QME | |
|---|---|---|---|---|---|---|---|---|
| Modelo 1 | Tempo | 0.948 | 18.21 | 0.944 | 155.033 | 156.95 | 51545.126 | 2868.27911315993 |
| Modelo 2 | Dose de íons | 0.032 | 517.641 | -0.048 | 196.035 | 197.952 | 875529.652 | 53647.9286674226 |
| Modelo 3 | Tempo + Dose de íons | 0.98 | 3 | 0.976 | 143.848 | 146.404 | 22225.01 | 1220.101393825 |
Tabela 4.2.18: Tabela da seleção de modelos e métricas
2.7.1.2 Coeficiente de Determinação Ajustado
Para evitar dificuldades na interpretação de $R^2$, alguns estatísticos preferem usar o $R_a^2$ ($R^2$ ajustado) , definido para uma equação com $p+1$ coeficientes como
$$R^2_a=1-\left(\frac{n-1}{n-(p+1)}\right)(1-R^2_p).$$
O $R_a^2$ não necessariamente aumenta com a adição de parâmetros no modelo. Na verdade se $s$ variáveis explicativas são incluídas no modelo (modelo com $p+s$ variáveis), o $R_a^2$ desse modelo excederá $R_a^2$ do modelo com p variáveis apenas se a estatística parcial F para testar a significância dos adicionais $s$ coeficientes passar de 1 (para mais detalhes ver Seber[1977]). Consequentemente, um critério para a seleção de um modelo ótimo é escolher o modelo que tem o $R_a^2$ máximo.
Exemplo 2.7.1.2
Para o exemplo na “Motivação 2” temos os seguintes coeficientes de determinação ajustado para cada um dos modelos possíveis, considerando os coeficientes de determinação múltipla no “Exemplo 2.7.1.1”.
- Modelo 1
$$R^2_a=1-\left(\frac{14-1}{14-2}\right)(1-0,948)=0,944.$$
- Modelo 2
$$R^2_a=1-\left(\frac{14-1}{14-2}\right)(1-0,032)=-0,04.$$
- Modelo 3
$$R^2_a=1-\left(\frac{14-1}{14-3}\right)(1-0,9798)=0,976.$$
Usando o software Action temos os seguintes resultados:
| Variáveis | R² Ajustado | CP | R² | AIC | BIC | PRESS | QME | |
|---|---|---|---|---|---|---|---|---|
| Modelo 1 | Tempo | 0.944 | 18.21 | 0.948 | 155.033 | 156.95 | 51545.126 | 2868.27911315993 |
| Modelo 2 | Dose de íons | -0.048 | 517.641 | 0.032 | 196.035 | 197.952 | 875529.652 | 53647.9286674226 |
| Modelo 3 | Tempo + Dose de íons | 0.976 | 3 | 0.98 | 143.848 | 146.404 | 22225.01 | 1220.101393825 |
Tabela 4.2.19: Tabela da seleção de modelos e métricas
2.7.1.3 Quadrado Médio dos Resíduos
O quadrado médio dos resíduos de um modelo de regressão é obtido por meio de
$$QME=SQE/(n-p-1), \quad \quad \quad \quad \quad \quad ~(2.7.1.3)$$
em que SQE é a soma dos quadrados dos resíduos. O QME sempre decresce conforme p aumenta. O quadrado médio do erro inicialmente decresce, estabiliza e eventualmente pode aumentar. Esse eventual aumento ocorre quando a redução do QME em adicionar um coeficiente para o modelo não é suficiente para compensar a perda nos graus de liberdade do denominador de (2.7.1.3).
O modelo que minimiza QME também maximizará $R_a^2$. Para entender isso, notamos que
$$R_a^2=1-\left(\frac{n-1}{n-p}\right)(1-R_p^2)$$
$$=1-\left(\frac{n-1}{n-p}\right)\left(\frac{SQE}{SQT}\right)$$
$$=1-\left(\frac{n-1}{SQT}\right)\left(\frac{SQE}{n-p}\right)$$
$$=1-\left(\frac{n-1}{SQT}\right)QME.$$
Assim, minimizar QME e maximar $R_a^2$ são equivalentes.
Exemplo 2.7.1.3
Utilizando o exemplo na “Motivação 2”, temos os quadrados médios dos erros de cada modelo possível:
- Modelo 1
$$QME=34.419/(14-1-1)=2.868.$$
- Modelo 2
$$QME=643.775/(14-1-1)=53.648.$$
- Modelo 3
$$QME=13.421/(14-2-1)=1.220.$$
Usando o Software Action temos os seguintes resultados:
| Variáveis | QME | CP | R² | R² Ajustado | AIC | BIC | PRESS | |
|---|---|---|---|---|---|---|---|---|
| Modelo 1 | Tempo | 2868.2791 | 18.21 | 0.948 | 0.944 | 155.033 | 156.95 | 51545.126 |
| Modelo 2 | Dose de íons | 53647.9287 | 517.641 | 0.032 | -0.048 | 196.035 | 197.952 | 875529.652 |
| Modelo 3 | Tempo + Dose de íons | 1220.1014 | 3 | 0.98 | 0.976 | 143.848 | 146.404 | 22225.01 |
Tabela 4.2.20: Tabela da seleção de modelos e métricas
2.7.1.4 Cp de Mallows
O critério $C_p$ de Mallows é baseado no conceito do erro quadrático médio (EQM) dos valores ajustados. O erro quadrático médio da previsão é
$$EQM=E(\hat{y}_i-E(y_i))^2=E(\hat{y}_i-E(\hat{y}_i)+E(E(\hat{y}_i)-E(y_i)))^2$$
$$=E(\hat{y}_i-E(\hat{y}_i))^2+(E(\hat{y}_i)-E(y_i))^2=\hbox{Var}(\hat{y}_i)+\hbox{vício}^2(\hat{y}_i),$$
em que $E(\hat{y}_i)-E(y_i)$ é o vicio. Assim, o EQM é a soma da variância de $\hat{y}_i$ e o vício ao quadrado. O EQM considerando os n valores amostrais é
$$\sum{E(\hat{y}_i-y_i)^2}= \sum{Var(\hat{y}_i)}+\sum{(E(\hat{y}_i)-E(y_i))^2}.$$
O critério $\Gamma_p$ é o erro quadrático médio dividido pela variância dos erros $\sigma^2$.
$$\Gamma_p=\left(\frac{1}{\sigma^2}\right)[\sum{Var(\hat{y}_i)}+\sum{(E(\hat{y}_i)-E(y_i))^2}], \quad \quad \quad \quad \quad ~(2.7.1.4)$$
em que $\sum{Var(\hat{y}_i)}=(p+1)\sigma^2$ e o valor esperado da soma dos quadrados dos erros é: $$E(SQE)=(E(\hat{y}_i)-E(y_i))^2+(n-(p+1))\sigma^2.$$
Substituindo esses valores em (2.7.1.4), obtemos
$$\Gamma_p=\left(\frac{1}{\sigma^2}\right)[E(SQE)-(n-(p+1))\sigma^2+(p+1)\sigma^2]$$
$$=\frac{E(SQE)}{\sigma^2}-n+2(p+1).$$
Como $\sigma^2$ é desconhecido, assumindo que o modelo que inclui todas as variáveis explicativas é tal que o QME é um estimador não viciado de $\sigma^2$ e substituindo E(SQE) pelo valor observado SQE, $\Gamma_p$ pode ser estimado por
$$C_p=\frac{SQE(p)}{QME}-n+2(p+1),$$
em que SQE(p) é a soma de quadrados dos erros do submodelo e QME é o quadrado médio do modelo com todas as variáveis explicativas.
Pode também ser mostrado que quando não há vício na estimativa do modelo com as p variáveis, $E(SQE)=(n-(p+1))\sigma^2$ e então,
$$E[C_p|Vicio=0]=\frac{(n-(p+1))\sigma^2}{\sigma^2}-n+2(p+1)=p+1,$$
em que $p+1$ é o número de parâmetros no modelo já que p é o número de variáveis explicativas mais o intercepto.
A estratégia usada para selecionar modelos com o critério $C_p$ é identificar modelos com $C_p$ próximo do número de parâmetros $(p+1)$.
Exemplo 2.7.1.4
Vamos calcular o $C_p$ para os todos os modelos possíveis do exemplo na “Motivação 2”.
- Modelo 1
A soma de quadrados dos erros do modelo apenas com Tempo como variável explicativa é 34.419. O quadrado médio do modelo completo é 1.220. Assim, o $C_p$ é dado por:
$$C_p=\frac{34.419}{1.220}-14+2\times 2=18,21.$$
- Modelo 2
A soma de quadrados dos erros do modelo apenas com Dose de íons como variável explicativa é 643.775 e o quadrado médio do modelo completo é 1.220. O $C_p$ do modelo 2 é:
$$C_p=\frac{643.775}{1.220}-14+2\times 2=517,6.$$
- Modelo 3
A soma de quadrados dos erros do modelo completo é 13.421 e o quadrado médio é 1.220.
$$C_p=\frac{13.421}{1.220}-14+2\times 3=3.$$
Usando o Software Action temos os seguintes resultados:
| Variáveis | CP | R² | R² Ajustado | AIC | BIC | PRESS | QME | |
|---|---|---|---|---|---|---|---|---|
| Modelo 1 | Tempo | 18.21 | 0.948 | 0.944 | 155.033 | 156.95 | 51545.126 | 2868.27911315993 |
| Modelo 2 | Dose de íons | 517.641 | 0.032 | -0.048 | 196.035 | 197.952 | 875529.652 | 53647.9286674226 |
| Modelo 3 | Tempo + Dose de íons | 3 | 0.98 | 0.976 | 143.848 | 146.404 | 22225.01 | 1220.101393825 |
Tabela 4.2.21: Tabela da seleção de modelos e métricas
2.7.1.5 AIC e BIC
O Critério de Informação de Akaike (AIC) é definido como
$$AIC_p=-2log(L_p)+2[(p+1)+1],$$
em que $L_p$ é a função de máxima verossimilhança do modelo e p é o número de variáveis explicativas consideradas no modelo.
O Critério de Informação Bayesiano (BIC) é definido como
$$BIC_p=-2log(L_p)+[(p+1)+1]log(n).$$
Tanto o AIC quanto o BIC aumentam conforme SQE aumenta. Além disso, ambos critérios penalizam modelos com muitas variáveis sendo que valores menores de AIC e BIC são preferíveis.
Como modelos com mais variáveis tendem a produzir menor SQE mas usam mais parâmetros, a melhor escolha é balancear o ajuste com a quantidade de variáveis.
Exemplo 2.7.1.5
Vamos calcular o AIC e BIC para os todos os modelos possíveis do exemplo na “Motivação 2”.
- Modelo 1
O $log(L_p)$ do modelo considerando apenas Tempo como variável explicativa é -74,51. Assim, o AIC e BIC é dado por
$$AIC_p=-2(-74,51)+2(3)=155,02.$$
$$BIC_p=-2(-74,51)+(3)log(14)=156,94.$$
- Modelo 2
O $log(L_p)$ do modelo considerando apenas Dose de íons como variável explicativa é -95,01. Assim, o AIC e BIC é dado por
$$AIC_p=-2(-95,01)+2(3)=196,02.$$
$$BIC_p=-2(-95,01)+(3)log(14)=197,94.$$
- Modelo 3
O $log(L_p)$ do modelo considerando Tempo e Dose de íons como variáveis explicativas é -67,92. Assim, o AIC e BIC é dado por
$$AIC_p=-2(-67,92)+2(4)=143,84.$$
$$BIC_p=-2(-67,92)+(4)log(14)=146,39.$$
Usando o software Action temos os seguintes resultados:
| Variáveis | CP | R² | R² Ajustado | AIC | BIC | PRESS | QME | |
|---|---|---|---|---|---|---|---|---|
| Modelo 1 | Tempo | 18.21 | 0.948 | 0.944 | 155.033 | 156.95 | 51545.126 | 2868.27911315993 |
| Modelo 2 | Dose de íons | 517.641 | 0.032 | -0.048 | 196.035 | 197.952 | 875529.652 | 53647.9286674226 |
| Modelo 3 | Tempo + Dose de íons | 3 | 0.98 | 0.976 | 143.848 | 146.404 | 22225.01 | 1220.101393825 |
Tabela 4.2.22: Tabela da seleção de modelos e métricas
2.7.1.6 Critério PRESS
O critério $PRESS_p$ (Prediction Error Sum of Squares) é definido por
$$PRESS_p=\sum_{i=1}^{n}(Y_i-\hat{Y}_{(i)})^2,$$
em que $\hat{Y}_{(i)}$ é o valor predito da regressão sem a i-ésima observação.
Podemos escrever o PRESS como
$$PRESS_p=\sum_{i=1}^{n}\left(\frac{Y_i-\hat{Y}_i}{1-h_{ii}}\right)^2,$$
em que o $h_{ii}$ é o i-ésimo valor da diagonal da matriz H.
No uso desse critério, modelo com menor $PRESS_p$ é preferível.
Exemplo 2.7.1.6
No exemplo na “Motivação 2” vamos calcular o critério PRESS para todos os modelos possíveis.
- Modelo 1 (apenas Tempo como variável explicativa):
Na Tabela 2.7.1.6.1 temos os valores da variável resposta Ganho, seus valores ajustados através do modelo 1 e também a diagonal da matriz H.
| Ganho (Y) | Valores Ajustados $(\hat{Y})$ | $h_{ii}$ |
|---|---|---|
| 1004 | 922,2 | 0,241 |
| 1636 | 1569,4 | 0,233 |
| 852 | 922,2 | 0,241 |
| 1506 | 1569,4 | 0,233 |
| 1272 | 1245,8 | 0,071 |
| 1270 | 1245,8 | 0,071 |
| 1269 | 1245,8 | 0,071 |
| 903 | 922,2 | 0,24 |
| 1555 | 1569,4 | 0,23 |
| 1260 | 1245,8 | 0,071 |
| 1146 | 1245,8 | 0,071 |
| 1276 | 1245,8 | 0,071 |
| 1225 | 1245,8 | 0,071 |
| 1321 | 1299,7 | 0,075 |
Tabela 4.2.23: Medidas para o cálculo do PRESS do modelo 1
$$PRESS_p=\sum_{i=1}^{14}\left(\frac{Y_i-\hat{Y}_i}{1-h_{ii}}\right)^2=51.545.$$
- Modelo 2 (apenas Dose de íons como variável explicativa):
Na Tabela 2.7.1.6.2 temos os valores da variável resposta Ganho, seus valores ajustados através do modelo 2 e também a diagonal da matriz H.
| Ganho (Y) | Valores Ajustados $(\hat{Y})$ | $h_{ii}$ |
|---|---|---|
| 1004 | 1301,683 | 0,197 |
| 1636 | 1301,683 | 0,197 |
| 852 | 1209,069 | 0,147 |
| 1506 | 1209,069 | 0,147 |
| 1272 | 1270,812 | 0,092 |
| 1270 | 1286,247 | 0,133 |
| 1269 | 1209,069 | 0,147 |
| 903 | 1255,376 | 0,073 |
| 1555 | 1255,376 | 0,073 |
| 1260 | 1301,683 | 0,197 |
| 1146 | 1193,634 | 0,216 |
| 1276 | 1255,376 | 0,073 |
| 1225 | 1190,546 | 0,233 |
| 1321 | 1255,376 | 0,073 |
Tabela 4.2.24: Medidas para o cálculo do PRESS do Modelo 2
$$PRESS_p=\sum_{i=1}^{14}\left(\frac{Y_i-\hat{Y}_i}{1-h_{ii}}\right)^2=875.529,6.$$
- Modelo 3 (modelo com Tempo e Dose de íons como variáveis explicativas):
Na Tabela 2.7.1.6.3 temos os valores da variável resposta Ganho, seus valores ajustados através do modelo 3 e também a diagonal da matriz H.
| Ganho (Y) | Valores Ajustados $(\hat{Y})$ | $h_{ii}$ |
|---|---|---|
| 1004 | 973,6 | 0,367 |
| 1636 | 1620,5 | 0,358 |
| 852 | 882,4 | 0,317 |
| 1506 | 1529,2 | 0,31 |
| 1272 | 1266,6 | 0,092 |
| 1270 | 1281,9 | 0,133 |
| 1269 | 1205,8 | 0,148 |
| 903 | 928 | 0,243 |
| 1555 | 1574,9 | 0,235 |
| 1260 | 1297,1 | 0,197 |
| 1146 | 1190,6 | 0,216 |
| 1276 | 1251,4 | 0,073 |
| 1225 | 1187,5 | 0,233 |
| 1321 | 1305,3 | 0,077 |
Tabela 2.4.25: Medidas para o cálculo do PRESS do Modelo 3
$$PRESS_p=\sum_{i=1}^{14}\left(\frac{Y_i-\hat{Y}_i}{1-h_{ii}}\right)^2=22.225.$$
Usando o software Action temos os seguintes resultados:
| Variáveis | PRESS | CP | R² | R² Ajustado | AIC | BIC | QME | |
|---|---|---|---|---|---|---|---|---|
| Modelo 1 | Tempo | 51545.126 | 18.21 | 0.948 | 0.944 | 155.033 | 156.95 | 2868.27911315993 |
| Modelo 2 | Dose de íons | 875529.652 | 517.641 | 0.032 | -0.048 | 196.035 | 197.952 | 53647.9286674226 |
| Modelo 3 | Tempo + Dose de íons | 22225.01 | 3 | 0.98 | 0.976 | 143.848 | 146.404 | 1220.101393825 |
Tabela 4.2.26: Tabela da seleção de modelos e métricas
2.7.2 Seleção Automática
Como a seleção de todas as regressões possíveis necessita de um considerável esforço computacional, outros métodos foram desenvolvidos para selecionar o melhor subconjunto de variáveis sequencialmente, adicionando ou removendo variáveis em cada passo.
O critério para a adição ou remoção de covariáveis é geralmente baseado na estatística F, comparando modelos com e sem as variáveis em questão. O AIC, assim como outros critérios, também podem ser utilizados na decisão de inserir e remover variáveis. Existem três procedimentos automáticos: (1) Método Forward, (2) Método Backward e (3) Método Stepwise.
2.7.2.1 Seleção Forward
Esse procedimento parte da suposição de que não há variável no modelo, apenas o intercepto. A ideia do método é adicionar uma variável de cada vez. A primeira variável selecionada é aquela com maior correlação com a resposta.
Procedimento:
- Ajustamos o modelo com a variável com maior correlação amostral com a variável resposta.
Supondo que essa variável seja $x_1$, calculamos a estatística F para testar se ela realmente é significativa para o modelo. A variável entra no modelo se a estatística F for maior do que o ponto crítico, chamado de $F_{in}$ ou F para entrada. Notemos que $F_{in}$ é calculado para um dado $\alpha$ crítico.
- Considerando que $x_1$ foi selecionado para o modelo, o próximo passo é encontrar uma variável com maior correlação com a resposta considerando a presença da primeira variável no modelo. Esta é chamada de correlação parcial e é a correlação dos resíduos do modelo $\hat{y}=\hat{\beta_0}+\hat{\beta_1}x_1$ com os resíduos do modelo $\hat{x_j}=\hat{\alpha_{0j}}+\hat{\alpha_{1j}}x_1$, j=2,3,…,p.
Vamos supor que a maior correlação parcial com y seja $x_2$. Isso implica que a maior estatística $F$ parcial é:
$$F=\frac{SQR(x_2|x_1)}{QME(x_1,x_2)}.$$
Se o valor da estatística é maior do que $F_{in}$, $x_2$ é selecionado para o modelo.
- O processo é repetido, ou seja, variável com maior correlação parcial com y é adicionada no modelo se sua estatística $F$ parcial for maior que $F_{in}$, até que não seja incluída mais nenhuma variável explicativa no modelo.
2.7.2.2 Seleção Backward
Enquanto o método Forward começa sem nenhuma variável no modelo e adiciona variáveis a cada passo, o método Backward faz o caminho oposto; incorpora inicialmente todas as variáveis e depois, por etapas, cada uma pode ser ou não eliminada.
A decisão de retirada da variável é tomada baseando-se em testes F parciais, que são calculados para cada variável como se ela fosse a última a entrar no modelo.
Procedimento:
- Para cada variável explicativa calcula-se a estatística $F$. Para a variável $x_k$, por exemplo,
$$F=\frac{SQR(x_k|x_1,…,x_{k-1})}{QME}.$$
O menor valor das estatísticas $F$ parciais calculadas é então comparado com o $F$ crítico, $F_{out}$, calculado para um dado valor $\alpha$ crítico. Se o menor valor encontrado for menor do que $F_{out}$, elimina-se do modelo a covariável responsável pelo menor valor da estatística $F$ parcial.
-
Ajusta-se novamente o modelo, agora com as $p-1$ variáveis. As estatísticas $F$ parciais são calculadas para esse modelo e o processo é repetido.
-
O algoritmo de eliminação termina quando a menor estatística $F$ parcial não for menor do que $F_{out}$.
2.7.2.3 Seleção Stepwise
Stepwise é uma modificação da seleção Forward em que cada passo todas as variáveis do modelo são previamente verificadas pelas suas estatísticas $F$ parciais. Uma variável adicionada no modelo no passo anterior pode ser redundante para o modelo por causa do seu relacionamento com as outras variáveis e se sua estatística $F$ parcial for menor que $F_{out}$, ela é removida do modelo.
Procedimento:
- Iniciamos com uma variável: aquela que tiver maior correlação com a variável resposta.
- A cada passo do forward, depois de incluir uma variável, aplica-se o backward para ver se será descartada alguma variável.
- Continuamos o processo até não incluir ou excluir nenhuma variável.
Assim, a regressão Stepwise requer dois valores de corte: $F_{in}$ e $F_{out}$. Alguns autores preferem escolher $F_{in}=F_{out}$ mas isso não é necessário. Se $F_{in}< F_{out}$: mais difícil remover que adicionar; se $F_{in}> F_{out}$: mais difícil adicionar que remover.
2.7.3 Algumas Considerações
A seleção de variáveis é um meio para se chegar a um modelo, mas não é a etapa final. O objetivo é construir um modelo que seja bom para obter predições ou que explique bem o relacionamento entre os dados.
Os métodos de seleção automática têm a vantagem de não necessitar de grande esforço computacional. Mas eles não indicam o melhor modelo respeitando algum critério (não retorna um conjunto de modelos em que o pesquisador tem o poder de escolha).
Já o método de todos os modelos possíveis identifica modelos que são melhores respeitando o critério que o pesquisador quiser.
É indicado, então, usar métodos passo a passo combinados com outros critérios.
Se por acaso existe um grande número de variáveis, é recomendado usar métodos de seleção automática para eliminar aquelas com efeitos insignificantes. E o conjunto reduzido de variáveis pode então ser investigado pelo método de todos os modelos possíveis.
A escolha do modelo final não é uma tarefa fácil. Além dos critérios formais, devemos responder às seguintes questões:
- O modelo faz sentido?
- O modelo é útil para o objetivo pretendido? Se, por exemplo, o custo da coleta dos dados de uma variável é exorbitante e impossível de ser obtido, isso resultará em um modelo sem utilidade.
- Todos os coeficientes são razoáveis, ou seja, os sinais e magnitude dos valores fazem sentido e os erros padrões são relativamente pequeno?
- A adequabilidade do modelo é satisfatória? Sem outliers, tem variância constante, normalidade e os dados são independentes?
Um princípio a ser levado em consideração é o “princípio da parcimônia”: modelos mais simples devem ser escolhidos aos mais complexos, desde que a qualidade do ajuste seja similar.
2.8 Análise de Resíduos na Regressão Linear Múltipla
Assim como na Regressão Linear Simples, as suposições do modelo descritas na Seção 2.1 precisam ser verificadas para que os resultados obtidos realmente sejam válidos. Para isso, é realizada a Análise dos Resíduos em que são avaliadas um conjunto de técnicas para investigar a adequabilidade do modelo. Na Seção 3 estão descritas algumas dessas técnicas.
Além disso, quando estamos lidando com mais de uma variável regressora, precisamos diagnosticar colinearidade ou multicolinearidade entre elas. Por isso, na Seção 3.6 temos técnicas para detectar variáveis correlacionadas que podem fazer com que os resultados sejam pouco confiáveis.