绘制分段增长曲线答案

【问题标题】：Plotting Piecewise growth curves绘制分段增长曲线
【发布时间】：2021-04-12 14:30:41
【问题描述】：

我正在尝试绘制类似于this first plot 的分段增长曲线。我使用了单独的斜率编码方案，并在时间 2 放置了一个断点

|时间 | 0 | 1 | 2 | 5 | 10 | 15 | 20|

|时间1 | 0 | 1 | 2 | 2 | 2 | 2 | 2 |

|时间2 | 0 | 0 | 0 | 1 | 2 | 3 | 4 |

我使用以下代码创建我的增长模型

m1 <- lmer(sdmtwr ~ time1 + time2 + (time1 | id) + (0 + time2 | id), data = SDMT, REML = FALSE)

我还在探索使用以下代码与 2 级分类预测器的交互

m2 <- lmer(sdmtwr ~ (time1 + time2)*edu + (time1 | id) + (0 + time2 | id), data = SDMT, REML = FALSE)

我尝试使用 ggplot2、sjPlot 和效果包创建绘图，但无济于事，由于编程经验有限，我不知所措。我只能为基线模型和交互模型分别绘制分段。

如果有人可以就相应的代码提供帮助，我将不胜感激！

编辑：这是 dput 摘要（编辑长度以显示 edu、time1 和 time2）

> dput(sdmt)
structure(list(id = c(3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 5L, 6L, 
6L, 6L, 28L, 28L, 28L, 28L, 28L, 28L, 28L, 62L, 62L, 62L, 62L, 
108L, 108L, 108L, 108L, 119L, 119L, 120L, 120L, 120L, 120L, 132L, 
132L, 132L, 132L, 132L, 148L, 148L, 148L, 148L, 148L, 148L, 175L, 
175L, 175L, 178L, 178L, 178L, 178L, 201L, 201L, 201L, 201L, 201L, 
201L, 201L, 253L, 253L, 253L, 253L, 327L, 327L, 327L, 327L, 336L, 
336L, 336L, 336L, 336L, 336L, 343L, 343L, 360L, 360L, 360L, 366L, 
366L, 366L), time = c(0L, 2L, 10L, 15L, 20L, 5L, 10L, 15L, 2L, 
2L, 15L, 20L, 0L, 1L, 2L, 5L, 10L, 15L, 20L, 5L, 10L, 15L, 20L, 
0L, 2L, 15L, 20L, 0L, 2L, 0L, 10L, 15L, 20L, 0L, 1L, 5L, 10L, 
20L, 1L, 2L, 5L, 10L, 15L, 20L, 0L, 1L, 2L, 0L, 1L, 2L, 5L, 0L, 
1L, 2L, 5L, 10L, 15L, 20L, 0L, 1L, 5L, 15L, 0L, 1L, 10L, 20L, 
0L, 1L, 5L, 10L, 15L, 20L, 0L, 10L, 1L, 5L, 10L, 0L, 10L, 15L
), sdmtwr = c(20L, 24L, 18L, 19L, 9L, 17L, 24L, 17L, 41L, 33L, 
27L, 29L, 31L, 29L, 26L, 29L, 32L, 20L, 19L, 40L, 42L, 46L, 38L, 
14L, 25L, 24L, 29L, 46L, 45L, 29L, 26L, 34L, 38L, 30L, 33L, 71L, 
52L, 51L, 29L, 33L, 50L, 55L, 40L, 39L, 32L, 34L, 35L, 28L, 37L, 
37L, 36L, 37L, 29L, 52L, 51L, 50L, 44L, 42L, 30L, 43L, 43L, 41L, 
33L, 46L, 49L, 38L, 52L, 50L, 48L, 49L, 49L, 50L, 40L, 39L, 18L, 
NA, 3L, 31L, 43L, 47L), time_seg1 = c(0, 2, 2, 2, 2, 2, 2, 2, 
2, 2, 2, 2, 0, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 0, 2, 2, 2, 0, 2, 
0, 2, 2, 2, 0, 1, 2, 2, 2, 1, 2, 2, 2, 2, 2, 0, 1, 2, 0, 1, 2, 
2, 0, 1, 2, 2, 2, 2, 2, 0, 1, 2, 2, 0, 1, 2, 2, 0, 1, 2, 2, 2, 
2, 0, 2, 1, 2, 2, 0, 2, 2), time_seg2 = c(0, 0, 2, 3, 4, 1, 2, 
3, 0, 0, 3, 4, 0, 0, 0, 1, 2, 3, 4, 1, 2, 3, 4, 0, 0, 3, 4, 0, 
0, 0, 2, 3, 4, 0, 0, 1, 2, 4, 0, 0, 1, 2, 3, 4, 0, 0, 0, 0, 0, 
0, 1, 0, 0, 0, 1, 2, 3, 4, 0, 0, 1, 3, 0, 0, 2, 4, 0, 0, 1, 2, 
3, 4, 0, 2, 0, 1, 2, 0, 2, 3), ed_dich = structure(c(2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, NA, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L), .Label = c("< HS", 
">= HS"), class = "factor")), row.names = c(NA, -80L), class = "data.frame")

【问题讨论】：

您能否通过将dput(SDMT) 的结果粘贴到代码块中来发布数据。我认为所有功能都不适合您的原因是这些功能不了解time1 和time2 之间的确定性关系。
当然，我创建了一个数据子集，仅显示 time1 和 time2 以及我的分类预测器 edu。
我们需要所有变量来运行模型，因此包括sdmtwr 和id 以及time 会很有用。如果你不想把所有东西都放完，你可以先tmp <- SDMT %>% select(id, time, sdmtwr, time1, time2, edu) %>% filter(id %in% unique(SDMT$id)[1:20]) 然后dput(tmp)。
我遇到了所有变量的空间问题，所以我使用了您提供的tmp 代码并粘贴了上面的摘要。

标签： r ggplot2 regression lme4

【解决方案1】：

我认为您想要的是分段线性样条曲线。您可以使用截断的幂基函数来执行此操作。在您的模型中，如果时间大于 2，您将包含时间和时间为 2 的函数，否则为 0。这使得分段线性函数在时间 = 2 处相遇。您可以在模型中执行以下操作：

library(lme4)
mod <- lmer(sdmtwr ~ time + I(ifelse(time > 2, time-2, 0)) + 
              (1 |id), data=tmp, REML=TRUE)

然后，您可以使用 ggeffects 包中的 ggpredict() 函数来生成绘图：

library(ggeffects)
g <- ggpredict(mod, "time")
plot(g)

注意：我无法让它在对时间变量产生随机影响的情况下运行，但如果有更多数据，也许你就能让它工作。

【讨论】：

这真的很有趣！我尝试使用几个包来实现样条函数，但它们都不起作用。这并不是我目前所需要的，但我将在数据上使用这种样条方法。我将如何解释 I(ifelse(time > 2, time-2, 0)) 系数？这是否代表 time2 的斜率？
@johnsanders 大多数时候，您使用 bs() 或 ns() 之类的样条从 splines 包实现样条 - 它们默认为三次样条。这个简单的分段线性模型很容易实现，没有其他开销。
是的，我尝试了 splines 和 segmented 包，我不知道它们默认为三次样条。感谢您的帮助。
@johnsanders 两个系数之和将是 2-T 时间的斜率，time 上的系数是 0-2 时间的斜率。