【问题标题】:Add age adjustment to geom_smooth为 geom_smooth 添加年龄调整
【发布时间】:2020-01-03 18:21:38
【问题描述】:

我需要在我添加到我的 ggscatter 图中的 geom_smooth 行中包含年龄调整。

我的数据看起来像~ table link

structure(list(Time = c(0L, 0L, 0L, 0L, 6L, 12L, 18L, 18L, 0L, 
12L, 18L, 6L), group = structure(c(1L, 1L, 2L, 2L, 1L, 3L, 3L, 
3L, 3L, 4L, 4L, 1L), .Label = c("A", "B", "C", "D"), class = "factor"), 
    Age = c(77, 70.2, 69.9, 65.7, 66.2, 66.7, 67.2, 67.7, 66.8, 
    67.8, 68.3, 68.8), Average = c(96L, 90L, 94L, 94L, 96L, 96L, 
    92L, 120L, 114L, 109L, 113L, 103L)), row.names = c(NA, 12L
), class = "data.frame")

我目前拥有的('Average' 值依赖于年龄..):

ggscatter(dtable, "Time","Average",conf.int = TRUE)+theme_bw()+
geom_smooth(aes(group=1),method='lm')+facet_wrap(~groups)

我想要的是这样的:

ggscatter(dtable, "Time","Average",conf.int = TRUE)+theme_bw()+
geom_smooth(aes(group=1),method='lm', adjust= ~age)+facet_wrap(~groups)

每组平均年龄调整

有什么建议吗?

【问题讨论】:

  • 您可以运行lm,然后使用预测函数将线性值附加到数据集上,然后绘制它们。如果你分享你的代码,我可以给你看。
  • 谢谢。我添加了一个示例表。我不确定应该在哪个部分添加预测函数以及如何确保它按组划分。我应该单独创建一个“公式”并将其插入 geom_smooth 吗?
  • 嗨@BarrMorgenstein,尝试使用dput(head(your_data, 20)),而不是照片链接,然后将该输出放在帖子中。这将有助于解决问题。
  • @MDEWITT - 添加:)

标签: r static-methods adjustment


【解决方案1】:

这就是我想你所追求的。

首先,我们需要拟合更复杂的模型,因为 ggplot 还没有用于多变量模型的功能

fit <- lm(Average ~ Time + group + Age, data = tdata)

然后我们可以使用 broom 包中的一些功能来添加预测和相关的标准错误。有了这些,我们可以使用 geom_line 和 geom_ribbon geoms 手动构建绘图

library(broom)
tdata %>% 
  bind_cols(augment(fit)) %>% 
  ggplot(aes(Time, Average))+
  geom_point()+
  geom_line(aes(x = Time, y = .fitted), size = 2, color = "blue")+
  geom_ribbon(aes(ymin = .fitted + .se.fit*2, ymax = .fitted - .se.fit*2), alpha = .2)+
  facet_wrap(~group)+
  theme_bw()

此外,如果您想查看汇总与非汇总估计值

fit_no_pool <- lm(Average ~ Time + group + Age, data = tdata)
fit_complete_pool <- lm(Average ~ Time + Age, data = tdata)

library(broom)
tdata %>% 
  bind_cols(augment(fit_no_pool) %>% setNames(sprintf("no_pool%s", names(.)))) %>% 
  bind_cols(augment(fit_complete_pool) %>% setNames(sprintf("pool%s", names(.)))) %>% 
  ggplot(aes(Time, Average))+
  geom_point()+
  # Non-Pooled Estimates
  geom_line(aes(x = Time, y = no_pool.fitted, color = "blue"), size = 2)+
  geom_ribbon(aes(ymin = no_pool.fitted + no_pool.se.fit*2, 
                  ymax = no_pool.fitted - no_pool.se.fit*2), alpha = .2)+
  # Pooled Estimates
  geom_line(aes(x = Time, y = pool.fitted, color = "orange"), size = 2)+
  geom_ribbon(aes(ymin = pool.fitted + pool.se.fit*2, 
                  ymax = pool.fitted - pool.se.fit*2), alpha = .2)+
  facet_wrap(~group)+
  scale_color_manual(name = "Regression", 
                       labels = c("Pooled", "Non-Pooled"), 
                     values = c("blue", "orange"))+
  theme_bw()

【讨论】:

  • 他想要进行年龄调整并按组分面 - 如果将按组分面,则可能值得将年龄和组都放入模型中的真实数据中
  • 谢谢!有用!如果可能,最后一个问题 - 在调整之前我使用简单的线性回归线并评估斜率 - 有没有一种方法可以让我在这里创建斜率或随时间减少百分比?
  • 您的意思是提取参数估计值吗?你可以通过summary(fit)broom::tidy(fit) 来做到这一点。
  • @MDEWITT - 我试图运行表格,包括+组而不是分面,以及表格中没有+组的子集 - 结果似乎不同..
【解决方案2】:

一种方法是将年龄作为模型中的附加预测变量来运行模型。然后使用predict 得到带有 CI 的预测值。附加到您的数据,然后使用 ggplot 进行绘图。我知道您想通过group 进行刻面,因此也值得将其放入您的模型中。只是一个想法。步骤将是相同的。

df <- structure(list(Time = c(0L, 0L, 0L, 0L, 6L, 12L, 18L, 18L, 0L, 
    12L, 18L, 6L), group = structure(c(1L, 1L, 2L, 2L, 1L, 3L, 3L, 
    3L, 3L, 4L, 4L, 1L), .Label = c("A", "B", "C", "D"), class = "factor"), 
    Age = c(77, 70.2, 69.9, 65.7, 66.2, 66.7, 67.2, 67.7, 66.8, 
    67.8, 68.3, 68.8), Average = c(96L, 90L, 94L, 94L, 96L, 96L, 
    92L, 120L, 114L, 109L, 113L, 103L)), row.names = c(NA, 12L
    ), class = "data.frame")


#model adjusted for age
mod <- lm(Average ~ Time + Age, data = df)
#get prediction with CIS
premod <- predict(mod, interval = "predict")
#append to data
df2 <- cbind(df,premod)
#add prediction to ggplot with scatter plot
ggplot(df2) + 
    geom_point(aes(x=Time,y=Average)) +
    geom_line(aes(x=Time, y = fit)) +
    geom_ribbon(aes(x = Time,ymin = lwr, ymax = upr), alpha = .1)+
    facet_wrap(~group)+
    theme_bw()

【讨论】:

    猜你喜欢
    • 2023-01-29
    • 1970-01-01
    • 1970-01-01
    • 2018-09-30
    • 2022-12-03
    • 2012-03-26
    • 1970-01-01
    • 1970-01-01
    • 2017-02-05
    相关资源
    最近更新 更多