The table present the different plot functions in the packages mgcv (Wood 2006, 2011) and itsadug for visualizing GAMM models.


Partial effect Sum of all effects Sum of “fixed” effects1
surface plot.gam() vis.gam()
pvisgam() fvisgam()
smooth plot.gam() plot_smooth()
group estimates plot.gam()2 plot_parametric()
random smooths get_random(), inspect_random()

1: include rm.ranef=TRUE to zero all random effects.

2: include all.terms=TRUE to visualize parametric terms.


Examples

library(itsadug)
library(mgcv)
data(simdat)

## Not run: 
# Model with random effect and interactions:
m1 <- bam(Y ~ Group + te(Time, Trial, by=Group)
    + s(Time, Subject, bs='fs', m=1, k=5),
    data=simdat)
# Simple model with smooth:
m2 <-  bam(Y ~ Group + s(Time, by=Group)
    + s(Subject, bs='re')
    + s(Subject, Time, bs='re'),
    data=simdat)

Summary model m1:

gamtabs(m1, type='html')
A. parametric coefficients Estimate Std. Error t-value p-value
(Intercept) 1.3765 0.5420 2.5398 0.0111
GroupAdults 2.7069 0.7782 3.4783 0.0005
B. smooth terms edf Ref.df F-value p-value
te(Time,Trial):GroupChildren 19.7349 21.9399 1651.5464 < 0.0001
te(Time,Trial):GroupAdults 19.6303 21.8781 1466.8418 < 0.0001
s(Time,Subject) 175.1741 178.0000 5803.4440 < 0.0001

Summary model m2:

gamtabs(m2, type='html')
A. parametric coefficients Estimate Std. Error t-value p-value
(Intercept) 2.0575 2.9497 0.6975 0.4855
GroupAdults 3.1264 4.1715 0.7495 0.4536
B. smooth terms edf Ref.df F-value p-value
s(Time):GroupChildren 8.6203 8.9584 1630.0297 < 0.0001
s(Time):GroupAdults 8.8372 8.9921 8798.1199 < 0.0001
s(Subject) 33.9912 34.0000 125335342.8798 < 0.0001
s(Subject,Time) 33.9878 34.0000 125318454.9121 < 0.0001

a. Surfaces

Plotting the partial effects of te(Time,Trial):GroupAdults and te(Time,Trial):GroupChildren.

pvisgam()

par(mfrow=c(1,2), cex=1.1)

pvisgam(m1, view=c("Time", "Trial"), select=1,
        main="Children", zlim=c(-12,12))
pvisgam(m1, view=c("Time", "Trial"), select=2,
        main="Adults", zlim=c(-12,12))

Notes:

  • Plots same data as plot(m1, select=1): partial effects plot, i.e., the plot does not include intercept or any other effects.

  • Make sure to set the zlim values the same when comparing surfaces

  • Use the argument cond to specify the value of other predictors in a more complex interaction. For example, for plotting a modelterm te(A,B,C) use something like pvisgam(model, view=c("A", "B"), select=1, cond=list(C=5)).

fvisgam()

Plotting the fitted effects of te(Time,Trial):GroupAdults and te(Time,Trial):GroupChildren.

par(mfrow=c(1,2), cex=1.1)

fvisgam(m1, view=c("Time", "Trial"), cond=list(Group="Children"),
        main="Children", zlim=c(-12,12), rm.ranef=TRUE)
fvisgam(m1, view=c("Time", "Trial"), cond=list(Group="Adults"),
        main="Adults", zlim=c(-12,12), rm.ranef=TRUE)