Econometria de Painel

Author

Igor Viveiros de Souza

Published

April 19, 2024

1 Uma rápida revisão: MQO e suas hipóteses

Regressão simples:

\(Y_i = \beta_0 + \beta_1 X_i + u_i\)

Hipóteses:

\(H_1\): Esperança condicional do erro é zero, isto é: \(E[u_i|X_i] = 0\),

\(H_2\): \((X_i,Y_i)\) são ensaios de uma distribuição conjunta,

\(H_3\): Valores extremos (outliers) são improváveis,

\(H_4\): não há multicolinearidade,

\(H_5\): \(COV(u_{i}u_{j}) = 0\) para \(i \neq j\),

\(H_6\): \(E(u_{i}^2|X_i) = \sigma^2_u\) \(\forall\) \(i\)

Sob estas hipóteses o \(MQO\) é um estimador não viciado, consistente e com distribuição normal em grandes amostras


Consistência: Propriedade assintótica (isto é para grandes amostras) que significa que à medida que a amostra cresce o estimador \(\hat{\beta}_{1,n}\) converge para o seu verdadeiro valor, \(\beta\). Isto é, \(\hat{\beta}_{1,n} \rightarrow \beta_1\) quando \(n \rightarrow \infty\)

Viés: Propriedade para pequenas e grandes amostras que diz que o estimador está centrado no verdadeiro valor do parâmetro, isto é, existe uma distribuição de probabilidade (geralmente Gaussiana) em torno do valor verdadeiro, ou seja: \(E[\hat{\beta}_{1,n}] = \beta_1\) para todo \(n\). Importante: Viés e consistência são conceitos distintos. Um estimador pode ser viesado e consistente ou também pode ser não viesado e inconsistente.

Eficiência: Um conceito associado à variância do estimador. Só faz sentido falarmos de eficiência quando estamos comparando dois estimadores quaisquer que sejam ambos não viciados e consistentes.


Observação:

  • Consistência é o requerimento mais importante que um método de estimação precisa possuir. Sem ele, mesmo com amostras absurdamente grandes não identificamos o verdadeiro valor do parâmetro de interesse.

  • Para que o \(MQO\) seja consistente, \(H_1\) tem que estar garantida. Como?

  • Sabemos que tudo aquilo que não está explícito num modelo de regressão é capturado pelo termo de erro \(u_i\).

  • A hipótese \(H1\) pede que o termo de erro e sua variável explicativa \(X_i\) sejam não correlacionadas (isto é \(COV(X_i,u_i)=0\)).

  • Em geral, o modelo de regressão simples é insuficiente para explicar diversos fenômenos, assim usamos um modelo de regressão múltipla:

\(Y_i = \beta_0 + \beta_1 X_{1i} + \beta_2 X_{2i} + \ldots + \beta_k X_{ik} + u_i\)


  • Modelo de regressão múltipla permite que variáveis importantes não sejam omitidas, isto é: controlamos por tudo aquilo que achamos relevante para explicar o a variável de interesse.

  • Ainda assim, em fenômenos sociais, em que diversas variáveis afetam o objeto de estudo, mesmo modelos super parametrizados podem sofrer com variáveis omitidas.

  • Regressões com dados em painel podem mitigar os efeitos deletérios da omissão de variáveis sem observá-las diretamente sob certas circunstâncias.

Se as variáveis omitidas são associadas a distintos níveis (interceptos) para cada unidade observada, técnicas de painel pode ajudar a idetificar efeitos corretamente. Por nível, entende-se informações que variem entre as unidades observadas (indivíduos, estados, empresas…) mas que sejam fixas no tempo.


2 Dados em painel versus dados em cross-section

Dados em cross-section: Uma amostra de \(n\) indivíduos em um ponto do tempo.

Dados em painel: Uma amostra de \(n\) indivíduos em múltiplos períodos.


  • Se os mesmos indivíduos aparecem em todos os períodos observados, temos um painel balanceado.

  • Se os mesmos indivíduos desaparecem em todos os períodos observados, temos um painel desbalanceado.

  • Se indivíduos distintos aparecem na amostra, temos um pseudo-painel.

Este curso abrangerá métricas para lidar com painéis balanceados.


Notação:

Como um dado em painel abrange \(n\) observações cross-sectionais (indivíduos) em \(T\) períodos, uma observação, em particular, é denotada com dois subscritos: \(i,t\)

\(Y_{i,t} = \beta_0 + \beta_1 X_{i,t} + u_{i,t}\)

Assim, o vetor \((Y_{i,t},X_{i,t},u_{i,t})\) é composto, respectivamente, pela variável dependente, variável explicativa e o termo de erro (choque aleatório) para o indivíduo \(i\) no período \(t\).


3 Exemplo: Impostos sobre cerveja e mortes no trânsito

  • Estima-se que há cerca de 40 mil fatalidades de trânsito nos EUA a cada ano,

  • Aproximadamente 25% envolvem motoristas embriagados,

  • Governo deseja reduzir este número, uma potencial medida: aumentar impostos sobre o consumo de álcool,

  • Dados disponíveis para 48 estados Norte-Americanos sobre mortes no trânsito e impostos sobre cerveja entre 1982 e 1988,

  • Pergunta de interesse: Qual o efeito do aumento dos impostos?

  • Vamos propor o seguinte modelo para esse fenômeno:

\(Fatalities_{i,t} = \beta_0 + \beta_1 beertax_{i,t} + u_{i,t}\)

Conjunto de dados disponível em: http://fmwww.bc.edu/ec-p/data/stockwatson/datasets.list.html, ou no pacote AER.

Code
rm(list=ls())
library(AER)
library(tidyverse)
library(ggplot2)
library(sjPlot)
library(stargazer)

data("Fatalities")
Fatalities$fat.tax <- 100000*Fatalities$fatal/Fatalities$pop
ggplot(data = Fatalities, aes(x = beertax, y = fat.tax )) +
   geom_point() +   
   geom_smooth(method = "lm", se = FALSE, color = "blue") +
   labs(title = "Fatalidades vs. impostos",
       x = "Impostos (%)",
       y = "Fatalidades - por 100k hab")

