尝试从模型矩阵中提取信息时,lmer出现不一致的自变量错误 [英] non-conformable arguments error from lmer when trying to extract information from the model matrix

查看:1260
本文介绍了尝试从模型矩阵中提取信息时,lmer出现不一致的自变量错误的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我有一些纵向数据,希望可以在指定的时间获得预测的均值.该模型包括2个项,它们的相互作用和时间变量的样条项.当我尝试获取预测的均值时,出现"mm%*%fixef(m4)的错误:参数不一致"

I have some longitudinal data from which I'd like to get the predicted means at specified times. The model includes 2 terms, their interaction and a spline term for the time variable. When I try to obtain the predicted means, I get "Error in mm %*% fixef(m4) : non-conformable arguments"

我使用了lmer的sleep数据集来说明我的问题.首先,我导入数据并为我的互动创建一个变量年龄"

I've used the sleep data set from lmer to illustrate my problem. First, I import the data and create a variable "age" for my interaction

sleep <- as.data.frame(sleepstudy)  #get the sleep data
# create fake variable for age with 3 levels
set.seed(1234567)
sleep$age <- as.factor(sample(1:3,length(sleep),rep=TRUE))

然后我运行我的lmer模型

Then I run my lmer model

library(lme4)
library(splines)
m4 <- lmer(Reaction ~ Days + ns(Days, df=4) + age + Days:age + (Days | Subject), sleep) 

最后,我创建获取预测均值所需的数据和矩阵

Finally, I create the data and matrix needed to obtain predicted means

#new data frame for predicted means
d <- c(0:9)  # make a vector of days = 0 to 9 to obtain predictions for each day
newdat <- as.data.frame(cbind(Days=d, age=rep(c(1:3),length(d))))
newdat$Days <- as.numeric(as.character(newdat$Days))
newdat$age <- as.factor(newdat$age)

# create a matrix 
mm<-model.matrix(~Days + ns(Days, df=4) + age + Days:age, newdat)  
newdat$pred<-mm%*%fixef(m4) 

这是我得到的错误: 毫米%*%fixef(m4)的错误:参数不一致

It's at this point that I get the error: Error in mm %*% fixef(m4) : non-conformable arguments

我可以使用预测来获取平均值

I can use predict to get the means

newdat$pred <- predict(m4, newdata=newdat, re.form=NA)

这很好用,但是我希望能够计算出一个置信区间,所以我需要一个合适的矩阵.

which works fine, but I want to be able to calculate a confidence interval, so I need a conformable matrix.

我在某处读到问题可能是lmer创建了别名(我找不到该帖子).此评论是关于不能将effect()用于类似任务的.我不太明白如何克服这个问题.而且,我记得帖子有点老了,希望别名问题可能不再相关.

I read somewhere that the problem may be that lmer creates aliases (I can't find that post). This comment was made with regards to not being able to use effect() for a similar task. I couldn't quite understand how to overcome this problem. Moreover, I recall that post was a little old and hoped the alias problem may no longer be relevant.

如果有人对我可能做错的事情有建议,我将很感激您的反馈.谢谢.

If anyone has a suggestion for what I may be doing wrong, I'd appreciate the feedback. Thanks.

推荐答案

这里有几件事.

  • 您需要删除列以使模型矩阵与实际拟合的固定效果矢量相称(即在删除共线列之后与实际用于拟合的模型矩阵相称)
  • 为引起进一步的困惑,您碰巧只采样了2和3岁(可能的{1,2,3}中)

我已经整理了一下代码...

I've cleaned up the code a little bit ...

library("lme4")
library("splines")
sleep <- sleepstudy  #get the sleep data
set.seed(1234567)
## next line happens to sample only 2 and 3 ...
sleep$age <- as.factor(sample(1:3,length(sleep),rep=TRUE))
length(levels(sleep$age))  ## 2

适合的型号:

m4 <- lmer(Reaction ~ Days + ns(Days, df=4) +
    age + Days:age + (Days | Subject), sleep)
## message; fixed-effect model matrix is 
##    rank deficient so dropping 1 column / coefficient

检查固定效果:

f1 <- fixef(m4)
length(f1)  ## 7
f2 <- fixef(m4,add.dropped=TRUE)
length(f2)  ## 8

我们可以使用固定效果的扩展版本(其中具有NA值),但是这只会通过在计算中传播NA值来弄乱我们.

We could use this extended version of the fixed effects (which has an NA value in it), but this would just mess us up by propagating NA values through the computation ...

检查模型矩阵:

X <- getME(m4,"X")
ncol(X)  ## 7
(which.dropped <- attr(getME(m4,"X"),"col.dropped"))
## ns(Days, df = 4)4 
##             6

用于预测均值的新数据框

New data frame for predicted means

d <- 0:9  
## best to use data.frame() directly, avoid cbind()
##   generate age based on *actual* levels in data
newdat <- data.frame(Days=d,
   age=factor(rep(levels(sleep$age),length(d))))

创建矩阵:

mm <- model.matrix(formula(m4,fixed.only=TRUE)[-2], newdat)
mm <- mm[,-which.dropped]   ## drop redundant columns
## newdat$pred <- mm%*%fixef(m4)    ## works now

由sianagh添加:用于获取置信区间并绘制数据的代码:

Added by sianagh: Code to obtain confidence intervals and plot the data:

predFun <- function(x) predict(x,newdata=newdat,re.form=NA)
newdat$pred <- predFun(m4)
bb <- bootMer(m4,
   FUN=predFun,
    nsim=200)  
## nb. this produces an error message on its first run, 
## but not on subsequent runs (using the development version of lme4)
bb_ci <- as.data.frame(t(apply(bb$t,2,quantile,c(0.025,0.975))))
names(bb_ci) <- c("lwr","upr")
newdat <- cbind(newdat,bb_ci)

情节:

plot(Reaction~Days,sleep)
with(newdat,
    matlines(Days,cbind(pred,lwr,upr),
            col=c("red","green","green"),
            lty=2,
            lwd=c(3,2,2)))

这篇关于尝试从模型矩阵中提取信息时,lmer出现不一致的自变量错误的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