Как предсказать и изобразить нелинейные переменные уклоны в lmer или glmer?

Question

Моя цель состоит в том, чтобы вычислить предсказанные значения из многоуровневой модели с переменным перехватом и переменным наклоном, используя lmer а также glmer функции lme4 пакет в R. Чтобы сделать это конкретным и понятным, я привожу здесь игрушечный пример с набором данных "mtcars":

Вот как я обычно создаю прогнозируемые значения из многоуровневой модели с переменным перехватом и переменным наклоном (этот код должен работать очень хорошо):

# loading in-built cars dataset
data(mtcars)

# the "gear" column will be the group-level factor, so we'll have cars nested 
# within "gear" type
mtcars$gear <- as.factor(mtcars$gear)

# fitting varying-slope, varying-intercept model
m <- lmer(mpg ~ 1 + wt + hp + (1 + wt|gear), data=mtcars)

# creating the prediction frame
newdata <- with(mtcars, expand.grid(wt=unique(wt),
                              gear=unique(gear),
                              hp=mean(hp)))

# calculating predictions
newdata$pred <- predict(m, newdata, re.form=~(1 + wt|gear))

# quick ggplot2 graph
p <- ggplot(newdata, aes(x=wt, y=pred, colour=gear))
p + geom_line() + ggtitle("Varying Slopes")

Вышеприведенный код R должен работать, но если я хочу создать и отобразить прогнозы на основе нелинейного переменного перехвата, переменного наклона, то это явно не сработает. Для простоты и воспроизводимости вот камень преткновения, использующий набор данных "mtcars":

# key question: how to create predictions if I want to examine a non-linear 
# varying slope?

# creating a squared term for a non-linear relationship
# NB: usually I use the `poly` function
mtcars$wtsq <- (mtcars$wt)^2

# fitting varying-slope, varying-intercept model with a non-linear trend
m <- lmer(mpg ~ 1 + wt + wtsq + hp + (1 + wt + wtsq|gear), data=mtcars)

# creating the prediction frame
newdata <- with(mtcars, expand.grid(wt=unique(wt),
                                wtsq=unique(wtsq),
                                gear=unique(gear),
                                hp=mean(hp)))

# calculating predictions
newdata$pred <- predict(m, newdata, re.form=~(1 + wt + wtsq|gear))

# quick ggplot2 graph 
# clearly not correct (see the graph below)
p <- ggplot(newdata, aes(x=wt, y=pred, colour=gear))
p + geom_line() + ggtitle("Varying Slopes")

Очевидно, что прогнозирующий фрейм настроен неправильно. Любые идеи о том, как создавать и отображать прогнозируемые значения при подгонке нелинейной многоуровневой модели с переменным наклоном и переменной наклоном в R? Спасибо!

Answer 1

Проблема в том, что когда вы используете expand.grid с обоими wt а также wt^2Вы создаете все возможные комбинации wt а также wt^2, Эта модификация вашего кода работает:

newdata <- with(mtcars, expand.grid(wt=unique(wt),
                                gear=unique(gear),
                                hp=mean(hp)))
newdata$wtsq <- newdata$wt^2

newdata$pred <- predict(m, newdata)

p <- ggplot(newdata, aes(x=wt, y=pred, colour=gear, group=gear))
p + geom_line() + ggtitle("Varying Slopes")