smooth.construct  R Documentation 
Smooth terms in a GAM formula are turned into smooth specification objects of
class xx.smooth.spec
during processing of the formula. Each of these objects is
converted to a smooth object using an appropriate smooth.construct
function. New smooth classes
can be added by writing a new smooth.construct
method function and a corresponding
Predict.matrix
method function (see example code below).
In practice, smooth.construct
is usually called via smooth.construct2
and the wrapper
function smoothCon
, in order to handle by
variables and
centering constraints (see the smoothCon
documentation if
you need to handle these things directly, for a user defined smooth class).
smooth.construct(object,data,knots)
smooth.construct2(object,data,knots)
object 
is a smooth specification object, generated by an
If

data 
For 
knots 
an optional data frame or list containing the knots relating to 
There are built in methods for objects with the following classes:
tp.smooth.spec
(thin plate regression splines: see tprs
);
ts.smooth.spec
(thin plate regression splines with shrinkagetozero);
cr.smooth.spec
(cubic regression splines: see cubic.regression.spline
;
cs.smooth.spec
(cubic regression splines with shrinkagetozero);
cc.smooth.spec
(cyclic cubic regression splines);
ps.smooth.spec
(Eilers and Marx (1986) style Psplines: see p.spline
);
cp.smooth.spec
(cyclic Psplines);
ad.smooth.spec
(adaptive smooths of 1 or 2 variables: see adaptive.smooth
);
re.smooth.spec
(simple random effect terms);
mrf.smooth.spec
(Markov random field smoothers for smoothing over discrete districts);
tensor.smooth.spec
(tensor product smooths).
There is an implicit assumption that the basis only depends on the knots and/or the set of unique covariate combinations; i.e. that the basis is the same whether generated from the full set of covariates, or just the unique combinations of covariates.
Plotting of smooths is handled by plot methods for smooth objects. A default mgcv.smooth
method
is used if there is no more specific method available. Plot methods can be added for specific smooth classes, see
source code for mgcv:::plot.sos.smooth
, mgcv:::plot.random.effect
, mgcv:::plot.mgcv.smooth
for example code.
The input argument object
, assigned a new class to indicate what type of smooth it is and with at least the
following items added:
X 
The model matrix from this term. This may have an 
S 
A list of positive semidefinite penalty matrices that apply to this term. The list will be empty if the term is to be left unpenalized. 
rank 
An array giving the ranks of the penalties. 
null.space.dim 
The dimension of the penalty null space (before centering). 
The following items may be added:
C 
The matrix defining any identifiability constraints on the term, for use when fitting. If this is 
Cp 
An optional matrix supplying alternative identifiability constraints for use when predicting. By default the fitting constrants are used. This option is useful when some sort of simple sparse constraint is required for fitting, but the usual sumtozero constraint is required for prediction so that, e.g. the CIs for model components are as narrow as possible. 
no.rescale 
if this is nonNULL then the penalty coefficient matrix of the smooth will not be
rescaled for enhaced numerical stability (rescaling is the default, because 
df 
the degrees of freedom associated with this term (when
unpenalized and unconstrained). If this is null then 
te.ok 

plot.me 
Set to 
side.constrain 
Set to 
L 
smooths may depend on fewer ‘underlying’ smoothing parameters than there are elements of

Usually the returned object will also include extra information required to define the basis, and used by
Predict.matrix
methods to make predictions using the basis. See
the Details
section for links to the information included for the built in smooth classes.
tensor.smooth
returned objects will additionally have each element of
the margin
list updated in the same way. tensor.smooths
also
have a list, XP
, containing reparameterization matrices for any 1D marginal terms
reparameterized in terms of function values. This list will have NULL
entries for marginal smooths that are not reparameterized, and is only long
enough to reach the last reparameterized marginal in the list.
User defined smooth objects should avoid having attributes names
"qrc"
or "nCons"
as these are used internally to provide
constraint free parameterizations.
Simon N. Wood simon.wood@rproject.org
Wood, S.N. (2003) Thin plate regression splines. J.R.Statist.Soc.B 65(1):95114
Wood, S.N. (2006) Low rank scale invariant tensor product smooths for generalized additive mixed models. Biometrics 62(4):10251036
The code given in the example is based on the smooths advocated in:
Ruppert, D., M.P. Wand and R.J. Carroll (2003) Semiparametric Regression. Cambridge University Press.
However if you want psplines, rather than splines with derivative based penalties, then the built in "ps" class is probably a marginally better bet. It's based on
Eilers, P.H.C. and B.D. Marx (1996) Flexible Smoothing with Bsplines and Penalties. Statistical Science, 11(2):89121
https://www.maths.ed.ac.uk/~swood34/
s
,get.var
, gamm
, gam
,
Predict.matrix
,
smoothCon
, PredictMat
## Adding a penalized truncated power basis class and methods
## as favoured by Ruppert, Wand and Carroll (2003)
## Semiparametric regression CUP. (No advantage to actually
## using this, since mgcv can happily handle nonidentity
## penalties.)
smooth.construct.tr.smooth.spec<function(object,data,knots) {
## a truncated power spline constructor method function
## object$p.order = null space dimension
m < object$p.order[1]
if (is.na(m)) m < 2 ## default
if (m<1) stop("silly m supplied")
if (object$bs.dim<0) object$bs.dim < 10 ## default
nk<object$bs.dimm1 ## number of knots
if (nk<=0) stop("k too small for m")
x < data[[object$term]] ## the data
x.shift < mean(x) # shift used to enhance stability
k < knots[[object$term]] ## will be NULL if none supplied
if (is.null(k)) # space knots through data
{ n<length(x)
k<quantile(x[2:(n1)],seq(0,1,length=nk+2))[2:(nk+1)]
}
if (length(k)!=nk) # right number of knots?
stop(paste("there should be ",nk," supplied knots"))
x < x  x.shift # basis stabilizing shift
k < k  x.shift # knots treated the same!
X<matrix(0,length(x),object$bs.dim)
for (i in 1:(m+1)) X[,i] < x^(i1)
for (i in 1:nk) X[,i+m+1]<(xk[i])^m*as.numeric(x>k[i])
object$X<X # the finished model matrix
if (!object$fixed) # create the penalty matrix
{ object$S[[1]]<diag(c(rep(0,m+1),rep(1,nk)))
}
object$rank<nk # penalty rank
object$null.space.dim < m+1 # dim. of unpenalized space
## store "tr" specific stuff ...
object$knots<k;object$m<m;object$x.shift < x.shift
object$df<ncol(object$X) # maximum DoF (if unconstrained)
class(object)<"tr.smooth" # Give object a class
object
}
Predict.matrix.tr.smooth<function(object,data) {
## prediction method function for the `tr' smooth class
x < data[[object$term]]
x < x  object$x.shift # stabilizing shift
m < object$m; # spline order (3=cubic)
k<object$knots # knot locations
nk<length(k) # number of knots
X<matrix(0,length(x),object$bs.dim)
for (i in 1:(m+1)) X[,i] < x^(i1)
for (i in 1:nk) X[,i+m+1] < (xk[i])^m*as.numeric(x>k[i])
X # return the prediction matrix
}
# an example, using the new class....
require(mgcv)
set.seed(100)
dat < gamSim(1,n=400,scale=2)
b<gam(y~s(x0,bs="tr",m=2)+s(x1,bs="ps",m=c(1,3))+
s(x2,bs="tr",m=3)+s(x3,bs="tr",m=2),data=dat)
plot(b,pages=1)
b<gamm(y~s(x0,bs="tr",m=2)+s(x1,bs="ps",m=c(1,3))+
s(x2,bs="tr",m=3)+s(x3,bs="tr",m=2),data=dat)
plot(b$gam,pages=1)
# another example using tensor products of the new class
dat < gamSim(2,n=400,scale=.1)$data
b < gam(y~te(x,z,bs=c("tr","tr"),m=c(2,2)),data=dat)
vis.gam(b)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.