if (FALSE) {
## tensor product `tesmd1' example
require(scam)
simu <- function(x,z) { -exp(4*x)/(1+exp(4*x))+2*sin(pi*z) }
xs <-seq(-1,3,length=30); zs <- seq(0,1,length=30)
pr <- data.frame(x=rep(xs,30),z=rep(zs,rep(30,30)))
truth <- matrix(simu(pr$x,pr$z),30,30)
set.seed(5)
n <- 500
x <- runif(n)*4-1
z <- runif(n)
f <- simu(x,z)
y <- f + rnorm(n)*.3
## fit model ...
b <- scam(y~s(x,z,bs="tesmd1",k=10))
summary(b)
old.par <- par(mfrow=c(2,2),mar=c(4,4,2,2))
plot(b,se=TRUE)
plot(b,pers=TRUE,theta = 30, phi = 40);title("tesmd1")
plot(y,b$fitted.values,xlab="Simulated data",ylab="Fitted data",pch=".",cex=3)
persp(xs,zs,truth,theta = 30, phi = 40);title("truth")
par(old.par)
vis.scam(b,theta = 30, phi = 40)
## example with cyclic cubic regression spline along the second covariate...
simu2 <- function(x,z) { -exp(4*x)/(1+exp(4*x))+sin(2*pi*z) }
xs <-seq(-1,3,length=30); zs <- seq(0,1,length=30)
pr <- data.frame(x=rep(xs,30),z=rep(zs,rep(30,30)))
truth <- matrix(simu2(pr$x,pr$z),30,30)
set.seed(2)
n <- 500
x <- runif(n)*4-1
z <- runif(n)
f <- simu2(x,z)
y <- f + rnorm(n)*.3
## fit model ...
b <- scam(y~s(x,z,bs="tesmd1",xt=list("cc"),k=10))
summary(b)
old.par <- par(mfrow=c(2,2),mar=c(4,4,2,2))
plot(b,se=TRUE)
plot(b,pers=TRUE,theta = 30, phi = 40);title("tesmd1, cyclic")
plot(y,b$fitted.values,xlab="Simulated data",ylab="Fitted data",pch=".",cex=3)
persp(xs,zs,truth,theta = 30, phi = 40);title("truth")
par(old.par)
vis.scam(b,theta = 30, phi = 40)
}
Run the code above in your browser using DataLab