qgam example

qgam-example.R

library(mgcv)
library(qgam)
library(tidyverse)

sim_dat <- tibble(time = 1:3000,
       fast_c = (time %% 200)/200,
       slow_c = (time %% 1150)/1150,
       y_true = cos(pi*time/3000) + 
         5 * sin(2*pi*fast_c) + 
         2 * sin(2*pi*slow_c) + 
         0.2 * sin(2*pi*slow_c) * sin(2*pi*fast_c),
       y = y_true)

plot(sim_dat$time, sim_dat$y, type = 'l')

# Add noise

sim_dat$y[2031:2310] <- sim_dat$y[2031:2310] + 5 * sin(seq(0, 6*pi, length.out = 280))

plot(sim_dat$time, sim_dat$y, type = 'l')

gam_ex <- gam(y ~ s(fast_c, bs = 'cc', k = 20) + 
                s(slow_c, bs = 'cc', k = 20) + 
                ti(fast_c, slow_c, bs = c('cc', 'cc')) +
                s(time), 
           data = sim_dat,
            rho = 0.95,
            method = "REML")

gratia::draw(gam_ex)

plot(sim_dat$time, sim_dat$y, type = 'l', main = "GAM")
lines(predict(gam_ex), col = 'red')

qgam_ex <- qgam(y ~ s(fast_c, bs = 'cc', k = 20) + 
                 s(slow_c, bs = 'cc', k = 20) + 
                 ti(fast_c, slow_c, bs = c('cc', 'cc')) +
                 s(time), 
                qu = 0.5,
               data = sim_dat,
               argGam = list(rho = 0.95,
               method = "REML"))

gratia::draw(qgam_ex)

plot(sim_dat$time, sim_dat$y, type = 'l', main = "qGAM")
lines(predict(qgam_ex), col = 'red')

plot(resid(gam_ex))
plot(qgam_ex$residuals)