Well it's not pretty, but I guess you could try something like this...
best,
simon
library(mgcv)
dat <- gamSim(1,n=200,dist="normal",scale=2)
b <- gam(y ~ s(x1)+s(x2),data=dat)
m <- 40
x1 <- seq(0,1,length=m)
x2 <- seq(0,1,length=m)
pd <- expand.grid(x1=x1,x2=x2)
fv <- predict(b,newdata=pd,se=TRUE)
zlim=range(dat$y)
pmat <- persp(x1,x2,matrix(fv$fit,m,m),theta=-120,zlim=zlim,zlab="y",col=NA)
par(new=TRUE)
persp(x1,x2,matrix(fv$fit-2*fv$se.fit,m,m),theta=-120,
zlim=zlim,col=NA,axes=FALSE,box=FALSE,border="grey")
par(new=TRUE)
persp(x1,x2,matrix(fv$fit,m,m),theta=-120,zlim=zlim,zlab="y",col=NA)
par(new=TRUE)
persp(x1,x2,matrix(fv$fit+2*fv$se.fit,m,m),theta=-120,
zlim=zlim,col=NA,axes=FALSE,box=FALSE,border="grey")
raw <- trans3d(dat$x1,dat$x2,dat$y,pmat)
fit <- trans3d(dat$x1,dat$x2,fitted(b),pmat)
points(raw$x,raw$y,pch=19,cex=.5)
lines(t(cbind(raw$x,fit$x,NA)),t(cbind(raw$y,fit$y,NA)))
On 09/11/2012 10:07 AM, Anna Zakrisson wrote:
> Hi,
>
> I have used the mgcv library to generate a simple additive model. I want to
> know how to plot the function on the raw data with confidence intervals whan
> I have TWO variables in the model. I get it to work with one variable but
> not with two. I am on the limit for what I understand in R, so be gentle. I
> have read the help file on predict.gam, but did not get any help out of it.
>
> #My model:
> Model <- gam(var ~ s(var1, k=4) + s(var2, k=4), data = mydata)
>
> # Plotting the model:
> par(mfrow=c(1,1))
> plot(Model ) #get a plot of the smoother
> plot(var ~ var1, data = mydata) #plotting the raw data
>
> ##predictions from M1
> pred.gam <- predict.gam(M1, se.fit=TRUE, type="response")
>
> #Plotting the function and confidence intervals.
> I <- order(mydata$Ncell)
> lines(mydata$Ncell[I], pred.gam$fit[I], lwd = 2)
> lines(mydata$Ncell[I], pred.gam$fit[I] + 2 * pred.gam$se.fit[I], lty = 2,
> lwd = 2)
> lines(mydata$Ncell[I], pred.gam$fit[I] - 2 * pred.gam$se.fit[I], lty = 2,
> lwd = 2)
>
> It all works fine if the model instead is:
> ModelSimple <- gam(var ~ s(var1, k=4), data = mydata)
>
> However, I get different results (for obvious reasons) if I predict from
> Model or ModelSimple.
> Please, help me!!!!
> Anna
>
> Anna Zakrisson Braeunlich
> PhD Student
>
> Department of Systems Ecology
> Stockholm University
> Svante Arrheniusv. 21A
> SE-106 91 Stockholm
>
> E-mail: a...@ecology.su.se
> Tel work: +46 (0)8 161103
> Mobile: +46-(0)700-525015
> Web site: http://www.ecology.su.se/staff/personal.asp?id=163
>
>> <((((º>`âEUR¢. . âEUR¢ `âEUR¢. .âEUR¢ `âEUR¢. . ><((((º>`âEUR¢. . âEUR¢
>> `âEUR¢. .âEUR¢
> `âEUR¢. .><((((º>
>
>
> [[alternative HTML version deleted]]
>
>
>
> ______________________________________________
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
