2011-08-19 20 views
10

Mam mały N duży panel T, który szacuję za pomocą plm (model liniowej regresji panelu), ze stałymi efektami.Czy istnieje funkcja przewidywania dla PLM w R?

Czy istnieje sposób na uzyskanie przewidywanych wartości dla nowego zestawu danych? (Chcę oszacować parametry w podzbiorze mojej próbki, a następnie użyć ich do obliczenia implikowanych wartości modelu dla całej próbki).

Dzięki!

+0

Wydaje się, że używasz 'lm' pod maską, więc czy próbowałeś wywoływać' predict.lm'? – James

+2

Podejrzewam, że autorzy wiedzą, że zwolnienie funkcji 'predict.plm' zachęciłoby ludzi, którzy nie rozumieją problemów statystycznych, do ślepego zastosowania jej, gdy założenia nie zostaną spełnione. IIRC, pakiet lme4 nie zapewnia również funkcji przewidywania, a autorzy tabeli twierdzą, że estymują zarówno składniki losowe, jak i stałe. –

+0

predict.lm nie działa. Przypuszczam, że istnieje sposób na wyodrębnienie współczynników i przechwyceń, ale wyobrażam sobie, że inni już napotkali ten problem –

Odpowiedz

7

Jest to (co najmniej) dwa sposoby w opakowaniu w celu wytworzenia oszacowania z PLM obiektów:

- fixef.plm: Wyodrębnianie stałymi

- pmodel.response: Funkcja celu wyekstrahowania model.response

Wydaje mi się, że autor (autorzy) nie są zainteresowani dostarczaniem szacunków dla "efektów losowych". Może to być kwestia "jeśli nie wiesz, jak to zrobić na własną rękę, to nie chcemy dać ci ostrego noża, żeby zbyt głęboko się skaleczyć".

2

napisałem funkcję o nazwie predict.out.plm, który może tworzyć prognozy dla oryginalnych danych z a dla manipulowane zbioru danych (z równych nazwy kolumn).

Parametr predict.out.plm oblicza a) przewidywany (dopasowany) wynik przekształconych danych i b) tworzy wynik zgodny z poziomem. Funkcja działa dla oszacowań Pierwszej Różnicy (FD) i Stałych Efektów (FE) przy użyciu plm. W przypadku FD tworzy on różny wynik w czasie, a dla FE tworzy wynik czasowy.

Funkcja jest w dużej mierze nietestowana i prawdopodobnie działa tylko z mocno zbalansowanymi ramkami danych.

Wszelkie sugestie i poprawki są bardzo mile widziane. Pomoc w opracowaniu małego pakietu R będzie bardzo cenna.

Funkcja predict.out.plm

