1. Home
  2. VSURFACE procedure

VSURFACE procedure

Fits a 2-dimensional spline surface using REML, and estimates its extreme point (D.B. Baird).


PRINT = string tokens What to print (description, model, components, effects, vcovariance, deviance, waldtests, extreme, confidence, monitoring); default desc, mode, comp, wald, extr
PLOT = string tokens What to plot (contour, surface); default * i.e. nothing
BASIS = string token Spline basis to use (thinplate, pspline, penalizedspline); default thin
KNOTS = scalar, variate or pointer Knots to be fitted in spline model, if a scalar, this is the total number of knots to be fitted; if a variate of length 2, this is the number of knots in the X1 and X2 directions; and if a pointer to 2 variates, these are the values for knots in the X1 and X2 directions; default 16
PENALTYMETHOD = string token Which tensor spline penalty to use (isotropic, semiconstrained, unconstrained); default unco
DEGREE = scalars Degree of polynomial used to form the underlying spline; default 1 for METHOD=penalizedspline and 3 for METHOD=pspline
DIFFORDER = scalars Differencing order for p-spline penalty; default 2
EXTREME = scalars Saves the estimated value of y at the extreme point
SEEXTREME = scalars Saves the standard error of the estimated value of y at the extreme point
TYPEEXTREME = string token Type of extreme to be identified (minimum, maximum); default maxi
PREDICTIONS = matrix or pointer Saves predictions
PMETHOD = string tokens Method of returning predictions (grid, list); default grid
NBOOT = scalars The number of bootstrap samples to estimate standard errors and confidence limits; default 100
NRETRIES = scalars Number of times to retry bootstrap sampling when the REML fit fails; default is the same value as NBOOT
SEED = scalars The seed used to initialize the randomization in the bootstrap sampling; default 0 continues an existing sequence or, if none, selects a seed automatically
CIPROBABILITY = scalar Probability level for confidence intervals for parameter estimates; default 0.95
COLOURS = text or variate Colours for the plots


Y = variates Y-variate to which the spline surface will be fitted
X1 = variates The first X-variate which defines the spline surface
X2 = variates The second X-variate which defines the spline surface
ESTIMATE = variates Estimated value of each x-variate at the extreme point
SE = variates Standard error of the estimated value of each x-variate at the extreme point
LEVELS = scalars, variates or pointers Number of values or values at which to evaluate each X for plots and predictions
TITLE = texts Title to use for graphs; default automatically made from the variate identifiers used for Y, X1 and X2
WINDOW = scalars Window number for the graphs; default 3
SCREEN = string tokens Whether to clear the screen before plotting or to continue plotting on the old screen (clear, keep); default clea
EXIT = scalars Exit code from the REML fit and location of extreme point


VSURFACE fits a spline surface defined by the X1 and X2 parameters to the Y variate, and estimates the extreme point within the region bounded by the values of x-variates. Parameters ESTIMATE and SE can save the estimated value of each x-variate, and their standard errors, at the extreme point. The y-value at the extreme point, and its standard error, can be saved by the EXTREME and SEEXTREME options. The TYPEEXTREME option specifies whether the extreme is a minimum or a maximum.

The BASIS option specifies whether to use thin-plate (the default), p-splines or penalized splines to construct the basis: p-splines or penalized splines are jointly known as tensor splines. Thin-plate splines are 2-dimensional cubic smoothing splines, and are formed using the THINPLATE procedure.

The positions of the knots used in the basis functions are specified by the KNOTS parameter. This can be if a scalar, specifying the total number of knots to be fitted; the procedure will then use equi-spaced knots divided proportionally to the number of distinct points in the two directions. Alternatively, it can be a variate of length 2 specifying the number of equi-spaced of knots in the X1 and X2 directions. Finally, it can be a pointer to 2 variates whose values are used for knots in the X1 and X2 directions.

The degree of polynomial used to form the underlying tensor spline basis functions is specified by the DEGREE option. This has a default of 3 for p-spline models, and 1 for penalized spline models. The DIFFORDER option specifies the differencing order to be used with p-spline models. This determines the strength of the penalty (for a given smoothness parameter). The default is to use second-order differencing. For a p-spline model, the underlying fixed polynomial in each dimension has degree d equal to DIFFORDER-1. For a penalized spline model, the underlying fixed polynomial in each dimension has degree d equal to the value specified by the DEGREE option. The tensor-spline basis is constructed via interactions of the one-dimensional spline bases, as detailed in the TENSORSPLINE procedure.

