# rm(list=ls())
options(OutDec = ",")
#==============================================================================
# Exemplo:
# Multicolinearidade (estimadores de minimos quadrados)
# Modelo:
# Emprego(y)=beta0+beta1*ano+beta2*preco+beta3*PIB+beta4*militar+epsilon
# preco = deflator do PIB.
#==============================================================================
dataxy <- read.table("../dados/d05_longley.txt",header=T)
#==============================================================================
# Correlacao entre as variaveis
#==============================================================================
print(cor(dataxy))
Employ Price GNP Armed Year Employ 1,0000000 0,9708985 0,9835516 0,4573074 0,9713295 Price 0,9708985 1,0000000 0,9915892 0,4647442 0,9911492 GNP 0,9835516 0,9915892 1,0000000 0,4464368 0,9952735 Armed 0,4573074 0,4647442 0,4464368 1,0000000 0,4172451 Year 0,9713295 0,9911492 0,9952735 0,4172451 1,0000000
#==============================================================================
# Ajuste com 15 observacoes
#==============================================================================
yfit1 <- lm(Employ ~ Year + Price + GNP + Armed,data=dataxy[1:15,])
print(summary(yfit1))
Call:
lm(formula = Employ ~ Year + Price + GNP + Armed, data = dataxy[1:15,
])
Residuals:
Min 1Q Median 3Q Max
-760,49 -281,76 -16,70 82,24 1182,69
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1,459e+06 7,142e+05 2,043 0,06825 .
Year -7,218e+02 3,700e+02 -1,951 0,07965 .
Price -1,811e+02 1,355e+02 -1,336 0,21101
GNP 9,107e-02 2,026e-02 4,495 0,00115 **
Armed -7,494e-02 2,611e-01 -0,287 0,77999
---
Signif. codes: 0 ‘***’ 0,001 ‘**’ 0,01 ‘*’ 0,05 ‘.’ 0,1 ‘ ’ 1
Residual standard error: 561,6 on 10 degrees of freedom
Multiple R-squared: 0,9798, Adjusted R-squared: 0,9717
F-statistic: 121 on 4 and 10 DF, p-value: 2,006e-08
#==============================================================================
# Ajuste com 16 observacoes
#==============================================================================
yfit2 <- lm(Employ ~ Year + Price + GNP + Armed,data=dataxy)
print(summary(yfit2))
Call:
lm(formula = Employ ~ Year + Price + GNP + Armed, data = dataxy)
Residuals:
Min 1Q Median 3Q Max
-905,8 -342,7 -107,6 216,8 1437,7
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1,169e+06 8,359e+05 1,399 0,18949
Year -5,765e+02 4,335e+02 -1,330 0,21049
Price -1,977e+01 1,389e+02 -0,142 0,88940
GNP 6,439e-02 1,995e-02 3,227 0,00805 **
Armed -1,015e-02 3,086e-01 -0,033 0,97436
---
Signif. codes: 0 ‘***’ 0,001 ‘**’ 0,01 ‘*’ 0,05 ‘.’ 0,1 ‘ ’ 1
Residual standard error: 667,3 on 11 degrees of freedom
Multiple R-squared: 0,9735, Adjusted R-squared: 0,9639
F-statistic: 101,1 on 4 and 11 DF, p-value: 1,346e-08
#==============================================================================
# Numero de condicao de XtX
#==============================================================================
X <- cbind(rep(1,16),as.matrix(dataxy[,2:5]))
XtX <- t(X)%*%X
for(i in 1:5) XtX[,i] <- XtX[,i]/sqrt(sum(XtX[,i]^2))
raizcar <- svd(XtX)$d
ncond <- sqrt(max(raizcar)/min(raizcar))
print(ncond) # Numero de condicao de XtX
[1] 5564519