Code
# A linha de regressão estimada encontra-se a seguir: 
lm(fat.tax~beertax, data = Fatalities) |> 
   tab_model(show.ci = F,
             dv.labels = "Fatalidades - por 100k hab")
  Fatalidades - por 100k hab
Predictors Estimates p
(Intercept) 18.53 <0.001
beertax 3.65 <0.001
Observations 336
R2 / R2 adjusted 0.093 / 0.091

Por que essa relação aparenta ser positiva?

  • A regressão acima usando os dados de 1982 à 1988 possue um viés de variável omitida. Mas qual?

  • É natural que os estados que tenham mais acidentes possuam maiores níveis de impostos. Por quê?

    • Impostos maiores são cobrados em mercados mais inelásticos;
    • Mercados são inelásticos quando, em geral, há um grande preferência dos consumidores naquele local pelo bem em questão;
  • A inelasticidade leva à níveis de consumo altos, logo:


Vamos revisitar a figura identificando os estados:

Code
ggplot(data = Fatalities, aes(x = beertax, y = fat.tax, color = state )) +
  geom_point() +
  labs(title = "Fatalidades vs. impostos",
       x = "Impostos (%)",
       y = "Fatalidades - por 100k hab") +
   geom_smooth(method = "lm", se = FALSE, color = "blue") 

Code
ggplot(data = Fatalities, aes(x = state, y = fat.tax, color = state)) + 
   geom_point() +
   scale_x_discrete(breaks ="") + 
   labs(title = "Fatalidades por estado: 1982-1988")

  • Percebemos que esses dados contêm clusters por estado e que, em média, estados com mais impostos ocupam níveis de fatalidades maiores,

  • Modelos que não consideram diferentes níveis produzem estimativas equivocadas dos efeitos de \(X\) em \(Y\)

  • Precisamos considerar então esses níveis em nossas regressões,

  • O nível pode ser visto como uma constante específica para cada estado, algo que não varia no tempo, assim, nosso modelo agora seria:

\(Fatalities_{i,t} = \beta_0 + \beta_1 beertax_{i,t} + \alpha_i + u_{i,t}\)

  • O termo \(\alpha_i\) representa o nível de consumo de cada estado, que está associado as preferências daquela localidade, que assumimos ser estável no tempo (pelo menos em uma janela curta),

  • Ponto principal aqui é que nosso MQO não considerou a existência de \(\alpha_i\) e como vimos, ele pode (e está!) correlacionado com a variável beertax, logo temos o endogeneidade, isto é: o erro de nossa regressão está correlacionado com a variável explicativa, violando \(H_1\),

  • É para lidar com situações como esta que a estrutura de dados de painel, com suas respectivas técnicas, existe.


4 Técnicas de painel: O estimador de primeira diferença

  • Se o problema identificado anteriomente é a existência de um nível, \(\alpha_i\), tudo que gostaríamos é remover (ou contornar) esse problema.

  • Uma potencial solução? Modelarmos variações no tempo. Ou seja, se fizermos a seguinte transformação nos dados:

\[\Delta Fatalities_{i,t} = Fatalities_{i,t} - Fatalities_{i,t-1} = \] \[(\beta_0 + \beta_1 beertax_{i,t} + \alpha_i + u_{i,t}) - (\beta_0 + \beta_1 beertax_{i,t-1} + \alpha_i + u_{i,t-1}) = \] \[\beta_1 \Delta beertax_{i,t} + \Delta u_{i,t}\]

  • Pensar em termos de taxa de variação, ou seja, uma estrutura antes e depois, teremos um modelo livre dos efeitos de nível (chamados de efeitos específicos) que mantêm as elasticidades originais,

  • Esta estratégia será chamada de estimador de primeira diferenças (FD),

  • Vamos visualizar o gráfico anterior utilizando as primeiras diferença.

Code
Fatalities <- Fatalities |>
   arrange(state,year) |>   # Organiza o banco em Estado e Ano
   group_by(state)  |>      # Declara que a transformação se dará intra grupo
   mutate(beertax.diff   = c(NA,diff(beertax)),         # Crio as variáveis
          fatalites.diff = c(NA,diff(100000*fatal/pop))
          ) 

ggplot(data = Fatalities, aes(x = beertax.diff, y =fatalites.diff)) +
  geom_point() +
  labs(title = "Fatalidades vs. impostos",
       x = "Impostos (%)",
       y = "Fatalidades - por 100k hab") +
   geom_smooth(method = "lm", se = FALSE, color = "blue") +
   stat_ellipse() + 
   stat_ellipse(level = .8, color = 3) + 
   stat_ellipse(level = .5, color = 2)  

  • Percebe-se um claro problema para o fenômeno com essa estrutura: a maior parte das varições de impostos está concentrada em torno de zero, o que é óbvio, pois mudanças de alíquotas, além de impopulares costumam necessitar de aprovação legislativa,

  • Neste cenário, é patente que teremos uma regressão não significativa! Não há muita variabilidade na variável explicativa,

Code
lm(fatalites.diff~beertax.diff, data = Fatalities) |> 
   tab_model(show.ci = F,
             dv.labels = "Fatalidades - por 100k hab")
  Fatalidades - por 100k hab
Predictors Estimates p
(Intercept) -0.03 0.792
beertax diff 0.14 0.962
Observations 288
R2 / R2 adjusted 0.000 / -0.003
  • Podemos tentar capturar mudanças mais longas usando a variação acumulada entre 1982 e 1988,

  • Vamos analisar os gráficos para 1982, 1988 e para a diferença, isto é 1988-1982:

Code
library(ggpubr) # Premite plotar vários gráficos numa mesma janela
# Recalculando as diferenças: 
Fatalities.subset <- Fatalities |> filter(year==1982|year==1988) |>
   arrange(state,year) |>   # Organiza o banco em Estado e Ano
   group_by(state)  |>      # Declara que a transformação se dará intra grupo
   mutate(beertax.diff   = c(NA,diff(beertax)),         # Crio as variáveis
          fatalites.diff = c(NA,diff(100000*fatal/pop))
          )

