4.3. Análise de Resíduos
Tanto na Regressão Linear Simples quanto na Regressão Múltipla, as suposições do modelo ajustado precisam ser validadas para que os resultados sejam confiáveis. Chamamos de Análise dos Resíduos um conjunto de técnicas utilizadas para investigar a adequabilidade de um modelo de regressão com base nos resíduos. Como visto anteriormente, o resíduo $(e_i)$ é dado pela diferença entre a variável resposta observada $(Y_i)$ e a variável resposta estimada $(\widehat{Y}_i)$, isto é
$$ e_i=Y_i-\widehat{Y_i}=Y_i-\widehat{\beta_0}-\widehat{\beta_1}x_{1i}-\dots-\widehat{\beta_p} x_{pi}\quad i=1,\dots,n.$$
A ideia básica da análise dos resíduos é que, se o modelo for apropriado, os resíduos devem refletir as propriedades impostas pelo termo de erro do modelo. Tais suposições são
$$Y = X \beta + \varepsilon,$$
em que $\varepsilon = (\varepsilon_1, \varepsilon_2, \cdots,\varepsilon_n)^\prime,$ com
i. $\varepsilon_i$ e $\varepsilon_j$ são independentes $(i\neq j)$;
ii. $Var(\varepsilon_i) = \sigma^2$ (constante);
iii. $\varepsilon_i \sim N(0,\sigma^2)$ (normalidade);
iv. Modelo é linear;
v. Não existir outliers (pontos atípicos) influentes.
Na Regressão Múltipla, além das suposições listadas acima, precisamos diagnosticar colinearidade e multicolinearidade entre as variáveis de entrada para que a relação existente entre elas não interfira nos resultados, causando inferências errôneas ou pouco confiáveis.
As técnicas utilizadas para verificar as suposições descritas acima podem ser informais (como gráficos) ou formais (como testes). As técnicas gráficas, por serem visuais, podem ser subjetivas e por isso técnicas formais são mais indicadas para a tomada de decisão. O ideal é combinar as técnicas disponíveis, tanto formais quanto informais, para o diagnóstico de problemas nas suposições do modelo.
Algumas técnicas gráficas para análise dos resíduos são:
-
Gráfico dos resíduos versus valores ajustados: verifica a homoscedasticidade do modelo, isto é, σ² constante.
-
Gráfico dos resíduos versus a ordem de coleta dos dados: avaliar a hipótese de independência dos dados.
-
Papel de probabilidade normal: verificar a normalidade dos dados. Ver detalhes em Papel de Probabilidade no conteúdo de Inferência.
-
Gráfico dos Resíduos Studentizados versus valores ajustados: verifica se existem outliers em Y.
-
Gráfico dos Resíduos Padronizados versus valores ajustados: verifica se existem outliers em Y.
-
Gráfico do Leverage (Diagonal da Matriz H): verifica se existem outliers em X.
Para a análise formal dos resíduos, podemos realizar os seguintes testes:
-
Testes de Normalidade em que detalhes estão contidos no conteúdo de Inferência.
-
Teste de Durbin-Watson para testar independência dos resíduos.
-
Teste de Breusch-Pagan e Goldfeld-Quandt para testar se os resíduos são homoscedásticos.
-
Teste de falta de ajuste para verificar se o modelo ajustado é realmente linear.
Para cada suposição do modelo, vamos mostrar detalhadamente as técnicas para diagnóstico.
3.1 Diagnóstico de Normalidade
A normalidade dos resíduos é uma suposição essencial para que os resultados do ajuste do modelo de regressão linear sejam confiáveis. Podemos verificar essa suposição por meio do gráfico de Papel de Probabilidade e por meio de testes tais como Shapiro-Wilk, Anderson-Darling e Kolmogorov-Smirnov. Para maiores detalhes, ver Testes de Normalidade no conteúdo de Inferência.
Exemplo 3.1.1
- Motivação 1:
Considerando o ajuste do modelo linear simples para os dados do exemplo na “Motivação 1”, vamos fazer o gráfico de Papel de Probabilidade e os testes de Shapiro-Wilk, Anderson-Darling e Kolmogorov-Smirnov para testar a normalidade dos resíduos.
Solução:
- Papel de Probabilidade
O vetor com os resíduos ordenados em ordem crescente, considerando o ajuste do modelo linear, é dado por
$$e=(-2,82; -2,66; -2,14; -1,14; -0,82; -0,82; -0,14; -0,14; -0,14; 0,02; 0,34; 0,34; 0,34; 0,34;$$ $$\quad 1,02; 1,02; 1,02; 1,18; 2,18; 3,02).$$
O vetor com os valores de $\Phi^{-1}(d_i)^\prime s$ é dado por
$$\Phi^{-1}(d)=(-1,82; -1,38; -1,12; -0,91; -0,74; -0,58; -0,44; -0,31; -0,19; -0,06; 0,06; 0,19; 0,31;$$ $$\quad 0,44; 0,58; 0,74; 0,91; 1,12; 1,38; 1,82),$$
em que $d_i=(i-0,3)/(n+0,4)$ para $i=1,\dots,n$ e $\Phi^{-1}(d_{(i)})$ é o quantil da distribuição normal padrão calculado em $d_{(i)}$. Neste exemplo, $n=20$. Assim, desenhando os pontos $(e_i,\Phi^{-1}(d_i))$, $i=1,\dots,20$, obtemos o gráfico de Papel de Probabilidade. Se a suposição de normalidade for adequada, esperamos um comportamento linear dos pontos.
Usando o software Action temos o seguinte resultado:
Figura 4.3.1: Gráfico de Papel de Probabilidade para os resíduos do modelo linear simples ajustado - Motivação 1.
Como os pontos seguem o comportamento da reta (não estão distantes dela), temos indícios de que os erros são normalmente distribuídos.
Testes de Normalidade
Em relação aos Testes de Normalidade, precisamos encontrar os valores de $F_{n}(e_i)$ e $F(e_i),$ $i=1,\dots,n.$ $F_{n}(e_{(i)})$ representa a função de distribuição acumulada empírica dos dados e é dada pela razão entre a posição $i$ e o valor total de resíduos, isto é,
$$F_{n}(e_{(i)})=i/n \hbox{(distribuição empírica acumulada)}.$$
$F(e_{(i)})$ representa a função de distribuição acumulada assumida para os dados e seu valor é obtido da tabela da distribuição normal padrão após transformarmos os dados pela relação
$$Z_{(i)}=\dfrac{e_{(i)}-\overline{e}}{s},$$
em que $\overline{e}$ é a média aritmética e s é o desvio padrão dos resíduos. Na Tabela 4.3.1 apresentamos os valores de $F_{n}(e_{(i)})$ e $F(e_{(i)})$ para os resíduos analisados.
| $e$ | $F(e)$ | $F_n(e)$ | $F_n(e)-F(e)$ | $F(e)-F_n(e-1)$ | $ln(F(e))$ | $ln(1-F(e))$ |
|---|---|---|---|---|---|---|
| -2,82 | 0,028 | 0,05 | 0,022 | 0,028 | -3,587 | -0,028 |
| -2,66 | 0,035 | 0,1 | 0,065 | -0,015 | -3,342 | -0,036 |
| -2,14 | 0,073 | 0,15 | 0,077 | -0,027 | -2,618 | -0,076 |
| -1,14 | 0,219 | 0,2 | -0,019 | 0,069 | -1,517 | -0,248 |
| -0,82 | 0,289 | 0,25 | -0,039 | 0,089 | -1,242 | -0,341 |
| -0,82 | 0,289 | 0,3 | 0,011 | 0,039 | -1,242 | -0,341 |
| -0,14 | 0,462 | 0,35 | -0,112 | 0,162 | -0,772 | -0,620 |
| -0,14 | 0,462 | 0,4 | -0,062 | 0,112 | -0,772 | -0,620 |
| -0,14 | 0,462 | 0,45 | -0,012 | 0,062 | -0,772 | -0,620 |
| 0,02 | 0,505 | 0,5 | -0,005 | 0,055 | -0,682 | -0,704 |
| 0,34 | 0,591 | 0,55 | -0,041 | 0,091 | -0,525 | -0,895 |
| 0,34 | 0,591 | 0,6 | 0,009 | 0,041 | -0,525 | -0,895 |
| 0,34 | 0,591 | 0,65 | 0,059 | -0,009 | -0,525 | -0,895 |
| 0,34 | 0,591 | 0,7 | 0,109 | -0,059 | -0,525 | -0,895 |
| 1,02 | 0,756 | 0,75 | -0,006 | 0,056 | -0,280 | -1,410 |
| 1,02 | 0,756 | 0,8 | 0,044 | 0,006 | -0,280 | -1,410 |
| 1,02 | 0,756 | 0,85 | 0,094 | -0,044 | -0,280 | -1,410 |
| 1,18 | 0,789 | 0,9 | 0,111 | -0,061 | -0,237 | -1,554 |
| 2,18 | 0,931 | 0,95 | 0,019 | 0,031 | -0,072 | -2,670 |
| 3,02 | 0,980 | 1 | 0,020 | 0,030 | -0,020 | -3,90 |
Tabela 4.3.1: Valores para cálculo das estatísticas de teste de Normalidade no exemplo da “Motivação1”.
- Kolmogorov-Smirnov
A estatística de teste para o teste de Kolmogorov-Smirnov é dada por
$$D_n=\max(0,11;0,162)=0,162.$$
Considerando $\alpha = 0,05$ e $n = 20$, encontramos pela tabela de valores críticos que o valor crítico neste caso é de aproximadamente $0,29$. Como $Dn = 0,162 < 0,29$, não temos evidências para rejeitar a hipótese de normalidade dos resíduos.
- Anderson-Darling
Utilizando a fórmula $(\star)$ no teste de Anderson-Darling no conteúdo de Inferência, temos que $$D=\sum_{i=1}^n[(2i-1)\ln(F(e_i))+(2(n-i)+1)\ln(1-F(e_i))]=-408,11.$$
Assim,
$$\displaystyle A^2=-n-\frac{D}{n}=-20-\frac{(-408,11)}{20}=0,405.$$
Para o cálculo do p-valor, precisamos encontrar a estatística de Anderson Darling modificada. Considerando μ e σ desconhecidos, temos que
$$A_m^2=A^2\times(1+(0,75/n)+(2,25/n^2))=0,405\times(1+(0,75/20)+(2,25/20^2))=0,405\times 1,043=0,422.$$
Desta forma, obtemos o p-valor aproximado analisamos a Tabela com quantis e valores da estatística de Anderson Darling. Como $A_m^2=0,422~<~0,56$ temos que p_valor é maior do que 15%. No R um valor aproximado do p-valor é encontrado por meio de interpolação.
Como p_valor é maior do que 15%, existe forte evidência pelo teste de Anderson-Darling de que os resíduos são normalmente distribuídos.
- Shapiro-Wilk
Para o cálculo da estatística teste de Shapiro-Wilk, precisamos dos valores contidos na Tabela 4.3.2, conforme apresentamos no teste de Shapiro-Wilk no conteúdo de Inferência.
| $i$ | $n-i+1$ | $a$ | $e_{n-i+1}$ | $e_i$ | $a(e_{n-i+1}-e_i)$ |
|---|---|---|---|---|---|
| 1 | 20 | 0,4734 | 3,02 | -2,82 | 2,76 |
| 2 | 19 | 0,3211 | 2,18 | -2,66 | 1,55 |
| 3 | 18 | 0,2565 | 1,18 | -2,14 | 0,85 |
| 4 | 17 | 0,2085 | 1,02 | -1,14 | 0,45 |
| 5 | 16 | 0,1686 | 1,02 | -0,82 | 0,31 |
| 6 | 15 | 0,1334 | 1,02 | -0,82 | 0,25 |
| 7 | 14 | 0,1013 | 0,34 | -0,14 | 0,05 |
| 8 | 13 | 0,0711 | 0,34 | -0,14 | 0,03 |
| 9 | 12 | 0,0422 | 0,34 | -0,14 | 0,02 |
| 10 | 11 | 0,014 | 0,34 | 0,02 |
Tabela 4.3.2: Medidas para o cálculo da estatística de Shapiro-Wilk
A estatística de teste é dada por
$$W=\dfrac{b^2}{\sum\limits_{i=1}^n(e_i-\bar{e})^2},$$
em que
$$b=\sum_{i=1}^{n/2}a_{n-i+1}\times (e_{n-i+1}-e_i).$$
Assim, segue que
$$W=\dfrac{(6,2839)^2}{41,16}=\dfrac{39,4872}{41,16}=0,959.$$
Como $W_{calc}=0,959 > W_{(0,05;20)}=0,905$, em que $W_{(0,05;20)}$ é obtido na tabela de valores críticos, dizemos que os resíduos são normalmente distribuídos com nível de significância de 5%.
Usando o software Action temos os seguintes resultados:
| Estatística | P-valor | |
|---|---|---|
| Anderson-Darling | 0.4053853 | 0.3202935 |
| Shapiro-Wilk | 0.9594731 | 0.5333777 |
| Kolmogorov-Smirnov | 0.1621102 | 0.1827284 |
Tabela 4.3.3: Testes de Normalidade
- Motivação 2:
Considerando agora os dados na “Motivação 2”, verificamos se os resíduos obtidos pelo ajuste do modelo de regressão linear múltipla segue distribuição normal utilizando o gráfico de Papel de Probabilidade e os testes de normalidade, como citados acima.
- Papel de Probabilidade
Para o gráfico de Papel de Probabilidade, devemos primeiramente ordenar os resíduos $e_i$, encontrar os valores de $\Phi^{-1}(d_i)$ e então apresentar no gráfico os pontos $(e_i,\Phi^{-1}(d_i))$, $i=1,\dots,n.$ Se a suposição de normalidade for adequada, esperamos um comportamento linear dos pontos no gráfico. O vetor com os resíduos ordenados, para os dados na “Motivação 2” é
$$e=(-44,58; -37,09; -30,36; -25,01; -23,23; -19,88; -11,87; 5,34; 15,48; 15,65; 24,56; 30,35;$$
$$37,46; 63,20).$$
Já o vetor com os valores de $\Phi^{-1}(d_i)$ é dado por
$$\Phi^{-1}(d)=(-1,66;-1,18;-0,89;-0,65;-0,45;-0,26;-0,09; 0,09; 0,26; 0,45; 0,65; 0,89; 1,18; 1,66),$$
em que $d_i=(i-0,3)/(n+0,4)$ para $i=1,\dots,n$ e $\Phi^{-1}(d_{(i)})$ é o quantil da distribuição normal padrão calculado no ponto $d_{(i)}$.
Usando o software Action temos o seguinte resultado:
Figura 4.3.2: Gráfico de Papel de Probabilidade para os resíduos do modelo linear simples ajustado - Motivação 2.
Como os pontos seguem o comportamento da reta (não estão distantes dela), temos indícios de que os erros são normalmente distribuídos.
Testes de Normalidade
Em relação aos testes de normalidade, precisamos obter os valores de $F_{n}(e_{(i)})$ e de $F(e_{(i)})$. O primeiro é calculado fazendo a razão entre a posição i e o valor total de resíduos, no caso $n=14$ (distribuição empírica acumulada). O segundo é encontrado na tabela da distribuição normal padrão, considerando a transformação dos dados dada pela relação
$$Z_{(i)}=\frac{e_{(i)}-\overline{e}}{s}$$
em que $\overline{e}$ é a média aritmética e s é o desvio padrão dos resíduos. Resumimos na tabela 3.1.3 os valores utilizados no cálculo das estatísticas de teste de normalidade.
| $e$ | $F(e_i)$ | $F_n(e)$ | $F_n(e_i)-F(e_i)$ | $F(e_i)-F_n(e_{i-1})$ | $ln(F(e_i))$ | $ln(1-F(e_i)$ |
|---|---|---|---|---|---|---|
| -44,58 | 0,0826 | 0,0714 | -0,0112 | 0,0826 | -2,493 | -0,086 |
| -37,08 | 0,1242 | 0,1428 | 0,0187 | 0,0528 | -2,086 | -0,132 |
| -30,36 | 0,1723 | 0,2143 | 0,0419 | 0,0295 | -1,758 | -0,189 |
| -25 | 0,2181 | 0,2857 | 0,0675 | 0,0038 | -1,5224 | -0,246 |
| -23,23 | 0,2348 | 0,3571 | 0,1223 | -0,0509 | -1,449 | -0,267 |
| -19,87 | 0,268 | 0,4286 | 0,1605 | -0,0891 | -1,316 | -0,312 |
| -11,87 | 0,3558 | 0,5 | 0,1441 | -0,0728 | -1,033 | -0,439 |
| 5,34 | 0,5660 | 0,5714 | 0,0054 | 0,066 | -0,569 | -0,835 |
| 15,47 | 0,685 | 0,6428 | -0,0421 | 0,1136 | -0,378 | -1,155 |
| 15,65 | 0,6869 | 0,7142 | 0,0274 | 0,0441 | -0,376 | -1,161 |
| 24,55 | 0,7776 | 0,7857 | 0,0081 | 0,0634 | -0,251 | -1,503 |
| 30,34 | 0,8275 | 0,8571 | 0,0296 | 0,0418 | -0,189 | -1,757 |
| 37,45 | 0,8781 | 0,9286 | 0,0504 | 0,021 | -0,129 | -2,105 |
| 63,2 | 0,9754 | 1 | 0,0246 | 0,0468 | -0,0249 | -3,705 |
Tabela 4.2.4: Valores para cálculo das estatísticas de teste de Normalidade
- Kolmogorov-Smirnov
Para o teste de Kolmogorov-Smirnov, temos que a estatística de teste é dada por
$$D_n=\max(0,1605;0,1136)=0,1605.$$
Considerando $\alpha=0,05$ e n = 14, encontramos pela tabela de valores críticos que o valor crítico é de aproximadamente 0,34. Como Dn = 0,1605 < 0,34, não temos evidências para rejeitar a hipótese de normalidade dos resíduos.
- Anderson-Darling
Utilizando a fórmula $(\star)$ no teste de Anderson-Darling no conteúdo de Inferência, temos que
$$D=\sum_{i=1}^n[(2i-1)\ln(F(e_i))+(2(n-i)+1)\ln(1-F(e_i))]= -200,3.$$
Assim, segue que $$\displaystyle A^2=-n-\frac{D}{n}=-14-\frac{(-200,3)}{14}=0,308.$$
Para o cálculo do p-valor, precisamos encontrar a estatística de Anderson Darling modificada. Considerando μ e σ desconhecidos, temos que a estatística modificada é dada por
$$A_m^2=A^2\times(1+(0,75/n)+(2,25/n^2))=0,308^2\times(1+(0,75/14)+(2,25/14^2))=0,308\times 1,065=0,328.$$
Novamente, para obter o P-valor aproximado analisamos a Tabela com quantis e valores da estatística de Anderson Darling. Como $A_m^2=0,328~<~0,56$ temos que p_valor é maior do que 15%. No R um valor aproximado do p-valor é obtido fazendo interpolação.
Como p_valor é maior do que 15%, existe forte evidência pelo teste de Anderson-Darling de que os resíduos são normalmente distribuídos.
- Shapiro Wilk
Para o cálculo da estatística de teste, precisamos dos valores contidos na Tabela 4.2.5.
| $i$ | $n-i+1$ | $a$ | $e_{n-i+1}$ | $e_i$ | $a(e_{n-i+1}-e_i)$ |
|---|---|---|---|---|---|
| 1 | 14 | 0,5251 | 63,201 | -44,5841 | 56,59796 |
| 2 | 13 | 0,3318 | 37,4588 | -37,0884 | 24,73476 |
| 3 | 12 | 0,246 | 30,3464 | -30,3643 | 14,93483 |
| 4 | 11 | 0,1802 | 24,5563 | -25,009 | 8,931667 |
| 5 | 10 | 0,124 | 15,6505 | -23,2338 | 4,821653 |
| 6 | 9 | 0,0727 | 15,4769 | -19,8784 | 2,57033 |
| 7 | 8 | 0,024 | 5,3414 | -11,8735 | 0,413158 |
Tabela 4.2.5: Medidas para cálculo da estatística de Shapiro-Wilk da Motivação 2
A estatística do teste de Shapiro-Wilk é dada por
$$W=\frac{b^2}{\sum\limits_{i=1}^n(e_i-\bar{e})^2},$$
em que
$$b=\sum_{i=1}^{n/2}a_{n-i+1}\times (e_{n-i+1}-e_i).$$
Assim, temos que
$$W=\frac{113,0043^2}{13.421,12}=0,9515,$$
Como $W_{calc}=0,9515 > W_{(0,05;14)}=0,874$, em que $W_{(0,05;14)}$ é dada pela tabela de valores críticos, não rejeitamos a suposição de normalidade. Portanto, concluímos que com um nível de significância de 5%, os resíduos são normalmente distribuídos.
Usando o software Action temos os seguintes resultados:
| Estatística | P-valor | |
|---|---|---|
| Anderson-Darling | 0.3076533977 | 0.5183067119 |
| Shapiro-Wilk | 0.9514879996 | 0.5840384217 |
| Kolmogorov-Smirnov | 0.160504667 | 0.4235051205 |
Tabela 4.2.6: Teste de Normalidade
3.2 Diagnóstico de Homocedasticidade
Homocedasticidade é o termo para designar variância constante dos erros experimentais $\varepsilon_i$ para observações distintas. Caso a suposição de homocedasticidade não seja válida, podemos listar alguns efeitos no ajuste do modelo:
- Os erros padrões dos estimadores, obtidos pelo Método dos Mínimos Quadrados, são incorretos e portanto a inferência estatística não é válida.
- Não podemos mais dizer que os Estimadores de Mínimos Quadrados são os melhores estimadores de mínima variância para β, embora ainda possam ser não viciados.
Vale ressaltar que a ausência de homoscedasticidade é chamada de heteroscedasticidade. Com isso, testamos as hipóteses:
$$ \begin{cases} H_0: \sigma^2_1=\sigma^2_2 = … =\sigma^2_n \cr H_1: \hbox{ pelo menos um dos } \sigma_i^{2^\prime} \hbox{s diferente, } \quad i=1,\ldots,n. \end{cases} $$
Dado um modelo de regressão linear $$y_i = \beta_0 + \beta_1 x_{1,i} + \cdots + \beta_p x_{p,i} + \varepsilon_i, \quad i=1, \cdots , n,$$ a hipótese de homocedasticidade implica que $$Var(\varepsilon \mid x_{1,i}, \cdots , x_{p,i}) = \sigma^2.$$ Neste caso, a variância do erro experimental é a mesma para todas as observações. Ao longo desta seção, vamos apresentar os diversos testes propostos na literatura. Quando a hipótese de homocedasticidade é violada temos que
$$Var(\varepsilon \mid x_{1,i}, \cdots , x_{p,i}) = \sigma^2_i, \quad i=1, \cdots , n$$
Neste caso, a vari6ancia de cada observação pode ser diferente. Visualmente, isto significa que se realizarmos o gráfico dos resíduos pelos valor ajustados, o comportamento não deve ser o mesmo ao longo dos valores ajustados.
3.2.1 Gráfico dos Resíduos versus Valores Ajustados
O gráfico dos resíduos versus valores ajustados (valores preditos) é uma das principais técnicas utilizadas para verificar as suposições dos resíduos. Além da detecção de heteroscedasticidade, esse gráfico pode indicar que não existe uma relação linear entra as variáveis explicativas com a variável resposta por meio de alguma tendência nos pontos. Por exemplo, se os pontos do gráfico formam uma parábola, é indicativo que termos de segundo grau sejam necessários.
Para o diagnóstico de heteroscedasticidade, tentamos encontrar alguma tendência no gráfico. Por isso, se os pontos estão aleatoriamente distribuídos em torno do 0, sem nenhum comportamento ou tendência, temos indícios de que a variância dos resíduos é homoscedástica. Já a presença de “funil” é um indicativo da presença de heteroscedasticidade.
Exemplo 3.2.1.1
Considerando o exemplo na “Motivação 1”, vamos fazer o gráfico dos resíduos versus valores ajustados.
Usando o software Action temos o seguinte resultado:
Figura 4.3.3: Gráfico dos resíduos versus valores ajustados (Motivação 1)
Exemplo 3.2.1.2
Considerando o exemplo na “Motivação 2”, vamos fazer o gráfico dos resíduos versus valores ajustados.
Usando o software Action temos o seguinte resultado:
Figura 4.3.4: Gráfico dos resíduos versus valores ajustados (Motivação 2)
3.2.2 Teste de Breusch-Pagan
Baseado no teste multiplicador de Lagrange, o teste de Breusch-Pagan é bastante utilizado para testar a hipótese nula de que as variâncias dos erros são iguais (homoscedasticidade) versus a hipótese alternativa de que as variâncias dos erros são uma função multiplicativa de uma ou mais variáveis, sendo que esta(s) variável(eis) pode(m) pertencer ou não ao modelo em questão. É indicado para grandes amostras e quando a suposição de normalidade nos erros é assumida.
A estatística de teste neste caso é obtida da seguinte maneira:
Inicialmente, ajustamos o modelo de regressão linear (simples ou múltiplo) e encontramos os resíduos $e=(e_1,\dots,e_n)$ e os valores ajustados $\widehat{y}=(\widehat{y_1},\dots,\widehat{y}_n)$. Em seguida, consideramos os resíduos ao quadrado e os padronizamos de modo que a média do vetor de resíduos padronizados, que denotaremos por u, seja 1. Esta padronização é feita dividindo cada resíduo ao quadrado pela SQE/n em que SQE é a Soma de Quadrados dos Resíduos do modelo ajustado e n é o número de observações. Desta forma, temos que cada resíduo padronizado é dado por
$$u_i=\dfrac{e_i^2}{SQE/n},~~i=1,\dots,n,$$
em que $SQE=\displaystyle\sum\limits_{i=1}^{n}{e_i^2}.$
Por fim, fazemos a regressão entre $u=(u_1,…,u_n)$ (variável resposta) e o vetor $\widehat{y}$ (variável explicativa) e obtemos a estatística $\chi^2_{BP}$ calculando a Soma de Quadrados da Regressão de u sobre $\widehat{y}$ e dividindo o valor encontrado por 2. A estatística usada no teste é a $\chi_{B P}^{2}$ studentizada que édada por
$$\chi_{B P \text { studentizada}}^{2}=\frac{\chi_{B P}^{2}}{\lambda}$$
em que $\lambda=\frac{\operatorname{Var}\left(\varepsilon^{2}\right)}{2 \operatorname{Var}(\varepsilon)^{2}}$.
Sob a hipótese nula, esta estatística tem distribuição qui-quadrada com 1 grau de liberdade. Resumidamente, se não existe heteroscedasticidade, é de se esperar que os resíduos ao quadrado não aumentem ou diminuam com o aumento do valor predito, $\widehat{y}$ e assim, a estatística de teste deveria ser insignificante.
Exemplo 3.2.2.1
Considerando o exemplo na “Motivação 1”, vamos realizar o teste de Breusch Pagan para testar a homocedasticidade dos resíduos.
A SQE (Soma de Quadrados do Erro) da regressão da variável Dureza (resposta) sobre a variável explicativa Temperatura é 41,16. Além disso, temos que n=20. Assim, dividindo os resíduos ao quadrado por $\dfrac{41,16}{20}=2,058$, encontramos os resíduos padronizados.
A SQR (Soma de Quadrados da Regressão) obtida da regressão de $u$ (vetor de resíduos padronizados) sobre o valor ajustado da regressão da variável Dureza sobre a variável Temperatura, $\widehat{y}$, é dada por 0,151. Logo, o valor da estatística de teste é dado por $\chi^2_{BP}=\dfrac{0,151}{2}=0,0784$.
Já o p_valor para este exemplo, isto é, $P(\chi^2> 0,0784)$ é $0,7794$. Considerando um nível de significância de 5%, não rejeitamos a hipótese de homoscedasticidade dos erros.
Usando o software Action temos os seguintes resultados:
Teste de Homocedasticidade - Breusch Pagan
| Estatística | GL | P-valor |
|---|---|---|
| 0.0784452527 | 1 | 0.7794155173 |
Tabela 4.3.7: Teste de Homocedasticidade - Breusch Pagan (Motivação 1)
Exemplo 3.2.2.2
No exemplo na “Motivação 2”, vamos calcular a estatística de Breusch Pagan para testar a homoscedasticidade dos resíduos.
A SQE (Soma de Quadrados do Erro) da regressão da variável Ganho (resposta) sobre as variáveis explicativas Tempo e Dose_de_íons é 13.421. Além disso, temos que n=14. Assim, dividindo os resíduos ao quadrado por $\dfrac{13.421}{14}=958,643$, encontramos os resíduos padronizados.
A SQR (Soma de Quadrados da Regressão) obtida da regressão de $u$ (vetor de resíduos padronizados) sobre o valor ajustado da regressão da variável Ganho sobre as variáveis Tempo e Dose_de_íons, $\widehat{y}$, é dada por 0,866. Logo, o valor da estatística de teste é dado por $\chi^2_{BP studentizada}=0,7993$.
Já o p_valor para este exemplo, isto é, $P(\chi^2 > 0,7993)$ é 0,3713. Considerando um nível de significância de 5%, não rejeitamos a hipótese de homoscedasticidade dos erros.
Usando o software Action temos os seguintes resultados:
| Estatística | GL | P-valor |
|---|---|---|
| 0.7992707163 | 1 | 0.3713115028 |
Tabela 4.3.8: Teste de Homocedasticidade - Breusch Pagan (Motivação 2)
3.2.3 Teste de Goldfeld-Quandt
O teste de Goldfeld-Quandt também é utilizado para testar a homoscedasticidade dos resíduos. Entre as limitações deste teste está a exigência de que a amostra seja relativamente grande.
O teste consiste inicialmente em ordenar as observações de acordo com a variável explicativa que se acredita a responsável pela heteroscedasticidade. Após isso, divide-se a amostra ordenada em 3 partes de tal forma que a parte do meio tenha aproximadamente 20% dos dados e que as partes 1 e 3 tenham quantidade de dados semelhantes. Então, ajusta-se um modelo de regressão com os dados da parte 1 (contendo os menores valores da variável explicativa utilizada na ordenação) e outro modelo de regressão com os dados da parte 3 (contendo os maiores valores da variável explicativa utilizada na ordenação). Por fim, testa-se a hipótese de que as variâncias dos erros em ambas regressões são iguais contra a hipótese de que a variância dos erros na parte 3 é maior do que a variância dos erros na parte 1, utilizando o teste F.
A estatística de teste neste caso é dada por
$$F_{GQ}=\frac{SQE^b/(n_{3}-(p+1))}{SQE^a/(n_{1}-(p+1))},$$
em que $SQE^a$ e $SQE^b$ são as somas de quadrados dos resíduos da regressão para o grupo inferior (parte 1) e para o grupo superior (parte 3), respectivamente, $n_{1}$ é o número de observações da parte 1 e $n_{3}$ é o número de observações da parte 3. Chamamos de d o número de observações omitidas (parte 2). Essa estatística tem distribuição $F_{(n_{3}-(p+1),n_{1}-(p+1))}$. Desta forma, considerando um nível de significância $\alpha=0,05$, rejeitamos a hipótese nula, ou seja, a hipótese de que as variâncias são iguais se $F_{GQ} > F_{\frac{\alpha}{2}}$.
Exemplo 3.2.3.1
Considerando o exemplo na “Motivação 1”, vamos realizar o teste de Goldfeld-Quandt para testar a homocedasticidade dos resíduos.
Temos apenas uma variável explicativa, Temperatura $(p=1)$. Desta forma, vamos ordenar os dados em relação a essa variável. Como $n=20$, descartamos $d=4$ observações, ou seja, 20% dos dados. Desta forma, a parte 1 possui $n_1=8$ observações e a parte 3 possui $n_3=8$ observações. A soma de quadrados dos resíduos do grupo superior (parte 3) é 13,2 enquanto que a soma de quadrados dos resíduos do grupo inferior (parte 1) é 7,867. Logo, a estatística de teste é
$$F_{GQ}=\frac{13,2/(8-2)}{7,867/(8-2)}=1,678.$$
$$\text{p-valor} = 2\min{\mathbb{P}[F \ > \ 1,678 | H_0];\mathbb{P}[F \ < \ 1,678 | H_0]}=0,545208844.$$
O p-valor do teste é 0,5452. Considerando um nível de confiança de 95%, não rejeitamos a hipótese de homoscedasticidade dos resíduos.
Usando o software Action temos os seguintes resultados:
| Variável | Estatística | GL1 | GL2 | P-valor |
|---|---|---|---|---|
| Temperatura | 1.67796610169501 | 6 | 6 | 0.54520884 |
Tabela 4.3.9: Teste de Homocedasticidade - Goldfeld Quandt (Motivação 1)
Exemplo 3.2.3.2
No exemplo na “Motivação 2”, vamos calcular a estatística de Goldfeld-Quandt para testar a homoscedasticidade dos resíduos.
Neste exemplo temos duas variáveis explicativas $(p=2)$. Desta forma, aplicamos o teste de Goldfeld-Quandt inicialmente ordenando os dados em relação a variável Tempo e posteriormente ordenando em relação a variável Dose_de_íons. Desta forma, conseguimos detectar qual das duas variáveis explicativas é a responsável pela heteroscedasticidade, caso exista.
Em ambos os casos (ordenação por Tempo ou Dose_de_íons), descartamos aproximadamente 20% das observações, isto é, 4 observações. Desta forma, temos que o número de observações restantes é 11. Quando o número de observações é impar, sempre deixamos a parte 3 com uma observação a mais do que a parte 1. Desta forma, temos $n_1=5$ e $n_3=6$.
Considerando a variável Tempo na ordenação dos dados, temos que a soma de quadrados dos resíduos do grupo inferior (parte 1) é 780,126 enquanto que a do grupo superior (parte 3) é 849,822. Logo, a estatística de teste fica
$$F_{GQ}=\frac{849,822/(6-3)}{780,126/(5-3)}=0,726.$$
$$\text{p-valor} = 2\min{\mathbb{P}[F \ > \ 0,726 | H_0];\mathbb{P}[F \ < \ 0,726 | H_0]}=0,752941559.$$
O p-valor do teste é 0,753. Logo, a um nível de confiança de 95%, não rejeitamos a hipótese de homoscedasticidade dos resíduos da Motivação 2.
Em relação a ordenação dos dados pela variável Dose_de_íons, temos que a soma de quadrados dos resíduos do grupo inferior (parte 1) é 2.592 enquanto que a do grupo superior (parte 3) é 8.837,965. Logo, a estatística de teste fica
$$F_{GQ}=\frac{8.837,965/(6-3)}{2.592/(5-3)}=2,286.$$
$$\text{p-valor} = 2\min{\mathbb{P}[F \ > \ 2,286 | H_0];\mathbb{P}[F \ < \ 2,286 | H_0]}=0,637566748.$$
O p-valor do teste é 0,638. Logo, a um nível de confiança de 95%, também não rejeitamos a hipótese de homoscedasticidade dos resíduos da Motivação 2.
Usando o software Action temos os seguintes resultados:
| Variável | Estatística | GL1 | GL2 | P-valor |
|---|---|---|---|---|
| Dose de ions | 2.27313910924677 | 3 | 2 | 6401498344179 |
| Tempo | 0.72622614253487 | 3 | 2 | 0.752941558959312 |
Tabela 4.3.9: Teste de Homocedasticidade - Goldfeld Quandt (Motivação 1)
3.3 Diagnóstico de Independência
Para verificar se os resíduos são independentes, podemos utilizar técnicas gráficas e testes. A seguir, temos o diagnóstico de independência por essas duas formas.
3.3.1 Gráfico dos Resíduos versus a Ordem de Coleta
Uma análise gráfica para verificar a hipótese de independência dos resíduos pode ser feita por meio do gráfico dos resíduos versus a ordem da coleta dos dados. Se ao avaliar o gráfico, percebemos uma tendência dos pontos, ou seja, se os pontos tiverem um comportamento que se repete em determinado ponto do gráfico, temos indícios de dependência dos resíduos.
Exemplo 3.3.1.1
Para os dados do exemplo na “Motivação 1”, fazer o gráfico dos resíduos do modelo ajustado versus a ordem da coleta e avaliar a suposição de independência dos resíduos.
Usando o software Action temos o seguinte resultado:
Figura 4.3.5: Gráfico dos resíduos versus a ordem da coleta considerando os dados do exemplo na Motivação 1.
Analisando o gráfico vemos que os pontos não parecem ter uma tendência e por isso temos indícios de independência dos erros.
Exemplo 3.3.1.2
Considerando os dados transformados do exemplo na “Motivação 2”, fazer o gráfico dos resíduos do modelo ajustado versus a ordem da coleta e avaliar a suposição de independência dos resíduos.
Usando o software Action temos o seguinte resultado:
Figura 4.3.6: Gráfico dos resíduos versus a ordem da coleta considerando os dados do exemplo na Motivação 2.
Novamente, verificamos pelo gráfico que os pontos não parecem ter uma tendência e por isso temos indícios de independência dos erros.
3.3.2 Teste de Durbin-Watson
O teste de Durbin-Watson é utilizado para detectar a presença de autocorrelação (dependência) nos resíduos de uma análise de regressão. Este teste é baseado na suposição de que os erros no modelo de regressão são gerados por um processo autoregressivo de primeira ordem, de acordo com $$\varepsilon_i=\rho\varepsilon_{i-1}+a_i,$$
em que $\varepsilon_i$ é o termo do erro do modelo na i-ésima observação, $a_i\overset{iid}{\sim}N(0,\sigma_{a}^2)$ e $\rho$ ($|\rho|< 1$) é o parâmetro de autocorrelação. Testamos a presença de autocorrelação por meio das hipóteses
$$\begin{cases} H_{0}: \rho=0. \cr H_{1}:\rho \neq 0. \end{cases}.$$
Sendo $e_i$ o resíduo associado à i-ésima observação, temos que a estatística do teste de Durbin-Watson é dada por
$$dw=\frac{\sum\limits_{i=2}^n(e_i-e_{i-1})^2}{\sum\limits_{i=1}^{n}e_i^2},$$
em que 0 ≤ $dw$ ≤ 4. A distribuição de $dw$ depende da matriz $X$. Entretanto, podemos tomar a decisão comparando o valor de $dw$ com os valores críticos $dL$ e $dU$ da Tabela de Durbin-Watson (Tabela 4.3.10). Assim,
se $0$ ≤ $dw$ < $dL$ então rejeitamos $H_0$ (dependência);
-
se $dL$ ≤ $dw$ ≤ $dU$ então o teste é inconclusivo;
-
se $dU$ < $dw$ < $4-dU$ então não rejeitamos $H_0$ (independência);
-
se $4-dU$ ≤ $dw$ ≤ $4-dL$ então o teste é inconclusivo;
-
se $4-dL$ < $dw$ ≤ $4$ então rejeitamos $H_0$ (dependência).
Quando $0~\leq~dw~<~d_L$ temos evidência de uma correlação positiva. Já quando $4-d_L ~<~dw~\leq~4$, a correlação é negativa. No caso em que não rejeitamos $H_0$, temos que não existe autocorrelação, ou seja, os resíduos são independentes. Podemos também tomar a decisão pelo p-valor.
Os valores críticos tabelados apresentados na Tabela 4.3.10 são geralmente utilizados para testar $\rho~=~0$ versus $\rho~>~0$. Desta forma, se para um determinado $\alpha$ utilizarmos os valores da Tabela 4.3.10 para testar $\rho~=~0$ versus $\rho~\neq~0$, o erro tipo I do teste em questão será $2\alpha$.
| Nível de significância | 1 | 1 | 2 | 2 | 3 | 3 | 4 | 4 | 5 | 5 | |
|---|---|---|---|---|---|---|---|---|---|---|---|
| $n$ | $d_L$ | $d_U$ | $d_L$ | $d_U$ | $d_L$ | $d_U$ | $d_L$ | $d_U$ | $d_L$ | $d_U$ | |
| 0,01 | 0,81 | 1,07 | 0,7 | 1,25 | 0,59 | 1,46 | 0,49 | 1,7 | 0,39 | 1,96 | |
| 15 | 0,025 | 0,95 | 1,23 | 0,83 | 1,4 | 0,71 | 1,61 | 0,59 | 1,84 | 0,48 | 2,09 |
| 0,05 | 1,08 | 1,36 | 0,95 | 1,54 | 0,82 | 1,75 | 0,69 | 1,97 | 0,56 | 2,21 | |
| 0,01 | 0,95 | 1,15 | 0,86 | 1,27 | 0,77 | 1,41 | 0,63 | 1,57 | 0,6 | 1,74 | |
| 20 | 0,025 | 1,08 | 1,28 | 0,99 | 1,41 | 0,89 | 1,55 | 0,79 | 1,7 | 0,7 | 1,87 |
| 0,05 | 1,2 | 1,41 | 1,1 | 1,54 | 1 | 1,68 | 0,9 | 1,83 | 0,79 | 1,99 | |
| 0,01 | 1,05 | 1,21 | 0,98 | 1,3 | 0,9 | 1,41 | 0,83 | 1,52 | 0,75 | 1,65 | |
| 25 | 0,025 | 1,13 | 1,34 | 1,1 | 1,43 | 1,02 | 1,54 | 0,94 | 1,65 | 0,86 | 1,77 |
| 0,05 | 1,2 | 1,45 | 1,21 | 1,55 | 1,12 | 1,66 | 1,04 | 1,77 | 0,95 | 1,89 | |
| 0,01 | 1,13 | 1,26 | 1,07 | 1,34 | 1,01 | 1,42 | 0,94 | 1,51 | 0,88 | 1,61 | |
| 30 | 0,025 | 1,25 | 1,38 | 1,18 | 1,46 | 1,12 | 1,54 | 1,05 | 1,63 | 0,98 | 1,73 |
| 0,05 | 1,35 | 1,49 | 1,28 | 1,57 | 1,21 | 1,65 | 1,14 | 1,74 | 1,07 | 1,83 | |
| 0,01 | 1,25 | 1,34 | 1,2 | 1,4 | 1,15 | 1,46 | 1,1 | 1,52 | 1,05 | 1,58 | |
| 40 | 0,025 | 1,35 | 1,45 | 1,3 | 1,51 | 1,25 | 1,57 | 1,2 | 1,63 | 1,15 | 1,69 |
| 0,05 | 1,44 | 1,54 | 1,39 | 1,6 | 1,34 | 1,66 | 1,29 | 1,72 | 1,23 | 1,79 | |
| 0,01 | 1,32 | 1,4 | 1,28 | 1,45 | 1,24 | 1,49 | 1,2 | 1,54 | 1,16 | 1,59 | |
| 50 | 0,025 | 1,42 | 1,5 | 1,38 | 1,54 | 1,34 | 1,59 | 1,3 | 1,64 | 1,26 | 1,69 |
| 0,05 | 1,5 | 1,59 | 1,46 | 1,63 | 1,42 | 1,67 | 1,38 | 1,72 | 1,34 | 1,7 |
Tabela 4.3.10: Valores críticos do teste de Durbin-Watson.
Exemplo 3.3.2.1
Utilizando o teste de Durbin Watson, verificar se os resíduos do modelo ajustado do exemplo na “Motivação 1” são independentes.
Para calcular a estatística de Durbin Watson, precisamos dos valores da Tabela 4.3.11.
| $e$ | $e_{t}-e_{t-1}$ | ${(e_{t}-e_{t-1})}^2$ | $e_t^2$ |
|---|---|---|---|
| -0,14 | - | - | 0,0196 |
| -0,14 | 0 | 0 | 0,0196 |
| -0,14 | 0 | 0 | 0,0196 |
| -1,14 | -1 | 1 | 1,2996 |
| -2,14 | -1 | 1 | 4,5796 |
| 3,02 | 5,16 | 26,62 | 9,1204 |
| 1,02 | -2 | 4 | 1,0404 |
| 0,02 | -1 | 1 | 0,0004 |
| 1,02 | 1 | 1 | 1,0404 |
| 1,02 | 0 | 0 | 1,0404 |
| 1,18 | 0,16 | 0,0256 | 1,3924 |
| -2,82 | -4 | 16 | 7,9524 |
| -0,82 | 2 | 4 | 0,6724 |
| 2,18 | 3 | 9 | 4,7524 |
| -0,82 | -3 | 9 | 0,6724 |
| 0,34 | 1,16 | 1,34 | 0,1156 |
| 0,34 | 0 | 0 | 0,1156 |
| 0,34 | 0 | 0 | 0,1156 |
| -2,66 | -3 | 9 | 7,0756 |
| 0,34 | 3 | 9 | 0,1156 |
| Soma: | 91,99 | 41,16 |
Tabela 4.3.11: Medidas para o cálculo da estatística de Durbin-Watson.
A estatística de teste é dada por $$dw=\frac{91,99}{41,16}=2,235.$$
Para encontrar os valores críticos $d_L$ e $d_U$ na Tabela 4.3.10, devemos considerar $\alpha=\alpha^*/2$, em que $\alpha^*$ é o nível de significância do teste. Para $\alpha^*=0,05$ temos que $\alpha=0,025$. Assim, com $n=20$ e $k=1$ variável explicativa, temos que $d_L=1,08$ e $d_U=1,28$. Como $d_U = 1,28<2,235\leq4-d_U=2,72$, não rejeitamos $H_0$ e portanto, podemos afirmar que com um nível de confiança de 95%, os resíduos são independentes.
Usando o software Action temos os seguintes resultados:
| Estatística | P-valor |
|---|---|
| 2.2351 | 0.7686 |
Tabela 4.3.12: Teste de Independência - Durbin-Watson (Motivação 1)
Exemplo 3.3.2.2
Analogamente ao exemplo 3.3.2.1, vamos testar pelo teste de Durbin Watson se os resíduos do modelo ajustado do exemplo na “Motivação 2” são independentes.
Para o cálculo da estatística, precisamos dos valores da Tabela 3.3.2.3.
| $e$ | $e_t-e_{t-1}$ | ${(e_t-e_{t-1})}^2$ | $e_t^2$ |
|---|---|---|---|
| 30,3464 | - | - | 920,9 |
| 15,4769 | -14,8695 | 221,1 | 239,53 |
| -30,3643 | -45,8412 | 2101,42 | 921,99 |
| -23,2338 | 7,1305 | 50,84 | 539,81 |
| 5,3414 | 28,5752 | 816,54 | 28,53 |
| -11,8735 | -17,2149 | 296,35 | 140,98 |
| 63,201 | 75,0744 | 5636,17 | 3994,36 |
| -25,009 | -88,2099 | 7780,99 | 625,45 |
| -19,8784 | 5,1305 | 26,32 | 395,15 |
| -37,0884 | -17,2099 | 296,18 | 1375,55 |
| -44,5841 | -7,4958 | 56,19 | 1987,75 |
| 24,5563 | 69,1405 | 4780,4 | 603,01 |
| 37,4588 | 12,9025 | 166,48 | 1403,16 |
| 15,6505 | -21,8083 | 475,6 | 244,94 |
| Soma | 22704,58 | 13421,11 |
Tabela 3.3.2.3: Medidas para o cálculo da estatística de Durbin-Watson
A estatística de teste é dada por $$dw=\frac{22704,58}{13421,11}=1,6917.$$
Analisando a Tabela 4.3.10 vemos que não há valores críticos $d_U$ e $d_L$ para $n=14$. Desta forma, consideramos os valores aproximados para n=15. Com um nível de significância $\alpha^{*}=0,05$ $(\alpha=0,025)$ e $k=2$ variáveis explicativas, temos que $d_L~\approx~0,83$ e $d_U~\approx~1,40$. Como $d_L~\approx~1,40~<~1,69~<~2,60~\approx4-d_L$, com um nível de significância de 5%, não rejeitamos $H_0$. Assim, dizemos que os resíduos são independentes.
Usando o software Action temos os seguintes resultados:
| Estatística | P-valor |
|---|---|
| 1.691706 | 0.709529 |
Tabela 4.3.13: Teste de Independência - Durbin-Watson (Motivação 2)
3.4 Pontos influentes e valores extremos
O valor extremo ou valor atípico (outlier) é uma observação com “alto” resíduo. Neste caso, temos uma observação cuja variável dependente é inesperada dado o valor da variável explicativa. Por outro lado, dizemos que uma observação é um ponto de alavanca (ou leverage alto) se a variável explicativa é extrema. Quando temos apenas uma variável explicativa, um valor extremo significa um valor baixo e/ou alto para a variável explicativa. Se tivermos um vetor de variáveis explicativas, um ponto extremo significa uma combinação de valores inesperados para as variáveis explicativa. Assim, vamos utilizar o resíduo para detectar valores extremos na resposta e o ponto de alavanca para detectarmos valores extremos ou inesperados para a variável explicativa.
Um ponto influente é uma observação que pode influenciar em qualquer parte da análise de regressão, como a previsão de uma resposta, o coeficiente angular ou o resultado de testes de hipóteses. A inclusão ou não desta observação na análise influencia em alguma parte da análise de regressão. Os valores extremos e os pontos de alavanca tem potencial para serem influentes, mas temos que investigar para avaliar o quanto eles são influentes. Se um valor extremo for influente, ele interfere sobre a função de regressão ajustada (a inclusão ou não do ponto modifica substancialmente os valores ajustados). Mas uma observação pode ser considerada um valor extremo e não ser um ponto influente. Da mesma forma, podemos ter pontos que influenciam na análise de regressão mas não são valores extremos ou pontos de alavanca.
3.4.1 Ponto de alavanca
A localização dos pontos no espaço das variáveis explicativas é importante para determinarmos propriedades do modelo de regressão. Pontos de alavanca ou valores extremos na variável explicativa ($X$) são detectados por meio da matriz chapéu $H$. A matriz hat é dada por
$$H=X(X^\prime X)X^\prime,$$
no qual $X$ é a matriz das variáveis explicativas ou a matriz de planejamento. Sabemos que a matriz chapéu $H$ determina a variâncias e as covariâncias da resposta estimada $\hat{y}$ e dos resíduos $e$,
$$Var(\hat{y})=\sigma^2H\quad \text{e} \quad Var(e)=\sigma^2(I-H).$$
Denotamos por $h_{ii}$ é o i-ésimo elemento da diagonal principal da matriz chapéu $H$, que é denominado ponto de alavanca (ou, leverage) da observação $i$. Este elemento pode ser intepretado como a alavanca exercida pela i-ésima observação $y_i$ na resposta prevista $\hat{y}_i$. O ponto de alavanca $h_{ii}$ pode ser calculado na forma
$$h_{ii}=x_i^{^\prime }(X^\prime X)^{-1}x_i,$$
em que $x_i^\prime $ é a i-ésima linha da matriz $X$. A diagonal da matriz chapéu $H$ é uma medida padronizada da distância da i-ésima observação para o centro do espaço definido pelas variáveis explicativa. Portanto, valor alto para $h_{ii}$ significa que o ponto está distante das outras observações em relação às variáveis explicativas.
Como $\sum_{i=1}^nh_{ii}=posto(X)=p+1$, a média dos $h_{ii}$ é $(p+1)/n$. É recomendado destacar as observações para as quais $h_{ii}> 2(p+1)/n$. Também podemos identificar pontos de alavanca altos através do box-plot do vetor de $h_{ii}$. É importante destacarmos que nem toda observação com $h_{ii}$ alto é influente na análise de regressão. Observe que a diagonal da matriz chapéu avalia a localização da observação no espaço das variáveis explicativas. Como consequência, podemos ter observações com ponto de alavanca alto, mas nada influenciarem nas estimativas dos parâmetros. Assim, é nteressante avaliarmos os pontos de alavanca em conjunto com os resíduos padronizados ou studentizados.
Exemplo 3.4.1.1
No exemplo da “Motivação 1”, vamos calcular os valores de $h_{ii}$ para verificar se alguma observação é um valor extremo em X.
Para o primeiro indivíduo, temos que $x_1^{^\prime }=[1~~220]$, em que o primeiro valor é referente ao intercepto e o segundo é o valor da temperatura observada do indivíduo.
A matriz $(X^\prime X)^{-1}$ é:
$$(X^\prime X)^{-1}~= \begin{bmatrix} 82,86 \quad -0,36 \cr -0,36 \quad 0,0016 \end{bmatrix} $$
Com isso,
$$h_{11}=[1 \quad 220] \begin{bmatrix} 82,86 \quad -0,36 \cr -0,36 \quad 0,0016 \cr \end{bmatrix} \begin{bmatrix} 1 \cr 220 \end{bmatrix} =0,14.$$
Para as demais observações calculamos os valores de $h_{ii}$ de forma análoga e obtemos os seguintes valores:
| $h_{ii}$ | valores |
|---|---|
| $h_1$ | 0,14 |
| $h_2$ | 0,14 |
| $h_3$ | 0,14 |
| $h_4$ | 0,14 |
| $h_5$ | 0,14 |
| $h_6$ | 0,06 |
| $h_7$ | 0,06 |
| $h_8$ | 0,06 |
| $h_9$ | 0,06 |
| $h_{10}$ | 0,06 |
| $h_{11}$ | 0,06 |
| $h_{12}$ | 0,06 |
| $h_{13}$ | 0,06 |
| $h_{14}$ | 0,06 |
| $h_{15}$ | 0,06 |
| $h_{16}$ | 0,14 |
| $h_{17}$ | 0,14 |
| $h_{18}$ | 0,14 |
| $h_{19}$ | 0,14 |
| $h_{20}$ | 0,14 |
Tabela 4.3.14: Valores da diagonal principal da matriz H da Motivação 1.
Usando o software Action temos o seguinte resultado:
Figura 4.3.7: Gráfico dos valores de $h_{ii}$ considerando os dados do exemplo na Motivação 1.
Observamos pela Tabela 4.3.14 e pela Figura 4.3.7 que nenhum $h_{ii}$ é maior do que $2(p+1)/n=2(2)/20= 0,2$. Por isso, temos que nenhuma observação é considerada outlier em $X$.
Exemplo 3.4.1.2
Para o exemplo na “Motivação 2” vamos calcular a diagonal principal da matriz H ($h_{ii}$) para verificar se alguma observação é um outlier em X. Consideramos os dados transformados.
Para o primeiro indivíduo, temos que $x_1^\prime=[1~~-1~~-1]$, em que o primeiro valor é referente ao intercepto, o segundo referente ao tempo e o terceiro é referente à dose de íons. Além disso,
$$(X^\prime X)^{-1}~= \begin{bmatrix} 0,0720 \quad -0,0019 \quad 0,0091 \cr -0,0019 \quad 0,1660 \quad 0,0004 \cr 0,0091 \quad 0,004 \quad 0,1429 \end{bmatrix} $$
Assim,
$$h_{11}=[1 \quad -1 \quad -1] \begin{bmatrix} 0,0720 \quad -0,0019 \quad 0,0091 \cr -0,0019 \quad 0,1660 \quad 0,0004 \cr 0,0091 \quad 0,004 \quad 0,1429 \end{bmatrix} \begin{bmatrix} 1 \cr -1 \cr -1 \end{bmatrix} =0,367$$
Calculando o $h_{ii}$ para as outras observações obtemos os seguintes valores:
| $h_{ii}$ | valores |
|---|---|
| $h_1$ | 0,36744 |
| $h_2$ | 0,35801 |
| $h_3$ | 0,316927 |
| $h_4$ | 0,310215 |
| $h_5$ | 0,092191 |
| $h_6$ | 0,133456 |
| $h_7$ | 0,147617 |
| $h_8$ | 0,242964 |
| $h_9$ | 0,234893 |
| $h_{10}$ | 0,19677 |
| $h_{11}$ | 0,216595 |
| $h_{12}$ | 0,072974 |
| $h_{13}$ | 0,233037 |
| $h_{14}$ | 0,07691 |
Tabela 4.3.15: Valores da diagonal principal da matriz H.
Usando o software Action temos o seguinte resultado:
Figura 4.3.8: Gráfico dos valores de $h_{ii}$ considerando os dados do exemplo na Motivação 2.
Observamos pela Tabela 4.3.15 e pela Figura 4.3.8 que nenhum $h_{ii}$ é maior do que $2(p+1)/n=2(3)/14=0,428$. Por isso, temos que nenhuma observação é considerada outlier em $X$.
3.4.2 Valores extremos na resposta: resíduos altos
Em análise de regressão, a resposta observada de alguns casos pode não corresponder ao modelo ajustado aos dados, neste caso, dizemos que esta resposta é um outlier ou valor extremo. Em regressão linear simples, o gráfico de dispersão pode nos mostrar quais observações não correspondem ao modelo ajustado. Entretanto, quando lidamos com modelos de regressão múltipla, devemos realizar uma análise mais cuidadosa. Os resíduos são definidos como $e_i=Y_i-\hat{Y}_i$, que corresponde a diferença entre o valor observado e o valor ajustado pelo modelo. Entretanto, para uma melhor detecção de outliers em Y, diversas formas “padronizadas” foram propostas.
3.4.2.1 Resíduos Normalizado
O resíduo normalizado, $d_i$, corresponde ao resíduo dividido pelo devio padrão estimado dos resíduos, $\sqrt{QME}$. $$d_i=\frac{e_i}{\sqrt{QME}}.$$
Se os erros têm distribuição normal, então aproximadamente 95% dos resíduos normalizado $(d_i)$ devem estar no intervalo de (-2,2). Resíduos fora desse intervalo podem indicar a presença de outliers.
3.4.2.2 Resíduos Padronizado
Existem inúmeras maneiras de se expressar o vetor de resíduos “e” que nos será útil. $$e=Y-\hat{Y}=Y-X\hat{\beta}=Y - HY=(I-H)Y.$$
A matriz de covariâncias dos resíduos é,
$$Cov[e]=Cov[(I-H)Y]=(I-H)Var(Y)(I-H)^\prime$$
$$=\sigma^2(I-H)(I-H)^\prime$$
$$=\sigma^2(I-H),$$
no qual $H=X(X^\prime X)^{-1} X^\prime$ é a matriz chapéu (hat) e $X$ a matriz do modelo.
Assim, definimos os resíduos padronizados por
$$r_i =\dfrac{e_i}{\sqrt{QME(1-h_{ii})}}, \qquad \qquad i=1,2,\ldots,n,$$
com $\hat{\sigma}^2=QME$ e $h_{ii}$ o $i$-ésimo elemento da matriz chapéu $H$.
Os resíduos padronizados tem variâncias constantes $Var(r_i)=1$ o que consequentemente torna muito prática a procura por outliers, que são observações distantes das demais. Em geral, as observações com resíduo padronizado fora do intervalo $-3\leq r_i \leq 3$ deve são considerados outliers.
3.4.2.3 Resíduo Studentizado
Um dos testes mais utilizados para detectar outlier é o teste da mudança na média. Suponha que a $i$-ésima resposta é um candidato a outlier. Para todas as outras respostas assumimos que $\mathbb{E} [Y \mid X=x_j] = x_j^\prime \beta$, mas para o caso $i$ a função média é dada por $\mathbb{E} [Y \mid X=x_i] = x_i^\prime \beta + \delta$. A resposta esperado para o caso $i$ é modificada pelo parâmetro $\delta$. Assim, podemos realizar o teste $H_0 : \delta = 0$ versus $H_1 : \delta \ne 0$ para verificar se a $i$-ésima respota é um outlier. Assumimos que a variância é constante, isto é, $Var(Y \mid X) = \sigma^2$.
Resposta com valores altos de resíduos são candidatos a outlier, porém, nem toda observação com resíduo alto é um outlier. Resíduos alto devem ocorrer conforme a distribuição de probabilidade dos dados (normal). Além disso, nem todo outlier é ruim, pois o outlier é detectado em função do modelo ajustado, que pode não ser adequado aos dados.
Suponha que a observação $i$ seja suspeita de ser um outlier. Para verificar esta suspeita procedemos da seguinte forma:
-
Tiramos a observação $i$ dos dados ficando com as $n-1$ observações restantes;
-
Com os dados reduzidos, estimamos os parâmetros $\beta$ e $\sigma^2$. Denotamos as estimativs por $\beta_{(i)}$ e $\sigma_{(i)}$ para evidenciar que retiramos a observação $i$ do conjunto de dados. Neste caso, o estimador para $\sigma^2_{(i)}$ tem $n-1-p-1$ graus de liberdade, para o modelo com o intercepto. Se denotarmos por $p^\prime$ o número de linhas da matriz $X$, obtemos que $p^\prime = p+1$ no caso com intercepto e $p^\prime = p$ no caso sem o intercepto. Neste caso, temos que $\hat{\sigma}^2_{(i)}$ tem $n-p^{\prime} -1$ graus de liberdade.
-
Para os dados sem a observação $i$, calculamos o valor ajustado $\hat{Y} _{i(i)} = x_i^\prime \beta_{(i)}$. Desde que a $i$-ésima observação não foi utilizada para estimar os parâmetros $\beta_{(i)}$ concluímos que $Y_i$ e $\hat{Y} _{i(i)}$ são independentes. Como consequência, obtemos que
$$Var [Y_i - \hat{Y} _{i(i)}] = \sigma^2 + \sigma^2 x_i^\prime (X_{(i)}^\prime X_{(i)})^{-1} x_i,$$
no qual $X_{(i)}$ é a matrix $X$ sem a linha $i$. A variância é estimada substituindo $\sigma^2$ por $\hat{\sigma}^2_{(i)}$ na expressão acima.
- Sabemos que $\mathbb{E} [Y_i - \hat{Y}_{i(i)}] = \delta$ que é zero sob $H_0$ e diferente de zero sob $H_1$. Assumindo distribuição normal para os resíduos, temos que sob $H_0: \delta = 0$, a estatística do teste
$$ t_i = \frac{Y_i - \hat{Y} _{i(i)}}{\hat{\sigma}\sqrt{1+x_i^\prime(X_{(i)}^\prime X_{(i)})^{-1} x_{i}}},$$
tem distribuição $t$-Student com $n-p^\prime -1$ graus de liberdade. Esta estatística pode ser escrita na forma
$$ t_i = r_i \sqrt{\frac{n-p^\prime -1}{n-p^\prime -r_{i}^2}} = \frac{e_i}{\hat{\sigma}_{(i)} \sqrt{1-h_{ii}}},$$
nos quais $e_i$ é o resíduo e $r_i$ o resíduo padronizado.
Na sequência, vamos mostrar que a fórmula acima é válida. Suponha que a matrix do modelo $X$ tenha colunas linearmente independentes. Denotamos por $X_{(i)}$ a matrix do modelo sem a linha $i$. Sabemos que
$$ (X_{(i)}^\prime X_{(i)})^{-1} = (X^{\prime} X)^{-1} - \frac{(X^\prime X)^{-1} x_{i}^\prime x_i (X^\prime X)^{-1}}{1-h_{ii}}.$$
Esta fórmula é conhecida desde Gauss em 1821 e tem uma longa história de aplicações na estatística. A partir desta fórmula, temos que
$$\hat{\beta} _{(i)} = \hat{\beta} - \frac{(X^\prime X)^{-1} x_{i}^\prime e_i}{1-h_{ii}}.$$
Assim, uma estimativa da variância é dada por
$$ \hat{\sigma}^2_{(i)} = \hat{\sigma}^2 \left( \frac{n-p^\prime -1}{n - p^\prime - r_i^2}\right)^{-1},$$
o que nos leva à fórmula do resíduo studentizado $t_i$. Observe que o resíduo padronizado $r_i$ e o resíduo studentizado $t_i$ carregam a mesma informação. Assim, ao considerarmos que o grau de liberdade $n - p^\prime -1$ é alto, sabemos que $t_i$ se comporta aproximadamente como a distribuição normal e com isso, $\mid t_i \mid \geq 3$ é considerado uma observação extrema (outlier).
Exemplo 3.4.2.1
Para os dados na “Motivação 1” vamos calcular os resíduos studentizados e padronizados.
Solução: Como vimos em “Estimação dos Parâmetros do Modelo Linear”, o QME do exemplo na motivação 1 é calculado no “Exemplo 1.2.2” e é dado por 2,286. Já a diagonal da matriz H foi calculada no “Diagnóstico de Outliers em X”. Assim, os resíduos são:
| Resíduos | Resíduos Studentizados | Resíduos Padronizados |
|---|---|---|
| -0.14 | -0.097 | -0.0998 |
| -0.14 | -0.097 | -0.0998 |
| -0.14 | -0.097 | -0.0998 |
| -1.14 | -0.8049 | -0.8129 |
| -2.14 | -1.5894 | -1.526 |
| 3.02 | 2.2898 | 2.0599 |
| 1.02 | 0.6854 | 0.6957 |
| 0.02 | 0.0133 | 0.0136 |
| 1.02 | 0.6854 | 0.6957 |
| 1.02 | 0.6854 | 0.6957 |
| 1.18 | 0.7966 | 0.8049 |
| -2.82 | -2.0972 | -1.9235 |
| -0.82 | -0.5483 | -0.5593 |
| 2.18 | 1.5429 | 1.4869 |
| -0.82 | -0.5483 | -0.5593 |
| 0.34 | 0.236 | 0.2425 |
| 0.34 | 0.236 | 0.2425 |
| 0.34 | 0.236 | 0.2425 |
Tabela 4.3.16: Resíduos da Motivação 1.
Usando o software Action temos os seguintes resultados:
Figura 4.3.9: Resíduos Padronizados e Studentizados vs Valores Ajustados (Motivação 1)
Exemplo 3.4.2.2
Usaremos novamente o exemplo da “Motivação 2” para fazer os resíduos padronizados e studentizados.
Solução: O QME neste caso foi calculado no “Exemplo 2.3.1” e é dado por $QME=1.220,1$. Além disso, a matriz H referente ao conjunto de dados é dada por
$$ H = \begin{bmatrix} \mathbf{0,37} & 0,03 & 0,14 & -0,19 & 0,13 & 0,16 & -0,03 & 0,26 & -0,08 & 0,20 & -0,06 & 0,09 & -0,07 & 0,06 \cr 0,03 & \mathbf{0,36} & -0,19 & 0,14 & 0,12 & 0,16 & -0,03 & -0,08 & 0,25 & 0,19 & -0,07 & 0,08 & -0,07 & 0,11 \cr 0,14 & -0,19 & \mathbf{0,32} & -0,02 & 0,03 & 0,01 & 0,15 & 0,23 & -0,11 & -0,02 & 0,18 & 0,06 & 0,18 & 0,04 \cr -0,19 & 0,14 & -0,02 & \mathbf{0,31} & 0,03 & 0,00 & 0,15 & -0,11 & 0,22 & -0,03 & 0,18 & 0,06 & 0,18 & 0,09 \cr 0,13 & 0,12 & 0,03 & 0,09 & \mathbf{0,09} & 0,11 & 0,03 & 0,08 & 0,08 & 0,12 & 0,02 & 0,08 & 0,01 & 0,08 \cr 0,16 & 0,16 & 0,01 & 0,00 & 0,11 & \mathbf{0,13} & 0,00 & 0,08 & 0,08 & 0,16 & -0,02 & 0,08 & -0,03 & 0,08 \cr -0,03 & -0,03 & 0,15 & 0,15 & 0,03 & 0,00 & \mathbf{0,15} & 0,06 & 0,06 & -0,03 & 0,18 & 0,06 & 0,18 & 0,06 \cr 0,26 & -0,08 & 0,23 & -0,11 & 0,08 & 0,08 & 0,06 & \mathbf{0,24} & -0,09 & 0,09 & 0,06 & 0,08 & 0,06 & 0,05 \cr -0,08 & 0,25 & -0,11 & 0,22 & 0,08 & 0,08 & 0,06 & -0,09 & \mathbf{0,24} & 0,08 & 0,06 & 0,07 & 0,05 & 0,10 \cr 0,20 & 0,19 & -0,02 & -0,03 & 0,12 & 0,16 & -0,03 & 0,09 & 0,08 & \mathbf{0,20} & -0,06 & 0,09 & -0,07 & 0,09 \cr -0,06 & -0,07 & 0,18 & 0,18 & 0,02 & -0,02 & 0,18 & 0,06 & 0,06 & -0,06 & \mathbf{0,22} & 0,06 & 0,23 & 0,06 \cr 0,09 & 0,08 & 0,06 & 0,06 & 0,08 & 0,08 & 0,06 & 0,08 & 0,07 & 0,09 & 0,06 & \mathbf{0,07} & 0,06 & 0,07 \cr -0,07 & -0,07 & 0,18 & 0,18 & 0,01 & -0,03 & 0,18 & 0,06 & 0,05 & -0,07 & 0,23 & 0,06 & \mathbf{0,23} & 0,06 \cr 0,06 & 0,11 & 0,04 & 0,09 & 0,08 & 0,08 & 0,06 & 0,05 & 0,10 & 0,09 & 0,06 & 0,07 & 0,06 & \mathbf{0,08} \end{bmatrix} $$
Desta forma, os resíduos são:
| Resíduos | Resíduos Studentizados | Resíduos Padronizados |
|---|---|---|
| 30.3464 | 1.103 | 1.0923 |
| 15.4769 | 0.5347 | 0.553 |
| -30.3643 | -1.0574 | -1.0518 |
| -23.2338 | -0.7869 | -0.8009 |
| 5.3414 | 0.1532 | 0.1605 |
| -11.8735 | -0.3503 | -0.3652 |
| 63.201 | 2.3162 | 1.9598 |
| -25.009 | -0.8099 | -0.8229 |
| -19.8784 | -0.6326 | -0.6506 |
| -37.0883 | -1.2094 | -1.1847 |
| -44.5841 | -1.5269 | -1.4421 |
| 24.5563 | 0.7137 | 0.7302 |
| 37.4588 | 1.2563 | 1.2245 |
| 15.6505 | 0.4491 | 0.4663 |
Tabela 4.3.18: Resíduos da Motivação 2.
Usando o software Action temos os seguintes resultados:
Figura 4.3.10: Resíduos Padronizados e Studentizados vs Valores Ajustados da Motivação 2
| Observação | t-Valor | P-valor | P-valor Bonferroni |
|---|---|---|---|
| 7 | 2.3162 | 0.0431 | 0.6027 |
Tabela 4.3.19: Outliers (pontos atípicos) da Motivação 2
3.4.3 Pontos Influentes
Um ponto é influente se sua exclusão do ajuste da regressão causa uma mudança substancial na análise de regressão, po exemplo, nos valores ajustados ou nas estimativas dos coeficientes do modelo. Por isso, técnicas foram desenvolvidas para identificar essas observações influentes.
3.4.3.1 DFFITS
A medida DFFITS mede a influência que a observação $i$ tem sobre seu próprio valor ajustado. Neste caso, medimos a influência da exclusão da i-ésima observação no seu valor previsto ou ajustado. Consideremos a medida
$$DFFITS_{(i)}=\frac{\hat{Y}_i-\hat{Y}_{i(i)}}{\sqrt{QME_{(i)}h_{ii}}},$$
isto é, a diferença dos valores preditos de $Y_i$ com e sem a observação $i$ (se $i$ está entre parênteses, significa que excluímos observação $i$), expressa em unidades de desvios padrão dos valores preditos de $Y_i$. Assim, essa técnica mede o quanto a exclusão da observação $i$ aumenta ou diminui seu valor predito. Dizemos que uma observação é um ponto é influente segundo a medida DFFITS se:
(1) $|DFFITS_{(i)}|~>~ 1$, para amostras pequenas ou médias.
(2) $|DFFITS_{(i)}|~>~ 2\sqrt{(p+1)/n}$, para amostras grandes, no qual $(p+1)$ é o número de parâmetros.
3.4.3.2 DFBETA
A medida DFBETA mede a influência da observação $i$ sobre o coeficiente de $X_j$. Esta medida de influência é deinida por $$DFBETA_{j(i)}=\frac{\hat{\beta}_j-\hat{\beta}_{j(i)}}{\sqrt{QME_{i}c_{jj}}},~~~~~~~~j=0,1,…,p,$$
em que $c_{jj}$ é o j-ésimo elemento da diagonal de $(X^\prime X)^{-1}$. Valor alto para a medida de influência DFBETA nos indica que a observação i influência na estimativa do coeficiente angular da variável explicativa $X_j$. São consideradas observações influentes aquelas que
(1) $|DFBETA|~>~1$, para amostras pequenas.
(2) $|DFBETA|~>~2/\sqrt{n}$, para amostras grandes.
3.4.3.3 Distância de Cook
A distância de Cook mede a influência da observação $i$ sobre todos $n$ valores ajustados $\hat{Y}_i$. Esta medida de influência é definida por:
$$D_i=\frac{e_i^2}{(p+1)QME}\frac{h_{ii}}{(1-h_{ii})^2}.$$
Percebemos que $D_i$ é grande quando ou o resíduo $e_i$ é grande, a leverage $h_{ii}$ é grande ou ambos. Destacamos as observações quando $D_i~>~1$.
Exemplo 3.4.3.1
Vamos verificar se as observações da “Motivação 1” são pontos influentes. Vale ressaltar que na análise dos resíduos estamos considerando o modelo $$\hbox{Dureza}=\beta_0+\beta_1\hbox{Temperatura}.$$
Solução:
| Observações | DFFIT | DF $\beta_0$ | DF $\beta_1$ | DCOOK |
|---|---|---|---|---|
| 1 | -0,039 | -0,032 | 0,031 | 8,11E-04 |
| 2 | -0,039 | -0,032 | 0,031 | 8,11E-04 |
| 3 | -0,039 | -0,032 | 0,031 | 8,11E-04 |
| 4 | -0,325 | -0,265 | 0,260 | 5,38E-02 |
| 5 | -0,641 | -0,523 | 0,514 | 1,90E-01 |
| 6 | 0,579 | 0,249 | -0,236 | 1,35E-01 |
| 7 | 0,173 | 0,075 | -0,071 | 1,54E-02 |
| 8 | 0,003 | 0,001 | -0,001 | 5,94E-06 |
| 9 | 0,173 | 0,075 | -0,071 | 1,54E-02 |
| 10 | 0,173 | 0,075 | -0,071 | 1,54E-02 |
| 11 | 0,201 | -0,078 | 0,082 | 2,07E-02 |
| 12 | -0,530 | 0,204 | -0,216 | 1,18E-01 |
| 13 | -0,139 | 0,053 | -0,057 | 9,98E-03 |
| 14 | 0,390 | -0,150 | 0,159 | 7,06E-02 |
| 15 | -0,139 | 0,053 | -0,057 | 9,98E-03 |
| 16 | 0,095 | -0,075 | 0,076 | 4,78E-03 |
| 17 | 0,095 | -0,075 | 0,076 | 4,78E-03 |
| 18 | 0,095 | -0,075 | 0,076 | 4,78E-03 |
| 19 | -0,831 | 0,654 | -0,667 | 2,93E-01 |
| 20 | 0,095 | -0,075 | 0,076 | 4,78E-03 |
Tabela 4.3.20: Medidas de Influência do exemplo na “Motivação 1”.
Usando o software Action temos os seguintes resultados:
Figura 4.3.11: Gráficos com os valores de DFFITS, D-COOK e DFBETAS considerando os dados da Motivação 1.
Pelos resultados da Tabela 4.3.20 e pela Figura 4.3.11 temos que nenhum DFFITS, D-COOK e DFBETAS é, em módulo, maior do que 1. Assim, temos que nenhuma observação do exemplo na “Motivação 1” é um ponto influente.
Exemplo 3.4.3.2
Usaremos novamente o exemplo da “Motivação 2” para verificar se as observações são pontos influentes. Vale a pena ressaltar que na análise dos resíduos estamos considerando o seguinte modelo
$\hbox{Ganho} = \beta_0 + \beta_1 \hbox{Dose de íons} + \beta_2 \hbox{Tempo}.$
Solução:
| Observações | DFFITS | DF $\beta_0$ | DF $\beta_1$ | DF $\beta_2$ | DCOOK |
|---|---|---|---|---|---|
| 1 | 0,841 | 0,745 | -0,492 | -0,573 | 0,231 |
| 2 | 0,399 | 0,055 | -0,236 | 0,268 | 0,057 |
| 3 | -0,720 | -0,018 | -0,352 | 0,526 | 0,171 |
| 4 | -0,528 | 0,418 | -0,262 | -0,382 | 0,096 |
| 5 | 0,049 | 0,022 | -0,023 | -0,001 | 0,001 |
| 6 | -0,137 | -0,084 | 0,094 | 0,002 | 0,007 |
| 7 | 0,964 | -0,535 | 0,692 | -0,010 | 0,222 |
| 8 | -0,459 | -0,257 | 0,037 | 0,384 | 0,072 |
| 9 | -0,351 | 0,130 | 0,027 | -0,291 | 0,043 |
| 10 | -0,599 | -0,417 | 0,478 | 0,008 | 0,115 |
| 11 | -0,803 | 0,518 | -0,657 | 0,007 | 0,192 |
| 12 | 0,200 | 0,036 | -0,029 | -0,004 | 0,014 |
| 13 | 0,693 | -0,456 | 0,577 | -0,005 | 0,152 |
| 14 | 0,130 | 0,005 | -0,018 | 0,029 | 0,006 |
Tabela 4.3.21: Medidas de Influência do exemplo na “Motivação 2”.
Usando o software Action temos os seguintes resultados:
Figura 4.3.12: Gráficos com os valores de DFFITS, D-COOK e DFBETAS considerando os dados da Motivação 2.
Pelos resultados da Tabela 4.3.21 e pela Figura 4.3.12 temos que nenhum DFFITS, D-COOK e DFBETAS é, em módulo, maior que 1. Assim, temos que nenhuma observação do exemplo na “Motivação 2” é um ponto influente.
3.5 Teste da Falta de Ajuste (Lack of Fit)
Neste módulo apresentamos um teste estatístico para avaliar a falta de ajuste em um modelo de regressão. O teste assume que a normalidade, independência e homocedasticiade da variância dos resíduos sejam válidos. Assim, após o ajuste, é importante verificar se o modelo linear é adequado. As Figuras 3.5.1 e 3.5.2 exemplificam modelos de regressão que ficaram bem ajustados e modelos com problema de ajuste.
Na Figura 4.3.13, o valor esperado de $Y$ dado $X$ ($E[Y\mid X]$) está alinhado com a média amostral de cada nível. Nesta situação, dizemos que a falta de ajuste de um modelo linear praticamente não existe, pois a variação da amostra em torno da reta de regressão é o erro aleatório devido a variação das observações das replicas.
Figura 4.3.13: Reta de regressão perfeitamente ajustada sem Falta de Ajuste
Na Figura 4.3.14 a reta de regressão não se ajusta perfeitamente aos dados e existe uma variância grande em torno do ajuste equivocado, neste caso observamos uma falta de ajuste. Observe que a reta não passa pela média de cada réplica.
Figura 4.3.14: Reta de regressão com Falta de Ajuste
Uma maneira formal de verificar o ajuste de um modelo linear é por meio do teste de falta de ajuste.
Esse teste requer medidas repetidas para um ou mais níveis de $X$. A seguir apresentamos o teste tanto para o modelo linear quanto para o múltiplo.
3.5.1 Análise da falta de ajuste nos modelos de Regressão Linear Simples
O modelo linear simples é:
$$Y_{ij}=\beta_0+\beta_1~x_i + \varepsilon \begin{cases}i=1,2,\ldots,m\cr j=1,2,\ldots,n_i \end{cases} \tag{3.5.1}$$
em que,
-
$Y_{ij}$ : representa a j-ésima observação para o i-ésimo valor da variável;
-
$\beta_0$ : onde a reta intercepta o eixo y;
-
$\beta_1$ : representa a inclinação da reta de regressão;
-
$x_i$: representa o i-ésimo valor da variável explicativa;
-
$\varepsilon_{ij}$: representa o erro aleatório associado à i-ésima e j-ésima observação;
-
$n_i$ representa o número de observações para o i-ésimo valor de x.
Supondo que temos $m$ diferentes valores da variável explicativa $(x_1,~x_2,~\ldots,~x_m)$ e que temos $n_i$ réplicas da variável resposta para cada valor da variável explicativa, ou seja,
$$\begin{cases}x_1~~\Longrightarrow~~Y_{11}~~Y_{12}~~Y_{13}~~\ldots~~Y_{1n_1} \cr x_2~~\Longrightarrow~~Y_{21}~~Y_{22}~~Y_{23}~~\ldots~~Y_{1n_2}\cr x_3~~\Longrightarrow~~Y_{31}~~Y_{32}~~Y_{33}~~\ldots~~Y_{1n_3} \cr \vdots~~~~~~\vdots~~~~~~\vdots~~~~~~\vdots~~~~~~\vdots~~~~~~\ddots~~~~~~\vdots \cr x_m~~\Longrightarrow~~Y_{m1}~~Y_{m2}~~Y_{m3}~~\ldots~~Y_{mn_m} \end{cases}$$
Vamos quebrar a soma de quadrados do erro em dois componentes.
Para entendermos a quebra, tomamos,
$$Y_{ij}-\hat{Y}_i= (Y_{ij}-\bar{Y}_{i.})+(\bar{Y}_{i.}-\hat{Y}_{i}).$$
A primeira parte reflete a variabilidade entre as observações da variável resposta para o mesmo $x_i$, enquanto que a segunda parte reflete os desvios das médias das observações (em cada $x_i$) para o modelo.
Assim, associamos a primeira parte ao erro puro (PE) e a segunda parte à falta de ajuste do modelo (lof).
$$SQE=\sum_{i=1}^{m}\sum_{j=1}^{n_i} (Y_{ij}-\hat{Y}_{i})^2$$
$$=\sum_{i=1}^{m}\sum_{j=1}^{n_i}(Y_{ij}-\bar{Y}_{i.}+\bar{Y}_{i.}-\hat{Y}_i)^2$$
$$=\sum_{i=1}^{m}\sum_{j=1}^{n_i}( (Y_{ij}-\bar{Y}_{i.})-(\hat{Y}_i - \bar{Y}_{i.}))^2$$
$$=\sum_{i=1}^{m}\sum_{j=1}^{n_i}(Y_{ij}-\bar{Y}_{i.})^2-2~\sum_{i=1}^{m}\sum_{j=1}^{n_i} (Y_{ij}-\bar{Y}_{i.})(\hat{Y}_{i}-\bar{Y}_{i.})+\sum_{i=1}^{m}\sum_{j=1}^{n_i}(\hat{Y}_{i}-\bar{Y}_{i.})^2$$
$$=\sum_{i=1}^{m}\sum_{j=1}^{n_i} (Y_{ij}-\bar{Y}_{i.})^2+\sum_{i=1}^{m}n_i(\hat{Y}_{i}-\bar{Y}_{i.})^2~,~~{pois}~~2~\sum_{i=1}^{m}\sum_{j=1}^{n_i} (Y_{ij}-\bar{Y}_{i.})(\hat{Y}_{i}-\bar{Y}_{i.})~=0$$
$$=SQ_{PE}+SQ_{LOF}.$$
Assim,
$$SQ_{PE}=\sum_{i=1}^{m}\sum_{j=1}^{n_i}(Y_{ij}-\bar{Y}_{i.})^2~~~{ e }~~~SQ_{LOF} =\sum_{i=1}^{m}n_i(\hat{Y}_{i}-\bar{Y}_{i.})^2.$$
O quadrado médio devido ao Falta de Ajuste é dado por
$$QM_{LOF}=\frac{SQ_{LOF}}{m-2}.$$
Assim o teste para a falta de ajuste é dada por
$$F_0=\frac{QM_{LOF}}{QM_{PE}}=\frac{\frac{(SQ_{LOF})}{(m-2)}}{\frac{SQ_{PE}}{(n-m)}}~~~~\sim F_{(m - 2; n - m)}$$
Temos que
$$E[QM_{LOF}] =\sigma^2 + \frac{\displaystyle\sum_{i=1}^{m}n_i [E(Y_i) - \beta_0 -\beta_1~x_i]^2}{m-2}. \tag{3.5.2}$$
Se a regressão linear for perfeita obtemos que $E(Y_i)=\beta_0+\beta_1~x_i$ e o segundo termos da equação (3.5.2) é zero. Neste caso,
$$E[QM_{LOF}]=\sigma^2$$
Se o modelo linear não for adequado, obtemos que $E(Y_i)\neq\beta_0+\beta_1~x_i$. Neste caso,
$$E[QM_{LOF}]> \sigma^2$$
Seja,
$$\begin{cases} H_0:E(Y_i)=\beta_0+\beta_1~x_i \hbox{ modelo linear adequado} \cr H_1: E(Y_i) \neq \beta_0 + \beta_1~x_i \hbox{ modelo linear inadequado} \end{cases}$$
Se $H_0$ é verdadeiro, obtemos que $F_0 \sim F_{(m - 2; n - m)}.$
Com isso, rejeitamos $H_0$ se $F_0> F_{(\alpha,m - 2; n - m)}.$
O P-valor é dado por
$$\text{P-valor}=P[F_{(m - 2; n - m)}> F_0].$$
A tabela ANOVA se resume em:
| Fonte | SQ | GL | QM | Estatística |
|---|---|---|---|---|
| Regressão | $SQR$ | 1 | $SQR$ | $\frac{SQR}{\frac{SQE}{n-2}}$ |
| Resíduo | $SQE$ | n-2 | $\frac{SQE}{n-2}$ | |
| Falta de Ajuste | $SQE - SQE(\text{puro})$ | m-2 | $\dfrac{SQE - SQE(\text{puro})}{m-2}$ | $\frac{\dfrac{SQE - SQE(\text{puro})}{m-2}}{\dfrac{SQE(\text{puro})}{n-m}}$ |
| Erro Puro | $SQE(\text{puro})$ | n-m | $\dfrac{SQE(\text{puro})}{n-m}$ | |
| Total | $SQT$ | n-1 |
Tabela 4.3.22: ANOVA para o teste de linearidade da regressão.
Exemplo 3.5.1
Vamos fazer a aplicação do Teste da Falta de Ajuste ao exemplo da “Motivação 1”. Nesse caso de regressão simples, temos quatro níveis $(m=4)$ da variável explicativa temperatura, com quatro replicações cada um.
Solução:
$$S_{xx}=\sum_{i=1}^{n}\sum_{j=0}^{n_i}(x_{ij}-\bar{x_{i.}})^2= 625$$
$$S_{yy}=\sum_{i=1}^{n}\sum_{j=0}^{n_i}(y_{ij}-\bar{y_{i.}})^2= 706,8$$
$$S_{xy}=\sum_{i=1}^{n}\sum_{j=0}^{n_i}(x_{ij}-\bar{x_{i.}})(y_{ij}-\bar{y_{i.}})=-645$$
Temos inicialmente $$(y_{ij}-\hat{y}_i)=(y_{ij}-\bar{y}_{i.})-(\bar{y}_{i.}-\hat{y}_{i}).$$
A partir daí podemos chegar na soma de quadrados separada para o Teste da Falta de Ajuste.
$$SQE=\sum_{i=1}^{m}\sum_{j=1}^{n_i}(y_{ij}-\bar{y}_{i.})^2 +\sum_{i=1}^{m}\sum_{j=1}^{n_i}(\hat{y}_{i}-\bar{y}_{i.})^2$$
$$=30,40+10,76=41,16.$$
$$QM_{LOF}=\frac{SQ_{LOF}}{m-2}=\frac{10,76}{2}=5,38.$$
$$QM_{PE}=\frac{SQ_{PE}}{n-m}=\frac{30,4}{16}=1,9.$$
Podemos então calcular a estatística F com base nos quadrados médios, sabendo que neste exemplo os valores de $m = 4$ e $n = 20.$
$$F_0=\frac{QM_{LOF}}{QM_{PE}}=\frac{5,38}{1,9}= 2,83~~~~\sim ~~~~F_{(2;16)}$$
Se $H_0$ é verdadeiro, obtemos que $F_0 \sim F_{(2;16)}.$ Com isso, rejeitamos $H_0$ se $F_0> F_{(2,83;2; 16)}.$
O p-valor é dado por
$$\text{P-valor}=P[F_{(2; 16)}> F_0]= 0,089.$$
| Fonte | SQ | GL | QM | Estatística | P-valor |
|---|---|---|---|---|---|
| Regressão | 665,64 | 1 | 665,64 | 350,3368 | 0 |
| Resíduo | 41,16 | 18 | 2,2867 | ||
| Falta de Ajuste | 10,76 | 2 | 5,38 | 2,8316 | 0,0885 |
| Erro Puro | 30,4 | 16 | 1,9 | ||
| Total | 706,8 | 19 |
Tabela 4.3.23: Teste de Falta de Ajuste (Motivação 1).
Portanto, não rejeitamos a hipótese de que o modelo linear é adequado.
3.5.2 Análise da falta de ajuste nos modelos na Regressão Linear Múltipla
No caso da regressão linear múltipla, nosso modelo é dado por
$$Y_{ij}=\beta_0+\beta_i~x_{ij}+\varepsilon_{ij} ~~\begin{cases}i=1,2,\ldots,m \cr j=1,2,\ldots,n_i \cr \end{cases} \tag{3.5.3}$$
em que,
-
$Y_{ij}$ : representa a j-ésima observação para o i-ésimo valor da variável;
-
$\beta_0$ : onde o hiperplano intercepta o eixo y;
-
$\beta_i$ : é o coeficiente da i-ésima variável regressora ($x_i$);
-
$x_{ij}$: o j-ésimo valor da variável explicativa $i$;
-
$\varepsilon_{ij}$: erro aleatório associado à j-ésima observação;
-
$n_i$ representa o número de observações para o i-ésimo valor de x.
Como temos $n_i$ observações para cada combinação i das variáveis explicativas, $x_i$, temos $n=\displaystyle\sum_{i=1}^{n}n_i$ observações.
O procedimento do teste envolve particionar a soma dos quadrados dos resíduos em dois componentes (da mesma forma que no caso linear simples),
$$SQ_E = SQ_{PE}+SQ_{LOF}.$$
em que $SQ_{PE}$ é a soma dos quadrados devido ao Erro Puro e $SQ_{LOF}$ é a soma dos quadrados devido à Falta de Ajuste.
Para desenvolver esta partição de $SQ_E$, note que o $(ij)$-ésimo resíduo é,
$$e_{ij}=Y_{ij}-\hat{Y}_i= (Y_{ij}-\bar{Y}_{i.})+(\bar{Y}_{i.}-\hat{Y}_{i}). \tag{3.5.4}$$
em que $\bar{Y}_i$ é a média de $n_i$ observações em $x_i$.
Com os dois lados elevados ao quadrado da equação (3.5.4) e somando através de i e j, temos:
$$\sum_{i=1}^m \sum_{j=1}^{n_i} (Y_{ij}-\hat{Y}_i)^2=\sum_{i=1}^m\sum_{j=1}^{n_i}(Y_{ij}-\bar{Y}_i)^2+\sum_{i=1}^{m}n_i(\bar{Y}_i-\hat{Y})^2. \tag{3.5.5}$$
O lado esquerdo da equação (3.5.5) é a soma dos quadrados dos resíduos tradicional. Os dois componentes do lado direito medem o puro erro e o Falta de Ajuste.
A soma dos quadrados do puro erro,
$$SQ_{PE}=\sum_{i=1}^m\sum_{j=1}^{n_i}(Y_{ij}-\bar{Y}_i)^2, \tag{3.5.6}$$
é obtido através da soma dos quadrados corrigidos das observações repetidas em cada nível $x_i$ e então somamos para todos os m níveis de x.
Como temos $(n_i - 1)$ graus de liberdade para a soma dos quadrados do puro erro em cada nível $x_i$, o total de graus de liberdade associados com a soma dos quadrados do puro erro é,$\displaystyle\sum_{i=1}^m (n_i -1)=n-m.$
A soma dos quadrados da Falta de Ajuste,
$$\sum_{i=1}^{m}n_i(\bar{Y}_i-\hat{Y})^2. \tag{3.5.7}$$
é uma soma dos desvios quadráticos ponderada entre a resposta média $\bar{Y}_i$ em cada nível e seu
correspondente valor ajustado.
Se o valor ajustado $\hat{Y}_i$ está próximo das correspondentes médias respostas $\bar{Y}_i$ então existe um grande indicativo de que o modelo linear é adequado. Se $\hat{Y}_i$ desvia muito de $\bar{Y}_i$, este por sua vez é um forte indício de que a regressão não é linear.
O valor esperado de $QM_{PE}$ é $\sigma^2$, e o valor esperado do $QM_{LOF}$ é dado por:
$$E(QM_{LOF})=\sigma^2+\frac{\sum_{i=1}^m n_i \left[ E(Y_i)-\beta_0-\sum_{j=1}^k \beta_j x_{ij} \right]^2}{m-p-1}.$$
em que $p$ é o número de variáveis explicativas da regressão.
O teste estatístico para o Falta de Ajuste é,
$$F_0=\frac{\frac{SQ_{LOF}}{(m-p-1)}}{\frac{SQ_{PE}}{(n-m)}}=\frac{QM_{LOF}}{QM_{PE}}.$$
Se a regressão linear for perfeita obtemos que $E(Y_i)=\beta_0+\beta_1~x_i$ e o segundo termos da equação é zero. Neste caso, $$E[QM_{LOF}]=\sigma^2.$$
Se o modelo linear não for adequado, obtemos que $E(Y_i)\neq\beta_0+\beta_1~x_i$. Neste caso, $$E[QM_{LOF}]> \sigma^2.$$
Seja,
$$\begin{cases} H_0:E(Y_i)=\beta_0+\beta_1~x_i \hbox{ modelo linear adequado} \cr H_1: E(Y_i) \neq \beta_0 + \beta_1~x_i \hbox{ modelo linear inadequado} \end{cases}$$
Se $H_0$ é verdadeiro, obtemos que $F_0 \sim F_{(m - p-1; n - m)}.$ Com isso, rejeitamos $H_0$ se $F_0 > F_{(\alpha,m - p-1; n - m)}$. O P-valor é dado por:
$$\text{P-valor}=P[F_{(m - p -1; n - m)}> F_0].$$
A Tabela ANOVA se resume em
| Fonte | SQ | GL | QM | Estatística |
|---|---|---|---|---|
| Regressão | $SQR$ | p | $QMR$ | $\frac{SQR}{s^2}$ |
| Resíduo | $SQE$ | n-p-1 | ||
| Falta de Ajuste | $SQE - SQE(puro)$ | m-p-1 | $\frac{SQE - SQE(puro)}{m-p-1}$ | $\frac{SQE - SQE(puro)}{s^2(m-p-1)}$ |
| Erro Puro | $SQE(puro)$ | n-m | $\frac{SQE(puro)}{n-m}$ | |
| Total | $SQT$ | n-1 |
Tabela 4.3.24: ANOVA para o teste de linearidade da regressão.
Exemplo 3.5.2
Vamos fazer a aplicação do Falta de Ajuste ao problema dado na “Motivação 2”. Como necessitamos de réplicas nos níveis, vamos considerar um novo conjunto de dados conforme a Tabela 4.3.25.
| Observação | Tempo | Dose | Ganho | Observação | Tempo | Dose | Ganho |
|---|---|---|---|---|---|---|---|
| 1 | 195 | 4 | 1004 | 16 | 225 | 4,7 | 1160 |
| 2 | 195 | 4 | 1010 | 17 | 225 | 4,3 | 1276 |
| 3 | 195 | 4,6 | 852 | 18 | 225 | 4,3 | 1270 |
| 4 | 195 | 4,6 | 849 | 19 | 225 | 4,72 | 1225 |
| 5 | 195 | 4,3 | 903 | 20 | 225 | 4,72 | 1220 |
| 6 | 195 | 4,3 | 920 | 21 | 230 | 4,3 | 1321 |
| 7 | 225 | 4,2 | 1272 | 22 | 230 | 4,3 | 1330 |
| 8 | 225 | 4,2 | 1285 | 23 | 230 | 4,5 | 1340 |
| 9 | 225 | 4,1 | 1270 | 24 | 230 | 4,5 | 1345 |
| 10 | 225 | 4,1 | 1280 | 25 | 255 | 4 | 1636 |
| 11 | 225 | 4,6 | 1269 | 26 | 255 | 4 | 1640 |
| 12 | 225 | 4,6 | 1200 | 27 | 255 | 4,6 | 1506 |
| 13 | 225 | 4 | 1260 | 28 | 255 | 4,6 | 1510 |
| 14 | 225 | 4 | 1250 | 29 | 255 | 4,3 | 1555 |
| 15 | 225 | 4,7 | 1146 | 30 | 255 | 4,3 | 1550 |
Tabela 4.3.25: Dados da “Motivação 2” com réplicas
Solução:
Temos inicialmente
$$ (y_{ij}-\hat{y}_i)=(y_{ij}-\bar{y}_{i.})-(\bar{y}_{i.}-\hat{y}_{i}).$$
A partir daí podemos chegar na soma de quadrados separada para o Falta de Ajuste.
$$SQE=\sum_{i=1}^m\sum_{j=1}^{n_i}(y_{ij}-\bar{y}_{i.})^2 +\sum_{i=1}^m \sum_{j=1}^{n_i}(\hat{y}_{ij}-\bar{y}_{i})^2$$
$$=2942+28587=31529.$$
$$QM_{LOF}=\frac{SQ_{LOF}}{m-p-1}=\frac{28587}{12}= 2382.$$
$$QM_{PE}=\frac{SQ_{PE}}{n-m}=\frac{2942}{15} =196.$$
Podemos então, calcular a estatística F, com base nos quadrados médios, sabendo que neste exemplo os valores de $m=15$ e $n=30$.
$$F_0=\frac{QM_{LOF}}{QM_{PE}} = \frac{2382}{196}= 12,15.$$
Se $H_0$ é verdadeiro, obtemos que $F_0 \sim F_{(12;15)}.$
Com isso, rejeitamos $H_0$ se $F_0> F_{(12; 15)}$.
O p-valor é dado por
$$\text{P-valor}=P[F_{(12; 15)}> F_0]= 0,000.$$
Temos a seguir a Tabela ANOVA com o Falta de Ajuste dos dados do exemplo da Motivação 2, com réplicas.
| Fonte | SQ | GL | QM | Estatística | P-valor |
|---|---|---|---|---|---|
| Regressão | 1307620 | 2 | 653810 | 559,89 | 0 |
| Resíduo | 31529 | 27 | 1168 | ||
| Falta de Ajuste | 28587 | 12 | 2382 | 12,15 | 0 |
| Erro Puro | 2942 | 15 | 196 | ||
| Total | 133914 | 29 |
Tabela 4.3.26: Tabela da ANOVA com a Falta de Ajuste com os dados da Tabela 4.3.25
Usando o software Action temos os seguintes resultados:
- Sem réplicas:
| Teste de Falta de Ajuste |
|---|
| Dados de entrada não contém réplicas. |
Tabela 4.3.27: Teste de Falta de Ajuste com os dados da Motivação 2 (sem réplicas)
- Com réplicas:
| GL | Soma de Quadrados | Quadrado Médio | Estat. F | P-valor | |
|---|---|---|---|---|---|
| Tempo | 1 | 1264886.4029802 | 1264886.4029802 | 6449.1149030262 | 0 |
| Dose | 1 | 42733.8932878853 | 42733.8932878853 | 217.8818488505 | 2e-10 |
| Resíduos | 27 | 31529.1703985826 | 1167.7470517994 | ||
| Falta de Ajuste | 12 | 28587.1703985826 | 2382.2641998819 | 12.1461464984 | 0.0000124497 |
| Erro Puro | 15 | 2942 | 196.1333333333 |
Tabela 4.3.27: Teste de Falta de Ajuste com os dados da Motivação 2 com réplicas.
3.6 Análise de Colinearidade e Multicolinearidade
Quando trabalhamos com mais de uma variável regressora, é muito importante verificar se essas variáveis explicativas são correlacionadas. Desta forma, se não houver nenhum relacionamento entre elas, dizemos que são ortogonais.
Na prática, é muito difícil que as variáveis de entrada sejam ortogonais e, felizmente, a falta de ortogonalidade não é séria. Mas se as variáveis forem muito correlacionadas, as inferências baseadas no modelo de regressão podem ser errôneas ou pouco confiáveis.
Por isso, é necessário verificar se as variáveis são altamente correlacionadas. Na literatura, os termos Colinearidade (Multicolinearidade) são utilizados para indicar a existência forte de correlação entre duas (ou mais) variáveis independentes. Entretanto, alguns autores designam de Colinearidade a existência de relação linear entre duas variável explicativa (matriz de correlação) e de Multicolinearidade a existência de relação linear entre uma variável explicativa e as demais.
Se a matriz $X’X$ é singular, isto é, algumas variáveis explicativas são combinações lineares de outras, temos Multicolinearidade e não há Estimadores de Mínimos Quadrados único para os parâmetros. Se $X’X$ é aproximadamente singular, temos Multicolinearidade aproximada.
3.6.1 Colinearidade
Algumas características no ajuste do modelo podem indicar problema de colinearidade. São elas:
-
$\hat{\beta}_i$ com sinal oposto ao esperado;
-
Grandes mudanças em $\hat{\beta}_i$ quando adicionamos ou excluímos variáveis ou observações;
-
$\hat{\beta}_i$ não significante para um $X_i$ teoricamente importante.
Para diagnosticar colinearidade, temos a seguinte opção:
- Verificar se a matriz de correlações das variáveis explicativas apresenta altas correlações. Se a correlação de duas variáveis for próxima de 1, indica problema.
3.6.1.1 Matriz de Correlação
A matriz de correlação apresenta a correlação das variáveis com todas as outras em consideração. Vamos supor que temos três variáveis: $X_1,X_2$ e $X_3.$ A matriz de correlação é dada por
$$r=\begin{bmatrix}r_{(x_1,x_1)}~~r_{(x_1,x_2)}~~r_{(x_1,x_2)}\cr r_{(x_2,x_1)}~~r_{(x_2,x_2)}~~r_{(x_2,x_3)}\cr r_{(x_3,x_1)}~~r_{(x_3,x_2)}~~r_{(x_3,x_3)} \end{bmatrix},$$
em que
$$r_{x_j,x_k}=\frac{\sum_{i=1}^n(x_{ji}-\bar{x}_j)(x_{ki}-\bar{x}_k)}{(n-1)s_{x_j}s_{x_k}}.$$
$s_{x_j}$ e $s_{x_k}$ são os desvios padrão de $X_j$ e $X_k$, respectivamente.
Se a correlação de duas variáveis é maior de 0,9, é indicativo de correlação forte entre elas.
Exemplo 3.6.1.1
No exemplo na “Motivação 2” vamos verificar se existe correlação entre as variáveis Tempo e Dose de íons.
Sejam $\bar{x_1}$ e $s_{x1}$ a média e desvio padrão da variável Tempo, respectivamente e sejam $\bar{x}_2$ e $s_{x2}$ a média e o desvio padrão da variável Dose de íons. A correlação é dada por
$$r_{x_1,x_2}=\frac{\sum_{i=1}^n(x_{1i}-\bar{x}_1)(x_{2i}-\bar{x}_2)}{(n-1)s_{x_1}s_{x_2}},$$
em que $\bar{x_1}$, $=225,36$, $\bar{x_2}=4,34$, $s_{x_1}=20,42$ e $s_{x_2}=0,26$.
$$r_{x_1,x_2}=\frac{\sum_{i=1}^n(x_{1i}-225,36)(x_{2i}-4,34)}{(14-1)(20,42)(0,26)}.$$
$$r_{x_1,x_2}= -0,00265.$$
A matriz de correlação é
$$r=\begin{bmatrix} 1~~~~~-0,002 \cr -0,002~~~~~~~1 \end{bmatrix}.$$
Assim, a correlação entre as variáveis Tempo e Dose de íons é baixa.
Usando o software Action temos os seguintes resultados:
Colinearidade
| Tempo | Dose de ions | |
|---|---|---|
| Tempo | 1 | -0.0026478463 |
| Dose.de.ions | -0.0026478463 | 1 |
Tabela 4.3.28: Matriz de colinearidade das variáveis Tempo e Dose
Figura 4.3.15: Gráfico da colinearidade
3.6.2 Multicolinearidade
A multicolinearidade é um problema no ajuste do modelo que pode causar impactos na estimativa dos parâmetros. Podemos diagnosticar Multicolinearidade por meio do VIF (Variance Inflation Factor).
3.6.2.1 VIF
Os elementos da diagonal principal de $(X^\prime X)^{-1}$ são também úteis para detectar multicolinearidade. O j-ésimo elemento da diagonal principal $(X^\prime X)^{-1}$, $C_{jj}$ pode ser escrito como $$C_{jj}=(1-R^2_j)^{-1},~~~~~~~~j=1,…,p.$$
em que $R_j^2$ é o coeficiente de determinação da regressão de $X_j$ sobre as outras variáveis explicativas.
$C_{jj}$ é chamado de fator de inflação da variância e outra notação usada é $VIF_j$. Assim, o $VIF_j$ é dado por
$$VIF_j=\frac{1}{1-R^2_j}.$$
Verificamos que $VIF_j$ mede o quanto a variância do coeficiente $\hat{\beta}_j$ é inflacionada por sua colinearidade.
Geralmente, o VIF é indicativo de problemas de multicolinearidade se VIF > 10.
Exemplo 3.6.2.1
Ainda considerando o exemplo na “Motivação 2”, vamos calcular o VIF das variáveis Dose e Tempo.
- Cálculo do VIF da variável Tempo:
O $R^2$ da variável Tempo considerando Dose de íons como variável explicativa é 0,000007. Assim, o VIF da variável Tempo é
$$VIF=\frac{1}{1-0,000007}=1,000007.$$
- Cálculo do VIF da variável Dose de íons
O $R^2$ da variável Dose de íons considerando Tempo como variável explicativa é 0,000007. Assim, o VIF da variável Dose de íons é
$$VIF=\frac{1}{1-0,000007}=1,000007.$$
Observe que os valores do VIF para as duas variáveis é o mesmo. Isso acontece pois o VIF mede a correlação da variável com todas as outras do modelo. Como temos apenas essas duas variáveis no modelo ajustado, obviamente o VIF das duas é o mesmo.
Como o VIF é menor do que 10, não temos o problema de multicolinearidade no exemplo da “Motivação 2”.
Usando o software Action temos os seguintes resultados:
| Tempo | Dose de ions | |
|---|---|---|
| VIF | 1 | 1 |
Tabela 4.3.29: Multicolinearidade das variáveis Tempo e Dose