I want to evaluate each component of the predictor from a GAM model separately using the option type="terms"
. As a sanity check, I compared the results to an evaluation of the total prediction using the option type="response"
.
It turns out that the results differ. Here is an example:
library(mgcv)
n<-200
sig <- 2
dat <- gamSim(1,n=n,scale=sig)
b<-gam(y~x0+s(I(x1^2))+s(x2)+offset(x3),da=dat)
nd <- data.frame(x0=c(.25,.5),x1=c(.25,.5),x2=c(.25,.5),x3=c(.25,.5))
a1 <- predict.gam(b,newdata=nd,type="response")
a2 <- rowSums(predict.gam(b,newdata=nd,type="terms")) + b$coefficients[1]
a1 - a2 # Should be zero!
# 1 2
# 0.25 0.50
Can anyone help me with this problem? Thank you very much for your help!
Your model:
y ~ x0 + s(I(x1^2)) + s(x2) + offset(x3)
has an offset term.
Offset will be considered by predict.gam
when type = "link"
or type = "response"
, but not considered when type = "terms"
.
a1 <- predict.gam(b, newdata=nd, type="response")
# 1 2
#11.178280 6.865068
a2 <- predict.gam(b, newdata=nd, type="terms")
# x0 s(I(x1^2)) s(x2)
#1 0.006878346 -1.8710120 5.6467813
#2 0.013756691 -0.6037635 -0.1905571
#attr(,"constant")
#(Intercept)
# 7.145632
So you have to add offset yourself:
a2 <- rowSums(a2) + b$coef[1] + nd$x3
# 1 2
#11.178280 6.865068
Now a1
and a2
are the same.
In case you wonder, I have documentation for you in ?predict.gam
:
type: ... When ‘type="terms"’ each component of the linear
predictor is returned seperately (possibly with standard
errors): this includes parametric model components, followed
by each smooth component, **but excludes any offset and any
intercept**.