r - heckman 模型的预测函数

Question

我使用sampleSelection包中的示例

## Greene( 2003 ): example 22.8, page 786
    data( Mroz87 )
    Mroz87$kids  <- ( Mroz87$kids5 + Mroz87$kids618 > 0 )
    # Two-step estimation
    test1  = heckit( lfp ~ age + I( age^2 ) + faminc + kids + educ,
                     wage ~ exper + I( exper^2 ) + educ + city, Mroz87 ) 
    # ML estimation
    test2 =  selection( lfp ~ age + I( age^2 ) + faminc + kids + educ,
                        wage ~ exper + I( exper^2 ) + educ + city, Mroz87 ) 
    pr2 <- predict(test2,Mroz87)
    pr1 <- predict(test1,Mroz87)

我的问题是预测功能不起作用。我收到此错误:

    Error in UseMethod("predict") : 
  no applicable method for 'predict' applied to an object of class "c('selection', 'maxLik', 'maxim', 'list')"

预测函数适用于许多模型，所以我想知道为什么赫克曼回归模型会出现错误。

------------更新------------我取得了一些进展，但我仍然需要你的帮助。我建立了一个原始的赫克曼模型进行比较:

data( Mroz87 )
Mroz87$kids  <- ( Mroz87$kids5 + Mroz87$kids618 > 0 )
test1  = heckit( lfp ~ age + I( age^2 ) + faminc + kids + educ,
                 wage ~ exper + I( exper^2 ) + educ + city, Mroz87[1:600,] ) 

之后我开始自己构建它。 Heckman 模型需要一个选择方程:

zi* = wi γ + ui
where zi =1 if zi* >0  and zi = 0 if zi* <=0
after you calculate yi = xi*beta +ei ONLY for the cases where zi*>0

我首先构建概率模型:

library(MASS)
#probit1 = probit(lfp ~ age + I( age^2 ) + faminc + kids + educ, Mroz87, x = TRUE, print.level = print.level -      1, iterlim = 30)
myprobit <- glm(lfp ~ age + I( age^2 ) + faminc + kids + educ, family = binomial(link = "probit"), 
                data = Mroz87[1:600,])
summary(myprobit)

该模型与 heckit 命令完全相同。

然后我构建了一个 lm 模型:

#get predictions for the variables (the data is not needed but I specify it anyway)
selectvar <- predict(myprobit,data = Mroz87[1:600,])
#bind the prediction to the table (I build a new one in my case)
newdata = cbind(Mroz87[1:600,],selectvar)
#Build an lm model for the subset where zi>0
lm1 = lm(wage ~ exper + I( exper^2 ) + educ + city , newdata, subset = selectvar > 0)
summary(lm1)

我现在的问题是 lm 模型与 heckit 创建的模型不太一样。我不知道为什么。有什么想法吗？

Answer 1

实现

这是 predict.selection 函数的实现 - 它产生 4 种不同类型的预测(解释 here ):

library(Formula)
library(sampleSelection)
predict.selection = function(objSelection, dfPred, 
                                 type = c('link', 'prob', 'cond', 'uncond')) {

  # construct the Formula object
  tempS = evalq(objSelection$call$selection)
  tempO = evalq(objSelection$call$outcome)

  FormHeck = as.Formula(paste0(tempO[2], '|', tempS[2], '~', tempO[3], '|', tempS[3]))

  # regressor matrix for the selection equation
  mXSelection = model.matrix(FormHeck, data = dfPred, rhs = 2)

  # regressor matrix for the outcome equation
  mXOutcome = model.matrix(FormHeck, data = dfPred, rhs = 1)

  # indices of the various parameters in selectionObject$estimate
  vIndexBetaS = objSelection$param$index$betaS
  vIndexBetaO = objSelection$param$index$betaO
  vIndexErr = objSelection$param$index$errTerms

  # get the estimates
  vBetaS = objSelection$estimate[vIndexBetaS]
  vBetaO = objSelection$estimate[vIndexBetaO]

  dLambda = objSelection$estimate[vIndexErr['rho']]*
    objSelection$estimate[vIndexErr['sigma']]

  # depending on the type of prediction requested, return
  # TODO allow the return of multiple prediction types
  pred = switch(type, 
         link = mXSelection %*% vBetaS,
         prob = pnorm(mXSelection %*% vBetaS),
         uncond = mXOutcome %*% vBetaO,
         cond = mXOutcome %*% vBetaO + 
           dnorm(temp <- mXSelection %*% vBetaS)/pnorm(temp) * dLambda)
  return(pred)
}

测试

假设您使用 MLE 估计以下 Heckman 样本选择模型:

data(Mroz87)

# define a new variable
Mroz87$kids  = (Mroz87$kids5 + Mroz87$kids618 > 0)

# create the estimation sample
Mroz87Est = Mroz87[1:600, ]

# create the hold out sample
Mroz87Holdout = Mroz87[601:nrow(Mroz87), ]

# estimate the model using MLE
heckML =  selection(selection = lfp ~ age + I(age^2) + faminc + kids + educ,
                    outcome = wage ~ exper + I(exper^2) + educ + city, data = Mroz87Est) 
summary(heckML)  

不同类型的预测计算如下:

vProb = predict(objSelection = heckML, dfPred = Mroz87Holdout, type = 'prob')
vLink = predict(objSelection = heckML, dfPred = Mroz87Holdout, type = 'link')
vCond = predict(objSelection = heckML, dfPred = Mroz87Holdout, type = 'cond')
vUncond = predict(objSelection = heckML, dfPred = Mroz87Holdout, type = 'uncond')

您可以在生成这些输出的平台上验证这些计算，例如 Stata .