g1 <- ggplot(data = Fatalities.subset |> filter(year==1982)
             , aes(x = beertax, y = 100000*fatal/pop)) +
        geom_point() +
        labs(title = "1982",
             x = "Impostos (%)",
             y = "Fatalidades - por 100k hab") +
        geom_smooth(method = "lm", se = FALSE, color = "blue") 

g2 <- ggplot(data = Fatalities.subset |> filter(year==1988)
             , aes(x = beertax, y = 100000*fatal/pop)) +
        geom_point() +
        labs(title = "1988",
             x = "Impostos (%)",
             y = "Fatalidades - por 100k hab") +
        geom_smooth(method = "lm", se = FALSE, color = "blue") 

g3 <- ggplot(data = Fatalities.subset |> filter(year==1988)
             , aes(x = beertax.diff, y =fatalites.diff)) +
        geom_point() +
        labs(title = "Diferença: 1988 - 1982",
             x = "Impostos (%)",
             y = "Fatalidades - por 100k hab") +
        geom_smooth(method = "lm", se = FALSE, color = "red") 

# Coloca os gráficos na mesma janela:
gga <- ggarrange( ggarrange(g1,g2,ncol=2), g3, nrow = 2)

# Coloca título comum:
annotate_figure(gga, top = text_grob("Fatalidades vs. impostos", 
               color = "red", face = "bold", size = 14))

As regressões para esses três períodos são:

Code
tab_model(
   lm(fat.tax~beertax, data = Fatalities.subset |> filter(year==1982)),
   lm(fat.tax~beertax, data = Fatalities.subset |> filter(year==1988)),
   lm(fatalites.diff~beertax.diff, data = Fatalities.subset |> filter(year==1988)),
   show.ci = F,
   dv.labels = c("1982","1988","1988-1982"))
  1982 1988 1988-1982
Predictors Estimates p Estimates p Estimates p
(Intercept) 20.10 <0.001 18.59 <0.001 -0.72 0.241
beertax 1.48 0.435 4.39 0.011
beertax diff -10.41 0.016
Observations 48 48 48
R2 / R2 adjusted 0.013 / -0.008 0.134 / 0.115 0.119 / 0.100
  • O estimador de FD foi a primeira proposta para lidar com heterogeneidade não observada fixa no tempo,

  • Entretanto, trata-se de uma alternativa pouco eficientes (seus desvios-padrão são maiores do que outras estratégias) produzindo, muitas vezes, estimativas não signidficativas, em especial quando a dimensão temporal é pequena.


5 Técnicas de painel: MQO com dummies (LSDV) o estimador within (FE)

  • Para o modelo \(Y_{i,t} = \beta_0 + \beta_1 X_{i,t} + \alpha_i + u_{i,t} \ \), \(\alpha_i\) é uma variável não observada que varia entre estados mas é fixa no tempo como, por exemplo, cultura local quanto dirigir embriagado,

  • Nosso parâmetro de interesse é \(\beta_1\). Se observássemos \(\alpha_i\) poderíamos colocá-lo diretamente na equação,

  • Como levar em consideração \(\alpha_i\) se não o observamos diretamente? Existem algumas formas de fazê-lo,

  • Podemos enxergar que \(\alpha_i\) como sendo um intercepto específico para cada unidade observada, que captura o impacto de todas as variáveis fixas no tempo em \(Y_{it}\), assim o chamamos de efeito fixo (efeito específico)


5.1 O estimador LSDV

Para estimar \(\alpha_i\) uma solução é a introdução de dummies de indivíduos, ou seja, criar um conjunto de \(n\) variáveis \(D_1,D_2,D_3,\ldots,D_n\), sendo:


\[ \begin{cases} D_1 = 1\; \text{ se } i=1 & \text{e} & D_1 = 0\; \text{ se } i \neq 1 \\ D_2 = 1\; \text{ se } i=1 & \text{e} & D_2 = 0\; \text{ se } i \neq 2 \\ &... \\ D_n = 1\; \text{ se } i=n & \text{e} & D_n = 0\; \text{ se } i \neq n \\ \end{cases} \]

O modelo de regressão passará a ser escrito como: \[Y_{i,t} = \beta_1 X_{i,t} + \alpha_1 D_1 + \alpha_2 D_2 + \ldots + \alpha_n D_n + u_{i,t}\] Ou alternativamente: \[Y_{i,t} = \beta_0 + \beta_1 X_{i,t} + \gamma_1 D_1 + \gamma_2 D_2 + \ldots + \gamma_n D_n + u_{i,t}\] em que \(\beta_0=\alpha_1\) e \(\gamma_i = \alpha_i - \beta_0\)

  • Lembre-se: uma dummy nada mais é do que a diferença de média entre os grupos de análise.

  • No nosso exemplo estudado, se colocarmos uma dummy para cada estado, estamos capturando tudo aquilo que é inerente àquele estado porém fixo no tempo! Ou seja: \(\alpha_i\).

Code
tab_model(
   (LSDV <- lm(fat.tax~beertax + state, data = Fatalities)), 
   dv.labels = "Fatalidades - tx. por 100k hab. - LSDV",
   show.ci = F
)
  Fatalidades - tx. por 100k hab. - LSDV
