Econometria Espacial em R

Aplicação à Convergência de Renda

O objetivo deste documento é prover um tutorial para as principais técnicas da Econometria Espacial. Dessa forma, os seguintes tópicos hão de ser abordados:

• Mapas com análise quantílica
• Matriz de Ponderação Espacial
• Autocorrelação espacial (Global e Local)
• Modelo de Defasagem Espacial (SAR) e Modelo de Erro Espacial (SEM)

1 Pacotes

install.packages(c("geobr","ggplot2", "dplyr","spdep","spatialreg", "normtest"), repos = "http://cran.us.r-project.org")

## package 'geobr' successfully unpacked and MD5 sums checked
## package 'ggplot2' successfully unpacked and MD5 sums checked
## package 'dplyr' successfully unpacked and MD5 sums checked
## package 'spdep' successfully unpacked and MD5 sums checked
## package 'spatialreg' successfully unpacked and MD5 sums checked
## package 'normtest' successfully unpacked and MD5 sums checked
## 
## The downloaded binary packages are in
##  C:\Users\Usuario\AppData\Local\Temp\Rtmpw9jbNp\downloaded_packages

library(spatialreg)
library(geobr)
library(ggplot2)
library(dplyr)
library(readr)
library(spdep)
library(classInt)
library(tseries)
library(lmtest)
library(MASS)
library(normtest)

2 Importação dos dados

Além de carregar os arquivos, faz-se manipulações em prol da organização do dataset principal joinedata. Atenção: ../Monografia 1/ refere-se ao Working Directory do Rstudio, de forma que o arquivo dataset4.csv precisa estar no diretório específico de cada computador.

#Leitura do mapa de MG
muni <- read_municipality(code_muni = "MG", year = 2001, showProgress = FALSE)

#Leitura dos dados dos municípios
data1 <- read.csv("../Monografia 1/dataset4.csv", sep = ";", header = TRUE, fileEncoding = "latin1")

data2 <- read.csv("../Monografia 1/ipeadata_anos_med.csv", sep = ";", header = TRUE, fileEncoding = "latin1")

data2 <- data2[-c(4:5)]

names(data1)[names(data1) == "Índice.de.Gini.2000"] <- "gini_2000"
names(data1)[names(data1) == "X..dos.ocupados.com.médio.completo...18.anos.ou.mais.2000"] <- "medcomp18_2002"
names(data1)[names(data1) == "X..dos.ocupados.com.superior.completo...18.anos.ou.mais.2000"] <- "supcomp18_2002"
names(data1)[names(data1) == "Esperança.de.vida.ao.nascer.2000"] <- "espvid_2002"
names(data1)[names(data1) == "Expectativa.de.anos.de.estudo.2000"] <- "expest_2002"

#Reunindo os datasets carregados
joinedata <- left_join(muni, data1, by = c("code_muni" = "Codigo"))

joinedata <- left_join(joinedata, data2, by = c("code_muni" = "Codigo"))

joinedata <- joinedata[,-c(20)]

3 Análise quantílica do PIB per capita municipal

Para o ano de 2002:

#2002
pibinho <- as.numeric(joinedata$PIB.per.capita.2002)

breaks_qt <- classIntervals(c(min(pibinho) - 1, pibinho), style = 'quantile',
                            n = 4)

joinedata <- mutate(joinedata, pib2002 = cut(PIB.per.capita.2002, breaks_qt$brks,
                                             dig.lab = 5,
                                             labels = c("Primeiro Quantil", "Segundo Quantil", "Terceiro Quantil","Quarto Quantil")))

no_axis <- theme(axis.title = element_blank(),
                 axis.text = element_blank(),
                 axis.ticks = element_blank(),
                 legend.title = element_blank(),
                 panel.grid = element_blank())

ggplot() +
  geom_sf(data = joinedata, aes(fill = pib2002), size = .15) +
  scale_fill_brewer(palette = "Oranges") +
  theme_minimal() +
  no_axis +
  labs(title = "Distribuição espacial do PIB per capita em 2002")

Para o ano de 2017:

pibinho_2 <- as.numeric(joinedata$PIB.per.capita.2017)

breaks_qt_2 <- classIntervals(c(min(pibinho_2) - 1, pibinho_2), style = 'quantile',
                              n = 4)

