plm: Verwenden von fixef (), um angepasste Werte für ein Zwei-Wege-Modell mit festen Effekten manuell zu berechnen

Bitte beachten Sie: Ich versuche, den Code zum Arbeiten zu bringenbeide Zeit & individuelle Fixeffekte und einunausgeglichen Datensatz. Der folgende Beispielcode funktioniert mit einem ausgeglichenen Datensatz.

Siehe bitte auch unten bearbeiten

Ich versuche, die angepassten Werte eines Modells mit festen Effekten (sowohl mit Einzel - als auch mit Zeiteffekten) manuell mit der zu berechnenplm Paket. Dies ist eher eine Übung, um zu bestätigen, dass ich die Mechanik des Modells und des Pakets versteheplm Objekt, aus den beiden verwandten Fragen (Hier undHier).

Von demplm Vignette (S. 2) ist das zugrunde liegende Modell:

y_it =Alpha + Beta_transponiert *x_it + (mu_i +Lambda_t +Epsilon_es)

Dabei ist mu_i die individuelle Komponente des Fehlerausdrucks (a.k.a. "individueller Effekt") und lambda_t ist der "Zeiteffekt".

Die festen Effekte können mit extrahiert werdenfixef() und ich dachte, ich könnte sie (zusammen mit den unabhängigen Variablen) verwenden, um die angepassten Werte für das Modell zu berechnen, indem ich (mit zwei unabhängigen Variablen) auf diese Weise verwende:

passen_it =Alpha + Beta_1 *x1 + Beta_2 *x2 + mu_i +Lambda_t

Hier scheitere ich - die Werte, die ich erhalte, liegen bei weitem nicht in der Nähe der angepassten Werte (die ich als Differenz zwischen den tatsächlichen Werten und den Residuen im Modellobjekt erhalte). Zum einen sehe ich nichtalpha irgendwo. Ich habe versucht, mit den festen Effekten zu spielen, die als Unterschiede zum ersten, zum Mittelwert usw. angezeigt wurden, ohne Erfolg.

Was mir fehlt Es könnte durchaus ein Missverständnis des Modells oder ein Fehler im Code sein, fürchte ich ... Vielen Dank im Voraus.

PS: Eine der verwandten Fragen deutet darauf hinpmodel.response() sollte mit meinem Problem zusammenhängen (und der Grund dafür ist keinplm.fit Funktion), aber seine Hilfeseite hilft mir nicht zu verstehen, was diese Funktion tatsächlich tut, und ich kann keine Beispiele finden, wie ich das Ergebnis interpretieren kann, das sie erzeugt.

Vielen Dank!

Beispielcode für das, was ich getan habe:

library(data.table); library(plm)

set.seed(100)
DT <- data.table(CJ(id=c("a","b","c","d"), time=c(1:10)))
DT[, x1:=rnorm(40)]
DT[, x2:=rnorm(40)]
DT[, y:=x1 + 2*x2 + rnorm(40)/10]
DT <- DT[!(id=="a" & time==4)] # just to make it an unbalanced panel
setkey(DT, id, time)

summary(plmFEit <- plm(data=DT, id=c("id","time"), formula=y ~ x1 + x2, model="within", effect="twoways"))

# Extract the fitted values from the plm object
FV <- data.table(plmFEit$model, residuals=as.numeric(plmFEit$residuals))
FV[, y := as.numeric(y)]
FV[, x1 := as.numeric(x1)]
FV[, x2 := as.numeric(x2)]

DT <- merge(x=DT, y=FV, by=c("y","x1","x2"), all=TRUE)
DT[, fitted.plm := as.numeric(y) - as.numeric(residuals)]

FEI <- data.table(as.matrix(fixef(object=plmFEit, effect="individual", type="level")), keep.rownames=TRUE) # as.matrix needed to preserve the names?
setnames(FEI, c("id","fei"))
setkey(FEI, id)
setkey(DT, id)
DT <- DT[FEI] # merge the fei into the data, each id gets a single number for every row

