mgcv：为“全局”分组平滑指定因子平滑交互

Question

我找到了Pedersen 等人的Hierarchical GAM 预印本（和GH repo ）。对于功能反应的组间变异建模非常有帮助，但我遇到了绊脚石。

我有一些时间序列数据（随着时间的推移计数），具有以下基本期望：

1. 计数和时间之间存在感兴趣的全局功能响应；
2. 全局函数在几个“固定”组 (2-5) 之间变化，可能具有类似的平滑参数；
3. 每年与该全局函数的“随机”偏差似乎存在系统趋势（即，功能响应的形状朝一致的方向漂移）；但
4. 年度偏差趋势对"固定"组的应用有所不同。

在这里，我尝试用一​​个人为的例子来说明，因为我担心我的描述不够充分。我为两个“固定”组/治疗模拟 20 个重复（年度）时间序列。在对照（trt == 0）中，年度变化很小。在治疗中（trt == 1），响应的形状存在系统趋势。

library(reshape)
library(dplyr)
library(ggplot2)
set.seed(2020)
n_yr <- 20
n_trt <- 2
n_x <- 10
dat <- tibble(yr = rep(seq(0, n_yr - 1), n_trt)) %>%
  expand.grid.df(tibble(x = rep(seq(0, 1, length.out = 10), n_trt),
                        tweak = c(rep(0, 15), rep(0.02, 5)),
                        trt = rep(0:1, each = n_x)), .) %>%
  mutate(e = rnorm(n_x * n_yr * n_trt, 0, 0.3),
         y = 3 + 9 * (x - tweak * yr) - 1 * (x - tweak * yr)^2 - 10 * (x - tweak * yr)^3 + e,
         fyr = factor(yr),
         ftrt = factor(trt),
         fty = factor(interaction(trt,yr)))
ggplot(dat, aes(x, y)) + geom_line(aes(group = yr, colour = yr)) +
    facet_grid(trt ~ .) + theme_bw() +
    scale_colour_gradient("Year", low = "#2166ac", high = "#b2182b")

捕获感兴趣的功能反应以及对各组之间不同的反应的年际变化进行建模的适当机制是什么？

我最好的猜测是为我的“全局”平滑建模因子平滑交互（使用处理作为by变量）。但是，如何对这些单独平滑的年度偏差进行建模呢？

来自单一全球响应的按年的简单因子平滑 ( fs; ?factor.smooth.interaction) 似乎不够：

library(mgcv)
m <- gam(y ~ s(x, by = ftrt, k = 5) + ftrt + s(x, fyr, bs = "fs", k = 5, m = 1),
         data = dat, method = "REML")
plot(m, pages = 1)

fs按治疗年份进行交互似乎更合适：

# Different smoothing parameters by treatment
m1 <- gam(y ~ s(x, by = ftrt, k = 5) + ftrt + s(x, fty, bs = "fs", k = 5, m = 1),
          data = dat, method = "REML")
p <- plot(m1, pages = 1)

下图更好地说明了治疗之间与治疗平滑度的年度偏差的差异，以及治疗组这些偏差形状的“趋势”：

df <- data.frame(p[[3]], 
                 fty = rep(levels(dat$fty), each = 100),
                 trt = rep(rep(0:1, each = 100), 20),
                 yr = rep(0:19, each = 200))
ggplot(df, aes(x, fit, group = fty)) +
  geom_line(aes(color = yr)) +
  scale_colour_gradient("Year", low = "#2166ac", high = "#b2182b") +
  facet_grid(trt ~ .) + theme_bw()

所以我的问题是：

1. 固定治疗组 ( ) 的一组基本因子平滑s(x, by = ftrt)与某种全局平滑 ( s(x)) 是否合适？
2. 是否有更好的方法来允许治疗组的平滑形状分别出现“随机”年度变化？
3. 有什么方法可以捕获某些组形状变化的“趋势”，即对平滑参数或类似参数的线性变化进行建模？