joinedata <- mutate(joinedata, pib2017 = cut(PIB.per.capita.2017, breaks_qt_2$brks,
                                             dig.lab = 5,
                                             labels = c("Primeiro Quantil", "Segundo Quantil", "Terceiro Quantil","Quarto Quantil")))

ggplot() +
  geom_sf(data = joinedata, aes(fill = pib2017), size = .15) +  
  scale_fill_brewer(palette = "Oranges") +
  theme_minimal() +
  no_axis +
  labs(title = "Distribuição espacial do PIB per capita em 2017")

Entre 2002 e 2017:

joinedata <- mutate(joinedata, txcrescm=log(PIB.per.capita.2017/PIB.per.capita.2002)/15)

txcrescm1 <- joinedata$txcrescm
breaks_qt_3 <- classIntervals(c(min(txcrescm1) - .001, txcrescm1), style = 'quantile',
                              n = 4)

joinedata <- mutate(joinedata, crescm = cut(txcrescm1, breaks_qt_3$brks,
                                            dig.lab = 4,
                                            labels = c("Primeiro Quantil", "Segundo Quantil", "Terceiro Quantil","Quarto Quantil")))

ggplot() +
  geom_sf(data = joinedata, aes(fill = crescm), size = .15) +
  scale_fill_brewer(palette = "Oranges") +
  theme_minimal() +
  no_axis +
  labs(title = "Distribuição espacial da Taxa de Crescimento Médio do PIB per capita", subtitle = "Período de 2002-17")

4 Criando uma Matriz de Pesos Espaciais

Para a aplicação das técnicas, é necessário definir uma matriz de ponderação espacial (W). Trata-se de um meio para representar a estrutura espacial dos dados. Conforme o procedimento de Baumont (2004), escolhe-se a matriz que gera o maior valor do I de Moran para os resíduos da estimação MQO. Tal matriz é a de convenção torre.

#Definindo uma matriz de convenção torre
mg.rook <- poly2nb(joinedata$geom, queen = FALSE)

#Tranformando a matriz em lista
rook.listw <- nb2listw(mg.rook)

#Mapa da Matriz de Pesos Espaciais
mg.coords <- cbind(joinedata$SEDE, joinedata$SEDE0)

plot(joinedata$geom, border = 'lightgrey')         
plot(mg.rook, mg.coords, add = TRUE, col = 'blue')

5 Análise Exploratória Espacial

5.1 I de Moran Global para PIB per capita

Para o ano de 2002:

#2002
logpib2002 <- log(joinedata$PIB.per.capita.2002)
moran.test(x = logpib2002, listw = rook.listw)

## 
##  Moran I test under randomisation
## 
## data:  logpib2002  
## weights: rook.listw    
## 
## Moran I statistic standard deviate = 22.36, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.4726340260     -0.0011737089      0.0004490232

Para o ano de 2017:

#2017
logpib2017 <- log(joinedata$PIB.per.capita.2017)
moran.test(x = logpib2017, listw = rook.listw)

## 
##  Moran I test under randomisation
## 
## data:  logpib2017  
## weights: rook.listw    
## 
## Moran I statistic standard deviate = 8.3156, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.1751235880     -0.0011737089      0.0004494782

5.2 I de Moran local (LISA) para PIB per capita

Para o ano de 2002:

###LISA 2002
locm <- localmoran(logpib2002, listw = rook.listw)

joinedata$spib02 <- scale(logpib2002) %>% as.vector()
joinedata$lag_spib02 <- lag.listw(rook.listw, joinedata$spib02)

joinedata$quad_sig <- NA
joinedata[(joinedata$spib02 >= 0 & joinedata$lag_spib02 >= 0) & (locm[, 5] <= 0.05), "quad_sig"] <- 1
joinedata[(joinedata$spib02 <= 0 & joinedata$lag_spib02 <= 0) & (locm[, 5] <= 0.05), "quad_sig"] <- 2
joinedata[(joinedata$spib02 >= 0 & joinedata$lag_spib02 <= 0) & (locm[, 5] <= 0.05), "quad_sig"] <- 3
joinedata[(joinedata$spib02 <= 0 & joinedata$lag_spib02 >= 0) & (locm[, 5] <= 0.05), "quad_sig"] <- 4
joinedata[(locm[, 5] >= 0.05), "quad_sig"] <- 5 

