Econometria: regressão com dados do Índice de Competitividade Global do WEF e o pacote data360r

Abstract

We analyse a regression model using the WEF Global Competitiveness Index (GCI) dataset and the World Bank’s package data360r.

Licença

This work is licensed under the Creative Commons Attribution-ShareAlike 4.0 International License. To view a copy of this license, visit http://creativecommons.org/licenses/by-sa/4.0/ or send a letter to Creative Commons, PO Box 1866, Mountain View, CA 94042, USA.

License: CC BY-SA 4.0

Citação

Sugestão de citação: FIGUEIREDO, Adriano Marcos Rodrigues. Econometria: regressão com dados do Índice de Competitividade Global do WEF e o pacote data360r. Campo Grande-MS,Brasil: RStudio/Rpubs, 2019. Disponível em http://rpubs.com/amrofi/competitiveness_data360r.

1 Introdução

Vamos examinar os metadados do conjunto de dados e identificar um conjunto de dados que queremos usar. Para este caso de uso, vamos nos concentrar no conjunto de dados do Índice de Competitividade Global (GCI) do WEF com o conjunto de dados id == 53. Um website importante destes dados é, por exemplo para o Brasil (BRA), https://tcdata360.worldbank.org/indicators/gci?country=BRA&indicator=632&viz=line_chart&years=2007,2017.

2 Dados

Nosso pacote básico é o data360r, que pode ser instalado a partir do CRAN (install.packages("data360r")) ou a partir do repositório github (devtools::install_github("mrpsonglao/data360r")). Nosso caso foi pelo github.

# install.packages('data360r') devtools::install_github('mrpsonglao/data360r')
library(data360r)

Faremos uma regressão rápida com o conjunto de dados do Índice de Competitividade Global (GCI) do WEF com o conjunto de dados id == 53. Primeiro puxarei os metadados do dataset desejado. Nesse caso, temos 101 observações com 7 variáveis.

df_usecase3_datasets <- get_metadata360(metadata_type = "datasets")

Esse dataframe contém a descrição das variáveis dentro desse conjunto de id == 53. Pegaremos uma variável específica, acompanhando o exemplo do #Case 3 em https://tcdata360.worldbank.org/tools/data360r. O indicador do WEF GCI é de id == 53 dentro do df_usecase3_datasets, ou com slug igual a global_competitiveness_index, na linha 52 do df_usecase3_datasets.

Portanto, escolherei para o último ano disponível (2017-2018). Existem informações desde 2007-2008. São vários indicadores e vários países, de modo que o usuário pode fazer um (View(df_usecase3_result)) caso deseje maiores detalhes.

library(tidyverse)
df_usecase3_result <- get_data360(dataset_id = c(53), output_type = "long") %>% filter(Period == 
    c("2017-2018"))

Para a regressão, manteremos apenas os indicadores GCI do WEF para Subindicator type == ‘Value’. Em seguida, reformulamos o quadro de dados resultante para que os indicadores sejam os nomes das colunas que estão sendo usados fazendo uso da função reshape::acast. Isso facilita a utilização dos dados do indicador nos modelos de regressão.

df2 <- filter(df_usecase3_result, df_usecase3_result$"Subindicator Type" == "Value", 
    !is.na(Observation))

knitr::kable(head(df2))

df3 <- as.data.frame(reshape2::acast(df2, df2$"Country ISO3" ~ df2$Indicator, value.var = "Observation"))

knitr::kable(head(df3))

Agora faremos a regressão com o dataset organizado (df3). A regressão nos indicadores “Inovação” (Innovation) e “Prontidão tecnológica” (Technological Readiness). Como o quadro de dados foi processado adequadamente, é fácil implementar a regressão nos indicadores WEF GCI 2017-2018. Antes de ajustar os dados a um modelo de regressão, primeiro vamos gerar um gráfico de dispersão para indicadores GCI WEF selecionados. Para simplificar, vamos nos concentrar nos dois indicadores citados.