The PENALTYMETHOD option controls the interaction between the one-dimensional spline bases. An unconstrained penalty (the default) allows a separate smoothing parameter for each term. In this case, the basis pointer has 2d+3 matrices, one for each term. With the semiconstrained penalty, the same smoothing parameter is imposed across the interaction of polynomials in the first dimension with random terms in the second, and for the interaction of random terms in the first dimension with polynomials in the second dimension. An isotropic penalty uses a single common penalty, and the terms are combined into a single matrix.

The PRINT option selects the output to be displayed:

description description of the data and spline basis to be fitted,
model description of model fitted,
components estimates of variance components and estimated parameters of covariance models,
effects estimates of the fixed and random effects,
vcovariance variance-covariance matrix of the estimated components,
deviance deviance of the fitted model (-2 × log-likelihood RL),
waldtests Wald tests for fixed terms,
extreme y and x-values of the extreme fitted value, with
confidence estimated confidence limits of the extreme y and x-values obtained from the bootstrap analysis, and
monitoring monitoring information at each iteration in the REML fitting and for each sample of the bootstrap analysis

The EXIT parameter saves a scalar containing the exit code from REML if the fit failed (-2, -1 or 1…8), or 9 if the extreme is on the boundary of the X1, X2 region (so the optimum may be outside the region), or 10 if the bootstrapping has not found NBOOT successful fits before NRETRIES failures (see below). EXIT will be 0 if an interior optimum has been found and any bootstrapping has been successful.

If standard errors or confidence limits are required, these are formed by bootstrapping the observations. The NBOOT option controls the number of bootstrap samples that are taken. If the REML fit for a sample fails, an extra sample will be taken until a total of NRETRIES samples have failed, in which case the procedure exits with parameter EXIT set to 10. The SEED option controls the randomization seed used for the bootstrapping, and the CIPROBABILITY controls the probability levels of the confidence limits. The value of NBOOT must be large enough that at least one sample falls outside the confidence limits on either side (i.e. NBOOT >= 2/(1 - CIPROBABILITY)).

The PREDICTIONS option can save predictions and fitted values from the fitted spline model. If option PMETHOD=list it saves both of these, while if PMETHOD=grid it saves just the grid of predictions. The LEVELS parameter specifies the values at which to form predictions. This can be a scalar giving the number of equi-spaced grid points between the minimum and maximum of each x-variate, or a variate of length 2 which contains the number of equi-spaced grid points in the X1 and X2 direction, or a pointer to two variates containing the grid points to be used for X1 and X2. The predictions are stored either in a matrix (the default if the structure type is not set) or in a pointer.

The PLOT option specifies which plots to display, with settings:
contour for a contour plot, and
surface for surface plot.

By default nothing is plotted. The COLOURS option specifies a text or variate to define the colours to use. (This is used as the setting of the PENFILL parameter of DCONTOUR and DSURFACE.) The default is a text containing the values ‘darkgreen‘ and ‘yellow‘. The TITLE, WINDOW and SCREEN parameters control the title, window and whether a new plot is started in similar manner to those used in DCONTOUR and DSURFACE. Note that if both surface and contour plots are produced, then SCREEN=keep will cause these to over-plot each other in the same window.



VSURFACE forms the spline basis functions using the THINPLATE or TENSORSPLINE procedures, and fits using REML. The extreme value from the fitted surface (over observations and grid points), is then found. Standard errors and confidence limits are formed by bootstrap resampling of the observations.

Action with RESTRICT

As in REML, either the y-variate or x-variates can be restricted to analyse a subset of the data. If more than one of Y, X1 or X2 are restricted, the restrictions must be consistent.

See also

Directive: REML.
Commands for: REML analysis of linear mixed models.


"Sinusoidal peak with diagonal trend on a 9 x 9 grid with 2 replicates"
VARIATE x1,x2,y; VALUES = !(9(1...9)2),!((1...9)18),!(162(*)); DECIMALS=2(0),2
CALC [SEED=5635] y = 10*SIN(x1/3)*SIN(x2/3) + (x1-4)*(x2-5)/5 + \

VSURFACE [PLOT=surface] y; X1=x1; X2=x2

VSURFACE [PRINT=extreme,confidence; PLOT=contour; EXTREME=max; \
          SEEXTREME=semax; PREDICTIONS=Predictions; NBOOT=20; NRETRIES=10; \
          CIPROB=0.9; SEED=345; COLOURS=!T('red','yellow')] \
          y; X1=x1; X2=x2; ESTIMATE=xpos; SE=sexpos; EXIT=exit; \
          TITLE='Thin-plate spline fit to sinusoidal peak with diagonal trend'
PRINT    max,semax
PRINT    xpos,sexpos




Updated on October 29, 2020

Was this article helpful?