Atividade2_Grupo5
Ana Júlia UbaldoRA RA 791689, Guilherme Ziegler Trintinalia RA 791865, Luiz Antonio Oliveira Ferreira RA: 790107, Fernando Lascane Arida RA: 791837, Thalles Roberto de Jesus Godoi RA 795298
08/04/2022
Carregando a base de dados do grupo:
## Lendo dados do Grupo
dados <-read_excel("..\\Dados\\Grupo5.xlsx",sheet = "dados", col_names = TRUE)
head(dados)## # A tibble: 6 x 13
## N UA lrenda renda sexo age age2 Defi Defc Demi Demc Desupi
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 12 6.03 415 0 29 841 0 0 1 0 0
## 2 2 166 7.50 1800 1 33 1089 0 1 0 0 0
## 3 3 259 6.52 678 0 22 484 0 0 1 0 0
## 4 4 372 6.62 750 1 27 729 0 1 0 0 0
## 5 5 435 5.63 280 1 16 256 1 0 0 0 0
## 6 6 536 6.52 679 1 34 1156 1 0 0 0 0
## # ... with 1 more variable: Desupc <dbl>Percebe-se pela análise da base de dados as seguintes variáveis: lrenda, renda, sexo, age, age2, Defi, Defc, Demi, Demc, Desupi e Desupc, das quais uma inspeção prévia nos permitiu verificar Irenda como o logaritmo natural de renda, age2 como o quadrado de age e as demais são variáveis do tipo binárias.
Dessa forma, é esperado que haja multicolineariade entre age e age2, dado que elas são a mesma variável com uma transformação matemática. O que justifica o aparecimento de age2 é o modelo teórico proposto que sustenta que a renda de um trabalhador pode ser melhor predita pela inserção dessa variável quadrada.
Plotando modelos de regressão para as variáveis explicativas:
Escolhemos três modelos de análise:
Linear - Linear: definido pela variável EqLinear Log - Linear: definido pela variável EqLogLinear Log - Log: definido pela variável EqLogLog
Ressalta-se que no caso do modelo EqLogLog as variáveis age e age2 foram transformadas. Todavia, é provavel que haja multicolinearidade perfeita nesta transformação.
options(scipen = 999)
# estimando o modelo linear
EqLinear <- lm(renda~sexo+age+age2+Defi+Defc+Demi + Demc + Desupi + Desupc, data = dados)
summary(EqLinear)##
## Call:
## lm(formula = renda ~ sexo + age + age2 + Defi + Defc + Demi +
## Demc + Desupi + Desupc, data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5180.0 -817.4 -261.9 356.6 12760.5
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4412.885 2762.337 -1.598 0.113
## sexo 700.854 554.543 1.264 0.209
## age 185.414 120.664 1.537 0.128
## age2 -1.645 1.430 -1.150 0.253
## Defi 51.823 1323.167 0.039 0.969
## Defc 677.164 1429.342 0.474 0.637
## Demi 1277.346 1857.425 0.688 0.493
## Demc 1845.055 1282.179 1.439 0.153
## Desupi 1992.718 1499.507 1.329 0.187
## Desupc 5902.875 1405.007 4.201 0.0000577 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2657 on 100 degrees of freedom
## Multiple R-squared: 0.3393, Adjusted R-squared: 0.2798
## F-statistic: 5.706 on 9 and 100 DF, p-value: 0.000002316# estimando o modelo loglinear
EqLogLinear <- lm(lrenda~sexo+age+age2+Defi+Defc+Demi + Demc + Desupi + Desupc, data = dados)
summary(EqLogLinear)##
## Call:
## lm(formula = lrenda ~ sexo + age + age2 + Defi + Defc + Demi +
## Demc + Desupi + Desupc, data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.39318 -0.37310 -0.03453 0.30482 1.96971
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.5885395 0.6333276 7.245 0.0000000000917 ***
## sexo 0.3617904 0.1271413 2.846 0.00538 **
## age 0.0852934 0.0276650 3.083 0.00265 **
## age2 -0.0008854 0.0003280 -2.700 0.00815 **
## Defi -0.0786331 0.3033656 -0.259 0.79601
## Defc 0.2301542 0.3277088 0.702 0.48412
## Demi 0.2219269 0.4258564 0.521 0.60343
## Demc 0.7456115 0.2939682 2.536 0.01275 *
## Desupi 0.9307101 0.3437955 2.707 0.00798 **
## Desupc 1.7710014 0.3221294 5.498 0.0000002951616 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6092 on 100 degrees of freedom
## Multiple R-squared: 0.514, Adjusted R-squared: 0.4702
## F-statistic: 11.75 on 9 and 100 DF, p-value: 0.000000000001986# estimando o modelo loglog
EqLogLog <- lm(lrenda~sexo+log(age)+log(age2)+Defi+Defc+Demi + Demc + Desupi + Desupc, data = dados)
summary(EqLogLog)##
## Call:
## lm(formula = lrenda ~ sexo + log(age) + log(age2) + Defi + Defc +
## Demi + Demc + Desupi + Desupc, data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.42078 -0.32959 -0.02206 0.32105 1.82451
##
## Coefficients: (1 not defined because of singularities)
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.5616 0.7644 5.968 0.0000000359 ***
## sexo 0.3988 0.1286 3.100 0.00250 **
## log(age) 0.5190 0.1857 2.794 0.00622 **
## log(age2) NA NA NA NA
## Defi -0.1487 0.3078 -0.483 0.63013
## Defc 0.1766 0.3333 0.530 0.59734
## Demi 0.1536 0.4330 0.355 0.72358
## Demc 0.7145 0.2993 2.387 0.01884 *
## Desupi 0.8979 0.3503 2.563 0.01184 *
## Desupc 1.7851 0.3283 5.437 0.0000003780 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6214 on 101 degrees of freedom
## Multiple R-squared: 0.4893, Adjusted R-squared: 0.4489
## F-statistic: 12.1 on 8 and 101 DF, p-value: 0.000000000005479# colocando todas as regressões estimadas juntas para análise:
screenreg(list(EqLinear, EqLogLinear, EqLogLog), digits = 4,
stars = c(0.01, 0.05, 0.1), ci.force = TRUE,
custom.model.names = c("Modelo Linear", "Modelo Log-Linear", "Modelo Log-Log"),
caption = "Multiple model types and single row", label = "tab:3",
include.adjrs = TRUE, include.bic = TRUE)##
## ===========================================================================
## Modelo Linear Modelo Log-Linear Modelo Log-Log
## ---------------------------------------------------------------------------
## (Intercept) -4412.8846 4.5885 * 4.5616 *
## [-9826.9650; 1001.1958] [ 3.3472; 5.8298] [ 3.0635; 6.0598]
## sexo 700.8542 0.3618 * 0.3988 *
## [ -386.0295; 1787.7378] [ 0.1126; 0.6110] [ 0.1467; 0.6508]
## age 185.4144 0.0853 *
## [ -51.0833; 421.9120] [ 0.0311; 0.1395]
## age2 -1.6448 -0.0009 *
## [ -4.4485; 1.1589] [-0.0015; -0.0002]
## Defi 51.8234 -0.0786 -0.1487
## [-2541.5353; 2645.1821] [-0.6732; 0.5160] [-0.7520; 0.4547]
## Defc 677.1639 0.2302 0.1766
## [-2124.2956; 3478.6234] [-0.4121; 0.8725] [-0.4766; 0.8299]
## Demi 1277.3459 0.2219 0.1536
## [-2363.1408; 4917.8327] [-0.6127; 1.0566] [-0.6951; 1.0022]
## Demc 1845.0555 0.7456 * 0.7145 *
## [ -667.9691; 4358.0800] [ 0.1694; 1.3218] [ 0.1278; 1.3011]
## Desupi 1992.7180 0.9307 * 0.8979 *
## [ -946.2612; 4931.6971] [ 0.2569; 1.6045] [ 0.2114; 1.5845]
## Desupc 5902.8754 * 1.7710 * 1.7851 *
## [ 3149.1114; 8656.6395] [ 1.1396; 2.4024] [ 1.1416; 2.4286]
## log(age) 0.5190 *
## [ 0.1550; 0.8831]
## ---------------------------------------------------------------------------
## R^2 0.3393 0.5140 0.4893
## Adj. R^2 0.2798 0.4702 0.4489
## Num. obs. 110 110 110
## ===========================================================================
## * Null hypothesis value outside the confidence interval.Como pode ser verificado na tabela que compara os três modelos, há uma perda de significância do modelo log-linear para o modelo log-log. De fato, apareceu um Na na variável age2, o que nos levou a deixar de lado este modelo e seguir trabalhando com os dois primeiros.
Teste Reset para forma funcional
A seguir, para os modelos linear e log-linear, performamos o teste reset.
#Teste Reset Para modelo Linear Linear
resettest(EqLinear, power = 2, type = "regressor", data = dados)##
## RESET test
##
## data: EqLinear
## RESET = 0.052486, df1 = 9, df2 = 91, p-value = 1#Teste Reset Para modelo Log Linear
resettest(EqLogLinear, power = 2, type = "regressor", data = dados)##
## RESET test
##
## data: EqLogLinear
## RESET = 0.14309, df1 = 9, df2 = 91, p-value = 0.9982Para ambos os modelos, p-valor > 0,05 de modo que rejeitamos a hipótese nula de que a forma funcional não esta adequada.
Avaliando a presença de Multicolinearidade.
Para avaliar a presença de Multicolinearidade, transformamos as variáveis do data frame em numéricas e calculamos Fiv e Tol.
# transformação das variáveis para numéricas
sexo<-as.numeric(dados$sexo)
age<-as.numeric(dados$age)
age2<-as.numeric(dados$age2)
Defi<-as.numeric(dados$Defi)
Defc<-as.numeric(dados$Defc)
Demi<-as.numeric(dados$Demi)
Demc<-as.numeric(dados$Demc)
Desupi<-as.numeric(dados$Desupi)
Desupc <-as.numeric(dados$Desupc)
# calculando o valor de Fiv e Tol
vi <- data.frame(sexo,age, age2, Defi, Defc, Demi, Demc, Desupi,Desupc)
pairs(vi)
fiv<-vif(vi)
FIV<-as.vector(fiv)
TOL<-as.vector(1/fiv)
fiv_tol<-round(cbind.data.frame(FIV, TOL),4)
fiv_tol## FIV TOL
## 1 1.1582 0.8634
## 2 37.6237 0.0266
## 3 37.2482 0.0268
## 4 4.5108 0.2217
## 5 3.3171 0.3015
## 6 1.8835 0.5309
## 7 6.0455 0.1654
## 8 2.6317 0.3800
## 9 3.4161 0.2927cor(as.matrix(dados[3:8]))## lrenda renda sexo age age2 Defi
## lrenda 1.00000000 0.84718406 0.04954888 0.11129888 0.06054600 -0.3427693
## renda 0.84718406 1.00000000 -0.02950961 0.15432003 0.12332561 -0.2114688
## sexo 0.04954888 -0.02950961 1.00000000 0.06760864 0.05081268 0.2004727
## age 0.11129888 0.15432003 0.06760864 1.00000000 0.98578677 0.2643811
## age2 0.06054600 0.12332561 0.05081268 0.98578677 1.00000000 0.2717951
## Defi -0.34276929 -0.21146882 0.20047271 0.26438108 0.27179508 1.0000000cor(as.matrix(vi))## sexo age age2 Defi Defc
## sexo 1.000000000 0.06760864 0.050812676 0.20047271 0.13277604
## age 0.067608644 1.00000000 0.985786768 0.26438108 -0.08261395
## age2 0.050812676 0.98578677 1.000000000 0.27179508 -0.07385213
## Defi 0.200472709 0.26438108 0.271795083 1.00000000 -0.18823055
## Defc 0.132776037 -0.08261395 -0.073852133 -0.18823055 1.00000000
## Demi -0.134693115 -0.18916212 -0.170785879 -0.09988065 -0.07111527
## Demc 0.006919711 -0.15868794 -0.158595062 -0.40408663 -0.28771069
## Desupi -0.088924847 -0.01189371 -0.010202578 -0.15348462 -0.10928143
## Desupc -0.237047450 0.01467205 -0.005192555 -0.19635081 -0.13980227
## Demi Demc Desupi Desupc
## sexo -0.13469311 0.006919711 -0.08892485 -0.237047450
## age -0.18916212 -0.158687942 -0.01189371 0.014672048
## age2 -0.17078588 -0.158595062 -0.01020258 -0.005192555
## Defi -0.09988065 -0.404086629 -0.15348462 -0.196350808
## Defc -0.07111527 -0.287710690 -0.10928143 -0.139802266
## Demi 1.00000000 -0.152667731 -0.05798793 -0.074183183
## Demc -0.15266773 1.000000000 -0.23460148 -0.300122524
## Desupi -0.05798793 -0.234601485 1.00000000 -0.113995831
## Desupc -0.07418318 -0.300122524 -0.11399583 1.000000000Como pode ser avaliado pelas tabelas de covariância e pelo valor de TOL há apenas multicolinearidade para os parâmetros age e age2, conforme era de se esperar.
Corrigir essa questão não é algo trivial, dado que a formula funcional esta adequada e ambos os parâmetros são significativos para o modelo Log linear, este será aquele que adotaremos daqui para frente.
Além disso, a partir de agora, uma nova equação será gerada, utiizando apenas os parâmetros significativos ja econtrados: sexo, age, age2, Def, Desupi, Desupc. Entende-se que essa mudança seja melhor para prever a variável dependente renda.
Escolhemos a variável EqlogLinearStar para representar este novo modelo que pode ser estimado abaixo:
EqlogLinearStar <- lm(lrenda~sexo+age+age2+Defi + Desupi + Desupc, data = dados)
summary(EqlogLinearStar)##
## Call:
## lm(formula = lrenda ~ sexo + age + age2 + Defi + Desupi + Desupc,
## data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.39287 -0.40214 -0.04482 0.36440 2.19461
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.1496600 0.5511373 9.344 0.00000000000000213 ***
## sexo 0.3499367 0.1309819 2.672 0.008775 **
## age 0.0859690 0.0286082 3.005 0.003335 **
## age2 -0.0009051 0.0003410 -2.654 0.009212 **
## Defi -0.6175125 0.1631134 -3.786 0.000258 ***
## Desupi 0.3806861 0.2284963 1.666 0.098742 .
## Desupc 1.2188520 0.1953044 6.241 0.00000000980656488 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6384 on 103 degrees of freedom
## Multiple R-squared: 0.4503, Adjusted R-squared: 0.4183
## F-statistic: 14.06 on 6 and 103 DF, p-value: 0.00000000001233Investigando Heterocedasticia
Uma vez escolhido o novo modelo funcional, seguimos para a investigação de heterocedasticia no modelo. Fizemos uma inspeção visual que sugere presença de heterocedasticia. É preciso tomar certo cuidado, pois com exceção de age e age2 as demais variáveis explicativas são binárias e não faria muito sentido uma inspeção visual.
Se elevarmos age ao quadrado, obtemos age2 de modo que podemos utilizar essa variável diretamente para a inspeção visual.
A seguir conduzimos o teste de Breusch-Pagan obtendo como resultado a ausência de heterocedasticia, com p-valor de 0,07279 > que 0,05. Em outras palavras ficamos com a hipótese nula de que os resíduos são homocedasticos. Todavia, fizemos a correção para critérios didáticos e não obtivemos melhoras para o método GLS.
resEqStarq <-(EqlogLinearStar$residuals^2)
lmrendaq <- (EqlogLinearStar$fitted.values^2)
ageq<- age^2
par(mfrow=c(1,2))
plot(lmrendaq, resEqStarq, xlab = "Quadrado da renda estimada",
ylab = "Resíduos quadráticos",
main="Estimativas para Reg. linear")
plot(ageq,resEqStarq, xlab = "Quadrado de idade estimada",
ylab = "Resíduos quadráticos",
main="Estimativas para Reg. linear")
bptest(EqlogLinearStar)##
## studentized Breusch-Pagan test
##
## data: EqlogLinearStar
## BP = 11.55, df = 6, p-value = 0.07279kable(tidy(bptest(EqlogLinearStar)),
caption="Breusch-Pagan heteroskedasticity test")| statistic | p.value | parameter | method |
|---|---|---|---|
| 11.55019 | 0.0727894 | 6 | studentized Breusch-Pagan test |
Corrigindo Heterocedasticia
O método escolhido foi o de Mínimos Quadrados Generalizados- GLS, uma vez que o tipo de variáveis binárias impossibilitaria utilizarmos o método de Mínimos Quadrados Ponderados (MQP) - WLS.
EqlogLinearStar <- lm(lrenda~sexo+age+age2+Defi + Desupi + Desupc, data = dados)
summary(EqlogLinearStar)##
## Call:
## lm(formula = lrenda ~ sexo + age + age2 + Defi + Desupi + Desupc,
## data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.39287 -0.40214 -0.04482 0.36440 2.19461
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.1496600 0.5511373 9.344 0.00000000000000213 ***
## sexo 0.3499367 0.1309819 2.672 0.008775 **
## age 0.0859690 0.0286082 3.005 0.003335 **
## age2 -0.0009051 0.0003410 -2.654 0.009212 **
## Defi -0.6175125 0.1631134 -3.786 0.000258 ***
## Desupi 0.3806861 0.2284963 1.666 0.098742 .
## Desupc 1.2188520 0.1953044 6.241 0.00000000980656488 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6384 on 103 degrees of freedom
## Multiple R-squared: 0.4503, Adjusted R-squared: 0.4183
## F-statistic: 14.06 on 6 and 103 DF, p-value: 0.00000000001233EqlogLinearStar_hc3<-coeftest(EqlogLinearStar, vcov = vcovHC(EqlogLinearStar, "HC3"))
class(EqlogLinearStar_hc3)## [1] "coeftest"library(texreg)
# Plotando os três modelos de regressão
screenreg(list(EqLogLinear,EqlogLinearStar, EqlogLinearStar_hc3),
digits = 4, stars = c(0.01, 0.05, 0.1),
custom.model.names = c("Modelo Log Linear", "Modelo Star","Star_hc3"),
caption = "Multiple model types and single row", label = "tab:2",
include.adjrs = TRUE, include.bic = TRUE)##
## =========================================================
## Modelo Log Linear Modelo Star Star_hc3
## ---------------------------------------------------------
## (Intercept) 4.5885 *** 5.1497 *** 5.1497 ***
## (0.6333) (0.5511) (0.5163)
## sexo 0.3618 *** 0.3499 *** 0.3499 **
## (0.1271) (0.1310) (0.1402)
## age 0.0853 *** 0.0860 *** 0.0860 ***
## (0.0277) (0.0286) (0.0286)
## age2 -0.0009 *** -0.0009 *** -0.0009 **
## (0.0003) (0.0003) (0.0004)
## Defi -0.0786 -0.6175 *** -0.6175 ***
## (0.3034) (0.1631) (0.1428)
## Defc 0.2302
## (0.3277)
## Demi 0.2219
## (0.4259)
## Demc 0.7456 **
## (0.2940)
## Desupi 0.9307 *** 0.3807 * 0.3807
## (0.3438) (0.2285) (0.2627)
## Desupc 1.7710 *** 1.2189 *** 1.2189 ***
## (0.3221) (0.1953) (0.2569)
## ---------------------------------------------------------
## R^2 0.5140 0.4503
## Adj. R^2 0.4702 0.4183
## Num. obs. 110 110
## =========================================================
## *** p < 0.01; ** p < 0.05; * p < 0.1Respondendo ao questionário
Qual é o melhor modelo estimado? Por quê? R: O modelo log-linear ajustado é o melhor. Pois ele estatisticamente significativo comm p-valor menor e variáveis com maior poder explicativo. Para este exercício, descartamos as variáveis não significativas, chegando ao modelo Log Linear Star como o melhor cujas variáveis são: sexo, age, age2, Defi, Desupi e Desupc
Interprete os parâmetros age e age2. Qual a relação entre as variáveis com a renda? R: Os parâmetros têm uma alta correlação entre si, pois age2 é age elevado ao quadrado. São os que apresentam tol menor que 0,1, e portanto, multicolinearidade. Entretanto, como são significativos, não foram adotadas medidas corretivas. Age tem um maior poder explicativo no modelo que consideramos o melhor, porém está abaixo de outros parâmetros, como desupc. Adotar age2 é condizente com a literatura, pois a medida que um indivíduo envelhece ele terá um papel negativo (inclinação negativa) para frear o aumento da renda. Desta forma, a renda aumenta com a idade até um certo ponto, quando passa a diminuir influenciada pelo parametro age2, num peso menor que age.
Interprete os parâmetros Desupc e Demc. Qual possui um maior impacto na renda do indivíduo? Explique. R: Para isso, teremos que olhar para os modelos linear e log-linear, pois o demc foi retirado no ajuste do log-linear. Em ambos os modelos, desupc mostra significância estatística maior que demc, isto é, temos o primeiro com valores menores que o segundo na coluna Pr(>|t|). Sendo Desupc uma variável binária com valor 1 se o chefe de família tem ensino superior completo - zero caso contrário. e Demc uma variável binária com valor 1 se o chefe de família tem ensino médio completo - zero caso contrário. Ensino superior tem maior impacto sobre a renda média de um indivíduo do que ensino médio, corroborando o que sugere a litera