#Mapeamento do LISA
breaks <- seq(1, 5, 1)

labels <- c("alto-alto", "baixo-baixo", "alto-baixo", "baixo-alto", "não signif.")

np <- findInterval(joinedata$quad_sig, breaks)

#Estabelecendo cores para cada tipo de cluster
colors <- c("red", "blue", "lightpink", "skyblue2", "white")
plot(joinedata$geom, col = colors[np])  #colors[np] coloca uma cor para cada município
mtext("I de Moran local 2002", cex = 1.5, side = 3, line = 1)
legend("topleft", legend = labels, fill = colors, bty = "n")

Para o ano de 2017:

locm_1 <- localmoran(logpib2017, listw = rook.listw)

joinedata$spib17 <- scale(logpib2017) %>% as.vector()
joinedata$lag_spib17 <- lag.listw(rook.listw, joinedata$spib17)

joinedata$quad_sig_1 <- NA
joinedata[(joinedata$spib17 >= 0 & joinedata$lag_spib17 >= 0) & (locm_1[, 5] <= 0.05), "quad_sig_1"] <- 1
joinedata[(joinedata$spib17 <= 0 & joinedata$lag_spib17 <= 0) & (locm_1[, 5] <= 0.05), "quad_sig_1"] <- 2
joinedata[(joinedata$spib17 >= 0 & joinedata$lag_spib17 <= 0) & (locm_1[, 5] <= 0.05), "quad_sig_1"] <- 3
joinedata[(joinedata$spib17 <= 0 & joinedata$lag_spib17 >= 0) & (locm_1[, 5] <= 0.05), "quad_sig_1"] <- 4
joinedata[(locm_1[, 5] >= 0.05), "quad_sig_1"] <- 5 

#Mapeamento do LISA
breaks_1 <- seq(1, 5, 1)

labels_1 <- c("alto-alto", "baixo-baixo", "alto-baixo", "baixo-alto", "não signif.")

np_1 <- findInterval(joinedata$quad_sig_1, breaks_1)

#Estabelecendo cores para cada tipo de cluster
colors_1 <- c("red", "blue", "lightpink", "skyblue2", "white")
plot(joinedata$geom, col = colors[np_1])  #colors[np] coloca uma cor para cada município
mtext("I de Moran local 2017", cex = 1.5, side = 3, line = 1)
legend("topleft", legend = labels, fill = colors, bty = "n")

6 Análise de convergência

6.1 Procedimento de estimação

Adota-se o procedimento sugerido por Florax, Folmer e Rey (2003):

1. Estimar o modelo clássico por Mínimos Quadrados Ordinários (MQO);
2. Testar a hipótese de autocorrelação espacial devido a uma defasagem espacial ou um erro espacial. Isso há de ser feito via estatísticas Multiplicador de Lagrange para ρ defasagem espacial (MLρ) e Multiplicador de Lagrange para o λ erro espacial (MLλ);
3. Se ambos os testes não são significativos, o modelo clássico é recomendado. Caso contrário, segue-se para o próximo passo;
4. Se ambos forem significativos, escolhe-se o modelo mais significante de acordo com as versões robustas desses testes: Multiplicador de Lagrange Robusto para ρ defasagem espacial (MLRρ) e Multiplicador de Lagrange Robusto para o λ erro espacial (MLRλ). Desse modo, caso MLRρ > MLRλ, escolhe-se o modelo com defasagem espacial como o mais apropriado. Caso contrário, escolhe-se o modelo de erro espacial como o mais apropriado.

6.2 Abordagem clássica

Os dados abaixo se referem aos indicadores contemplados na modelagem econométrica.

#Os dados foram extraídos das seguintes fontes: IMRS(FJP), PNUD e IPEA
pib_growth <- log(joinedata$PIB.per.capita.2017/joinedata$PIB.per.capita.2002)/15

ln_pib02 <- log(joinedata$PIB.per.capita.2002)

gini_00 <- joinedata$gini_2000

anos_med_est <- joinedata$anos_med_2000

dat <- data.frame(pib_growth, ln_pib02, anos_med_est, gini_00)

#Visualizando as primeiras linhas 
head(dat)