load("df3.Rdata")
# o arquivo está em meu github em https://github.com/amrofi/datasets

Como parece que o site ou o pacote mudaram de nomes das variáveis, para preservar o algoritmo eu salvei o arquivo no diretorio de projeto.

library(ggplot2)
qplot(x = `12th pillar Innovation`, y = `9th pillar Technological readiness`, main = "Gráfico de dispersão entre Inovação x Prontidão tecnológica,
      2017-2018, países do GCI WEF", 
    xlab = "Inovação - pilar 12", ylab = "Prontidão Tecnológica - pilar 9", data = df3)

O gráfico de dispersão sugere que a inovação aumenta quadraticamente com a prontidão tecnológica. Vamos ajustar um modelo de regressão quadrática a esses pontos de dados. Com base nos resultados resumidos, o modelo é um bom ajuste. Também geramos os gráficos suplementares do modelo para ver se os resultados fazem sentido.

3 Modelo quadrático

mod_usecase3_quad <- lm(df3$`12th pillar Innovation` ~ poly(df3$`9th pillar Technological readiness`, 
    2) + df3$`1st pillar Institutions` + df3$`2nd pillar Infrastructure`, data = df3)
summary(mod_usecase3_quad)

Call:
lm(formula = df3$`12th pillar Innovation` ~ poly(df3$`9th pillar Technological readiness`, 
    2) + df3$`1st pillar Institutions` + df3$`2nd pillar Infrastructure`, 
    data = df3)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.83633 -0.25348 -0.00226  0.24333  1.06116 

Coefficients:
                                                   Estimate Std. Error t value
(Intercept)                                         1.30264    0.30309   4.298
poly(df3$`9th pillar Technological readiness`, 2)1  2.83413    0.94459   3.000
poly(df3$`9th pillar Technological readiness`, 2)2  2.21608    0.41620   5.325
df3$`1st pillar Institutions`                       0.38144    0.06498   5.870
df3$`2nd pillar Infrastructure`                     0.17087    0.07643   2.236
                                                   Pr(>|t|)    
(Intercept)                                        3.34e-05 ***
poly(df3$`9th pillar Technological readiness`, 2)1  0.00323 ** 
poly(df3$`9th pillar Technological readiness`, 2)2 4.27e-07 ***
df3$`1st pillar Institutions`                      3.38e-08 ***
df3$`2nd pillar Infrastructure`                     0.02706 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.3668 on 131 degrees of freedom
Multiple R-squared:  0.8213,    Adjusted R-squared:  0.8158 
F-statistic: 150.5 on 4 and 131 DF,  p-value: < 2.2e-16

mod_usecase3_quad$AIC <- AIC(mod_usecase3_quad)
mod_usecase3_quad$BIC <- BIC(mod_usecase3_quad)
par(mfrow = c(2, 2))
plot(mod_usecase3_quad)

4 Modelo log-log

mod_usecase3_loglog <- lm(log(df3$`12th pillar Innovation`) ~ log(df3$`9th pillar Technological readiness`) + 
    log(df3$`1st pillar Institutions`) + log(df3$`2nd pillar Infrastructure`), data = df3)
summary(mod_usecase3_loglog)

Call:
lm(formula = log(df3$`12th pillar Innovation`) ~ log(df3$`9th pillar Technological readiness`) + 
    log(df3$`1st pillar Institutions`) + log(df3$`2nd pillar Infrastructure`), 
    data = df3)

Residuals:
      Min        1Q    Median        3Q       Max 
-0.260358 -0.069793 -0.000115  0.076542  0.269855 

Coefficients:
                                               Estimate Std. Error t value
(Intercept)                                    0.009349   0.064182   0.146
log(df3$`9th pillar Technological readiness`)  0.278406   0.077962   3.571
log(df3$`1st pillar Institutions`)             0.638497   0.067009   9.529
log(df3$`2nd pillar Infrastructure`)          -0.028972   0.082934  -0.349
                                              Pr(>|t|)    
