Broom.mixed exp 模型预测
Broom.mixed exp model predictions
我想寻求一些帮助来绘制我的模型的预测值以及由 lmer() 的估计生成的方程。
所以,我拥有的数据是不同老鼠在不同日子里的质量体积。每只老鼠测量该体积的时间点不同。
那么,我使用的模型是:
m1 <- lmer(lVolume ~ Country*Day + (1|Rat))
我这样做是因为我对 exp(fitted)
值感兴趣,然后获得该模型的指数方法而不是使用非线性混合效应模型(目前)
为了绘制此模型的预测,我的尝试是:
m1%>%
augment() %>%
clean_names() %>%
ggplot(data = .,
mapping = aes(x = day,
y = exp(l_volume),
group = rat)) +
geom_point(alpha = 0.5) +
geom_line(alpha = 0.5) +
geom_point(aes(y = exp(fitted)),
color = "red") +
geom_line(aes(y = exp(fitted)),
color = "red") +
expand_limits(x = 0 , y = 0)
这里我绘制了更多的老鼠,但是,如您所见,(0,0) 与 lmer 的预测相差太远。我想知道如何绘制我的模型生成的预测以查看 (0,200) 的点。我从这里尝试通过创建一个新的数据框来提示,然后再次使用预测(m1,newdata = new_df)进行绘图,但是我不知道如何创建这个数据框,因为我有 20 只老鼠而且我不知道如何将其扩展到 predict()。
我的尝试:
pframe <- data.frame(Day=seq(0, 200, length.out=101))
pframe$continuous_outcome <- predict(m1, newdata = pframe, level = 0)
ggplot(data, aes(Day,lVolume)) +
geom_point() +
geom_line(data=pframe)
但是我得到一个错误:
Error in eval(predvars, data, env) : object 'Rat' not found
而且,还有一种方法可以绘制你从每次估计中生成的方程式,即,从每只老鼠你有一组固定和随机的估计量,我如何绘制方程式(红色曲线) lmer 正在为每只老鼠生成?
原来 predict
比 broom.mixed::augment
更容易使用。
构建预测
(Rat/Country/Days 0-150 的所有组合(一天达到 200 导致一些极端预测超出了垂直比例)
library(tidyverse)
dc <- distinct(dplyr::select(dat1, Rat, Country))
pframe <- (with(dat1,
expand_grid(Rat = unique(Rat),
Day = 0:150))
%>% full_join(dc, by = "Rat")
%>% mutate(lVolume = predict(m1, newdata = .))
)
将数据和预测合并到一个数据框中(您不必这样做,但它使图例变得容易)
comb <- dplyr::bind_rows(list(data = dat1, model = pframe),
.id = "type")
剧情:
ggplot(comb, aes(Day, exp(lVolume), colour = type)) +
geom_point(alpha = 0.2) +
geom_line(aes(group = interaction(type, Rat))) +
scale_colour_manual(values = c("black", "red"))
重建数据:
dat0 <- list(
list("rat1", vol=c(78,304,352,690,952,1250), days = c(89,110,117,124,131,138), country = "Chile"),
list("rat2", vol=c(202,440,520,870,1380), days = c(75,89,96,103,110), country = "Chile"),
list("rat3", vol=c(186,370,620,850,1150), days = c(75,89,96,103,110), country = "Chile"),
list("rat4", vol=c(92,250,430,450,510,850,1000,1200), days = c(47,61,75,82,89,97,103,110), country = "England"),
list("rat5", vol=c(110,510,710,1200), days = c(47,61,75,82), country = "England"),
list("rat6", vol=c(115,380,480,540,560,850,1150,1350), days = c(47,61,75,82,89,97,103,110), country = "England"))
dat1 <- purrr::map_dfr(dat0,
~ data.frame(Rat = .[[1]],
lVolume = log(.$vol), Day = .$days,
Country = .$country))
m1 <- lmer(lVolume ~ Country*Day + (1|Rat), data = dat1)
我想寻求一些帮助来绘制我的模型的预测值以及由 lmer() 的估计生成的方程。
所以,我拥有的数据是不同老鼠在不同日子里的质量体积。每只老鼠测量该体积的时间点不同。
那么,我使用的模型是:
m1 <- lmer(lVolume ~ Country*Day + (1|Rat))
我这样做是因为我对 exp(fitted)
值感兴趣,然后获得该模型的指数方法而不是使用非线性混合效应模型(目前)
为了绘制此模型的预测,我的尝试是:
m1%>%
augment() %>%
clean_names() %>%
ggplot(data = .,
mapping = aes(x = day,
y = exp(l_volume),
group = rat)) +
geom_point(alpha = 0.5) +
geom_line(alpha = 0.5) +
geom_point(aes(y = exp(fitted)),
color = "red") +
geom_line(aes(y = exp(fitted)),
color = "red") +
expand_limits(x = 0 , y = 0)
这里我绘制了更多的老鼠,但是,如您所见,(0,0) 与 lmer 的预测相差太远。我想知道如何绘制我的模型生成的预测以查看 (0,200) 的点。我从这里尝试通过创建一个新的数据框来提示,然后再次使用预测(m1,newdata = new_df)进行绘图,但是我不知道如何创建这个数据框,因为我有 20 只老鼠而且我不知道如何将其扩展到 predict()。
我的尝试:
pframe <- data.frame(Day=seq(0, 200, length.out=101))
pframe$continuous_outcome <- predict(m1, newdata = pframe, level = 0)
ggplot(data, aes(Day,lVolume)) +
geom_point() +
geom_line(data=pframe)
但是我得到一个错误:
Error in eval(predvars, data, env) : object 'Rat' not found
而且,还有一种方法可以绘制你从每次估计中生成的方程式,即,从每只老鼠你有一组固定和随机的估计量,我如何绘制方程式(红色曲线) lmer 正在为每只老鼠生成?
原来 predict
比 broom.mixed::augment
更容易使用。
构建预测
(Rat/Country/Days 0-150 的所有组合(一天达到 200 导致一些极端预测超出了垂直比例)
library(tidyverse)
dc <- distinct(dplyr::select(dat1, Rat, Country))
pframe <- (with(dat1,
expand_grid(Rat = unique(Rat),
Day = 0:150))
%>% full_join(dc, by = "Rat")
%>% mutate(lVolume = predict(m1, newdata = .))
)
将数据和预测合并到一个数据框中(您不必这样做,但它使图例变得容易)
comb <- dplyr::bind_rows(list(data = dat1, model = pframe),
.id = "type")
剧情:
ggplot(comb, aes(Day, exp(lVolume), colour = type)) +
geom_point(alpha = 0.2) +
geom_line(aes(group = interaction(type, Rat))) +
scale_colour_manual(values = c("black", "red"))
重建数据:
dat0 <- list(
list("rat1", vol=c(78,304,352,690,952,1250), days = c(89,110,117,124,131,138), country = "Chile"),
list("rat2", vol=c(202,440,520,870,1380), days = c(75,89,96,103,110), country = "Chile"),
list("rat3", vol=c(186,370,620,850,1150), days = c(75,89,96,103,110), country = "Chile"),
list("rat4", vol=c(92,250,430,450,510,850,1000,1200), days = c(47,61,75,82,89,97,103,110), country = "England"),
list("rat5", vol=c(110,510,710,1200), days = c(47,61,75,82), country = "England"),
list("rat6", vol=c(115,380,480,540,560,850,1150,1350), days = c(47,61,75,82,89,97,103,110), country = "England"))
dat1 <- purrr::map_dfr(dat0,
~ data.frame(Rat = .[[1]],
lVolume = log(.$vol), Day = .$days,
Country = .$country))
m1 <- lmer(lVolume ~ Country*Day + (1|Rat), data = dat1)