Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

Chapter Contents
Chapter Introduction
NAG Toolbox

# NAG Toolbox: nag_opt_lsq_uncon_covariance (e04yc)

## Purpose

nag_opt_lsq_uncon_covariance (e04yc) returns estimates of elements of the variance-covariance matrix of the estimated regression coefficients for a nonlinear least squares problem. The estimates are derived from the Jacobian of the function $f\left(x\right)$ at the solution.
This function may be used following any one of the nonlinear least squares functions nag_opt_lsq_uncon_mod_func_comp (e04fc), nag_opt_lsq_uncon_mod_func_easy (e04fy), nag_opt_lsq_uncon_quasi_deriv_comp (e04gb), nag_opt_lsq_uncon_mod_deriv_comp (e04gd), nag_opt_lsq_uncon_quasi_deriv_easy (e04gy), nag_opt_lsq_uncon_mod_deriv_easy (e04gz), nag_opt_lsq_uncon_mod_deriv2_comp (e04he) or nag_opt_lsq_uncon_mod_deriv2_easy (e04hy).

## Syntax

[v, cj, ifail] = e04yc(job, m, fsumsq, s, v, 'n', n)
[v, cj, ifail] = nag_opt_lsq_uncon_covariance(job, m, fsumsq, s, v, 'n', n)

## Description

nag_opt_lsq_uncon_covariance (e04yc) is intended for use when the nonlinear least squares function, $F\left(x\right)={f}^{\mathrm{T}}\left(x\right)f\left(x\right)$, represents the goodness-of-fit of a nonlinear model to observed data. The function assumes that the Hessian of $F\left(x\right)$, at the solution, can be adequately approximated by $2{J}^{\mathrm{T}}J$, where $J$ is the Jacobian of $f\left(x\right)$ at the solution. The estimated variance-covariance matrix $C$ is then given by
 $C=σ2JTJ-1, JTJ nonsingular,$
where ${\sigma }^{2}$ is the estimated variance of the residual at the solution, $\stackrel{-}{x}$, given by
 $σ2=Fx- m-n ,$
$m$ being the number of observations and $n$ the number of variables.
The diagonal elements of $C$ are estimates of the variances of the estimated regression coefficients. See the E04 Chapter Introduction, Bard (1974) and Wolberg (1967) for further information on the use of $C$.
When ${J}^{\mathrm{T}}J$ is singular then $C$ is taken to be
 $C=σ2 JTJ †,$
where ${\left({J}^{\mathrm{T}}J\right)}^{†}$ is the pseudo-inverse of ${J}^{\mathrm{T}}J$, and
 $σ2 = F x- m-k , k=rank ​J$
but in this case the argument ifail is returned as nonzero as a warning to you that $J$ has linear dependencies in its columns. The assumed rank of $J$ can be obtained from ifail.
The function can be used to find either the diagonal elements of $C$, or the elements of the $j$th column of $C$, or the whole of $C$.
nag_opt_lsq_uncon_covariance (e04yc) must be preceded by one of the nonlinear least squares functions mentioned in Purpose, and requires the arguments fsumsq, s and v to be supplied by those functions (e.g., see nag_opt_lsq_uncon_mod_func_comp (e04fc)). fsumsq is the residual sum of squares $F\left(\stackrel{-}{x}\right)$ and s and v contain the singular values and right singular vectors respectively in the singular value decomposition of $J$. s and v are returned directly by the comprehensive functions nag_opt_lsq_uncon_mod_func_comp (e04fc), nag_opt_lsq_uncon_quasi_deriv_comp (e04gb), nag_opt_lsq_uncon_mod_deriv_comp (e04gd) and nag_opt_lsq_uncon_mod_deriv2_comp (e04he), but are returned as part of the workspace argument w (from one of the easy-to-use functions). In the case of nag_opt_lsq_uncon_mod_func_easy (e04fy), s starts at ${\mathbf{w}}\left(\mathit{ns}\right)$, where
 $ns=6× n+2× m+m× n+1+ max1, n × n-1 / 2$
and in the cases of the remaining easy-to-use functions, s starts at ${\mathbf{w}}\left(\mathit{ns}\right)$, where
 $ns=7× n+2× m+m× n+n× n+1 / 2+1+ max1, n × n-1 / 2 .$
The argument v starts immediately following the elements of s, so that v starts at ${\mathbf{w}}\left(\mathit{nv}\right)$, where
 $nv=ns+n .$
