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

NAG Toolbox: nag_nearest_correlation_h_weight (g02aj)

Purpose

nag_nearest_correlation_h_weight (g02aj) computes the nearest correlation matrix, using element-wise weighting in the Frobenius norm and optionally with bounds on the eigenvalues, to a given square, input matrix.

Syntax

[g, h, x, iter, norm_p, ifail] = g02aj(g, alpha, h, 'n', n, 'errtol', errtol, 'maxit', maxit)
[g, h, x, iter, norm_p, ifail] = nag_nearest_correlation_h_weight(g, alpha, h, 'n', n, 'errtol', errtol, 'maxit', maxit)

Description

nag_nearest_correlation_h_weight (g02aj) finds the nearest correlation matrix, X$X$, to an approximate correlation matrix, G$G$, using element-wise weighting, this minimizes H(GX)F ${‖H\circ \left(G-X\right)‖}_{F}$, where C = AB$C=A\circ B$ denotes the matrix C$C$ with elements Cij = Aij × Bij${C}_{ij}={A}_{ij}×{B}_{ij}$.
You can optionally specify a lower bound on the eigenvalues, α$\alpha$, of the computed correlation matrix, forcing the matrix to be strictly positive definite, if 0 < α < 1$0<\alpha <1$.
Zero elements in H$H$ should be used when you wish to put no emphasis on the corresponding element of G$G$. The algorithm scales H$H$ so that the maximum element is 1$1$. It is this scaled matrix that is used in computing the norm above and for the stopping criteria described in Section [Accuracy].
Note that if the elements in H$H$ vary by several orders of magnitude from one another the algorithm may fail to converge.

References

Borsdorf R and Higham N J (2010) A preconditioned (Newton) algorithm for the nearest correlation matrix IMA Journal of Numerical Analysis 30(1) 94–107
Jiang K, Sun D and Toh K-C (To appear) An inexact accelerated proximal gradient method for large scale linearly constrained convex SDP
Qi H and Sun D (2006) A quadratically convergent Newton method for computing the nearest correlation matrix SIAM J. Matrix AnalAppl 29(2) 360–385

Parameters

Compulsory Input Parameters

1:     g(ldg,n) – double array
ldg, the first dimension of the array, must satisfy the constraint ldgn$\mathit{ldg}\ge {\mathbf{n}}$.
G$G$, the initial matrix.
2:     alpha – double scalar
The value of α$\alpha$.
If alpha < 0.0${\mathbf{alpha}}<0.0$, 0.0$0.0$ is used.
Constraint: alpha < 1.0${\mathbf{alpha}}<1.0$.
3:     h(ldh,n) – double array
ldh, the first dimension of the array, must satisfy the constraint ldhn$\mathit{ldh}\ge {\mathbf{n}}$.
The matrix of weights H$H$.
Constraint: h(i,j)0.0${\mathbf{h}}\left(\mathit{i},\mathit{j}\right)\ge 0.0$, for all i$i$ and j = 1,2,,n$j=1,2,\dots ,n$, ij$i\ne j$.

Optional Input Parameters

1:     n – int64int32nag_int scalar
Default: The first dimension of the arrays h, g and the second dimension of the arrays h, g. (An error is raised if these dimensions are not equal.)
The order of the matrix G$G$.
Constraint: n > 0${\mathbf{n}}>0$.
2:     errtol – double scalar
The termination tolerance for the iteration. If errtol0.0${\mathbf{errtol}}\le 0.0$ then n × sqrt(machine precision) is used. See Section [Accuracy] for further details.
Default: 0.0$0.0$
3:     maxit – int64int32nag_int scalar
Specifies the maximum number of iterations to be used.
If maxit0${\mathbf{maxit}}\le 0$, 200$200$ is used.
Default: 0$0$

ldg ldh ldx

Output Parameters

1:     g(ldg,n) – double array
ldgn$\mathit{ldg}\ge {\mathbf{n}}$.
A symmetric matrix (1/2)(G + GT)$\frac{1}{2}\left(G+{G}^{\mathrm{T}}\right)$ with the diagonal set to I$I$.
2:     h(ldh,n) – double array
ldhn$\mathit{ldh}\ge {\mathbf{n}}$.
A symmetric matrix (1/2) (H + HT) $\frac{1}{2}\left(H+{H}^{\mathrm{T}}\right)$ with its diagonal elements set to zero and the remaining elements scaled so that the maximum element is 1.0$1.0$.
3:     x(ldx,n) – double array
ldxn$\mathit{ldx}\ge {\mathbf{n}}$.
Contains the nearest correlation matrix.
4:     iter – int64int32nag_int scalar
The number of iterations taken.
5:     norm_p – double scalar
The value of H(GX)F${‖H\circ \left(G-X\right)‖}_{F}$ after the final iteration.
6:     ifail – int64int32nag_int scalar
${\mathrm{ifail}}={\mathbf{0}}$ unless the function detects an error (see [Error Indicators and Warnings]).

Error Indicators and 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.