(Intercept)                                   0.884404    
log(df3$`9th pillar Technological readiness`) 0.000497 ***
log(df3$`1st pillar Institutions`)             < 2e-16 ***
log(df3$`2nd pillar Infrastructure`)          0.727394    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.1087 on 132 degrees of freedom
Multiple R-squared:  0.7688,    Adjusted R-squared:  0.7636 
F-statistic: 146.3 on 3 and 132 DF,  p-value: < 2.2e-16

mod_usecase3_loglog$AIC <- AIC(mod_usecase3_loglog)
mod_usecase3_loglog$BIC <- BIC(mod_usecase3_loglog)
par(mfrow = c(2, 2))
plot(mod_usecase3_loglog)

5 Resumo

suppressMessages(library(stargazer))
stargazer(list(mod_usecase3_quad, mod_usecase3_loglog), title = "Título: Resultados das Regressões Quadrática e Log-Log", 
    align = TRUE, type = "text", style = "all", keep.stat = c("AIC", "BIC", "rsq", 
        "adj.rsq", "n"), column.labels = c("Quadrática", "Log-Log"), dep.var.labels.include = FALSE)

Título: Resultados das Regressões Quadrática e Log-Log
======================================================================
                                              Dependent variable:     
                                          ----------------------------
                                            Quadrática      Log-Log   
                                               (1)            (2)     
----------------------------------------------------------------------
`9th pillar Technological readiness`, 2)1    2.834***                 
                                             (0.945)                  
                                            t = 3.000                 
                                            p = 0.004                 
`9th pillar Technological readiness`, 2)2    2.216***                 
                                             (0.416)                  
                                            t = 5.325                 
                                           p = 0.00000                
`1st pillar Institutions`                    0.381***                 
                                             (0.065)                  
                                            t = 5.870                 
                                           p = 0.00000                
`2nd pillar Infrastructure`                  0.171**                  
                                             (0.076)                  
                                            t = 2.236                 
                                            p = 0.028                 
`9th pillar Technological readiness`)                      0.278***   
                                                            (0.078)   
                                                           t = 3.571  
                                                          p = 0.0005  
`1st pillar Institutions`)                                 0.638***   
                                                            (0.067)   
                                                           t = 9.529  
                                                           p = 0.000  
`2nd pillar Infrastructure`)                                -0.029    
                                                            (0.083)   
                                                          t = -0.349  
                                                           p = 0.728  
Constant                                     1.303***        0.009    
                                             (0.303)        (0.064)   
                                            t = 4.298      t = 0.146  
                                           p = 0.00004     p = 0.885  
----------------------------------------------------------------------
Observations                                   136            136     
R2                                            0.821          0.769    
Adjusted R2                                   0.816          0.764    
Akaike Inf. Crit.                            120.050       -211.667   
Bayesian Inf. Crit.                          137.526       -197.103   
======================================================================
Note:                                      *p<0.1; **p<0.05; ***p<0.01

6 Teste de Especificação do modelo quadrático: RESET e outliers

Os modelos rejeitam H0 e acusam erro de especificação. Vamos checar eventuais outliers.

library(lmtest)
resettest(mod_usecase3_quad, power = 2:3, type = "fitted")

    RESET test

data:  mod_usecase3_quad
RESET = 4.2149, df1 = 2, df2 = 129, p-value = 0.01686

resettest(mod_usecase3_loglog, power = 2:3, type = "fitted")

    RESET test

data:  mod_usecase3_loglog
RESET = 9.7877, df1 = 2, df2 = 130, p-value = 0.0001098

7 Investigação de outliers das regressões

O teste para outliers foi o Bonferonni do pacote car. Conforme resultados abaixo, foi identificado outlier da regressão da observação 58, correspondente a ISR = Israel.

library(car)
outlierTest(mod_usecase3_quad)

