Question

J'ai un ensemble de données de panel dans R (temps et section transversale) et je voudrais calculer des erreurs standard qui sont regroupées en deux dimensions, car mes résidus sont corrélés dans les deux sens.En cherchant sur Google, j'ai trouvé http://thetarzan.wordpress.com / 2011/06/11 / clustered-standard-errors-in-r / qui fournit une fonction pour ce faire.Cela semble un peu ad hoc alors je voulais savoir s'il y a un package qui a été testé et est-ce que cela?

Je sais que sandwich fait des erreurs standard HAC, mais il ne fait pas de double clustering (c'est-à-dire le long de deux dimensions).

Était-ce utile?

La solution

Le paquet rms de Frank Harrell (qui s'appelait auparavant Design) a une fonction que j'utilise souvent lors du clustering: robcov.

Voir cette partie de ?robcov, par exemple.

cluster: a variable indicating groupings. ‘cluster’ may be any type of
      vector (factor, character, integer).  NAs are not allowed.
      Unique values of ‘cluster’ indicate possibly correlated
      groupings of observations. Note the data used in the fit and
      stored in ‘fit$x’ and ‘fit$y’ may have had observations
      containing missing values deleted. It is assumed that if any
      NAs were removed during the original model fitting, an
      ‘naresid’ function exists to restore NAs so that the rows of
      the score matrix coincide with ‘cluster’. If ‘cluster’ is
      omitted, it defaults to the integers 1,2,...,n to obtain the
      "sandwich" robust covariance matrix estimate.

Autres conseils

Pour les régressions de panel, le package plm peut estimer les SE groupés selon deux dimensions.

Utilisation de M. Résultats de référence de Petersen :

require(foreign)
require(plm)
require(lmtest)
test <- read.dta("http://www.kellogg.northwestern.edu/faculty/petersen/htm/papers/se/test_data.dta")

##Double-clustering formula (Thompson, 2011)
vcovDC <- function(x, ...){
    vcovHC(x, cluster="group", ...) + vcovHC(x, cluster="time", ...) - 
        vcovHC(x, method="white1", ...)
}

fpm <- plm(y ~ x, test, model='pooling', index=c('firmid', 'year'))

Pour que vous puissiez maintenant obtenir des SE en cluster:

##Clustered by *group*
> coeftest(fpm, vcov=function(x) vcovHC(x, cluster="group", type="HC1"))

t test of coefficients:

            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 0.029680   0.066952  0.4433   0.6576    
x           1.034833   0.050550 20.4714   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

##Clustered by *time*
> coeftest(fpm, vcov=function(x) vcovHC(x, cluster="time", type="HC1"))

t test of coefficients:

            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 0.029680   0.022189  1.3376   0.1811    
x           1.034833   0.031679 32.6666   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

##Clustered by *group* and *time*
> coeftest(fpm, vcov=function(x) vcovDC(x, type="HC1"))

t test of coefficients:

            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 0.029680   0.064580  0.4596   0.6458    
x           1.034833   0.052465 19.7243   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Pour plus de détails, voir:


Cependant, ce qui précède ne fonctionne que si vos données peuvent être forcées à un pdata.frame. Il échouera si vous avez "duplicate couples (time-id)". Dans ce cas, vous pouvez toujours regrouper, mais uniquement selon une dimension.

Faites croire à plm que vous disposez d'un ensemble de données de panel approprié en spécifiant un seul un index:

fpm.tr <- plm(y ~ x, test, model='pooling', index=c('firmid'))

Pour que vous puissiez maintenant obtenir des SE en cluster:

##Clustered by *group*
> coeftest(fpm.tr, vcov=function(x) vcovHC(x, cluster="group", type="HC1"))

t test of coefficients:

            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 0.029680   0.066952  0.4433   0.6576    
x           1.034833   0.050550 20.4714   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Vous pouvez également utiliser cette solution de contournement pour regrouper par une dimension supérieure ou à un niveau supérieur (par exemple, industry ou country). Cependant, dans ce cas, vous ne pourrez pas utiliser le group (ou time) effects, qui est la principale limite de l'approche.


Le package multiwayvcov est une autre approche qui fonctionne à la fois pour le panel et d'autres types de données. Il permet le double clustering, mais aussi le clustering à des dimensions plus élevées. Selon le site Web des packages, il s'agit d'une amélioration par rapport au code d'Arai:

  • Traitement transparent des observations supprimées en raison de leur absence
  • Clustering intégralement multidimensionnel (ou n-way, ou n-dimensionnel ou multidimensionnel)

Utilisation des données Petersen et du cluster.vcov():

library("lmtest")
library("multiwayvcov")

data(petersen)
m1 <- lm(y ~ x, data = petersen)

coeftest(m1, vcov=function(x) cluster.vcov(x, petersen[ , c("firmid", "year")]))
## 
## t test of coefficients:
## 
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 0.029680   0.065066  0.4561   0.6483    
## x           1.034833   0.053561 19.3206   <2e-16 ***
## ---
## Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

La fonction d'Arai peut être utilisée pour regrouper les erreurs standard.Il a une autre version pour le clustering en plusieurs dimensions:

mcl <- function(dat,fm, cluster1, cluster2){
          attach(dat, warn.conflicts = F)
          library(sandwich);library(lmtest)
          cluster12 = paste(cluster1,cluster2, sep="")
          M1  <- length(unique(cluster1))
          M2  <- length(unique(cluster2))   
          M12 <- length(unique(cluster12))
          N   <- length(cluster1)          
          K   <- fm$rank             
          dfc1  <- (M1/(M1-1))*((N-1)/(N-K))  
          dfc2  <- (M2/(M2-1))*((N-1)/(N-K))  
          dfc12 <- (M12/(M12-1))*((N-1)/(N-K))  
          u1j   <- apply(estfun(fm), 2, function(x) tapply(x, cluster1,  sum)) 
          u2j   <- apply(estfun(fm), 2, function(x) tapply(x, cluster2,  sum)) 
          u12j  <- apply(estfun(fm), 2, function(x) tapply(x, cluster12, sum)) 
          vc1   <-  dfc1*sandwich(fm, meat=crossprod(u1j)/N )
          vc2   <-  dfc2*sandwich(fm, meat=crossprod(u2j)/N )
          vc12  <- dfc12*sandwich(fm, meat=crossprod(u12j)/N)
          vcovMCL <- vc1 + vc2 - vc12
          coeftest(fm, vcovMCL)}

Pour des références et des exemples d'utilisation, voir:

Licencié sous: CC-BY-SA avec attribution
Non affilié à StackOverflow
scroll top