ifail = 1${\mathbf{ifail}}=1$
Constraint: n > 0${\mathbf{n}}>0$.
ifail = 2${\mathbf{ifail}}=2$
Constraint: ldgn$\mathit{ldg}\ge {\mathbf{n}}$.
ifail = 3${\mathbf{ifail}}=3$
Constraint: ldhn$\mathit{ldh}\ge {\mathbf{n}}$.
ifail = 4${\mathbf{ifail}}=4$
Constraint: ldxn$\mathit{ldx}\ge {\mathbf{n}}$.
ifail = 5${\mathbf{ifail}}=5$
Constraint: alpha < 1.0${\mathbf{alpha}}<1.0$.
ifail = 6${\mathbf{ifail}}=6$
On entry, one or more of the off-diagonal elements of h were negative.
ifail = 7${\mathbf{ifail}}=7$
Routine fails to converge in _$_$ iterations.
Increase maxit or check the call to the function.
W ifail = 8${\mathbf{ifail}}=8$
ifail = 999${\mathbf{ifail}}=-999$
Dynamic memory allocation failed.

Accuracy

The returned accuracy is controlled by errtol and limited by machine precision. If ei${e}_{i}$ is the value of norm_p at the i$i$th iteration, that is
 ei = ‖H ∘ (G − X)‖F , $ei = ‖H∘ (G-X) ‖ F ,$
where H$H$ has been scaled as described above, then the algorithm terminates when:
 (|ei − ei − 1|)/( 1 + max (ei,ei − 1) ) ≤ errtol . $|ei-ei-1| 1+ max (ei,ei-1) ≤ errtol .$

Arrays are internally allocated by nag_nearest_correlation_h_weight (g02aj). The total size of these arrays is 15 × n + 5 × n × n + max (2 × n × n + 6 × n + 1,120 + 9 × n)$15×{\mathbf{n}}+5×{\mathbf{n}}×{\mathbf{n}}+\mathrm{max}\phantom{\rule{0.125em}{0ex}}\left(2×{\mathbf{n}}×{\mathbf{n}}+6×{\mathbf{n}}+1,120+9×{\mathbf{n}}\right)$ double elements and 5 × n + 3$5×{\mathbf{n}}+3$ integer elements. All allocated memory is freed before return of nag_nearest_correlation_h_weight (g02aj).

Example

```function nag_nearest_correlation_h_weight_example
g = [ 2, -1,  0,  0;
-1,  2, -1,  0;
0, -1,  2, -1;
0,  0, -1,  2];
h = [ 0.0, 10.0,  0.0,  0.0;
10.0,  0.0,  1.5,  1.5;
0.0,  1.5,  0.0,  0.0;
0.0,  1.5,  0.0,  0.0];
alpha = 0.04;
[g, h, x, iter, norm_p, ifail] = nag_nearest_correlation_h_weight(g, alpha, h);

if (ifail == 0)
fprintf('\nReturned H Matrix\n');
disp(h);
fprintf('Nearest Correlation Matrix\n');
disp(x);
fprintf('Number of iterations taken:     %d\n', iter);
fprintf('Norm value: %26.4f\n', norm_p);
fprintf('Alpha:      %26.4f\n', alpha);

[~, w, info] = nag_lapack_dsyev('n', 'u', x);
fprintf('\nEigenvalues of X\n');
disp(w');
end
```
```

Returned H Matrix
0    1.0000         0         0
1.0000         0    0.1500    0.1500
0    0.1500         0         0
0    0.1500         0         0

Nearest Correlation Matrix
1.0000   -0.9229    0.7734    0.0026
-0.9229    1.0000   -0.7843   -0.0000
0.7734   -0.7843    1.0000   -0.0615
0.0026   -0.0000   -0.0615    1.0000

Number of iterations taken:     66
Norm value:                     0.1183
Alpha:                          0.0400

Eigenvalues of X
0.0769    0.2637    1.0031    2.6563

```
```function g02aj_example
g = [ 2, -1,  0,  0;
-1,  2, -1,  0;
0, -1,  2, -1;
0,  0, -1,  2];
h = [ 0.0, 10.0,  0.0,  0.0;
10.0,  0.0,  1.5,  1.5;
0.0,  1.5,  0.0,  0.0;
0.0,  1.5,  0.0,  0.0];
alpha = 0.04;
[g, h, x, iter, norm_p, ifail] = g02aj(g, alpha, h);

if (ifail == 0)
fprintf('\nReturned H Matrix\n');
disp(h);
fprintf('Nearest Correlation Matrix\n');
disp(x);
fprintf('Number of iterations taken:     %d\n', iter);
fprintf('Norm value: %26.4f\n', norm_p);
fprintf('Alpha:      %26.4f\n', alpha);

[~, w, info] = f08fa('n', 'u', x);
fprintf('\nEigenvalues of X\n');
disp(w');
end
```
```

Returned H Matrix
0    1.0000         0         0
1.0000         0    0.1500    0.1500
0    0.1500         0         0
0    0.1500         0         0

Nearest Correlation Matrix
1.0000   -0.9229    0.7734    0.0026
-0.9229    1.0000   -0.7843   -0.0000
0.7734   -0.7843    1.0000   -0.0615
0.0026   -0.0000   -0.0615    1.0000

Number of iterations taken:     66
Norm value:                     0.1183
Alpha:                          0.0400

Eigenvalues of X
0.0769    0.2637    1.0031    2.6563

```