predict.out.plm<-function(
    estimate, 
    formula, 
    data, 
    model="fd", 
    pname="y", 
    pindex=NULL, 
    levelconstr=T 
){ 
    # estimate=e.fe 
    # formula=f 
    # data=d 
    # model="within" 
    # pname="y" 
    # pindex=NULL 
    # levelconstr=T 
    #get index of panel data 
    if (is.null(pindex) && class(data)[1]=="pdata.frame") { 
    pindex<-names(attributes(data)$index) 
    } else { 
    pindex<-names(data)[1:2] 
    } 
    if (class(data)[1]!="pdata.frame") { 
    data<-pdata.frame(data) 
    } 
    #model frame 
    mf<-model.frame(formula,data=data) 
    #model matrix - transformed data 
    mn<-model.matrix(formula,mf,model) 

    #define variable names 
    y.t.hat<-paste0(pname,".t.hat") 
    y.l.hat<-paste0(pname,".l.hat") 
    y.l<-names(mf)[1] 

    #transformed data of explanatory variables 
    #exclude variables that were droped in estimation 
    n<-names(estimate$aliased[estimate$aliased==F]) 
    i<-match(n,colnames(mn)) 
    X<-mn[,i] 

    #predict transformed outcome with X * beta 
    # p<- X %*% coef(estimate) 
    p<-crossprod(t(X),coef(estimate)) 
    colnames(p)<-y.t.hat 

    if (levelconstr==T){ 
    #old dataset with original outcome 
    od<-data.frame(
     attributes(mf)$index, 
     data.frame(mf)[,1] 
    ) 
    rownames(od)<-rownames(mf) #preserve row names from model.frame 
    names(od)[3]<-y.l 

    #merge old dataset with prediciton 
    nd<-merge(
     od, 
     p, 
     by="row.names", 
     all.x=T, 
     sort=F 
    ) 
    nd$Row.names<-as.integer(nd$Row.names) 
    nd<-nd[order(nd$Row.names),] 

    #construct predicted level outcome for FD estiamtions 
    if (model=="fd"){ 
     #first observation from real data 
     i<-which(is.na(nd[,y.t.hat])) 
     nd[i,y.l.hat]<-NA 
     nd[i,y.l.hat]<-nd[i,y.l] 
     #fill values over all years 
     ylist<-unique(nd[,pindex[2]])[-1] 
     ylist<-as.integer(as.character(ylist)) 
     for (y in ylist){ 
     nd[nd[,pindex[2]]==y,y.l.hat]<- 
      nd[nd[,pindex[2]]==(y-1),y.l.hat] + 
      nd[nd[,pindex[2]]==y,y.t.hat] 
     } 
    } 
    if (model=="within"){ 
     #group means of outcome 
     gm<-aggregate(nd[, pname], list(nd[,pindex[1]]), mean) 
     gl<-aggregate(nd[, pname], list(nd[,pindex[1]]), length) 
     nd<-cbind(nd,groupmeans=rep(gm$x,gl$x)) 
     #predicted values + group means 
     nd[,y.l.hat]<-nd[,y.t.hat] + nd[,"groupmeans"] 
    } 
    if (model!="fd" && model!="within") { 
     stop('funciton works only for FD and FE estimations') 
    } 
    } 
    #results 
    results<-p 
    if (levelconstr==T){ 
    results<-list(results,nd) 
    names(results)<-c("p","df") 
    } 
    return(results) 
} 

Testowanie funkcji:

##packages 
library(plm) 

##test dataframe 
#data structure 
N<-4 
G<-2 
M<-5 
d<-data.frame(
    id=rep(1:N,each=M), 
    year=rep(1:M,N)+2000, 
    gid=rep(1:G,each=M*2) 
) 
#explanatory variable 
d[,"x"]=runif(N*M,0,1) 
#outcome 
d[,"y"] = 2 * d[,"x"] + runif(N*M,0,1) 
#panel data frame 
d<-pdata.frame(d,index=c("id","year")) 

##new data frame for out of sample prediction 
dn<-d 
dn$x<-rnorm(nrow(dn),0,2) 

##estimate 
#formula 
f<- pFormula(y ~ x + factor(year)) 
#fixed effects or first difffernce estimation 
e<-plm(f,data=d,model="within",index=c("id","year")) 
e<-plm(f,data=d,model="fd",index=c("id","year")) 
summary(e) 

##fitted values of estimation 
#transformed outcome prediction 
predict(e) 
c(pmodel.response(e)-residuals(e)) 
predict.out.plm(e,f,d,"fd")$p 
# "level" outcome prediciton 
predict.out.plm(e,f,d,"fd")$df$y.l.hat 
#both 
predict.out.plm(e,f,d,"fd") 

##out of sampel prediciton 
predict(e,newdata=d) 
predict(e,newdata=dn) 
# Error in crossprod(beta, t(X)) : non-conformable arguments 
# if plm omits variables specified in the formula (e.g. one year in factor(year)) 
# it tries to multiply two matrices with different length of columns than regressors 
# the new funciton avoids this and therefore is able to do out of sample predicitons 
predict.out.plm(e,f,dn,"fd")