2012-05-21 18 views
6

मेरी समस्या यह है: मुझे NA मिलता है जहां मुझे मजबूत मानक त्रुटियों की गणना में कुछ मूल्य प्राप्त करना चाहिए।पैनल डेटा रिग्रेशन: मजबूत मानक त्रुटियां

मैं क्लस्टर-मजबूत मानक त्रुटियों के साथ एक निश्चित प्रभाव पैनल प्रतिगमन करने की कोशिश कर रहा हूं। इसके लिए, मैं Arai (2011) का पालन करता हूं जो पी पर है। 3 Stock/ Watson (2006) का अनुसरण करता है (बाद में Econometrica में प्रकाशित, जिनके पास पहुंच है)। मैं नीचे की पूर्वाग्रह के खिलाफ (M/(M-1)*(N-1)/(N-K) द्वारा स्वतंत्रता की डिग्री को सही करना चाहता हूं क्योंकि क्लस्टर की संख्या सीमित है और मेरे पास असंतुलित डेटा है।

इसी तरह की समस्याएं स्टैक ओवरव्लो और [3] पर क्रॉसविलिडेटेड पर [1, 2] से पहले पोस्ट की गई हैं।

अराइ (और 1 कड़ी में जवाब) कार्यों के लिए निम्न कोड (मैं कुछ आगे टिप्पणी साथ नीचे मेरी डेटा उपलब्ध कराने) का उपयोग करता है:

gcenter <- function(df1,group) { 
    variables <- paste(
     rep("C", ncol(df1)), colnames(df1), sep=".") 
    copydf <- df1 
    for (i in 1:ncol(df1)) { 
     copydf[,i] <- df1[,i] - ave(df1[,i], group,FUN=mean)} 
    colnames(copydf) <- variables 
    return(cbind(df1,copydf))} 

# 1-way adjusting for clusters 
clx <- function(fm, dfcw, cluster){ 
    # R-codes (www.r-project.org) for computing 
    # clustered-standard errors. Mahmood Arai, Jan 26, 2008. 
    # The arguments of the function are: 
    # fitted model, cluster1 and cluster2 
    # You need to install libraries `sandwich' and `lmtest' 
    # reweighting the var-cov matrix for the within model 
    library(sandwich);library(lmtest) 
    M <- length(unique(cluster)) 
    N <- length(cluster)   
    K <- fm$rank       
    dfc <- (M/(M-1))*((N-1)/(N-K)) 
    uj <- apply(estfun(fm),2, function(x) tapply(x, cluster, sum)); 
    vcovCL <- dfc*sandwich(fm, meat=crossprod(uj)/N)*dfcw 
    coeftest(fm, vcovCL) } 

, जहां gcenter मतलब से विचलन की गणना करता है (निश्चित प्रभाव)। मैं फिर अपने क्लस्टर वैरिएबल (मैंने अपना डेटा 'डेटा' नाम दिया है) DS_CODE के साथ प्रतिक्रिया जारी रखी है।

centerdata <- gcenter(data, data$DS_CODE) 
datalm <- lm(C.L1.retE1M ~ C.MCAP_SEC + C.Impact_change + C.Mom + C.BM + C.PD + C.CashGen + C.NITA + C.PE + C.PEdummy + factor(DS_CODE), data=centerdata) 
M <- length(unique(data$DS_CODE)) 
dfcw <- datalm$df/(datalm$df - (M-1)) 

और मैं अपने regressors के लिए कुछ मान केवल शुरुआत में मिलता है,

clx(datalm, dfcw, data$DS_CODE) 

हालांकि गणना करना चाहते हैं जब मैं uj गणना करने के लिए (उपरोक्त सूत्र clx देखें) विचरण के लिए चाहते हैं, तो, फिर शून्य के बहुत सारे। यदि यह इनपुट uj भिन्नता के लिए उपयोग किया जाता है, केवल NAs परिणाम।

मेरे डेटा

अपने डेटा विशेष संरचना का हो सकता है और मैं इस समस्या को समझ नहीं सकता के बाद से, मैं एक link हॉटमेल से के रूप में पूरे बात पोस्ट। इसका कारण यह है कि अन्य डेटा (अरि (2011) से लिया गया) के साथ मेरी समस्या नहीं होती है। गड़बड़ी के लिए अग्रिम में खेद है, लेकिन अगर आप इसे फिर से देख सकते हैं तो मैं बहुत आभारी हूं। फ़ाइल एक 5 एमबी .txt फ़ाइल है जिसमें पूरी तरह से डेटा है।

+0

अराइ के कागज अब आपके लिंक के नीचे मौजूद हैं। क्या आप वास्तविक लिंक प्रदान कर सकते हैं? – MERose

उत्तर

2

कुछ समय के आसपास खेलने के बाद, यह मेरे लिए काम करता है और मुझे देता है:

      Estimate Std. Error t value Pr(>|t|)  
(Intercept)   4.5099e-16 5.2381e-16 0.8610 0.389254  
C.MCAP_SEC   -5.9769e-07 1.2677e-07 -4.7149 2.425e-06 *** 
C.Impact_change  -5.3908e-04 7.5601e-05 -7.1306 1.014e-12 *** 
C.Mom     3.7560e-04 3.3378e-03 0.1125 0.910406  
C.BM     -1.6438e-04 1.7368e-05 -9.4645 < 2.2e-16 *** 
C.PD     6.2153e-02 3.8766e-02 1.6033 0.108885  
C.CashGen    -2.7876e-04 1.4031e-02 -0.0199 0.984149  
C.NITA    -8.1792e-02 3.2153e-02 -2.5438 0.010969 * 
C.PE     -6.6170e-06 4.0138e-06 -1.6485 0.099248 . 
C.PEdummy    1.3143e-02 4.8864e-03 2.6897 0.007154 ** 
factor(DS_CODE)130324 -5.2497e-16 5.2683e-16 -0.9965 0.319028  
factor(DS_CODE)130409 -4.0276e-16 5.2384e-16 -0.7689 0.441986  
factor(DS_CODE)130775 -4.4113e-16 5.2424e-16 -0.8415 0.400089 
... 

यह सवाल क्यों यह आप के लिए नहीं है के साथ हमें छोड़ देता है। मुझे लगता है कि यह आपके डेटा के प्रारूप के साथ कुछ करने के लिए है। क्या सबकुछ संख्यात्मक है? मैं स्तंभ वर्गों परिवर्तित और यह मेरे लिए इस तरह दिखता है:

str(dat) 
'data.frame': 48251 obs. of 12 variables: 
$ DS_CODE  : chr "902172" "902172" "902172" "902172" ... 
$ DNEW   : num 2e+05 2e+05 2e+05 2e+05 2e+05 ... 
$ MCAP_SEC  : num 78122 71421 81907 80010 82462 ... 
$ NITA   : num 0.135 0.135 0.135 0.135 0.135 ... 
$ CashGen  : num 0.198 0.198 0.198 0.198 0.198 ... 
$ BM   : num 0.1074 0.1108 0.097 0.0968 0.0899 ... 
$ PE   : num 57 55.3 63.1 63.2 68 ... 
$ PEdummy  : num 0 0 0 0 0 0 0 0 0 0 ... 
$ L1.retE1M : num -0.72492 0.13177 0.00122 0.07214 -0.07332 ... 
$ Mom   : num 0 0 0 0 0 ... 
$ PD   : num 5.41e-54 1.51e-66 3.16e-80 2.87e-79 4.39e-89 ... 
$ Impact_change: num 0 -10.59 -10.43 0.7 -6.97 ... 

क्या आप के लिए str(data) वापसी करता है?

+0

आपके प्रयास और आपके उत्तर के लिए बहुत बहुत धन्यवाद! 'DS_CODE' और * int * के लिए' DNEW' के लिए मेरा 'str (डेटा) 'रिटर्न * फैक्टर *। अन्य सभी परिणाम एक ही हैं .... लेकिन: यह अजीब चीज है: अगर मैं * कम * डेटा सेट का उपयोग करता हूं तो यह अब काम करता है (मैंने आपको केवल मेरे अन्य वैरिएबल और आर पंक्ति संख्याओं के बिना छोटा डेटा सेट दिया है)। बड़े सेट के साथ, मुझे * यूजे * की गणना में 'एनएएस' की 1 सिंगल पंक्ति मिलती है। यदि मैं पंक्ति संख्याओं के बिना अपना पूरा डेटा सेट निर्यात करता हूं ('row.names = FALSE'), इसे फिर से आयात करें और प्रतिगमन करें, यह बड़े डेटा सेट के साथ काम करता है। मुझे नहीं पता क्यों ... – Jan

+0

खुशी है कि यह अब काम करता है। –

0

plm पैकेज पैनल प्रतिगमन के लिए क्लस्टर एसई का अनुमान लगा सकता है। मूल डेटा अब उपलब्ध नहीं है, इसलिए यहां डमी डेटा का उपयोग करके एक उदाहरण दिया गया है।

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

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

##Arellano clustered by *group* SEs 
> coeftest(fpm, vcov=function(x) vcovHC(x, cluster="group", type="HC0")) 

t test of coefficients: 

      Estimate Std. Error t value Pr(>|t|)  
(Intercept) 0.029680 0.066939 0.4434 0.6575  
x   1.034833 0.050540 20.4755 <2e-16 *** 
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

आप lm मॉडल (plm के बजाय) का उपयोग कर रहे हैं, तो multiwayvcov पैकेज मदद मिल सकती है।

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

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

> coeftest(m1, vcov=function(x) cluster.vcov(x, petersen[ , c("firmid")], 
    df_correction=FALSE)) 

t test of coefficients: 

      Estimate Std. Error t value Pr(>|t|)  
(Intercept) 0.029680 0.066939 0.4434 0.6575  
x   1.034833 0.050540 20.4755 <2e-16 *** 
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

अधिक जानकारी के लिए देखें:

यह भी देखें:

संबंधित मुद्दे