Predictors Estimates p
(Intercept) 34.78 <0.001
beertax -6.56 0.001
state [az] -5.68 0.034
state [ar] -6.55 0.003
state [ca] -15.09 <0.001
state [co] -14.84 <0.001
state [ct] -18.62 <0.001
state [de] -13.08 <0.001
state [fl] -2.68 0.055
state [ga] 5.25 0.005
state [id] -6.69 0.010
state [il] -19.62 <0.001
state [in] -14.62 <0.001
state [ia] -15.44 <0.001
state [ks] -12.23 <0.001
state [ky] -12.18 <0.001
state [la] -8.47 <0.001
state [me] -11.08 <0.001
state [md] -17.06 <0.001
state [ma] -21.10 <0.001
state [mi] -14.85 <0.001
state [mn] -18.97 <0.001
state [ms] -0.29 0.845
state [mo] -12.96 <0.001
state [mt] -3.60 0.173
state [ne] -15.22 <0.001
state [nv] -6.01 0.037
state [nh] -12.54 <0.001
state [nj] -21.06 <0.001
state [nm] 4.26 0.095
state [ny] -21.87 <0.001
state [nc] -2.90 0.016
state [nd] -16.23 <0.001
state [oh] -16.74 <0.001
state [ok] -5.45 0.001
state [or] -11.68 <0.001
state [pa] -17.67 <0.001
state [ri] -22.65 <0.001
state [sc] 5.57 <0.001
state [sd] -10.04 <0.001
state [tn] -8.76 0.001
state [tx] -9.17 <0.001
state [ut] -11.64 <0.001
state [vt] -9.66 <0.001
state [va] -12.90 <0.001
state [wa] -16.60 <0.001
state [wv] -8.97 <0.001
state [wi] -17.59 <0.001
state [wy] -2.29 0.466
Observations 336
R2 / R2 adjusted 0.905 / 0.889

  • É natural que esta estratégia apresente um R² alto, afinal colocamos 47 dummies na equação portanto, não se empolgue…

  • No modelo acima um procedimento que entrega estimativas idênticas é o chamado estimador whitin,

  • Para aplicar o estimador within faz-se um tratamento prévio dos dados removendo as médias específicas de cada indivíduo, isto é desviando \(Y_{it}\) e \(X_{it}\) de suas médias,

  • Após a remoção das médias específicas aplica-se o MQO nos dados transformados,

  • Este é um procedimento especialmente interessante quando o número de unidades modeladas \(n\) for grande demais face a dimensão temporal \(T\).


5.2 O estimador within

  • Essência do estimador within e fazer uma regressão com dados desviados da média,

  • Calcula-se, para todas as variáveis do modelo, as médias de cada indívuo, isto é:

\[\overline{Y}_i = \beta_1 \overline{X}_{i} + \alpha_i + \overline{u}_{i}\], em que: \(\overline{Y}_i = \frac{\sum_{t=1}^T Y_{i,t}}{T}\) e \(\overline{X}_i = \frac{\sum_{t=1}^T X_{i,t}}{T}\). Depois, as subtraímos dos dados originais:

\[Y_{i,t} - \overline{Y}_i = (\beta_1 X_{i,t} + \alpha_i + u_{i,t}) - (\beta_1 \overline{X}_{i} + \alpha_i + \overline{u}_{i})\] rearranjando, tem-se:

\[ \underbrace{(Y_{i,t} - \overline{Y}_i)}_{\tilde{y}_{i,t}} = \beta_1 \underbrace{(X_{i,t} - \overline{X}_i)}_{\tilde{x}_{i,t}} + \underbrace{(u_{i,t} - \overline{u}_i)}_{\tilde{u}_{i,t}}\]

  • Agora, estima-se a equação: \(\tilde{y}_{i,t} = \beta_1 \tilde{x}_{i,t} + \tilde{u}_{i,t}\) por MQO.

  • Onde foi para o efeito específico? Ele foi removido quando subtraímos da média,

  • Essa transformação é idêntica a colocar uma dummie para cada estado,

  • Para visualizarmos isso vamos comparar as duas estratégias usando o pacote plm que contêm funções diversas para estimativas de dados em painel.


  • O primeiro passo é declarar o dado como um painel a fim de que o R entenda quais dimensões são o indivíduo e o tempo.
Code
library(plm) # Biblioteca com as funcionalidades para painel
# Converta os dados para o formato de painel
panel_data <- pdata.frame(Fatalities, index = c("state", "year"))

class(panel_data)
[1] "pdata.frame" "data.frame" 
Code
# Os dados ficarão organizados por state e year:
head(panel_data[,1:9],14)
        state year spirits unemp   income   emppop   beertax baptist  mormon
al-1982    al 1982    1.37  14.4 10544.15 50.69204 1.5393795 30.3557 0.32829
al-1983    al 1983    1.36  13.7 10732.80 52.14703 1.7889907 30.3336 0.34341
al-1984    al 1984    1.32  11.1 11108.79 54.16809 1.7142856 30.3115 0.35924
al-1985    al 1985    1.28   8.9 11332.63 55.27114 1.6525424 30.2895 0.37579
al-1986    al 1986    1.23   9.8 11661.51 56.51450 1.6099070 30.2674 0.39311
al-1987    al 1987    1.18   7.8 11944.00 57.50988 1.5599999 30.2453 0.41123
al-1988    al 1988    1.17   7.2 12368.62 56.83453 1.5014436 30.2233 0.43018
az-1982    az 1982    1.97   9.9 12309.07 56.89330 0.2147971  3.9589 4.91910
az-1983    az 1983    1.90   9.1 12693.81 57.55363 0.2064220  3.8901 4.83107
az-1984    az 1984    2.14   5.0 13265.93 60.37902 0.2967033  3.8226 4.74461
az-1985    az 1985    1.86   6.5 13726.70 58.64853 0.3813559  3.7562 4.65971
az-1986    az 1986    1.78   6.9 14107.33 60.28018 0.3715170  3.6910 4.57632
az-1987    az 1987    1.72   6.2 14241.00 60.21506 0.3600000  3.6269 4.49442
az-1988    az 1988    1.68   6.3 14408.08 60.49767 0.3464870  3.5640 4.41399
  • Notem que as linhas agora são nomeadas com a sigla e o ano, ou seja, agora o R entende que temos um painel, vamos estimar o modelo within, também conhecido como efeitos fixos (EF) com o comando a seguir:

Code
EF <- plm(fat.tax~beertax,
          data = panel_data,
          model = "within")

tab_model(EF, LSDV,
          show.ci = F,
          digits = 4,
          dv.labels = c("EF","LSDV")
)
  EF LSDV
