Page 66 - Compact Numerical Methods For Computers

Handling larger problems 55
(n + 1) by (n + g) which initially has all elements in the first n rows equal to zero.
Each of the m observations or rows of (A, b) can be loaded in succession into the
(n + 1)th row of the working array. The Givens’ rotations will suitably fill up the
workspace and create the triangular matrix R. The same workspace suffices for the
orthogonalisation. Note that the elements of the vectors d are automatically left
2
in the last g elements of row (n + 1) of the working array when the first n have
been reduced to zero. Since there are only (m – n) components in each d , but m
2
rows to process, at least n of the values left in these positions will be zero. These
will not necessarily be the first n values.
T
A further feature of this method is that the residual sum of squares, r r, is
equal to the sum of squared terms This can be shown quite easily since

T T
r r = (b – Ax) (b – Ax)
T T T T T T
= b b – b Ax – x A b + x A Ax. (4.21)
By using the normal equations (2.22) the last two terms of this expression cancel
leaving
T
T
T
r r = b b – b Ax. (4.22)
If least-squares problems with large numbers of observations are being solved
via the normal equations, expression (4.22) is commonly used to compute the
T T T
residual sum of squares by accumulating b b, A A and A b with a single pass
through the data. In this case, however, (4.22) almost always involves the
subtraction of nearly equal numbers. For instance, when it is possible to approxi-
T
mate b very closely with Ax, then nearly all the digits in b b will be cancelled by
T T
those in b Ax, leaving a value for r r with very few correct digits.
For the method using rotations, on the other hand, we have
(4.23)
and

(4.24)

by equation (4.12). Hence, by substitution of (4.23) and (4.24) into (4.22) we
obtain
(4.25)
The cancellation is now accomplished theoretically with the residual sum of
squares computed as a sum of positive terms, avoiding the digit cancellation.
The result (4.25) is derived on the assumption that (4.14) holds. In the
rank-deficient case, as shown by k zero or ‘small’ singular values, the vector f in
equation (4.15) can be decomposed so that

(4.26)

where f is of order (n – k) and f of order k. Now equation (4.24) will have the
1 2
form
(4.27)