\[E_i\sim N\left(0,\sigma^2\right)\iff y_i\sim N\left(\beta_0+\beta_1 x_i,\sigma^2\right)\]
$
O que é um modelo estatístico?
Um modelo probabilístico que estabelece relações matemáticas entre algumas variáveis estatísticas
Uma forma de introduzir suposições sobre a observação de dados amostrais
O contexto
\(y\) – uma variável resposta quantitativa
\({\bf x}=(x_1,\ldots,x_k)'\) – um qualquer número de preditores, covariáveis ou variáveis explicativas
Nota
Diferentes tipos de modelos (modelos de regressão, ANOVA, . . .) são necessários para modelar diferentes tipos de variáveis respostas e para incluir diferentes tipos de covariáveis
Algumas questões importantes
Qual o grau de associação entre essas covariáveis e a variável resposta?
Qual a precisão com que podemos predizer valores não observados de \(y\), futuros ou passados?
Princípio fundamental da modelação
Variação observada em \(y\)
=
Variação explicável + Variação aleatória
Forma funcional
\[y = g(\mathbf{x}, \boldsymbol{\beta}) + E\]
onde \(E\) é chamado o erro aleatório
O modelo de regressão linear é caraterizado por:
\[g(\mathbf{x}, \boldsymbol{\beta})=\mathbf{x}' \boldsymbol{\beta}\]
com \(E\sim N(0, \sigma^2)\), não correlacionados \(\forall \mathbf{x}\)
Nota
O modelo é linear nos parâmetros!
Os dados
Um conjunto de pontos observáveis \((x_i,y_i)\), \(i=1,\ldots,n\)
Modelo de Regressão Linear Simples (MRLS) \[y_i=\beta_0+\beta_1x_i+E_i,\, i=1,\ldots,n\] \(\beta_0,\, \beta_1\) – parâmetros da curva de regressão
\(E_i\) – erro aleatório associado a \(y_i=(y\mid x=x_i)\)
Pressupostos usuais do MRLS
\(Var[E_i]=\sigma^2\) \(\iff Var[y_i]=\sigma^2\)
\(E_i\)’s não correlacionadas \(\iff\) \(y_i\)’s não correlacionadas
\(E_i\sim\) normal
\[E_i\sim N\left(0,\sigma^2\right)\iff y_i\sim N\left(\beta_0+\beta_1 x_i,\sigma^2\right)\]
Para que servem os pressupostos?
. Tornam o modelo interpretável
\(E[y\mid x]=\beta_0+\beta_1 x\) — curva de regressão
\(\beta_0=E[y\mid x=0]\)
\(\beta_1=E[y\mid x=x_0+1]-E[y\mid x=x_0]\)
. Tornam simples a inferência estatística
Pelo método da máxima verossimilhança:
\[\left\lbrace\begin{array}{lcl} \hat {\beta }_1 &=& \frac{\sum_{i = 1}^n {x_i y_i } - n\bar {x}\bar {y}}{\sum_{i = 1}^n {x_i^2 } - n\bar {x}^2}\\ \hat{\beta }_0 &=& \bar{y} - \hat{\beta }_1 \bar{x}\\ \hat{\sigma}_{MV}^2 & = & \frac{\sum_{i=1}^{n}{\left(y_i-(\hat{\beta}_0+\hat{\beta}_1 x_i)\right)^2}}{n} \end{array}\right.\]
Nota
\(E[\hat{\beta}_i]=\beta_i\), \(i=1,2\)
\(E[\hat{\sigma}_{MV}^2]=\frac{n-2}{n}\sigma^2\implies \hat{\sigma}^2=\frac{n}{n-2}\hat{\sigma}_{MV}^2\) é um estimador centrado de \(\sigma^2\) que passaremos a usar
Mais ainda, temos:
\[\hat{\beta}_1 \sim N\left(\beta_1, \frac{\sigma^2}{Sxx}\right)\] \[\hat{\beta}_0 \sim N\left(\beta_0, \sigma^2 \left(\frac{1}{n} + \frac{\bar{x}^2}{Sxx} \right) \right)\]
em que \(Sxx = \sum_{i=1}^n (x_i-\bar{x})^2\)
e
\[\frac{(n-2)\hat{\sigma}^2}{\sigma^2}\sim \chi_{(n-2)}^2\]
\(\hat{\sigma}^2\) é independente de \(\hat{\beta}_0\) e \(\hat{\beta}_1\)
Os resultados anteriores conduzem a:
\[\frac{\hat{\beta}_1-\beta_1}{\hat{\sigma}\sqrt{\frac{1}{Sxx}}} = \frac{\hat{\beta}_1-\beta_1}{se(\hat{\beta}_1)} \sim t_{(n-2)}\]
e
\[\frac{\hat{\beta}_0-\beta_0}{\hat{\sigma}\sqrt{\frac{1}{n} + \frac{\bar{x}^2}{Sxx}}} = \frac{\hat{\beta}_0-\beta_0}{se(\hat{\beta}_0)} \sim t_{(n-2)}\]
que permitem construir intervalos de confiança ou testar hipóteses sobre \(\beta_1\) e \(\beta_0\).
Para que um modelo ajustado possa ser útil é necessário avaliar:
a adequação dos pressupostos
É possível encontrar algo que ponha em causa algum dos pressupostos?
a qualidade do ajustamento aos dados observados
O modelo permite explicar grande parte da variação observada em \(y\)?
a capacidade preditiva do modelo
É possível obter predições precisas?
Os dados
Um conjunto de pontos observáveis \((\mathbf{x}_i,y_i)\), \(i=1,\ldots,n\) com \(\mathbf{x}_i=(x_{i,1}, \ldots,x_{i,k})\).
Vamos, por agora, considerar que todas as \(k\) covariáveis são numéricas (VN).
Modelo de Regressão Linear Múltipla (MRLM) \[y_i = \beta_0 + \sum_{j=1}^{k}\beta_j x_{i,j} + E_i,\,i=1,\ldots,n\]
\(E_i\) — erro aleatório associado a \(y_i=(y\mid \mathbf{x}=\mathbf{x}_i)\) com \(E_i\sim N(0, \sigma^2)\), não correlacionadas.
O modelo anterior é apenas o caso particular mais simples em que todas as covariáveis são usadas diretamente.
Regressão de 1ª ordem
\[y_i = \beta_0 + \sum_{j=1}^{k}\beta_j x_{i,j} + E_i,\,i=1,\ldots,n\]
\[\begin{pmatrix}y_1\\ \vdots\\ y_n\end{pmatrix}=\begin{pmatrix} 1 & x_{1,1} & \cdots & x_{1,k}\\ \vdots & \vdots & \ddots & \vdots\\ 1 & x_{n,1} & \cdots & x_{n,k} \end{pmatrix} \begin{pmatrix}\beta_0\\ \vdots\\ \beta_{k}\end{pmatrix}+ \begin{pmatrix}E_1\\ \vdots\\ E_n\end{pmatrix}\]
Frequentemente é útil incluir termos de maior ordem:
| 1ª ordem | potências-(1) | \(x_1,\ldots, x_{k}\) |
| 2ª ordem | potências-(2) | \(x_1^2,\ldots, x_{k}^2\) |
| interações-(1,1) | \(x_1x_2,\ldots, x_1x_{k}, x_2x_3,\ldots, x_2x_{k},\ldots\) | |
| 3ª ordem | potências-(3) | \(x_1^3,\ldots, x_{k}^3\) |
| interações-(1,2) | \(x_1x_2^2,\ldots, x_1x_{k}^2, x_2x_3^2,\ldots, x_2x_{k}^2,\ldots\) | |
| interações-(2,1) | \(x_1^2x_2,\ldots, x_1^2x_{k}, x_2^2x_3,\ldots, x_2^2x_{k},\ldots\) | |
| . . . | . . . | . . . |
ou mesmo outras transformações das covariáveis.
O MRLM pode ser escrito em forma matricial:
\[\underset{n\times 1}{\mathbf{y}} = \underset{n\times p}{\mathbf{X}}\quad \underset{p\times 1}{\boldsymbol{\beta}}+\underset{n\times 1}{\mathbf{E}}\]
\[\begin{pmatrix}y_1\\ \vdots\\ y_n\end{pmatrix}=\begin{pmatrix} 1 & x_{1,1} & \cdots & x_{1,p-1}\\ \vdots & \vdots & \ddots & \vdots\\ 1 & x_{n,1} & \cdots & x_{n,p-1} \end{pmatrix} \begin{pmatrix}\beta_0\\ \vdots\\ \beta_{p-1}\end{pmatrix}+ \begin{pmatrix}E_1\\ \vdots\\ E_n\end{pmatrix}\]
\(\mathbf{X}\) — a matriz de delineamento
Daqui em diante:
\(y\) pode representar uma transformação da variável resposta original
\(\mathbf{X}\) pode ser a matriz de delineamento de qualquer transformação das covariáveis \((\mathbb{R}^k\rightarrow \mathbb{R}^{p-1})\)
cada coluna de \(\mathbf{X}\) fica associada a um termo do modelo
Regressão polinomial
\[y_i = \sum_{j=0}^{p-1}\beta_j x_{i}^j + E_i,\,i=1,\ldots,n\]
\[\begin{pmatrix}y_1\\ \vdots\\ y_n\end{pmatrix}=\begin{pmatrix} 1 & x_{1} & x_{1}^2 & \cdots & x_1^{p-1}\\ \vdots & \vdots & \vdots & \ddots & \vdots\\ 1 & x_{n} & x_{n}^2 & \cdots & x_n^{p-1} \end{pmatrix} \begin{pmatrix}\beta_0\\ \vdots\\ \beta_{p-1}\end{pmatrix}+ \begin{pmatrix}E_1\\ \vdots\\ E_n\end{pmatrix}\]
Nota
A interpretação usual de \(\beta_i\), \(i=1,\ldots, p-1\), como sendo a variação em \(E[y\mid \mathbf{x}]\) devida ao aumento de uma unidade em \(x_i\), mantendo todas as outras covariáveis fixas, isto é \[E[y\mid x_1,\ldots, x_i=x+1,\ldots, x_{k}]-E[y\mid x_1,\ldots, x_i=x,\ldots, x_{k}]\] só é válida para o modelo de 1ª ordem.
Para qualquer outro modelo de regressão, todas as interpretações decorrem da sua superfície de resposta particular:
\[E[y\mid \mathbf{x}] = \mathbf{x}\boldsymbol{\beta}\]
Notas
A distribuição normal multivariada
\(\mathbf{y}\sim N_m(\boldsymbol{\mu}, \boldsymbol{\Sigma}),\,\mathbf{y}\in\mathbb{R}^m\)
\[f_{\mathbf{Y}}(\mathbf{y}) = (2\pi)^{-\frac{n}{2}}|\boldsymbol{\Sigma}|^{-\frac{1}{2}}\exp\left\{-\frac{1}{2}(\mathbf{y}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\mathbf{y}-\boldsymbol{\mu})\right\}\]
\(E[\mathbf{y}]=\boldsymbol{\mu}\)
\(Cov[\mathbf{y}]=\boldsymbol{\Sigma}\)
Estimação pontual
Maximizando \[\log \mathcal{L}(\boldsymbol{\beta},\sigma^2)=-\frac{n}{2}\log \sigma^2-\frac{1}{2\sigma^2}(\mathbf{y}-\mathbf{X}\boldsymbol{\beta})'(\mathbf{y}-\mathbf{X}\boldsymbol{\beta}) + C\] obtêm-se os estimadores de MV:
\(\hat{\boldsymbol{\beta}}=(\mathbf{X}'\mathbf{X})^{-1}\mathbf{X}'\mathbf{y}\)
\(\hat{\sigma}^2_{MV}=\frac{(\mathbf{y}-\mathbf{X}\hat{\boldsymbol{\beta}})'(\mathbf{y}-\mathbf{X}\hat{\boldsymbol{\beta}})}{n}\)
\(\hat{\mathbf{y}} = \hat{E}[\mathbf{y}\mid \mathbf{x}] = \mathbf{X}\hat{\boldsymbol{\beta}}\) — valores ajustados
Usando a expressão do estimador \(\hat{\boldsymbol{\beta}}\) temos:
\[\hat{\mathbf{y}}=\mathbf{X}(\mathbf{X}'\mathbf{X})^{-1}\mathbf{X}'\mathbf{y}= \mathbf{H}\mathbf{y}\]
\(\mathbf{H}_{n\times n}\) é chamada a matriz de projeção (ou matriz chapéu)
\(\mathbf{e} = \mathbf{y} - \hat{\mathbf{y}} = \mathbf{y} - \mathbf{H}\mathbf{y} = (\mathbf{I}_n-\mathbf{H})\mathbf{y}\) — resíduos
\(\mathbf{e}'\mathbf{e}=(\mathbf{y}-\mathbf{X}\hat{\boldsymbol{\beta}})'(\mathbf{y}-\mathbf{X}\hat{\boldsymbol{\beta}}) = SSE\) — soma dos quadrados dos resíduos
Sob os pressupostos do modelo:
\(\hat{\boldsymbol{\beta}} \sim N_p\left(\boldsymbol{\beta}, \sigma^2(\mathbf{X}'\mathbf{X})^{-1}\right) \implies \hat{\beta}_k \sim N\left(\beta_k, \sigma^2(\mathbf{X}'\mathbf{X})^{-1}_{k,k}\right)\)
\(\frac{SSE}{\sigma^2}=\frac{n\hat{\sigma}^2_{MV}}{\sigma^2}\sim \chi^2_{(n-p)}\)
\(\hat{\boldsymbol{\beta}}\) e \(\hat{\sigma}^2_{MV}\) são independentes
De 2. temos que \(E[SSE]=\sigma^2(n-p)\) o que mostra que \[E[\hat{\sigma}^2_{MV}]=E\left[\frac{SSE}{n}\right]=\frac{n-p}{n}\sigma^2\]
Assim, passaremos a usar o estimador \[\hat{\sigma}^2=\frac{SSE}{n-p}=MSE\]
Em particular, os resultados anteriores implicam
\[\frac{\hat{\boldsymbol{\beta}}_k-\boldsymbol{\beta}_k}{se(\hat{\boldsymbol{\beta}}_k)}\sim t_{(n-p)}\] com \(se(\hat{\boldsymbol{\beta}}_k)=\hat{\sigma}\sqrt{(\mathbf{X}'\mathbf{X})^{-1}_{kk}}\)
Estes resultados distribucionais podem ser usados para calcular intervalos de confiança ou testar hipóteses sobre cada um dos parâmetros de regressão individualmente.
Exemplo
Para analisar a procura de um certo produto foram registados os níveis médios de educação e de rendimento e ainda uma medida do grau de urbanização em cada uma de 9 regiões.
demand education income urbanization
1 167 11.2 31.9 42.2
2 174 10.6 13.2 48.6
3 161 10.6 28.7 42.6
4 162 10.4 26.1 39.0
5 141 9.3 30.1 34.7
6 175 10.8 8.5 44.5
7 164 10.7 24.3 39.1
8 174 10.0 18.6 40.1
9 186 12.0 20.4 45.9
Call:
lm(formula = demand ~ education, data = product)
Residuals:
Min 1Q Median 3Q Max
-9.06 -5.98 -2.17 5.22 15.53
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 28.9 43.7 0.66 0.530
education 13.0 4.1 3.17 0.016 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 8.66 on 7 degrees of freedom
Multiple R-squared: 0.59, Adjusted R-squared: 0.531
F-statistic: 10.1 on 1 and 7 DF, p-value: 0.0157Inferências para a resposta média
O estimador de \(E[\mathbf{y}\mid \mathbf{x}_0]=\mathbf{x}'_0\boldsymbol{\beta}\) é naturalmente \(\mathbf{x}'_0\hat{\boldsymbol{\beta}}\), e:
\[\frac{\mathbf{x}'_0\hat{\boldsymbol{\beta}}-\mathbf{x}'_0\boldsymbol{\beta}}{se(\mathbf{x}'_0\hat{\boldsymbol{\beta}})}\sim t_{(n-p)}\]
com \(se(\mathbf{x}'_0\hat{\boldsymbol{\beta}})=\hat{\sigma}\sqrt{\mathbf{x}'_0(\mathbf{X}'\mathbf{X})^{-1}\mathbf{x}_0}\)
Um intervalo de confiança com nivel de confiança \(\gamma\) é dado por:
\[\mathbf{x}'_0\hat{\boldsymbol{\beta}}\pm F_{t_{(n-p)}}^{-1}\left(\frac{1+\gamma}{2}\right)\times se(\mathbf{x}'_0\hat{\boldsymbol{\beta}})\]
Cuidado com as extrapolações!
Predição de uma nova observação
O objetivo é predizer uma nova observação sob a validade do modelo, isto é
\[y_0\sim N(\mathbf{x}'_0\boldsymbol{\beta}, \sigma^2)\]
independente dos dados observados.
O estimador do seu valor esperado, \(\hat{y}_0=\mathbf{x}'_0\hat{\boldsymbol{\beta}}\), é agora chamado um preditor.
Temos ainda
\[\frac{\hat{y}_0-y_0}{se(\hat{y}_0-y_0)}\sim t_{(n-p)}\]
com \(se(\hat{y}_0-y_0)=\hat{\sigma}\sqrt{1+\mathbf{x}'_0(\mathbf{X}'\mathbf{X})^{-1}\mathbf{x}_0}\)
Um intervalo de predição ao nível de confiança \(\gamma\) é dado por
\[\hat{y}_0\pm F_{t_{(n-p)}}^{-1}\left(\frac{1+\gamma}{2}\right)\times se(\hat{y}_0-y_0)\]
Cuidado com as extrapolações!
Exemplo – banda de confiança
Exemplo – bandas de confiança e predição
A significância de um modelo de regressão
Consideremos as hipóteses
\[H_0: \forall i>0 : \beta_i = 0\] contra \[H_1: \exists i>0 : \beta_i \neq 0\]
Nota
\(H_0\) significa que nenhum dos termos contribui para explicar a variação de \(y\) sob o modelo ou
\[H_0: E[y\mid \mathbf{x}] = \beta_0\]
Sob \(H_0\), \(\hat{\beta}_0=\bar{y}\).
Pode mostrar-se que, para qualquer modelo de regressão:
\[\sum\left(y_i-\bar{y}\right)^2=\sum\left(y_i-\hat{y}_i\right)^2+\sum\left(\hat{y}_i-\bar{y}\right)^2\]
\[\iff SST=SSE+SSR\]
Sob \(H_0\),
\[F=\frac{\frac{SSR}{p-1}}{\frac{SSE}{n-p}} = \frac{MSR}{MSE}\sim \mathcal{F}_{(p-1, n-p)}\]
e \(H_0\) pode ser rejeitada se \(F_o > F^{-1}_{\mathcal{F}_{(p-1, n-p)}}(1-\alpha)\).
Relacionado com o teste anterior, já conhecíamos o coeficiente de determinação:
\[R^2=\frac{SSR}{SST}=1-\frac{SSE}{SST}\]
Notas
\(F=\frac{R^2}{1-R^2} \frac{n-p}{p-1}\)
\(SST\) é fixo e não depende do modelo
\(SSR\) não pode diminuir se se incluirem mais termos no modelo
\(R^2\) nunca diminui com a inclusão de termos!
Por causa desta inflação do \(R^2\), é preferível usar o:
Coeficiente de determinação ajustado
\[R^2_a=1-\frac{MSE}{MST}=1-\frac{(n-1) SSE}{(n-p) SST}\]
Notas
perde-se a interpretação mas \(R^2_a\) penaliza o aumento da complexidade do modelo
\(R^2_a\) pode ser negativo, ao contrário de \(R^2\)
Call:
lm(formula = demand ~ ., data = product)
Residuals:
1 2 3 4 5 6 7 8 9
0.876 -0.972 -1.492 0.762 -4.612 -5.462 -2.129 11.653 1.377
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 60.014 36.191 1.66 0.158
education 10.718 4.530 2.37 0.064 .
income -0.751 0.395 -1.90 0.116
urbanization 0.240 1.012 0.24 0.822
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 6.29 on 5 degrees of freedom
Multiple R-squared: 0.845, Adjusted R-squared: 0.753
F-statistic: 9.11 on 3 and 5 DF, p-value: 0.0181
Call:
lm(formula = demand ~ education + income, data = product)
Residuals:
Min 1Q Median 3Q Max
-5.873 -2.730 0.065 1.111 11.481
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 63.021 31.114 2.03 0.089 .
education 11.517 2.777 4.15 0.006 **
income -0.816 0.261 -3.12 0.021 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 5.77 on 6 degrees of freedom
Multiple R-squared: 0.844, Adjusted R-squared: 0.791
F-statistic: 16.2 on 2 and 6 DF, p-value: 0.00383Comparação de modelos encaixados
Um modelo R (modelo reduzido) diz-se encaixado no modelo A (modelo alargado) se o conjunto de termos de R é um subconjunto ou um subespaço do conjunto de termos de A.
Para comparar modelos encaixados
\[H_0: \text{Modelo R}\ \ \text{contra}\ \ H_1: \text{Modelo A}\] podemos usar a estatística de teste
\[F = \frac{\frac{SSE(R) - SSE(A)}{gl_R-gl_A}}{\frac{SSE(A)}{gl_A}}\stackrel{H_0}{\sim}\mathcal{F}_{(gl_R-gl_A,\ gl_A)}\]
e \(H_0\) pode ser rejeitada se \(F_o > F^{-1}_{\mathcal{F}_{\left(gl_R-gl_A,\ gl_A\right)}}(1-\alpha)\)
Comparação de 2 modelos encaixados
Testes F sequenciais
Somas de quadrados extra
Cada modelo de regressão \(M\) conduz a uma decomposição particular da SST: \[SST = SSR(M) + SSE(M)\]
Para o modelo nulo temos SSE(\(\emptyset\)) = SST. A inclusão de termos irá aumentar a SSR enquanto a SSE diminui.
Para identificar termos significantes é útil analisar as suas contribuições nesta transferência de variação.
Para isso, definem-se as somas de quadrados extra como:
\[SSR(A) - SSR(R) = SSE(R) - SSE(A)\]
Por exemplo, as somas de quadrados parciais ou sequenciais:
\[SSR(x_{j+1}\mid x_1,\ldots,x_{j}) = SSR(x_1,\ldots,x_{j+1}) - SSR(x_1,\ldots,x_{j})\]
As somas de quadrados extra permitem definir os coeficientes de determinação parciais.
Por exemplo,
\[R^2_{3\mid 1,2}=\frac{SSR(x_3\mid x_1, x_2)}{SSE(x_1, x_2)}\]
é o coeficiente de determinação parcial entre \(y\) e \(x_3\), dado que \(x_1\) e \(x_2\) estão já incluídas no modelo, que é interpretado como a proporção da variação de \(y\) não explicada pelo modelo reduzido que é explicada pela inclusão de \(x_3\).
Analysis of Variance Table
Response: demand
Df Sum Sq Mean Sq F value Pr(>F)
education 1 754 754 19.06 0.0072 **
income 1 325 325 8.21 0.0352 *
urbanization 1 2 2 0.06 0.8221
Residuals 5 198 40
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Analysis of Variance Table
Response: demand
Df Sum Sq Mean Sq F value Pr(>F)
urbanization 1 822 822 20.78 0.0061 **
education 1 116 116 2.93 0.1475
income 1 143 143 3.62 0.1157
Residuals 5 198 40
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Call:
lm(formula = demand ~ education + income + urbanization, data = product)
Residuals:
1 2 3 4 5 6 7 8 9
0.876 -0.972 -1.492 0.762 -4.612 -5.462 -2.129 11.653 1.377
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 60.014 36.191 1.66 0.158
education 10.718 4.530 2.37 0.064 .
income -0.751 0.395 -1.90 0.116
urbanization 0.240 1.012 0.24 0.822
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 6.29 on 5 degrees of freedom
Multiple R-squared: 0.845, Adjusted R-squared: 0.753
F-statistic: 9.11 on 3 and 5 DF, p-value: 0.0181
Call:
lm(formula = demand ~ urbanization + education + income, data = product)
Residuals:
1 2 3 4 5 6 7 8 9
0.876 -0.972 -1.492 0.762 -4.612 -5.462 -2.129 11.653 1.377
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 60.014 36.191 1.66 0.158
urbanization 0.240 1.012 0.24 0.822
education 10.718 4.530 2.37 0.064 .
income -0.751 0.395 -1.90 0.116
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 6.29 on 5 degrees of freedom
Multiple R-squared: 0.845, Adjusted R-squared: 0.753
F-statistic: 9.11 on 3 and 5 DF, p-value: 0.0181Covariáveis categorizadas
Quando incluímos covariáveis categorizadas é útil analisar a superfície de resposta:
\[E[y\mid \mathbf{x}] = \mathbf{X} \boldsymbol{\beta}\]
Covariáveis binárias (2 níveis)
Suponhamos que \(p=3\) e que \(x_2\) é binária (0/1). A superfície de resposta do modelo de 1ª ordem é:
\[\begin{align*} E[y\mid \mathbf{x}] & = \beta_0 + \beta_1 x_1 + \beta_2 x_2 =\\[0.25cm] & = \begin{cases}\beta_0+\beta_1 x_1, & x_2=0\\[0.25cm](\beta_0 + \beta_2)+\beta_1 x_1, & x_2=1\end{cases} \end{align*}\]
\(\leadsto\) estamos a ajustar 2 superfícies de resposta paralelas!
\(E[y\mid x_1, x_2=1] - E[y\mid x_1, x_2=0] = \beta_{2}\)
Com \(q\) covariáveis binárias teríamos \(2^q\) superfícies de resposta paralelas.
Exemplo
Exemplo
Exemplo
Exemplo
Para remover o paralelismo é necessário incluir interações:
\[\begin{align*} E[y\mid \mathbf{x}] & = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \beta_3 x_1 x_2 =\\[0.25cm] & = \begin{cases}\beta_0+\beta_1 x_1, & x_2=0\\[0.25cm](\beta_0 + \beta_2)+(\beta_1 + \beta_3) x_1, & x_2=1\end{cases} \end{align*}\]
\(E[y\mid x_1, x_2=1] - E[y\mid x_1, x_2=0] = \beta_{2} + \beta_3 x_1\)
Exemplo
Mais de 2 níveis
Suponhamos agora que \(x_2\) é uma covariável categorizada ou fator com 3 níveis. Codificando esses níveis como 0, 1 ou 2 temos:
\[\begin{align*} E[y\mid \mathbf{x}] & = \beta_0 + \beta_1 x_1 + \beta_2 x_2 =\\[0.25cm] & = \begin{cases}\beta_0+\beta_1 x_1, & x_2=0\\[0.25cm](\beta_0 + \beta_2)+\beta_1 x_1, & x_2=1\\[0.25cm](\beta_0 + 2\beta_2)+\beta_1 x_1, & x_2=2\end{cases} \end{align*}\]
Além de paralelas, as superfícies de resposta estão igualmente espaçadas e ordenadas. Isto é, em geral, indesejável.
Poderíamos associar uma variável binária a cada nível do fator, \(x_2\rightsquigarrow (x_{2,0}, x_{2,1}, x_{2,2})\), mas isso não conduziria a uma matriz de delineamento de característica completa.
Em vez disso usam-se apenas 2 variáveis indicadoras ou variáveis mudas:
| Nível | \(x_{2,1}\) | \(x_{2,2}\) |
|---|---|---|
| 0 | 0 | 0 |
| 1 | 1 | 0 |
| 2 | 0 | 1 |
Isto é chamada a codificação de nível de referência em que o nível codificado com \(\mathbf{x}=\mathbf{0}\) é o nível de referência.
\[\begin{align*} E[y\mid \mathbf{x}] & = \beta_0 + \beta_1 x_1 + \beta_{2,1} x_{2,1} + \beta_{2,2} x_{2,2} =\\[0.25cm] & = \begin{cases}\beta_0+\beta_1 x_1, & x_2=0\\[0.25cm](\beta_0 + \beta_{2,1})+\beta_1 x_1, & x_2=1\\[0.25cm](\beta_0 + \beta_{2,2})+\beta_1 x_1, & x_2=2\end{cases} \end{align*}\]
\(E[y\mid x_1, x_2=i] - E[y\mid x_1, x_2=0] = \beta_{2,i},\ i=1,2\)
Também aqui, a inclusão de interações pode remover o paralelismo.
\[\begin{align*} E[Y\mid \mathbf{x}] & = \beta_0 + \beta_1 x_1 + \beta_{2,1} x_{2,1} + \beta_{2,2} x_{2,2} + \\[0.25cm] & + \beta_{3} x_1 x_{2,1} + \beta_4 x_1 x_{2,2} \\[0.25cm] & = \begin{cases}\beta_0+\beta_1 x_1, & x_2=0\\[0.25cm](\beta_0 + \beta_{2,1})+(\beta_1 + \beta_3) x_1, & x_2=1\\[0.25cm](\beta_0 + \beta_{2,2}) + (\beta_1 + \beta_4) x_1, & x_2=2\end{cases} \end{align*}\]
\(E[y\mid x_1, x_2=i] - E[y\mid x_1, x_2=0] = \beta_{2,i} + \beta_{i+2}x_1,\ i=1,2\)
Notas
Outra aplicação de variáveis categorizadas
Regressão por troços
Seja \(x_1\) uma variável numérica e \(x_2=I_{[c, +\infty[}(x_1)\).
Considere
\[\begin{align*} E[Y\mid \mathbf{x}] & = \beta_0 + \beta_1 x_1 + \beta_2 (x_1-c)x_2 =\\[0.25cm] & = \begin{cases}\beta_0+\beta_1 x_1, & x_1<c\\[0.25cm](\beta_0 - c \beta_2)+(\beta_1 + \beta_2) x_1, & x_1\geq c\end{cases} \end{align*}\]
\(\leadsto\) modelo de regressão por troços contínuo
Considere
\[\begin{align*} E[Y\mid \mathbf{x}] & = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \beta_3 (x_1-c)x_2 =\\[0.25cm] & = \begin{cases}\beta_0+\beta_1 x_1, & x_1<c\\[0.25cm](\beta_0 + \beta_2 - c \beta_3)+(\beta_1 + \beta_3) x_1, & x_1\geq c\end{cases} \end{align*}\]
\(\leadsto\) modelo de regressão por troços descontínuo
Potenciais problemas
O preditor linear não é adequado
A estrutura do erro aleatório não é adequada
Variância não constante (heterocedasticidade)
Correlação nos erros
Não-normalidade
\(\leadsto\) inclusão de novos termos no modelo, transformações de \(y\) ou das covariáveis, tentar outros modelos.
\(\leadsto\) requer mais investigação: há uma explicação? É um erro? Ou uma fraude? Considerar a remoção com cautela.
\(\leadsto\) Estandardizar as covariáveis ou usar outros modelos de regressão. Considerar a remoção de covariáveis com cautela.
Muitos dos meios de diagnóstico são baseados nos resíduos \(\mathbf{e}=\mathbf{y}-\hat{\mathbf{y}}\):
Métodos gráficos
resíduos em função de \(\hat{\mathbf{y}}\) e \(\mathbf{x}_i\)
Histogramas, box-plots, QQ-plots, . . .
Medidas numéricas, testes de hipóteses, . . .
Estratégia usual
Começar com alguns gráficos e usar testes apenas para clarificar indicações suspeitas ou pouco claras dadas pelos gráficos.
Resíduos
Sob os pressupostos do modelo, um gráfico dos resíduos deve ter o aspeto de uma nuvem simétrica de pontos em torno de zero com variação constante.
No entanto,
\[\mathbf{e}=(\mathbf{I} - \mathbf{H})\mathbf{y}\sim N_n(\mathbf{0}, \sigma^2(\mathbf{I} - \mathbf{H}))\]
\(\implies e_i\sim N(0, \sigma^2(1-h_{ii}))\)
\(\therefore\) os resíduos são v.a.’s correlacionadas com diferentes variâncias
Resíduos Studentizados
\[r_i=\frac{e_i}{se(e_i)}=\frac{e_i}{\sqrt{MSE (1-h_{ii})}},\ i=1,\ldots, n\]
são chamados os resíduos internamente Studentizados ou estandardizados.
Nota: rstandard() no R
Seja \[d_i=y_i-\hat{y}_{-i}\] em que \(\hat{y}_{-i}\) é o valor ajustado após a remoção da iª observação (deleted residual) e \(r_i^*\) o correspondente resíduo Studentizado.
Aparentemente, o cálculo destes resíduos deverá ser dispendioso!
Felizmente, tem-se \(d_i=\frac{e_i}{1-h_{ii}}\) e,
\[\begin{align*} r_i^* & = e_i\sqrt{\frac{n-p-1}{SSE (1-h_{ii})-e_i^2}}=\\[0.25cm] & = r_i \sqrt{\frac{n-p-1}{n-p- r_i^2}}\sim t_{(n-p-1)},\ i=1,\ldots, n \end{align*}\]
são chamados os resíduos externamente Studentizados.
Nota: rstudent() no R
Exemplo
Num gráfico de resíduos devemos olhar para:
a existência de algum padrão ou um desvio sistemático de zero (A, B1 and B2)
assimetria (B3)
valores extremos (C)
alguma dependência em alguma covariável relativamente aos pontos anteriores
Exemplos
Exemplos
Outros gráficos
Observações atípicas e alavancagem (leverage)
Observações influentes: alavancagem \(\times\) resíduo
Lembremos que \(Var[e_i]=\sigma^2 (1- h_{ii})\).
\(\leadsto\) um valor elevado de \(h_{ii}\) torna a \(Var[e_i]\) pequena e, consequentemente, a superfície de regressão estimada será mais fortemente “atraída” por \((y_i, \mathbf{x}_i)\).
Chamamos a \(h_{ii}\) a alavancagem da iª observação.
A distância de Cook é definida como
\[D_i=\frac{(\hat{\mathbf{y}} - \hat{\mathbf{y}}_{-i})'(\hat{\mathbf{y}} - \hat{\mathbf{y}}_{-i})}{p MSE}=\frac{1}{p} (r^*_i)^2\frac{h_{ii}}{1-h_{ii}}\]
e é usada como uma medida agregada da influência de cada observação.
Nota Usualmente os valores de \(D_i\) são comparados com os quantis de uma distribuição \(\mathcal{F}_{(p, n-p)}\).
leverage std.residual cooks.d
High leverage, Low residual 0.4096 -0.163 0.00926
Low leverage, High residual 0.0478 3.616 0.32819
High leverage, High residual 0.4096 -3.916 5.31837