Skip to content

plotting gam curve against predictors

3 messages · Sonya Ku, Suresh Krishna, Gavin Simpson

Sonya Ku
# 
Hi 
I am just beginning to learn R and have fitted several GAM to my
species presense/absence data.

I have used plot(x,y) using fitted.values as a y variable against
predictors. However, it is hard to see general relationships where
there is wide spread in predicted values for any x.

So I'd like to 
1) plot general curves, not individual fitted values against predictors
2) have some indication of confidence interval 95% but this is optional

would greatly appreciate any tips

Sonya.
1 day later
Suresh Krishna
#
Sonya Ku wrote:
What exactly do you mean by "general curves" or "general relationships" ?

-s.
Gavin Simpson
#
Suresh Krishna wrote:
I think Sonya meant plotting the response curve of the species over the 
full environmental gradient in question rather then plotting a line that 
goes through the fitted values. If you have a small number of samples or 
they are clumped along that gradient, just plotting a line through the 
fitted points can result in a strange response curve.

The simplest way is to create a new set of predictors variable(s) that 
are evenly spread across the range of your environmental data. Then use 
the predict method to predict the abundance of each species for each 
value of the predictor you created.

e.g. assuming you used mgcv for GAMs (not stated), using dummy code:

# fit the model
gam.mod <- gam(spp ~ s(pH), family = binomial(link = "logit"))
# create some new predictor data aross range of your measured values
# Here we create 100 evenly spaced values
new.dat <- data.frame(pH = seq(min(pH), max(pH), length = 100))
# predict abundance for these new values
gam.pred <- predict(gam.mod, newdata = new.dat, type = "response")
# plot the response curve
plot(new.dat, gam.pred, type = "l")

if you want the standard errors as well then:

# note se.fit = TRUE, gam.pred contains to items in a list $fit and $se.fit
gam.pred <- predict(gam.mod, newdata = new.dat, type = "response",
	se.fit = TRUE)
# so we have to use pam.pred$fit for the predicted values
plot(new.dat, gam.pred$fit, type = "l")
# and +/- the se.fit for the standard errors
lines(new.dat, gam.pred$fit + gam.pred$se.fit, lty = "dotted")
lines(new.dat, gam.pred$fit, type = "l", lty = "dotted")

See ?predict.gam for more information.

I assume the procedure will be similar if you use package gam instead, 
the details of the predict methods may differ though.

HTH

Gav