##   pib_growth ln_pib02 anos_med_est gini_00
## 1 0.12041798 8.190410          4.4    0.50
## 2 0.09907125 8.314511          4.8    0.58
## 3 0.09022412 8.153669          4.0    0.51
## 4 0.09038039 7.807095          3.7    0.55
## 5 0.10067212 7.670079          3.0    0.55
## 6 0.09674629 7.645924          2.6    0.64

Vamos criar funções para os cálculos de velocidade de convergência (% ao ano) e meia-vida (anos). Note-se que as funções foram criadas conforme a abordagem de Barro e Sala-i-Martin (1992).

velconv <- function(betaconv) {
  vel <- -log(1 + 10 * betaconv)/10
  return(vel*100)
}

#meia-vida

halflife <- function(betaconv) {
  tau <- -log(2)/log(1 + betaconv)
  return(tau)
}

Convergência absoluta:

#Regressão simples
reg1 <- lm(pib_growth~ln_pib02, data = dat)
summary(reg1)

## 
## Call:
## lm(formula = pib_growth ~ ln_pib02, data = dat)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.080399 -0.021047 -0.006103  0.016379  0.190626 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.408057   0.016706   24.43   <2e-16 ***
## ln_pib02    -0.037969   0.002043  -18.59   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03403 on 851 degrees of freedom
## Multiple R-squared:  0.2888, Adjusted R-squared:  0.2879 
## F-statistic: 345.5 on 1 and 851 DF,  p-value: < 2.2e-16

#Visualizando a reg1
ggplot(dat, aes(x = ln_pib02, y = pib_growth)) +
  geom_point() +
  stat_smooth(method = "lm", col = "red") +
  ylab("Taxa de Crescimento 2002-2017")

Faz-se o diagnóstico da regressão para verificar se existe heterocedasticidade. Ademais, investiga-se a normalidade dos resíduos.

sresid <- studres(reg1)
jarque.bera.test(sresid)

## 
##  Jarque Bera Test
## 
## data:  sresid
## X-squared = 927.92, df = 2, p-value < 2.2e-16

plot(reg1, 2) 

De acordo com a análise gráfica e o teste de Jarque-Bera, os resíduos não seguem uma distribuição normal.

bptest(reg1) 

## 
##  studentized Breusch-Pagan test
## 
## data:  reg1
## BP = 1.2054, df = 1, p-value = 0.2722

dat$resi <- reg1$residuals
ggplot(data = dat, aes(y = resi, x = ln_pib02)) +
  geom_point(col = 'blue') + geom_abline(slope = 0) 

Tanto a análise gráfica quanto o teste de Breusch-Pagan indicam que os resíduos são homocedásticos. A seguir são feitos os cálculos para velocidade e meia-vida.

format(velconv(-0.03796852), digits = 4)

## [1] "4.775"

format(halflife(-0.03796852), digits = 5)

## [1] "17.907"

Convergência condicional:

reg2 <- lm(pib_growth~ln_pib02 + anos_med_est + gini_00,
           data = dat)
summary(reg2)

## 
## Call:
## lm(formula = pib_growth ~ ln_pib02 + anos_med_est + gini_00, 
##     data = dat)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.080225 -0.020815 -0.006505  0.015677  0.192197 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   0.434379   0.021917  19.820   <2e-16 ***
## ln_pib02     -0.042154   0.002785 -15.134   <2e-16 ***
## anos_med_est  0.003276   0.001531   2.141   0.0326 *  
## gini_00      -0.010741   0.018239  -0.589   0.5561    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03397 on 849 degrees of freedom
## Multiple R-squared:  0.293,  Adjusted R-squared:  0.2905 
## F-statistic: 117.3 on 3 and 849 DF,  p-value: < 2.2e-16

• Diagnóstico dos resíduos

sresid2 <- studres(reg2)
jarque.bera.test(sresid2) 

## 
##  Jarque Bera Test
## 
## data:  sresid2
## X-squared = 978.52, df = 2, p-value < 2.2e-16

Os resíduos não seguem uma distribuição normal.

bptest(reg2)

## 
##  studentized Breusch-Pagan test
## 
## data:  reg2
## BP = 4.5795, df = 3, p-value = 0.2053

Não há heterocedasticidade.

• Cálculos de velocidade e meia-vida

format(velconv(-0.042154466), digits = 5)

## [1] "5.4739"

format(halflife(-0.042154466), digits = 5)