No Studentized residuals with Bonferroni p < 0.05
Largest |rstudent|:
   rstudent unadjusted p-value Bonferroni p
58  3.05004          0.0027737      0.37722

outlierTest(mod_usecase3_loglog)

No Studentized residuals with Bonferroni p < 0.05
Largest |rstudent|:
   rstudent unadjusted p-value Bonferroni p
58 2.565249           0.011436           NA

8 Testes para Média zero e normalidade dos resíduos

Embora as observações 58 e 118 apareçam em destaque no qqplot, eles ainda estão dentro da banda de confiança. Os testes de normalidade Jarque-Bera não indicam violação do pressuposto de normalidade dos resíduos.

# obtendo residuos do modelo da tendencia cúbica
u.hat.quad <- resid(mod_usecase3_quad)
u.hat.loglog <- resid(mod_usecase3_loglog)
# teste de media zero dos residuos
t.test(u.hat.quad)

    One Sample t-test

data:  u.hat.quad
t = -3.4254e-17, df = 135, p-value = 1
alternative hypothesis: true mean is not equal to 0
95 percent confidence interval:
 -0.06127377  0.06127377
sample estimates:
    mean of x 
-1.061282e-18 

t.test(u.hat.loglog)

    One Sample t-test

data:  u.hat.loglog
t = -1.4306e-16, df = 135, p-value = 1
alternative hypothesis: true mean is not equal to 0
95 percent confidence interval:
 -0.01823166  0.01823166
sample estimates:
   mean of x 
-1.31878e-18 

# testes dos residuos da regressao cubica
library(car)
## normalidade dos residuos no R
qqPlot(mod_usecase3_quad)  #car

[1]  58 118

qqPlot(mod_usecase3_loglog)

[1]  58 118

# teste de jarque-bera para normalidade

(JB.regquad <- tseries::jarque.bera.test(u.hat.quad))

    Jarque Bera Test

data:  u.hat.quad
X-squared = 0.87481, df = 2, p-value = 0.6457

(JB.regloglog <- tseries::jarque.bera.test(u.hat.loglog))

    Jarque Bera Test

data:  u.hat.loglog
X-squared = 0.64692, df = 2, p-value = 0.7236

9 Teste de White para heterocedasticidade do modelo quadrático

Usarei o pacote lmtest para realizar o teste. A hipótese nula é H0: resíduos homocedásticos. A rejeição de H0 leva à adoção de modelo robusto a heterocedasticidade dos resíduos.

# mod_usecase3_quad <- lm(df3$`12th pillar Innovation` ~ poly(df3$`9th pillar
# Technological readiness`, 2)+ df3$`1st pillar Institutions`+ df3$`2nd pillar
# Infrastructure`)
library(lmtest)
mquad <- mod_usecase3_quad
bptest(mquad, ~poly(df3$`9th pillar Technological readiness`, 2) + df3$`1st pillar Institutions` + 
    df3$`2nd pillar Infrastructure` + I(df3$`9th pillar Technological readiness`^2) + 
    I(df3$`9th pillar Technological readiness`^4) + I(df3$`1st pillar Institutions`^2) + 
    I(df3$`2nd pillar Infrastructure`^2))

    studentized Breusch-Pagan test

data:  mquad
BP = 12.75, df = 7, p-value = 0.07845

A correção de White para o modelo robusto será:

vcov.white0 <- car::hccm(mod_usecase3_quad, type = c("hc1"))
coeftest(mod_usecase3_quad, vcov.white0)

t test of coefficients:

                                                   Estimate Std. Error t value
(Intercept)                                        1.302639   0.325046  4.0076
poly(df3$`9th pillar Technological readiness`, 2)1 2.834125   1.051113  2.6963
poly(df3$`9th pillar Technological readiness`, 2)2 2.216079   0.430958  5.1422
df3$`1st pillar Institutions`                      0.381438   0.062005  6.1517
df3$`2nd pillar Infrastructure`                    0.170869   0.083825  2.0384
                                                    Pr(>|t|)    
