Regressão Linear
Matrizes
Camila Bezerra
09/10/2020
1 Regressão Linear em Forma Matricial
Modelos de regressão linear podem ser expressos na forma matricial. A VD é representada por uma matriz nx1, onde n é o número de observações. As VIs são uma matriz de dimensão nxk, k é o número de VIs + o intercepto, diante disso, a primeira coluna da matriz X é formada por 1 porque se refere à constante. Os betas são expressos numa matriz kx1 e o termo de erro numa matriz nx1.
Regressão em forma matricial
Neste tutorial vamos estimar uma regressão linear com cálculo matricial no R. Mas antes vamos passar rapidamente sobre as notações matriciais.
\[\hat{y} = \boldsymbol X\hat{\beta} + \hat{u}\]
O objetivo da regressão linear (MQO, OLS) é diminuir os resíduos. Os resíduos são a diferença entre y e o y estimado. A seguir vamos demonstrar como chegamos à equação para encontrar o beta que minimiza a soma dos e^2. Acesse essas anotações para ver isso mais detalhadamente.
A notação abaixo representa a soma dos e^2, como os erros são a diferença entre o y real e o y estimado, substituímos o termo de erro por (y-Xb) na equação a seguir:
\[\sum{e^2} = e'e = (y- \boldsymbol X \beta)'(y- \boldsymbol X \beta)\]
Ou:
\[\sum{e^2} = e'e = y'y -2\beta' \boldsymbol X'y + \beta' \boldsymbol X' \boldsymbol X \beta \]
Para encontrarmos o beta que diminui a soma dos erros ao quadrado precisamos fazer a derivada da última equação. Como resultado da derivação temos a seguinte equação:
\[-2 \boldsymbol X' y + 2\boldsymbol X' \boldsymbol X \beta =0 \]
Agora vamos isolar o beta para chegar na equação final. Que é o inverso da X transposta multiplicada por X, multiplicada pela transposta de X e y.
\[\hat{\beta} = (X'X)^{-1}X'y\]
1.0.1 Estimando os betas
A seguir aplicaremos esssa equação no R. Primeiro criaremos as matrizes com os dados da lista da aula 4.
#voto do partido do incumbente no 1 turno
y<-matrix(c(4.7, 54.3, 53.1, 23.2, 48.6, 46.9, 41.6, 1.2))
#PIB per capita(%)
X1<-c(1.373, 3.628, -1.184, 1.710, 2.835, 6.524, -0.352, 0.328)
# Presidente busca a reeleicao?
#dummy reeleicao:
##0: nao busca
##1: busca
X2<-c(0,1,0,1,0,1,0,0)
# Matriz X
# c(rep(1, 8): para criarmos a coluna com 1s da constante
X<-matrix(c(rep(1, 8),X1,X2),nrow=8)
X## [,1] [,2] [,3]
## [1,] 1 1.373 0
## [2,] 1 3.628 1
## [3,] 1 -1.184 0
## [4,] 1 1.710 1
## [5,] 1 2.835 0
## [6,] 1 6.524 1
## [7,] 1 -0.352 0
## [8,] 1 0.328 0Agora resolveremos a equação para encontrarmos o beta. Lembre-se que o beta é uma matriz kx1, desse modo, a primeira linha se refere ao beta_0 (intercepto), a segunda ao beta_1 e a terceira ao beta_2.
\[\hat{\beta} = (X'X)^{-1}X'y\]
#solve: retorna o inverso
#t: retorna a matriz transpposta
beta <- solve(t(X) %*% X) %*% t(X) %*% y
print(beta)## [,1]
## [1,] 29.048521
## [2,] 1.319131
## [3,] 7.202301Quando rodamos a regressão com a função lm() no R chegamos ao mesmo resultado.
modelo1 <- lm(y ~ X1 + X2)
print(modelo1)##
## Call:
## lm(formula = y ~ X1 + X2)
##
## Coefficients:
## (Intercept) X1 X2
## 29.049 1.319 7.2021.0.2 Multicolinearidade Perfeita
Para demonstrar como a multicolinearidade perfeita é um problema que nos impede de estimar os coeficientes da regressão, criaremos uma nova variável dummy que codifica se o Presidente não está buscando reeleição e depois vamos inserir a dummy anterior no modelo também.
# Presidente nao busca a reeleicao?
#dummy reeleicao:
##0: busca
##1: nao busca
X3<-c(1,0,1,0,1,0,1,1)
Xd <-matrix(c(rep(1,8), X1, X2, X3), nrow =8)Quando tentamos estimar os betas com as duas variáveis sobre a reeleição retorna um erro. Isto ocorre porque as duas dummies geram uma relação linear perfeita com a coluna da constante na matriz X. Visto que, a soma das duas dummies sempre dá 1, e a coluna da constante é formada por 1. Isto gera uma matriz singular e ela não é inversível, pois tem determinante nulo. Em vista disso, só manteremos a primeira dummy.
print(Xd)## [,1] [,2] [,3] [,4]
## [1,] 1 1.373 0 1
## [2,] 1 3.628 1 0
## [3,] 1 -1.184 0 1
## [4,] 1 1.710 1 0
## [5,] 1 2.835 0 1
## [6,] 1 6.524 1 0
## [7,] 1 -0.352 0 1
## [8,] 1 0.328 0 1solve(t(Xd)%*% Xd) %*% t(Xd) %*% yQuando rodamos a regressão com a função lm() o coeficiente da variável 3 retorna como NA. Justamente porque a matriz é singular e não é possível estimar os betas.
lm(y~X1+X2 + X3)##
## Call:
## lm(formula = y ~ X1 + X2 + X3)
##
## Coefficients:
## (Intercept) X1 X2 X3
## 29.049 1.319 7.202 NA1.0.3 Estimando os resíduos
Como visto anteriormente os resíduos são a diferença entre o y e o y estimado. Abaixo estamos estimando os resíduos através de matrizes e com a função lm() do R. Alcançamos os mesmos resultados das duas formas.
\[\hat{u} = y - \boldsymbol X \hat{\beta} \]
# estimando os resíduos
u_hat <-y - X%*%beta
print(u_hat)## [,1]
## [1,] -26.159688
## [2,] 13.263370
## [3,] 25.613330
## [4,] -15.306536
## [5,] 15.811742
## [6,] 2.043166
## [7,] 13.015813
## [8,] -28.281196# salvando os residuos do modelo 1 em um objeto
u_hat2 <- modelo1$residuals
# imprimindo os resultados
print(u_hat2)## 1 2 3 4 5 6 7
## -26.159688 13.263370 25.613330 -15.306536 15.811742 2.043166 13.015813
## 8
## -28.2811961.0.4 Estimando o y (y_hat)
Agora vamos estimar o y com cálculo matricial. Também podemos somar os fitted.values (X * beta) e os resíduos (y – y_hat) do modelo 1 para estimarmos o y.
\[\hat{y} = \boldsymbol X\hat{\beta} + \hat{u}\]
y_hat <- X %*% beta + u_hat
print(y_hat)## [,1]
## [1,] 4.7
## [2,] 54.3
## [3,] 53.1
## [4,] 23.2
## [5,] 48.6
## [6,] 46.9
## [7,] 41.6
## [8,] 1.2#modelo1$fitted.values = X * beta
#modelo1$residuals = resíduos
y_hat2 <- modelo1$fitted.values + modelo1$residuals
print(y_hat2)## 1 2 3 4 5 6 7 8
## 4.7 54.3 53.1 23.2 48.6 46.9 41.6 1.2