## [1] "16.094"

6.3 Testes de Dependência Espacial

Realiza-se o teste de Moran para os resíduos da regressão MQO. Em seguida são feitos testes de Multiplicadores de Lagrange para identificar como a autocorrelação espacial se manifesta. As versões robustas desses testes também são aplicadas.

lm.morantest(reg1, listw = rook.listw)

## 
##  Global Moran I for regression residuals
## 
## data:  
## model: lm(formula = pib_growth ~ ln_pib02, data = dat)
## weights: rook.listw
## 
## Moran I statistic standard deviate = 2.5864, p-value = 0.004849
## alternative hypothesis: greater
## sample estimates:
## Observed Moran I      Expectation         Variance 
##     0.0530817578    -0.0017304748     0.0004491225

Há dependência espacial nos resíduos da regressão MQO.

lm.LMtests(reg1, listw = rook.listw, test = "all")

## 
##  Lagrange multiplier diagnostics for spatial dependence
## 
## data:  
## model: lm(formula = pib_growth ~ ln_pib02, data = dat)
## weights: rook.listw
## 
## LMerr = 6.2095, df = 1, p-value = 0.01271
## 
## 
##  Lagrange multiplier diagnostics for spatial dependence
## 
## data:  
## model: lm(formula = pib_growth ~ ln_pib02, data = dat)
## weights: rook.listw
## 
## LMlag = 0.7834, df = 1, p-value = 0.3761
## 
## 
##  Lagrange multiplier diagnostics for spatial dependence
## 
## data:  
## model: lm(formula = pib_growth ~ ln_pib02, data = dat)
## weights: rook.listw
## 
## RLMerr = 11.866, df = 1, p-value = 0.0005717
## 
## 
##  Lagrange multiplier diagnostics for spatial dependence
## 
## data:  
## model: lm(formula = pib_growth ~ ln_pib02, data = dat)
## weights: rook.listw
## 
## RLMlag = 6.4397, df = 1, p-value = 0.01116
## 
## 
##  Lagrange multiplier diagnostics for spatial dependence
## 
## data:  
## model: lm(formula = pib_growth ~ ln_pib02, data = dat)
## weights: rook.listw
## 
## SARMA = 12.649, df = 2, p-value = 0.001792

Há uma diferença em relação aos resultados das versões dos testes LM. As versões não robustas mostram que há uma dependência espacial quando a alternativa é o modelo de erro espacial, porém, não há evidência de dependência espacial quando a alternativa é o modelo de defasagem espacial. No que concerne às versões robustas, ambos os testes são significativos ao nível de 5% e, além disso, o modelo com melhor ajuste é o modelo de erro espacial.

6.4 Abordagem Espacial

As hipóteses de convergência são testadas com o modelo SAR.

Convergência absoluta:

reg3 <- lagsarlm(pib_growth~ln_pib02, data = dat, listw = rook.listw)
summary(reg3)

## 
## Call:spatialreg::lagsarlm(formula = formula, data = data, listw = listw, 
##     na.action = na.action, Durbin = Durbin, type = type, method = method, 
##     quiet = quiet, zero.policy = zero.policy, interval = interval, 
##     tol.solve = tol.solve, trs = trs, control = control)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0794172 -0.0212852 -0.0059807  0.0164471  0.1904992 
## 
## Type: lag 
## Coefficients: (asymptotic standard errors) 
##               Estimate Std. Error z value  Pr(>|z|)
## (Intercept)  0.3983434  0.0203749  19.551 < 2.2e-16
## ln_pib02    -0.0373192  0.0022182 -16.824 < 2.2e-16
## 
## Rho: 0.04537, LR test value: 0.84508, p-value: 0.35795
## Asymptotic standard error: 0.046953
##     z-value: 0.96628, p-value: 0.3339
## Wald statistic: 0.9337, p-value: 0.3339
## 
## Log likelihood: 1674.611 for lag model
## ML residual variance (sigma squared): 0.0011539, (sigma: 0.033968)
## Number of observations: 853 
## Number of parameters estimated: 4 
## AIC: -3341.2, (AIC for lm: -3342.4)
## LM test for residual autocorrelation
## test value: 12.116, p-value: 0.00049986

• Cálculo do Pseudo \(r^{2}\)