Predictors Estimates p Estimates p
beertax -6.5587 0.001 -6.5587 0.001
(Intercept) 34.7763 <0.001
state [az] -5.6773 0.034
state [ar] -6.5495 0.003
state [ca] -15.0947 <0.001
state [co] -14.8428 <0.001
state [ct] -18.6226 <0.001
state [de] -13.0760 <0.001
state [fl] -2.6813 0.055
state [ga] 5.2460 0.005
state [id] -6.6902 0.010
state [il] -19.6162 <0.001
state [in] -14.6154 <0.001
state [ia] -15.4393 <0.001
state [ks] -12.2322 <0.001
state [ky] -12.1752 <0.001
state [la] -8.4712 <0.001
state [me] -11.0795 <0.001
state [md] -17.0644 <0.001
state [ma] -21.0975 <0.001
state [mi] -14.8453 <0.001
state [mn] -18.9721 <0.001
state [ms] -0.2908 0.845
state [mo] -12.9626 <0.001
state [mt] -3.6039 0.173
state [ne] -15.2218 <0.001
state [nv] -6.0077 0.037
state [nh] -12.5445 <0.001
state [nj] -21.0575 <0.001
state [nm] 4.2637 0.095
state [ny] -21.8667 <0.001
state [nc] -2.9047 0.016
state [nd] -16.2344 <0.001
state [oh] -16.7442 <0.001
state [ok] -5.4506 0.001
state [or] -11.6800 <0.001
state [pa] -17.6747 <0.001
state [ri] -22.6505 <0.001
state [sc] 5.5717 <0.001
state [sd] -10.0372 <0.001
state [tn] -8.7566 0.001
state [tx] -9.1747 <0.001
state [ut] -11.6395 <0.001
state [vt] -9.6604 <0.001
state [va] -12.9018 <0.001
state [wa] -16.5952 <0.001
state [wv] -8.9675 <0.001
state [wi] -17.5927 <0.001
state [wy] -2.2850 0.466
Observations 336 336
R2 / R2 adjusted 0.041 / -0.120 0.905 / 0.889

5.3 Efeitos fixos de tempo

  • Assim como unidade distintas podem possuir sua própria média \(\alpha_i\), o fenômeno sob estudo também pode sofrer variações específicas em um determinado período de tempo,

  • Estamos lidando com choques que seriam comuns a todas unidades mas que variam de um ano para o outro,

  • Podemos adicionar esses efeitos de tempo nessas regressões a fim de capturar variáveis omitidas que são fixas nas cross-sections mas variáveis no tempo,

  • Típicos exemplos de efeitos de tempo: condições macroeconômicas, mudanças de legislação. O modelo passaria a ser:

\[Y_{i,t} = \beta_1 X_{i,t} + \alpha_i + \lambda_t + u_{i,t}\]

  • Da mesma forma que fizemos anteriormente, podemos capturar tal efeito criando-se T variáveis dummy correspondentes à cada período de tempo:

\[ \begin{cases} B_1 = 1\; \text{ se } t=1 & \text{e} & B_1 = 0\; \text{ se } t \neq 1 \\ B_2 = 1\; \text{ se } t=2 & \text{e} & B_2 = 0\; \text{ se } t \neq 2 \\ &... \\ B_T = 1\; \text{ se } t=T & \text{e} & B_T = 0\; \text{ se } t \neq T \\ \end{cases} \] * O modelo passaria a ser:

\[Y_{i,t} = \beta_1 X_{i,t} + \alpha_1 D_1 + \alpha_2 D_2 + \ldots + \alpha_n D_n + \lambda_{t_1}B_1 + \ldots + \lambda_{T}B_T + u_{i,t}\]

Vamos voltar ao nosso exemplo para verficar se existem tendências para as mortes no trânsito a cada ano:

Code
gplots::plotmeans(fat.tax ~ year, data = panel_data)

Code
ggplot(data = Fatalities, aes(x = beertax, y = fat.tax)) + 
         geom_point() +
         labs(title = "Fatalidades") +
         geom_smooth(aes(x = beertax, y = fat.tax,group=year,color=year), 
                     method = "lm", se = FALSE)

O modelo que leva em consideração os efeitos fixos de indivíduo e tempo é o conhecido como two ways fixed effects. Ele pode ser obtido tanto através de dummies como através de através do cálculo de desvios de médias individuais e temporais.

As duas estratégias vão produzir estimativas idênticas para as elasticidades:

Code
tab_model(
   TWFE    <- plm(fat.tax ~ beertax, data = panel_data, model = "within",
                  effect = "twoways"),
   FE.t    <- plm(fat.tax ~ beertax + year, data = panel_data, model = "within"),                      
   LSDV.it <- lm(fat.tax  ~ beertax + state + year, data = Fatalities),
   dv.labels = c("TWFE","FE.t","LSDV.it"),
   digits = 4,
   show.ci = F
)
  TWFE FE.t LSDV.it