FET <- data.table(as.matrix(fixef(object=plmFEit, effect="time", type="level")), keep.rownames=TRUE) # as.matrix needed to preserve the names?
setnames(FET, c("time","fet"))
FET[, time := as.integer(time)] # fixef returns time as character
setkey(FET, time)
setkey(DT, time)
DT <- DT[FET] # merge the fet into the data, each time gets a single number for every row

# calculate the fitted values (called calc to distinguish from those from plm)
DT[, fitted.calc := as.numeric(coef(plmFEit)[1] * x1 + coef(plmFEit)[2]*x2 + fei + fet)]
DT[, diff := as.numeric(fitted.plm - fitted.calc)]

all.equal(DT$fitted.plm, DT$fitted.calc)

Meine Sitzung ist wie folgt:

R version 3.2.2 (2015-08-14)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 8 x64 (build 9200)

locale:
[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United States.1252    LC_MONETARY=English_United States.1252 LC_NUMERIC=C                          
[5] LC_TIME=English_United States.1252    

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] plm_1.4-0           Formula_1.2-1       RJSONIO_1.3-0       jsonlite_0.9.17     readxl_0.1.0.9000   data.table_1.9.7    bit64_0.9-5         bit_1.1-12          RevoUtilsMath_3.2.2

loaded via a namespace (and not attached):
 [1] bdsmatrix_1.3-2  Rcpp_0.12.1      lattice_0.20-33  zoo_1.7-12       MASS_7.3-44      grid_3.2.2       chron_2.3-47     nlme_3.1-122     curl_0.9.3       rstudioapi_0.3.1 sandwich_2.3-4  
[12] tools_3.2.2  

Bearbeiten: (22.02.2015) Da dies auf Interesse gestoßen ist, werde ich versuchen, dies weiter zu klären. Ich habe versucht, ein "Fixed Effects" - Modell (a.k.a. "in" - oder "Least Squares - Dummy - Variablen") wie dasplm paket vignette nennt es auf S.3, oberster Absatz) - gleiche Steigung (en), verschiedene Abschnitte.

Dies ist dasselbe wie das Ausführen einer normalen OLS-Regression nach dem Hinzufügen von Dummys fürtime undid. Mit dem folgenden Code kann ich die angepassten Werte aus dem duplizierenplm Paket mit Basislm(). Bei den Dummies ist es eindeutig, dass die ersten Elemente von id und time die Gruppe sind, mit der verglichen werden soll. Was ich immer noch nicht kann, ist, wie ich die Einrichtungen desplm Paket, um das gleiche zu tun, kann ich leicht mit erreichenlm().

# fit the same with lm() and match the fitted values to those from plm()
lmF <- lm(data = DT, formula = y ~ x1 + x2 + factor(time) + factor(id))
time.lm <- coef(lmF)[grep(x = names(coef(lmF)), pattern = "time", fixed = TRUE)]
time.lm <- c(0, unname(time.lm)) # no need for names, the position index corresponds to time

id.lm <- coef(lmF)[grep(x = names(coef(lmF)), pattern = "id", fixed = TRUE)]
id.lm <- c(0, unname(id.lm))
names(id.lm) <- c("a","b","c","d") # set names so that individual values can be looked up below when generating the fit

DT[, by=list(id, time), fitted.lm := coef(lmF)[["(Intercept)"]]  +  coef(lmF)[["x1"]] * x1  +  coef(lmF)[["x2"]] * x2  +  time.lm[[time]]  +  id.lm[[id]]]
all.equal(DT$fitted.plm, DT$fitted.lm)

Hoffe, dies ist nützlich für andere, die interessiert sein könnten. Das Problem könnte etwas darüber sein, wieplm undfixef mit dem fehlenden Wert umgehen, den ich absichtlich geschaffen habe. Ich habe versucht, mit dem zu spielentype= Parameter vonfixef aber ohne Wirkung.

Antworten auf die Frage(4)

Ihre Antwort auf die Frage