1. Home
2. VSURFACE procedure

# VSURFACE procedure

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

### Options

 `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

### Parameters

 `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

### Description

`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.

Options: `PRINT`, `PLOT`, `BASIS`, `KNOTS`, `PENALTYMETHOD`, `DEGREE`, `DIFFORDER`, `EXTREME`, `SEEXTREME`, `TYPEEXTREME`, `PREDICTIONS`, `PMETHOD`, `NBOOT`, `NRETRIES`, `SEED`, `CIPROBABILITY`, `COLOURS`.
Parameters: `Y`, `X1`, `X2`, `ESTIMATE`, `SE`, `LEVELS`, `TITLE`, `WINDOW`, `SCREEN`, `EXIT`.

### Method

`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.

Directive: `REML`.
Procedures: `RQUADRATIC`, `THINPLATE`, `TENSORSPLINE`.
Commands for: REML analysis of linear mixed models.

### Example

```CAPTION 'VSURFACE example'; STYLE=major
"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 + \
GRNORMAL(162;0;0.2)

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
PRINT    [RLWIDTH=4] Predictions; FIELDWIDTH=6; DECIMALS=2```

Updated on October 29, 2020