Predictors Estimates p Estimates p Estimates p
beertax -6.3998 0.001 -6.3998 0.001 -6.3998 0.001
year [1983] -0.7990 0.038 -0.7990 0.038
year [1984] -0.7242 0.060 -0.7242 0.060
year [1985] -1.2398 0.001 -1.2398 0.001
year [1986] -0.3786 0.327 -0.3786 0.327
year [1987] -0.5090 0.193 -0.5090 0.193
year [1988] -0.5180 0.192 -0.5180 0.192
(Intercept) 35.1137 <0.001
state [az] -5.4686 0.050
state [ar] -6.3853 0.005
state [ca] -14.8519 <0.001
state [co] -14.6153 <0.001
state [ct] -18.4013 <0.001
state [de] -12.8426 <0.001
state [fl] -2.6005 0.068
state [ga] 5.1162 0.007
state [id] -6.4896 0.016
state [il] -19.3850 <0.001
state [in] -14.4014 <0.001
state [ia] -15.2428 <0.001
state [ks] -12.0431 <0.001
state [ky] -11.9479 <0.001
state [la] -8.3366 <0.001
state [me] -10.9424 <0.001
state [md] -16.8407 <0.001
state [ma] -20.8802 <0.001
state [mi] -14.6650 <0.001
state [mn] -18.7649 <0.001
state [ms] -0.1991 0.896
state [mo] -12.7540 <0.001
state [mt] -3.3977 0.218
state [ne] -15.0291 <0.001
state [nv] -5.7816 0.054
state [nh] -12.3893 <0.001
state [nj] -20.8122 <0.001
state [nm] 4.4611 0.093
state [ny] -21.6288 <0.001
state [nc] -2.8507 0.019
state [nd] -16.0376 <0.001
state [oh] -16.5474 <0.001
state [ok] -5.3361 0.002
state [or] -11.4540 <0.001
state [pa] -17.4574 <0.001
state [ri] -22.4173 <0.001
state [sc] 5.5358 <0.001
state [sd] -9.8820 <0.001
state [tn] -8.5467 0.002
state [tx] -8.9855 0.001
state [ut] -11.4972 <0.001
state [vt] -9.5038 <0.001
state [va] -12.7520 <0.001
state [wa] -16.3713 <0.001
state [wv] -8.7774 0.001
state [wi] -17.3593 <0.001
state [wy] -2.0346 0.534
Observations 336 336 336
R2 / R2 adjusted 0.036 / -0.149 0.080 / -0.096 0.909 / 0.891

6 Técnicas de painel: o modelo de efeitos aleatórios

  • O estimador de efeitos aleatórios (efeitos específicos aleatório), ou EA, foi delineado para uma situação bem específica: quando o efeito individual, \(\alpha_i\), não está correlacionado com nenhuma das variáveis explicativas \(X_{k,it}\),

  • No nosso exemplo discutido até aqui, suponha que taxação sobre a cerveja seja determinada de maneira aleatória e desconsiderando as preferências (elasticidades) dos consumidores locais. Sim, isso parece um tanto irracional por parte dos governos, por isso apenas suponha…

  • Neste caso os níveis ainda podem existir porém não estão mais correlacionados com a variável de interesse \(X_{i,t}\), isto é: \(COV(\alpha_i,X_{i,t}) = 0\),

  • O processo gerador de dados agora seria dado por:

Neste cenário, para o modelo: \(Y_{i,t} = \beta_1 X_{i,t} + \alpha_i + \lambda_t + u_{i,t}\), se pudermos assumir que \(COV(\alpha_i,X_{i,t}) = 0\) então o estimador de efeitos aleatórios se torna a melhor escolha!

Importante: Na situação descrita acima o estimador Within continua sendo consistente porém menos eficiente do que o EA. No entanto, se não pudermos garantir que \(COV(\alpha_i,X_{i,t}) = 0\) o EA será inconsistente e a escolha natural será o Within.


  • Para implentar o EA deve-se calcular a seguinte informação: \(\eta = 1 - \left(\frac{\sigma_{u}^{2}}{T\sigma_{\alpha}^{2} + \sigma_{u}^{2} } \right)^{1/2}\) onde \(\sigma_{u}^{2} = VAR(u_{it})\) e \(\sigma_{\alpha}^{2} = VAR(\alpha_i)\).

  • Posteriormente, multiplica-se \(\eta\) pela equação de média individual: \(\eta\bar{Y}_i = \eta\beta_1 \bar{X}_{i} + \eta(\alpha_i + \bar{u}_{i})\)

  • E aplica-se MQO na seguinte equação:

\[Y_{i,t} - \eta\bar{Y}_i = (\beta_1 X_{i,t} + \alpha_i + u_{i,t}) - \eta(\beta_1 \bar{X}_{i} + \alpha_i + \bar{u}_{i})\]

  • Como \(\sigma_{u}^{2}\) e \(\sigma_{\alpha}^{2}\) são desconhecidos a priori, implementar este estimador requer um primeiro passo para estimá-los que é o calculo de um Pooled (MQO normal)

  • Uma curiosidade sobre EA é que se \(\eta = 1\) teremos um estiador idêntico ao Within e, caso \(\eta = 0\), situação em que só seria possível quando \(\sigma_{\alpha}^{2} = 0\), o EA colapsa para o Pooled.


7 Qual modelo escolher? Testes para identificar a natureza (ou existência) do efeito específico

  • Em geral, o pesquisador não sabe definir a natureza do termo \(\alpha_i\). Isto é: ele existe? Se existe é correlacionado com a matriz de variáveis explicativas?

  • É importante ressaltar que caso \(\sigma_{\alpha}^{2} = 0\), situação em que deveríamos usar um pooled (MQO), o estimador within continua sendo uma opção válida (consistente) porém menos eficiente. O contrário não é verdadeiro,

  • Diante disto, temos uma séries de testes para averiguar tais hipóteses. O primeiro deles é um teste-F que compara p modelo Wihtin com o Pooled. Sua estatística de teste é \(F_{(n,n(T-1)-k)}\) onde \(n\) é número de unidades na cross-section sob análise e \(k\) o número de variáveis colocadas no modelo. O teste é tal que:

\[H_0: \text{efeitos específicos não significativos, isto é: } \alpha_i = 0 \forall \ i\]

\[H_1: \text{efeitos específicos significativos, isto é: } \exists \ \alpha_i \neq 0 \ \forall \ i\]

  • Caso a existência de efeitos específicos seja aceita, isto é: no teste anterior, rejeita-se \(H_0\), devemos verificar se tais efeitos estão correlacionados a matriz dos \(k\) regressores inclusos no modelo.

  • Ou seja se devemos optar pelo EA ou o Within. Para tanto, realizamos o teste de Hausman. A estrutura do teste é tal que:

\[H_0: \text{efeitos específicos não correlacionados com regressores (EA)} \]

\[H_1: \text{efeitos específicos correlacionados com regressores (Within)}\]


  • Vamos implementar, em nosso exemplo, essa estratégia de testes a fim de verificar qual modelo mais adequado,

  • A primeira etapa é, poratanto, verificar a existência de efeitos específicos, \(\alpha_i\), significativos, ou seja: vamos decidir entre um estimador de painel ou um MQO (pooled):