1-(reg3$SSE/(var(dat$pib_growth)*(length(dat$pib_growth)-1)))

## [1] 0.2897397

• Diagnóstico dos resíduos

sresid3 <- reg3$residuals
jb.norm.test(sresid3)

## 
##  Jarque-Bera test for normality
## 
## data:  sresid3
## JB = 848.78, p-value < 2.2e-16

Os resíduos não seguem uma distribuição normal.

bptest.sarlm(reg3)

## 
##  studentized Breusch-Pagan test
## 
## data:  
## BP = 1.1892, df = 1, p-value = 0.2755

Não há heterocedasticidade.

• Cálculos de velocidade e meia-vida

format(velconv(-0.03731924), digits = 5)

## [1] "4.6712"

format(halflife(-0.03731924), digits = 5)

## [1] "18.225"

Convergência condicional:

reg4 <- lagsarlm(pib_growth~ln_pib02 + anos_med_est + gini_00,
                 data = dat, listw = rook.listw)
summary(reg4)

## 
## Call:spatialreg::lagsarlm(formula = formula, data = data, listw = listw, 
##     na.action = na.action, Durbin = Durbin, type = type, method = method, 
##     quiet = quiet, zero.policy = zero.policy, interval = interval, 
##     tol.solve = tol.solve, trs = trs, control = control)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0806223 -0.0208830 -0.0065462  0.0161757  0.1921122 
## 
## Type: lag 
## Coefficients: (asymptotic standard errors) 
##                Estimate Std. Error  z value Pr(>|z|)
## (Intercept)   0.4237165  0.0238887  17.7371  < 2e-16
## ln_pib02     -0.0415501  0.0028480 -14.5894  < 2e-16
## anos_med_est  0.0034779  0.0015327   2.2692  0.02326
## gini_00      -0.0127193  0.0182302  -0.6977  0.48536
## 
## Rho: 0.06125, LR test value: 1.5189, p-value: 0.21778
## Asymptotic standard error: 0.047278
##     z-value: 1.2955, p-value: 0.19514
## Wald statistic: 1.6784, p-value: 0.19514
## 
## Log likelihood: 1677.507 for lag model
## ML residual variance (sigma squared): 0.0011457, (sigma: 0.033848)
## Number of observations: 853 
## Number of parameters estimated: 6 
## AIC: -3343, (AIC for lm: -3343.5)
## LM test for residual autocorrelation
## test value: 7.9686, p-value: 0.0047597

• Cálculo do Pseudo \(r^{2}\)

1-(reg4$SSE/(var(dat$pib_growth)*(length(dat$pib_growth)-1)))

## [1] 0.2947616

• Diagnóstico dos resíduos

sresid4 <- reg4$residuals
jb.norm.test(sresid4)

## 
##  Jarque-Bera test for normality
## 
## data:  sresid4
## JB = 895.44, p-value < 2.2e-16

Os resíduos não seguem uma distribuição normal.

bptest.sarlm(reg4)

## 
##  studentized Breusch-Pagan test
## 
## data:  
## BP = 4.5764, df = 3, p-value = 0.2056

Não há heterocedasticidade.

• Cáculos de velocidade e meia-vida

format(velconv(-0.041550063), digits = 4)

## [1] "5.37"

format(halflife(-0.041550063), digits = 5)

## [1] "16.333"

As hipóteses de convergência são testadas com o modelo SEM.

Convergência absoluta:

reg5 <- errorsarlm(pib_growth~ln_pib02, data = dat, listw = rook.listw)
summary(reg5)

## 
## Call:spatialreg::errorsarlm(formula = formula, data = data, listw = listw, 
##     na.action = na.action, Durbin = Durbin, etype = etype, method = method, 
##     quiet = quiet, zero.policy = zero.policy, interval = interval, 
##     tol.solve = tol.solve, trs = trs, control = control)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.078681 -0.021427 -0.006243  0.016390  0.187796 
## 
## Type: error 
## Coefficients: (asymptotic standard errors) 
##               Estimate Std. Error z value  Pr(>|z|)
## (Intercept)  0.4161555  0.0177186  23.487 < 2.2e-16
## ln_pib02    -0.0389591  0.0021663 -17.984 < 2.2e-16
## 
## Lambda: 0.14163, LR test value: 6.4906, p-value: 0.010844
## Asymptotic standard error: 0.052211
##     z-value: 2.7126, p-value: 0.0066759
## Wald statistic: 7.3581, p-value: 0.0066759
## 
## Log likelihood: 1677.433 for error model
## ML residual variance (sigma squared): 0.0011425, (sigma: 0.0338)
## Number of observations: 853 
## Number of parameters estimated: 4 
## AIC: -3346.9, (AIC for lm: -3342.4)