(Intercept)                                        0.0001024 ***
poly(df3$`9th pillar Technological readiness`, 2)1 0.0079324 ** 
poly(df3$`9th pillar Technological readiness`, 2)2 9.656e-07 ***
df3$`1st pillar Institutions`                      8.663e-09 ***
df3$`2nd pillar Infrastructure`                    0.0435222 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Ou de outra forma,

library(sandwich)
modrob <- coeftest(mod_usecase3_quad, vcov = vcovHC(mod_usecase3_quad, type = "HC1"))

suppressMessages(library(stargazer))
stargazer(list(mod_usecase3_quad, mod_usecase3_loglog, modrob), title = "Título: Resultados das Regressões Quadrática e Log-Log", 
    align = TRUE, type = "text", style = "all", keep.stat = c("AIC", "BIC", "rsq", 
        "adj.rsq", "n"), column.labels = c("Quadrática", "Log-Log", "Quad.Robusto"), 
    dep.var.labels.include = FALSE, model.numbers = FALSE, model.names = FALSE)

Título: Resultados das Regressões Quadrática e Log-Log
=============================================================================
                                                  Dependent variable:        
                                          -----------------------------------
                                          Quadrática   Log-Log   Quad.Robusto
-----------------------------------------------------------------------------
`9th pillar Technological readiness`, 2)1  2.834***                2.834***  
                                            (0.945)                (1.051)   
                                           t = 3.000              t = 2.696  
                                           p = 0.004              p = 0.008  
`9th pillar Technological readiness`, 2)2  2.216***                2.216***  
                                            (0.416)                (0.431)   
                                           t = 5.325              t = 5.142  
                                          p = 0.00000            p = 0.00000 
`1st pillar Institutions`                  0.381***                0.381***  
                                            (0.065)                (0.062)   
                                           t = 5.870              t = 6.152  
                                          p = 0.00000             p = 0.000  
`2nd pillar Infrastructure`                 0.171**                0.171**   
                                            (0.076)                (0.084)   
                                           t = 2.236              t = 2.038  
                                           p = 0.028              p = 0.044  
`9th pillar Technological readiness`)                  0.278***              
                                                       (0.078)               
                                                      t = 3.571              
                                                      p = 0.0005             
`1st pillar Institutions`)                             0.638***              
                                                       (0.067)               
                                                      t = 9.529              
                                                      p = 0.000              
`2nd pillar Infrastructure`)                            -0.029               
                                                       (0.083)               
                                                      t = -0.349             
                                                      p = 0.728              
Constant                                   1.303***     0.009      1.303***  
                                            (0.303)    (0.064)     (0.325)   
                                           t = 4.298  t = 0.146   t = 4.008  
                                          p = 0.00004 p = 0.885   p = 0.0002 
-----------------------------------------------------------------------------
Observations                                  136        136                 
R2                                           0.821      0.769                
Adjusted R2                                  0.816      0.764                
Akaike Inf. Crit.                           120.050    -211.667              
Bayesian Inf. Crit.                         137.526    -197.103              
=============================================================================
Note:                                             *p<0.1; **p<0.05; ***p<0.01

10 Teste de multicolinearidade do modelo quadrático

Não houve evidência de multicolinearidade prejudicial ao modelo.

regquad.vif <- vif(mod_usecase3_quad)
regquad.vif

                                                      GVIF Df GVIF^(1/(2*Df))
poly(df3$`9th pillar Technological readiness`, 2) 8.506647  2        1.707810
df3$`1st pillar Institutions`                     3.252166  1        1.803376
df3$`2nd pillar Infrastructure`                   8.503978  1        2.916158

Referências

Ma. Regina Paz S. Onglao (2018). data360r: Wrapper for TC/Govdata360 API. R package version 1.0.2. Disponível em: https://github.com/mrpsonglao/data360r