Code
PO <- plm(fat.tax~beertax,
          data = panel_data,
          model = "pooling")
pFtest(EF,PO)

    F test for individual effects

data:  fat.tax ~ beertax
F = 52.179, df1 = 47, df2 = 287, p-value < 2.2e-16
alternative hypothesis: significant effects
  • No teste acima, rejeitamos \(H_0\), ou seja, devemos considerar a presença de efeitos específicos,

  • Agora, vamos checar se esses efeitos são correlacionados com a variável beertax, isto é: \(COV(beertax_i,\alpha) = 0\):

Code
EA <- plm(fat.tax~beertax,
          data = panel_data,
          model = "random")
phtest(EF,EA)

    Hausman Test

data:  fat.tax ~ beertax
chisq = 18.353, df = 1, p-value = 1.835e-05
alternative hypothesis: one model is inconsistent
  • Rejeitamos \(H_0\), isso indica que \(COV(beertax_i,\alpha) \neq 0\), ou seja: Os impostos sobre cerveja estão correlacionados com nível (determinado pelas preferências locais) e, portanto, devemos usar o EF

  • Importante: nos dois teste acima, o modelo EF deve vir primeiro!

  • Também podemos testar a necessidade de efeitos fixos de tempo. Para isso usaremos o modelo EF (só com efeitos fixos de indivíduo) contra o TWFE (que caputura também os efeitos de tempo):

Code
pFtest(TWFE,EF)

    F test for twoways effects

data:  fat.tax ~ beertax
F = 2.0117, df1 = 6, df2 = 281, p-value = 0.0642
alternative hypothesis: significant effects
Code
plmtest(EF, c("time"), type=("bp"))

    Lagrange Multiplier Test - time effects (Breusch-Pagan)

data:  fat.tax ~ beertax
chisq = 1.8545, df = 1, p-value = 0.1733
alternative hypothesis: significant effects
  • Os testes acima indicam que controlar para efeitos fixos de tempo trazem poucos ganhos

  • A tabela a seguir resume a circunstâncias em que cada estimador é a melhor escolha

\(\alpha_i = 0\) \(\alpha_i \neq 0\) com \(COV(\alpha_i,X_{it})=0\) \(\alpha_i \neq 0\) com \(COV(\alpha_i,X_{it}) \neq 0\)
Consistência Eficiência Consistência Eficiência Consistência Eficiência
Pooled X X X - - -
Efeito Aleatório X - X X - -
Within X - X - X X
FD X - X - X -

8 Diagnósticos pós estimação:

Tomemos o modelo mais completo: \(Y_{i,t} = \beta_1 X_{i,t} + \alpha_i + \lambda_t + u_{i,t}\)

  • As hipóteses que garantem a consistência, ausência de vício e eficiência são:

    • \(H_1\): \(E[u_{i,t}|X_{i,1},\ldots, X_{i,T}, \alpha_i, \lambda_t] = 0\),

    • \(H_2\): \(\left(X_{i,1}, \ldots, X_{i,T}; Y_{i,1}, \ldots, Y_{i,T} \right)\) são independentes nas cross-sections,

    • \(H_3\): grandes outliers são improváveis,

    • \(H_4\): não há multicolinearidade,

    • \(H_5\): \(COV\left(u_{i,t}u_{i,s}|X_{i,1},\ldots, X_{i,T}, \alpha_i, \lambda_t\right) = 0\) para \(t\neq s\),

    • \(H_6\): \(E\left(u_{i,t}^2|X_{i,1},\ldots, X_{i,T}, \alpha_i, \lambda_t\right) = \sigma^2_u\).

  • \(H_1\) é a hipótese mais importante, é a que garante consistência,

  • \(H_3\) garante ausência de viés,

  • \(H_2\) (indivíduos observados são independentes) é facilmente preenchida mas, se violada, afeta eficiência,

  • \(H_4\) afeta eficiência,

  • \(H_5\) (ausência de correlação temporal intra indivíduo) também afeta eficiência e é facilmente violada: é razoável assumir que os erros são correlacionado no tempo,

  • \(H_6\) (hipótese de homocedasticidade) também afeta eficiência e é facilmente violada.

Temos, então que \(H_2, \ H_5\) e \(H_6\) quando violadas afetam eficiência.

INEFICIÊNCIA –> Desvios-padrão inadequados –> testes de significância comprometidos


Temos testes para cada uma dessas circunstâncias e caso elas estejam presentes precisaremos corrigir os desvios padrão a fim de realizarmos uma inferência correta!

Comecemos com um teste para \(H_2\): teste de independência cross-sectional. Vamos avaliar os dois modelos de efeitos fixos, isto é: com e sem efeitos fixos de tempo.

Temos implementados no pacote plm, dois testes: O teste Lagrange-Multiplier e o teste do Pesaran:

Code
pcdtest(EF, test = c("lm"))

    Breusch-Pagan LM test for cross-sectional dependence in panels

data:  fat.tax ~ beertax
chisq = 1899.2, df = 1128, p-value < 2.2e-16
alternative hypothesis: cross-sectional dependence
Code
pcdtest(TWFE, test = c("lm"))

    Breusch-Pagan LM test for cross-sectional dependence in panels

data:  fat.tax ~ beertax
chisq = 1987.5, df = 1128, p-value < 2.2e-16
alternative hypothesis: cross-sectional dependence
Code
pcdtest(EF, test = c("cd"))

    Pesaran CD test for cross-sectional dependence in panels

data:  fat.tax ~ beertax
z = 5.436, p-value = 5.449e-08
alternative hypothesis: cross-sectional dependence
Code
pcdtest(TWFE, test = c("cd"))

    Pesaran CD test for cross-sectional dependence in panels

data:  fat.tax ~ beertax
z = 1.1956, p-value = 0.2318
alternative hypothesis: cross-sectional dependence
  • É interessante notar que o TWEF não aponta, pelo teste do Pesaran, depencdência cross-sectional.

  • Um teste de heterocedasticidade tradiocional pode ser aplicado, no caso em tela, podemos usar o Breusch-Pagan:

