我有一個小N大T面板,我估計通過plm(面板線性迴歸模型),具有固定的影響。R中是否有PLM的預測函數?
有什麼辦法可以獲得新數據集的預測值嗎? (我想 估計我樣本子集上的參數,然後使用這些來計算整個樣本的模型隱含值)。
謝謝!
我有一個小N大T面板,我估計通過plm(面板線性迴歸模型),具有固定的影響。R中是否有PLM的預測函數?
有什麼辦法可以獲得新數據集的預測值嗎? (我想 估計我樣本子集上的參數,然後使用這些來計算整個樣本的模型隱含值)。
謝謝!
有(至少)兩種方法中的包,以產生從PLM對象估計:
- fixef.plm:提取固定效應
- pmodel.response:一個函數來提取模型。響應
在我看來,作者(s)沒有興趣提供「隨機效應」的估計。這可能是一個問題,「如果你不知道如何自己做,那麼我們不希望給你一把鋒利的刀把自己切得太深。」
我寫了一個叫predict.out.plm
功能,可以創建爲原數據和一個操縱數據集的預測(以相等的列名)。
predict.out.plm
計算a)轉換數據的預測(擬合)結果和b)構建根據水平結果的結果。該函數適用於使用plm
的一階差分(FD)估算和固定效應(FE)估算。對於FD來說,它隨着時間的推移創造了不同的結果,對於FE來說它創造了時間貶低的結果。
該功能很大程度上未經測試,可能只適用於強平衡的數據幀。
任何建議和更正都非常受歡迎。幫助開發一個小R包將非常感激。
功能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)
}
測試的功能:
##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")
貌似有一個新的軟件包,爲多種車型,包括PLM做樣品中的預測值
https://cran.r-project.org/web/packages/prediction/prediction.pdf
它似乎在引擎蓋下使用'lm',所以你嘗試調用'predict.lm'嗎? – James
我懷疑作者知道釋放一個'predict.plm'函數會鼓勵那些不瞭解統計問題的人盲目地應用它,當這些假設沒有得到滿足時。 IIRC,lme4軟件包也沒有提供預測功能,而且plm的作者指出他們正在評估隨機和固定組件。 –
predict.lm不起作用。我想有一種方法可以提取係數和攔截,但我想其他人已經遇到過這個問題 –