0

我正在学习在 R 中实现健壮的 glms,但是当我有一个由于共线性而删除了某些列的模型时,我无法弄清楚为什么我无法让 glmrob 从我的回归模型中预测值。特别是当我使用 predict 函数从 glmrob 预测值时,它总是为所有值提供 NA 。使用 glm 从相同的数据和模型预测值时,我没有观察到这一点。我使用什么数据似乎并不重要——只要拟合模型中有一个 NA 系数(并且 NA 不是系数向量中的最后一个系数),预测就不起作用。

这种行为适用于我尝试过的所有数据集和模型,其中由于共线性而删除了内部列。我包含了一个假数据集,其中从模型中删除了两列,这在系数列表中给出了两个 NA。glm 和 glmrob 都给出几乎相同的系数,但 predict 仅适用于 glm 模型。所以我的问题是:对于会阻止我的 glmrob 模型生成预测值的稳健回归,我不了解什么?

library(robustbase)

#Make fake data with two categorial predictors
df <- data.frame("category" = rep(c("A","B","C"),each=6))
df$location <- rep(1:6,each=3)
val <- rep(c(500,50,5000),each=6)+rep(c(50,100,25,200,100,1),each=3)
df$value <- rpois(NROW(df),val)

#note that predict works if we omit the newdata parameter. However I need the newdata param
#so I use the original dataframe here as a stand-in.  
mod <- glm(val ~ category + as.factor(location), data=df, family=poisson)
predict(mod, newdata=df) # works fine

mod <- glmrob(val ~ category + as.factor(location), data=df, family=poisson)
predict(mod, newdata=df) #predicts NA for all values
4

1 回答 1

1

我一直在研究这个并得出结论,问题不在于我对鲁棒回归的理解,而在于鲁棒基础包中的错误。predict.lmrob 函数在预测之前没有从模型中正确选择必要的系数。它需要选择前 x 个非 NA 系数(其中 x=模型矩阵的秩)。相反,它只选择前 x 个系数而不检查它们是否为 NA。这就解释了为什么这个问题只出现在 NA 不是系数向量中的最后一个系数的模型中。

为了解决这个问题,我使用以下方法复制了 predict.lmrob 源:

getAnywhere(predict.lmrob)

并创建了我自己的替换功能。在这个函数中,我对代码进行了一次修改:

...
p <- object$rank
if (is.null(p)) {
    df <- Inf
    p <- sum(!is.na(coef(object)))
    #piv <- seq_len(p) # old code
    piv <- which(!is.na(coef(object))) # new code
}
else {
    p1 <- seq_len(p)
    piv <- if (p) 
        qr(object)$pivot[p1]
}
...

我已经使用此更改运行了数百个数据集,并且效果很好。

于 2017-05-24T00:40:54.503 回答