Code
bptest(EF, studentize=F)

    Breusch-Pagan test

data:  EF
BP = 17.184, df = 1, p-value = 3.393e-05
Code
bptest(TWFE, studentize=F)

    Breusch-Pagan test

data:  TWFE
BP = 17.184, df = 1, p-value = 3.393e-05
  • Os dois modelos apresentam heterocedasticidade.
Code
Fatalities$Pred.EF   <- fitted(EF)
Fatalities$Pred.TWFE <- fitted(TWFE)
Fatalities$res.EF    <- resid(EF)
Fatalities$res.TWFE  <- resid(TWFE)

g1 <- ggplot(data = Fatalities, aes(x = Pred.EF, y = res.EF^2)) + 
         geom_point() +
         labs(title = "Resíduos ao quadrado - EF") +
         geom_smooth(method = "loess", se = FALSE)

g2 <- ggplot(data = Fatalities, aes(x = Pred.TWFE, y = res.TWFE^2)) + 
         geom_point() +
         labs(title = "Resíduos ao quadrado - TWFE") +
         geom_smooth(method = "loess", se = FALSE)

ggarrange(g1,g2)

  • Um teste de auto correlação serial, mostra que o problema está presente
Code
pbgtest(EF)

    Breusch-Godfrey/Wooldridge test for serial correlation in panel models

data:  fat.tax ~ beertax
chisq = 50.668, df = 7, p-value = 1.068e-08
alternative hypothesis: serial correlation in idiosyncratic errors
Code
pbgtest(TWFE)

    Breusch-Godfrey/Wooldridge test for serial correlation in panel models

data:  fat.tax ~ beertax
chisq = 54.188, df = 7, p-value = 2.159e-09
alternative hypothesis: serial correlation in idiosyncratic errors
  • A função –vcovHC- possui três estimadores de covariância consistentes com a presença de heterocedasticidade:

    • “white1” - para heterocedasticidade geral, mas sem correlação serial. Recomendado para efeitos aleatórios,
    • “white2” - é o “white1” restrito a uma variância comum dentro dos grupos. Recomendado para efeitos aleatórios,
    • “arellano” - heteroscedasticidade e correlação serial. Recomendado para efeitos efeitos fixos,
  • As seguintes opções se aplicam:

    • HC0 - heterocedasticidade consistente (O estimador de white padrão),
    • HC1,HC2, HC3 - Recomendado para amostras pequenas. HC3 dá menos peso a observações observações influentes (outliers),
    • HC4 - amostras pequenas com observações influentes,
  • Vamos ver, a seguir, como realizar as correções:

Code
coeftest(EF, vcovHC) # O default é HC0

t test of coefficients:

        Estimate Std. Error t value Pr(>|t|)  
beertax  -6.5587     2.8837 -2.2744  0.02368 *
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Code
coeftest(EF, vcovHC(EF, type = "HC3"))

t test of coefficients:

        Estimate Std. Error t value Pr(>|t|)  
beertax  -6.5587     2.9779 -2.2025  0.02843 *
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Code
coeftest(EF, vcovHC(EF, method ="arellano"))

t test of coefficients:

        Estimate Std. Error t value Pr(>|t|)  
beertax  -6.5587     2.8837 -2.2744  0.02368 *
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Code
coeftest(EF, vcovHC(EF, method ="arellano",type = "HC3"))

t test of coefficients:

        Estimate Std. Error t value Pr(>|t|)  
beertax  -6.5587     2.9779 -2.2025  0.02843 *
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  • Um truque para comparar os diferentes erros-padrão usando o comando “tab_model”, é criar um novo objeto do tipo “model” mas substituir a matrix de variância pela desejada
Code
ARL0  <- vcovHC(EF, method = "arellano") # O default é HC0
ARL3  <- vcovHC(EF, method = "arellano", type = "HC3") 
WHT0  <- vcovHC(EF, method ="white1")
WHT3  <- vcovHC(EF, method ="white1", type = "HC3")

EF.ARL0 <- EF; EF.ARL0$vcov <- ARL0
EF.ARL3 <- EF; EF.ARL3$vcov <- ARL3
EF.WHT0 <- EF; EF.WHT0$vcov <- WHT0
EF.WHT3 <- EF; EF.WHT3$vcov <- WHT3


tab_model(EF.ARL0, EF.ARL3, EF.WHT0, EF.WHT3,
          dv.labels = c("Arellano HC0","Arellano HC1","White HC0", "White HC3"), 
          show.se = T,
          show.ci = F,
          digits = 4
)
  Arellano HC0 Arellano HC1 White HC0 White HC3
Predictors Estimates std. Error p Estimates std. Error p Estimates std. Error p Estimates std. Error p
beertax -6.5587 2.8837 0.024 -6.5587 2.9779 0.028 -6.5587 1.8787 0.001 -6.5587 1.9512 0.001
Observations 336 336 336 336
R2 / R2 adjusted 0.041 / -0.120 0.041 / -0.120 0.041 / -0.120 0.041 / -0.120
  • Os estimadores de matriz robusta para objetos do tipo “plm” já criam, automaticamente, clusters na variável de unidade de análise (grupo)
Code
rbind(
coeftest(EF, vcov = vcovHC(EF, method = "arellano", type = "HC0")),
coeftest(EF, vcov = vcovHC(EF, method = "arellano", type = "HC0", cluster = "group")),
coeftest(EF, vcov = vcovHC(EF, method = "arellano", type = "HC0", cluster = "time"))
) |> `row.names<-`(c("beertax (defaut)", "beertax (group)", "beertax (time)"))
                  Estimate Std. Error   t value     Pr(>|t|)
beertax (defaut) -6.558737  2.8836811 -2.274432 2.367803e-02
beertax (group)  -6.558737  2.8836811 -2.274432 2.367803e-02
beertax (time)   -6.558737  0.9457327 -6.935085 2.691410e-11