• Cálculo do Pseudo \(r^{2}\)

1-(reg5$SSE/(var(dat$pib_growth)*(length(dat$pib_growth)-1)))

## [1] 0.2967513

• Diagnóstico dos resíduos

sresid5 <- reg5$residuals
jb.norm.test(sresid5)

## 
##  Jarque-Bera test for normality
## 
## data:  sresid5
## JB = 822.62, p-value < 2.2e-16

Os resíduos não seguem uma distribuição normal

bptest.sarlm(reg5)

## 
##  studentized Breusch-Pagan test
## 
## data:  
## BP = 1.0088, df = 1, p-value = 0.3152

Não há heterocedasticidade

• Cálculos de velocidade e meia-vida

format(velconv(-0.0389591), digits = 4)

## [1] "4.936"

format(halflife(-0.0389591), digits = 5)

## [1] "17.443"

Convergência condicional:

reg6 <- errorsarlm(pib_growth~ln_pib02 + anos_med_est + gini_00,
                 data = dat, listw = rook.listw)
summary(reg6)

## 
## Call:spatialreg::errorsarlm(formula = formula, data = data, listw = listw, 
##     na.action = na.action, Durbin = Durbin, etype = etype, method = method, 
##     quiet = quiet, zero.policy = zero.policy, interval = interval, 
##     tol.solve = tol.solve, trs = trs, control = control)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0784408 -0.0207665 -0.0064811  0.0162547  0.1895434 
## 
## Type: error 
## Coefficients: (asymptotic standard errors) 
##                Estimate Std. Error  z value Pr(>|z|)
## (Intercept)   0.4402667  0.0224090  19.6469  < 2e-16
## ln_pib02     -0.0427327  0.0028497 -14.9957  < 2e-16
## anos_med_est  0.0031817  0.0015837   2.0090  0.04454
## gini_00      -0.0121071  0.0183595  -0.6594  0.50961
## 
## Lambda: 0.13401, LR test value: 5.8511, p-value: 0.015567
## Asymptotic standard error: 0.052387
##     z-value: 2.5581, p-value: 0.010524
## Wald statistic: 6.5439, p-value: 0.010524
## 
## Log likelihood: 1679.673 for error model
## ML residual variance (sigma squared): 0.0011369, (sigma: 0.033718)
## Number of observations: 853 
## Number of parameters estimated: 6 
## AIC: -3347.3, (AIC for lm: -3343.5)

• Cálculo do Pseudo \(r^{2}\)

1-(reg6$SSE/(var(dat$pib_growth)*(length(dat$pib_growth)-1)))

## [1] 0.3001609

• Diagnóstico dos resíduos

sresid6 <- reg6$residuals
jb.norm.test(sresid6)

## 
##  Jarque-Bera test for normality
## 
## data:  sresid6
## JB = 873.98, p-value < 2.2e-16

Os resíduos não seguem uma distribuição normal

bptest.sarlm(reg6)

## 
##  studentized Breusch-Pagan test
## 
## data:  
## BP = 4.9392, df = 3, p-value = 0.1763

Não há heterocedasticidade

• Cálculos de velocidade e meia-vida

format(velconv(-0.0427351910), digits = 5)

## [1] "5.5748"

format(halflife(-0.0427351910), digits = 6)

## [1] "15.8705"

7 Referências

• Barro, Robert J., and Xavier Sala-i-Martin. 1992. Convergence. Journal of Political Economy 100(2): 223-251.

• BurkeyAcademy

• Create a LISA map Exercise

• Florax, Raymond., Folmer, Henk., and Sergio Rey. 2003. Specification searches in spatial econometrics: the relevance of Hendry’s methodology. Regional Science and Urban Economics 33: 557-579.

• Mendez C. (2020). Spatial regression analysis in R. R Studio/Rpubs. Disponível em https://rpubs.com/quarcs-lab/tutorial-spatial-regression