The R Project SVN R-packages

Rev

Rev 6283 | Blame | Compare with Previous | Last modification | View Log | Download | RSS feed

\name{smooth.construct.ds.smooth.spec}
\alias{smooth.construct.ds.smooth.spec}
\alias{Predict.matrix.duchon.spline}
\alias{Duchon.spline}
%- Also NEED an `\alias' for EACH other topic documented here.
\title{Low rank Duchon 1977 splines}

\description{Thin plate spline smoothers are a special case of the isotropic splines discussed in Duchon (1977). A subset of this more 
general class can be invoked by terms like \code{s(x,z,bs="ds",m=c(1,.5)} in a \code{\link{gam}} model formula.
In the notation of Duchon (1977) m is given by \code{m[1]} (default value 2), while s is given by \code{m[2]} (default value 0).

Duchon's (1977) construction generalizes the usual thin plate spline penalty as follows. The usual TPS penalty is given by the integral 
of the squared Euclidian norm of a vector of mixed partial mth order derivatives of the function w.r.t. its arguments. Duchon re-expresses 
this penalty in the Fourier domain, and then weights the squared norm in the integral by the Euclidean norm of the fourier frequencies, 
raised to the power 2s. s is a user selected constant taking integer values divided by 2. If d is the number of arguments of the smooth, 
then it is required that -d/2 < s < d/2. To obtain continuous functions we further require that m + s > d/2. If s=0 then the usual thin plate 
spline is recovered. 

The construction is amenable to exactly the low rank approximation method given in Wood (2003) to thin plate splines, with similar 
optimality properties, so this approach to low rank smoothing is used here. For large datasets the same subsampling approach as is used in the 
\code{\link{tprs}} case is employed here to reduce computational costs. 

These smoothers allow the use of lower orders of derivative in the penalty than conventional thin plate splines, 
while still yielding continuous functions. For example, we can set m = 1 and s = d/2 - .5 in order to use first derivative penalization 
for any d (which has the advantage that the dimension of the null space of unpenalized functions is only d+1). 

}

\usage{
\method{smooth.construct}{ds.smooth.spec}(object, data, knots)
\method{Predict.matrix}{duchon.spline}(object, data)
}

\arguments{ 
\item{object}{a smooth specification object, usually generated by a term \code{s(...,bs="ds",...)}.}

\item{data}{a list containing just the data (including any \code{by} variable) required by this term, 
            with names corresponding to \code{object$term} (and \code{object$by}). The \code{by} variable 
            is the last element.} 

\item{knots}{a list containing any knots supplied for basis setup --- in same order and with same names as \code{data}. 
             Can be \code{NULL}}

}

\value{ An object of class \code{"duchon.spline"}. In addition to the usual elements of a 
smooth class documented under \code{\link{smooth.construct}}, this object will contain:
\item{shift}{A record of the shift applied to each covariate in order to center it around zero and 
avoid any co-linearity problems that might otehrwise occur in the penalty null space basis of the term. }
\item{Xu}{A matrix of the unique covariate combinations for this smooth (the basis is constructed by first stripping 
out duplicate locations).}
\item{UZ}{The matrix mapping the smoother parameters back to the parameters of a full Duchon spline.}
\item{null.space.dimension}{The dimension of the space of functions that have zero wiggliness according to the 
wiggliness penalty for this term.}
}

\details{  The default basis dimension for this class is \code{k=M+k.def} where \code{M} is the null space dimension 
(dimension of unpenalized function space) and \code{k.def} is 10 for dimension 1, 30 for dimension 2 and 100 for higher dimensions. 
This is essentially arbitrary, and should be checked, but as with all penalized regression smoothers, results are statistically 
insensitive to the exact choise, provided it is not so small that it forces oversmoothing (the smoother's 
degrees of freedom are controlled primarily by its smoothing parameter).

The constructor is not normally called directly, but is rather used internally by \code{\link{gam}}. 
To use for basis setup it is recommended to use \code{\link{smooth.construct2}}.  

For these classes the specification \code{object} will contain
information on how to handle large datasets in their \code{xt} field. The default is to randomly
subsample 2000 `knots' from which to produce a reduced rank eigen approximation to the full basis, 
if the number of unique predictor variable combinations in excess of 2000. The default can be
modified via the \code{xt} argument to \code{\link{s}}. This is supplied as a
list with elements \code{max.knots} and \code{seed} containing a number
to use in place of 2000, and the random number seed to use (either can be
missing). Note that the random sampling will not effect the state of R's RNG. 

For these bases \code{knots} has two uses. Firstly, as mentioned already, for large datasets 
the calculation of the \code{tp} basis can be time-consuming. The user can retain most of the advantages of the 
approach by supplying  a reduced set of covariate values from which to obtain the basis - 
typically the number of covariate values used will be substantially 
smaller than the number of data, and substantially larger than the basis dimension, \code{k}. This approach is 
the one taken automatically if the number of unique covariate values (combinations) exceeds \code{max.knots}.
The second possibility is to avoid the eigen-decomposition used to find the spline basis altogether and simply use 
the basis implied by the chosen knots: this will happen if the number of knots supplied matches the 
basis dimension, \code{k}. For a given basis dimension the second option is 
faster, but gives poorer results (and the user must be quite careful in choosing knot locations).


}
\seealso{\code{\link{Spherical.Spline}}}

\references{
Duchon, J. (1977) Splines minimizing rotation-invariant semi-norms in Solobev spaces. in W. Shemp and K. Zeller (eds) 
Construction theory of functions of several variables, 85-100, Springer, Berlin.

Wood, S.N. (2003) Thin plate regression splines. J.R.Statist.Soc.B 65(1):95-114
}

\author{ Simon N. Wood \email{simon.wood@r-project.org}}

\examples{
require(mgcv)
eg <- gamSim(2,n=200,scale=.05)
attach(eg)
op <- par(mfrow=c(2,2),mar=c(4,4,1,1))
b0 <- gam(y~s(x,z,bs="ds",m=c(2,0),k=50),data=data)  ## tps
b <- gam(y~s(x,z,bs="ds",m=c(1,.5),k=50),data=data)  ## first deriv penalty
b1 <- gam(y~s(x,z,bs="ds",m=c(2,.5),k=50),data=data) ## modified 2nd deriv

persp(truth$x,truth$z,truth$f,theta=30) ## truth
vis.gam(b0,theta=30)
vis.gam(b,theta=30)
vis.gam(b1,theta=30)

detach(eg)

}
\keyword{models} \keyword{regression}%-- one or more ..