Calculates singular value decompositions of matrices.
|Printed output required (
||Matrices to be decomposed|
||Left-hand matrix of each decomposition|
||Singular values (middle) matrix|
||Right-hand matrix of each decomposition|
Suppose that we have a rectangular matrix A with m rows and n columns, and that p is the minimum of m and n. The singular value decomposition can be defined as
m An = mUp p Sp pVn′
The diagonal matrix S contains the p singular values of A, ordered such that
s1 ≥ s2 ≥ … ≥ sp ≥ 0
The matrices U and V contain the left and right singular vectors of A, and are orthonormal:
U′U = V′V = Ip
The smaller of U and V will be orthogonal. So, if A has more rows than columns, m>n, p=n and VV′=Ip.
The least-squares approximation of rank r to A can be formed as
Ar = Ur Sr Vr′
where Ur and Vr are the first r columns of U and V, and Sr contains the first r singular values of A (Eckart & Young 1936).
INMATRIX parameter specifies the matrices to be decomposed. The algorithm uses Householder transformations to reduce A to bi-diagonal form, followed by a QR algorithm to find the singular values of the bi-diagonal matrix (Golub & Reinsch 1971). The other parameters allow you to save the component parts of the decomposition:
RIGHT for U, S and V respectively.
Genstat will decide how many columns and singular values r to store, and will store that number for any of the components that you specify. If none of the matrices in the
RIGHT lists has been declared in advance, the full number of singular values (r=p) is stored; otherwise Genstat sets r to the maximum number of columns contained in any of the matrices. If r<p, the first r singular values will be saved, along with the corresponding columns of singular vectors.
One practical application of the singular value decomposition is to form generalized inverses of matrices. If you use the singular value decomposition you obtain the Moore-Penrose generalized inverse, sometimes called the pseudo-inverse, and this is the method used by the
Eckart, C. & Young, G. (1936). The approximation of one matrix by another of lower rank. Psychometrika, 1, 211-218.
Golub, G.H. & Reinsch, C. (1971). Singular value decomposition and least squares solutions. Numerische Mathematik, 14, 403-420.
Commands for: Calculations and manipulation,
" Genstat example SVD-1: Singular value decomposition of matrices. The SVD directive calculates singular value decompostitions of matrices using Householder transformations to reduce the input matrix to bi-diagonal form, followed by a QR algorithm to find the singular values of the bi-diagonal matrix. For this example, for matrix h, we wish to find u, v and s such that: h = u s v (16x7) (16x7) (7x7) (7x7) We need to declare and read the values of matrix h but u,s and v will be declared and calculated by the SVD directive. " MATRIX [ROWS=16; COLUMNS=7] h READ h 16.5 24.8 106 147 1112 905 494 4.2 13.3 122 90 982 669 954 11.6 24.7 340 242 808 609 645 18.1 34.2 184 293 1668 901 602 6.9 41.5 173 191 1534 1368 780 13.0 35.7 477 220 1566 1183 788 2.5 8.8 68 103 1017 724 468 3.6 12.7 42 28 1457 1102 637 16.8 26.6 289 186 1509 787 697 10.8 43.2 255 226 1494 955 765 9.7 51.8 286 355 1902 1386 862 10.3 39.7 266 283 1056 1036 776 9.4 19.4 522 267 1674 1392 848 5.0 23.0 157 144 1530 1281 488 5.1 22.9 85 148 1206 756 483 12.5 27.6 524 217 1496 1003 739 : PRINT h " Carry out the decomposition, printing the singular values and saving the component parts into u, s and v. Print the components one after the other. " SVD [PRINT=s] h; LEFT=u; SINGULAR=s; RIGHT=v PRINT [SERIAL=yes] s,u,v " To save a subset of the total number of singular values, the diagonal matrix s, must be declared explicitly. From the previous output, it is apparent that singular values 6 and 7 may be ommitted. To select only the first five singular values, declare ss to be a 5x5 diagonal matrix before carrying out the decomposition. Since diagonal matrices are square matrices, Genstat needs only the number of rows to be declared. " DIAGONALMATRIX [ROWS=5] ss " Carry out the decomposition and print the results. " SVD h; LEFT=us; SINGULAR=ss; RIGHT=vs PRINT [SERIAL=y] ss,us,vs