For all the easy-to-use functions the argument ldv must be supplied as n. Thus a call to nag_opt_lsq_uncon_covariance (e04yc) following nag_opt_lsq_uncon_mod_func_easy (e04fy) can be illustrated as
``` .
.
.
[x, fsumsq, user, ifail] = e04fy(m, lsfun1, x);
.
.
.
ns = 6*n + 2*m + m*n + 1 + max((1,(n*(n-1))/2);
nv = ns + n;
[v, cj, ifail] = e04yc(job, m, fsumsq, w(ns:nv-1), w(nv:numel(w)));```
where the arguments m, n, fsumsq and the $\left(n+{n}^{2}\right)$ elements ${\mathbf{w}}\left(\mathit{ns}\right),{\mathbf{w}}\left(\mathit{ns}+1\right),\dots ,{\mathbf{w}}\left(\mathit{nv}+{{\mathbf{n}}}^{2}-1\right)$ must not be altered between the calls to nag_opt_lsq_uncon_mod_func_easy (e04fy) and nag_opt_lsq_uncon_covariance (e04yc). The above illustration also holds for a call to nag_opt_lsq_uncon_covariance (e04yc) following a call to one of nag_opt_lsq_uncon_quasi_deriv_easy (e04gy), nag_opt_lsq_uncon_mod_deriv_easy (e04gz) or nag_opt_lsq_uncon_mod_deriv2_easy (e04hy), except that $\mathit{ns}$ must be computed as
 $ns=7× n+2× m+m× n+ n× n+1 /2+1+ max 1, n × n-1 / 2 .$

## References

Bard Y (1974) Nonlinear Parameter Estimation Academic Press
Wolberg J R (1967) Prediction Analysis Van Nostrand

## Parameters

### Compulsory Input Parameters

1:     $\mathrm{job}$int64int32nag_int scalar
Which elements of $C$ are returned as follows:
${\mathbf{job}}=-1$
The $n$ by $n$ symmetric matrix $C$ is returned.
${\mathbf{job}}=0$
The diagonal elements of $C$ are returned.
${\mathbf{job}}>0$
The elements of column job of $C$ are returned.
Constraint: $-1\le {\mathbf{job}}\le {\mathbf{n}}$.
2:     $\mathrm{m}$int64int32nag_int scalar
The number $m$ of observations (residuals ${f}_{i}\left(x\right)$).
Constraint: ${\mathbf{m}}\ge {\mathbf{n}}$.
3:     $\mathrm{fsumsq}$ – double scalar
The sum of squares of the residuals, $F\left(\stackrel{-}{x}\right)$, at the solution $\stackrel{-}{x}$, as returned by the nonlinear least squares function.
Constraint: ${\mathbf{fsumsq}}\ge 0.0$.
4:     $\mathrm{s}\left({\mathbf{n}}\right)$ – double array
The $n$ singular values of the Jacobian as returned by the nonlinear least squares function. See Description for information on supplying s following one of the easy-to-use functions.
5:     $\mathrm{v}\left(\mathit{ldv},{\mathbf{n}}\right)$ – double array
ldv, the first dimension of the array, must satisfy the constraint if ${\mathbf{job}}=-1$, $\mathit{ldv}\ge {\mathbf{n}}$.
The $n$ by $n$ right-hand orthogonal matrix (the right singular vectors) of $J$ as returned by the nonlinear least squares function. See Description for information on supplying v following one of the easy-to-use functions.

### Optional Input Parameters

1:     $\mathrm{n}$int64int32nag_int scalar
Default: the dimension of the array s and the second dimension of the array v. (An error is raised if these dimensions are not equal.)
The number $n$ of variables $\left({x}_{j}\right)$.
Constraint: $1\le {\mathbf{n}}\le {\mathbf{m}}$.

### Output Parameters

1:     $\mathrm{v}\left(\mathit{ldv},{\mathbf{n}}\right)$ – double array
If ${\mathbf{job}}\ge 0$, v is unchanged.
If ${\mathbf{job}}=-1$, the leading $n$ by $n$ part of v stores the $n$ by $n$ matrix $C$. When nag_opt_lsq_uncon_covariance (e04yc) is called with ${\mathbf{job}}=-1$ following an easy-to-use function this means that $C$ is returned, column by column, in the ${n}^{2}$ elements of w given by ${\mathbf{w}}\left(\mathit{nv}\right),{\mathbf{w}}\left(\mathit{nv}+1\right),\dots ,{\mathbf{w}}\left(\mathit{nv}+{{\mathbf{n}}}^{2}-1\right)$. (See Description for the definition of $\mathit{nv}$.)
2:     $\mathrm{cj}\left({\mathbf{n}}\right)$ – double array
If ${\mathbf{job}}=0$, cj returns the $n$ diagonal elements of $C$.
If ${\mathbf{job}}=j>0$, cj returns the $n$ elements of the $j$th column of $C$.
If ${\mathbf{job}}=-1$, cj is not referenced.
3:     $\mathrm{ifail}$int64int32nag_int scalar
${\mathbf{ifail}}={\mathbf{0}}$ unless the function detects an error (see Error Indicators and Warnings).

## Error Indicators and Warnings

Note: nag_opt_lsq_uncon_covariance (e04yc) may return useful information for one or more of the following detected errors or warnings.
Errors or warnings detected by the function:

Cases prefixed with W are classified as warnings and do not generate an error of type NAG:error_n. See nag_issue_warnings.

${\mathbf{ifail}}=1$
 On entry, ${\mathbf{job}}<-1$, or ${\mathbf{job}}>{\mathbf{n}}$, or ${\mathbf{n}}<1$, or ${\mathbf{m}}<{\mathbf{n}}$, or ${\mathbf{fsumsq}}<0.0$, or $\mathit{ldv}<{\mathbf{n}}$.
