From: <arthur.gilmour_at_DPI.NSW.GOV.AU>

Date: Wed, 7 Nov 2007 12:30:30 +1100

Date: Wed, 7 Nov 2007 12:30:30 +1100

Dear Kevin,

It will be very hard to do the prediction without using the spl()

coefficients used by ASReml.

These are returned in the .res file with standalone ASReml though I can't

see them in the returned asreml() object..

The coding is based on Verbyla et al Applied Stats 1999 48:269-311 and

described in the following tex code (from the Old ASReml reference manaul)

\section{Cubic\index{Cubic, see splines}\label{DET:SPLINE}

splines\index{splines}}

\special{ pdf:outline 2 << /Title (Cubic Splines) /Dest [@thispage /FitH

@ypos]>>}

A recent advance in \AAR\ is the implementation of

cubic smoothing splines as described by Verbyla {\it et al.}\footnote{

Verbyla, A.P., Cullis. B.R., Kenward, M.G. and Welham, S.J. (1999).

Smoothing splines in the analysis of designed experiments

and longitudinal data. {\it Applied Statistics} {\bf 48:} 269-311. }

A cubic spline is a smooth curve over an interval $[a,b]$

formed by linking segments of cubic polynomials at {\sl knot}

points such that the whole curve and it's first and second differentials

are continuous over the interval. See Green \& Silverman\footnote{

Green, P.J. and Silverman, B.W. (1994), Nonparametric Regression and

Generalized Linear models. Chapman \& Hall, 182pp}

for a general treatment of splines.

Verbyla {\em et al.} show that a cubic spline can be fitted as a

linear model of the form

$$\bmy=\X\bmb + \Z_s\bmu_s + \bmeta$$

where $\X\bmb$ fits a linear regression

(intercept and slope), $\Z_s=\Delta(\Delta^T\Delta)^{-1},$ both

$\bmu\;\sim$ N$(0, \sigma^2_s\G)$ and $\bmeta\;\sim$ N$(0, \sigma^2\I)$

are vectors of random effects, and

$\G$ and $\Delta$ are functions of the first differences between the

$N$ ordered

knot points. \AAR\ transforms this model to use $\Z=k^{-1.5}\Z_s\G^{1/2}$

so that the fitted effects are

$\tilde{\bmu}=\G^{-1/2}\tilde{\bmu}_sk^{1.5}$ where

$k$ is the averge distance between knots (range/$(N-1)$)

and $k$ standardizes $\Z$ so the variance component is not poorly scaled.

In \AAR, $\G^{1/2}$ is formed as a cholesky factor (unlike the

Genstat and Splus implementations which use an eigen factorization).

For $h_i=x_{i+1}-x_i$ being the differences between five ordered knot

points,

$$\Delta^T =

{\left( \begin{array}{ccccc}

1/h_1 & -( 1/h_1 + 1/h_2) & 1/h_2 & 0 & 0 \\

0 & 1/h_2 & -( 1/h_2 + 1/h_3) & 1/h_3 & 0 \\

0 & 0 & 1/h_3 & -( 1/h_3 + 1/h_4) & 1/h_4

\end{array}\right) }

$$

and

$$\G =

{\left( \begin{array}{ccc}

{h_1 +h_2\over 3} & {h_2\over 6} &0 \\*[0.15cm]

{h_2\over 6} & {h_2 +h_3\over 3} & {h_3\over 6} \\*[0.15cm]

0 & {h_3\over 6} & {h_3 +h_4\over 3}

\end{array}\right). }

$$

\AAR\ performs the scaling by actually using $h_i=(x_{i+1}-x_i)/k$ in

these expressions

where $k$ is $(x_N-x_1)/(N-1)$.

Note that

\AAR\

returns $\tilde{\bmu}$ (rather than $\tilde{\bmu}_s$) in the solution

vector.

It returns the knot points ($\bmx$) and the

fitted spline curvature ($\Z\tilde{\bmu}=\Z_s\tilde{\bmu}_s$) values

above mini plots of the curvature against the knot points

% and $\bmx$ (correpsonding to $\tilde{\bmu}$ and $\Z$) is presented

%above the spline plots

in the {\tt .asr} file. $\Z$ and $\bmx$ are

also returned in the {\tt .res} file. A {\tt .pin} file can be

used to predict the full spline $\X\hat{\bmb}+\Z\tilde{\bmu}$ (=

$\X\hat{\bmb}+\Z_s\tilde{\bmu}_s).$

The interpolated spline at $x$ between two knot points

$x_l$ and $x_h$ is

\begin{eqnarray*}

\frac{ (x-x_l)y_h + (x_h-x)y_l}{(x_h-x_l)}

- \hfill\hbox{\hspace*{6cm}}\\\frac{(x-x_l)(x_h-x)}{ 6}\left[

\left(1+\frac{x-x_l}{x_h-x_l}\right)d_h +

\left(1+\frac{x_h-x}{x_h-x_l}\right)d_l\right]

\end{eqnarray*}

where $\bmy= \X\hat{\bmb}+\Z\tilde{\bmu}$ is thus the vector of predicted

points,

and for

$$\bmdel = \G^{-1}\tilde{\bmu}_s

= \G^{-1}\Delta^T\Delta(\Delta^T\Delta)^{-1}\tilde{\bmu}_s

= \G^{-1}\Delta^T\Z_s\tilde{\bmu}_s,$$

$\bmd= (0\ \ \bmdel'\ \ 0)'$ is vector of second derivatives at the knot

points.

Extract from .res file showing spl coeffs at knot pouints, and then

interpolated at predict points

spl(x,7) Z matrix for specified KNOT points: Scale factor 0.0061

118.00 1.3072 1.4651 0.9457 0.4209 0.1926

484.00 -0.9158 -0.3924 -0.0087 0.0500 0.0352

664.00 -0.8386 -1.3059 -0.4781 -0.1323 -0.0421

1004.00 -0.3283 -0.9049 -1.3647 -0.4768 -0.1883

1231.00 0.0124 -0.1593 -0.5303 -0.7068 -0.2859

1372.00 0.2240 0.3038 0.1852 -0.1233 -0.3466

1582.00 0.5391 0.9936 1.2509 0.9684 0.6351

spl(x,7) has 5 levels. Z matrix standardized for scale.

118.00000 1.30719 1.46505 0.94568 0.42087 0.19262

150.00000 1.04362 1.31451 0.85826 0.38953 0.17864

200.00000 0.63940 1.07798 0.72211 0.34043 0.15682

etc

Arthur

*><><><><><><><><><><><><><><><><><><><><><><><>
*

Received on Sun Nov 07 2007 - 12:30:30 EST

This webpage is part of the ASReml-l discussion list archives 2004-2010. More information on ASReml can be found at the VSN website. This discussion list is now deprecated - please use the VSN forum for discussion on ASReml. (These online archives were generated using the hypermail package.)