W  ${\mathbf{ifail}}=2$
The singular values are all zero, so that at the solution the Jacobian matrix $J$ has rank $0$.
W  ${\mathbf{ifail}}>2$
At the solution the Jacobian matrix contains linear, or near linear, dependencies amongst its columns. In this case the required elements of $C$ have still been computed based upon $J$ having an assumed rank given by ${\mathbf{ifail}}-2$. The rank is computed by regarding singular values $SV\left(j\right)$ that are not larger than $10\epsilon ×SV\left(1\right)$ as zero, where $\epsilon$ is the machine precision (see nag_machine_precision (x02aj)). If you expect near linear dependencies at the solution and are happy with this tolerance in determining rank you should call nag_opt_lsq_uncon_covariance (e04yc) with ${\mathbf{ifail}}={\mathbf{1}}$ in order to prevent termination. It is then essential to test the value of ifail on exit from nag_opt_lsq_uncon_covariance (e04yc).
$\mathbf{\text{Overflow}}$
If overflow occurs then either an element of $C$ is very large, or the singular values or singular vectors have been incorrectly supplied.
${\mathbf{ifail}}=-99$
${\mathbf{ifail}}=-399$
Your licence key may have expired or may not have been installed correctly.
${\mathbf{ifail}}=-999$
Dynamic memory allocation failed.

## Accuracy

The computed elements of $C$ will be the exact covariances corresponding to a closely neighbouring Jacobian matrix $J$.

When ${\mathbf{job}}=-1$ the time taken by nag_opt_lsq_uncon_covariance (e04yc) is approximately proportional to ${n}^{3}$. When ${\mathbf{job}}\ge 0$ the time taken by the function is approximately proportional to ${n}^{2}$.

## Example

This example estimates the variance-covariance matrix $C$ for the least squares estimates of ${x}_{1}$, ${x}_{2}$ and ${x}_{3}$ in the model
 $y=x1+t1x2t2+x3t3$
using the $15$ sets of data given in the following table:
 $y t1 t2 t3 0.14 1.0 15.0 1.0 0.18 2.0 14.0 2.0 0.22 3.0 13.0 3.0 0.25 4.0 12.0 4.0 0.29 5.0 11.0 5.0 0.32 6.0 10.0 6.0 0.35 7.0 9.0 7.0 0.39 8.0 8.0 8.0 0.37 9.0 7.0 7.0 0.58 10.0 6.0 6.0 0.73 11.0 5.0 5.0 0.96 12.0 4.0 4.0 1.34 13.0 3.0 3.0 2.10 14.0 2.0 2.0 4.39 15.0 1.0 1.0$
The program uses $\left(0.5,1.0,1.5\right)$ as the initial guess at the position of the minimum and computes the least squares solution using nag_opt_lsq_uncon_mod_func_easy (e04fy). See the function document nag_opt_lsq_uncon_mod_func_easy (e04fy) for further information.
```function e04yc_example

fprintf('e04yc example results\n\n');

global y t;

% Fit data
n = 3;
m = int64(15);
y = [0.14, 0.18, 0.22, 0.25, 0.29, 0.32, 0.35, 0.39, 0.37, ...
0.58, 0.73, 0.96, 1.34, 2.10, 4.39];

t = [1.0, 15.0, 1.0;
2.0, 14.0, 2.0;
3.0, 13.0, 3.0;
4.0, 12.0, 4.0;
5.0, 11.0, 5.0;
6.0, 10.0, 6.0;
7.0,  9.0, 7.0;
8.0,  8.0, 8.0;
9.0,  7.0, 7.0;
10.0, 6.0, 6.0;
11.0, 5.0, 5.0;
12.0, 4.0, 4.0;
13.0, 3.0, 3.0;
14.0, 2.0, 2.0;
15.0, 1.0, 1.0];

% Initial guess
x = [0.5;  1;  1.5];

% Minimize
[x, fsumsq, fvec, fjac, s, v, niter, nf, user, ifail] = ...
e04fc( ...
m, @lsqfun, @lsqmon, x);

fprintf('The sum of squares is %12.4f\n', fsumsq)
fprintf('at the point\n');
fprintf('%12.4f', x);
fprintf('\n\n');

job = int64(0);

[~, cj, ifail] = e04yc(...
job, m, fsumsq, s, v);

fprintf('Estimated variances of the sample regression coefficients are:\n');
fprintf('%12.4f', cj);
fprintf('\n');

function [iflag,fvecc, user] = lsqfun(iflag, m, n, xc, user)

global y t;

fvecc = zeros(m,1);
for i=1:double(m)
fvecc(i) = xc(1) + t(i,1)/(xc(2)*t(i,2)+xc(3)*t(i,3)) - y(i);
end

function [user] = lsqmon(m, n, xc, fvecc, fjacc, ljc, ...
```
```e04yc example results

The sum of squares is       0.0082
at the point
0.0824      1.1330      2.3437

Estimated variances of the sample regression coefficients are:
0.0002      0.0948      0.0878
```