A new algorithm for solving subset regression problems is described. The algorithm performs a QR decomposition with a new column-pivoting strategy, which permits subset selection directly from the originally defined regression parameters. This, in combination with a number of extensions of the new technique, makes the method a very flexible tool for analyzing subset regression problems in which the parameters have a physical meaning

Verhaegen, M. H.

English

NASA Technical Reports Server

NASA T-..knir?al Memorandum 89442 
~ _ _  ~~~ 
On the Reliable and Flexible 
Solution of Practical Subset 
Regression Problems 
M. Verhaegen 
(NASA-TE-89442) CN T E E  B E l I A E L E  AND N87-24 156 
IIEXIBLE S G L U l l C b  C I  EHACIJCAL S U E S E T  
EEGBESSJCI  P b C B L E L S  (IASA) 7 F A v a i l :  
E93S fiC A02/rrP A01 CSCL 128 Unclas 
81/66 o o a o i w  
March 1987 
National Aeronautics and 
Space Administration 
https://ntrs.nasa.gov/search.jsp?R=19870014723 2020-03-20T11:00:13+00:00Z
NASA Technical Memorandum 89442 
On the Reliable and Flexible 
Solution of Practical Subset 
Regression Problems 
M. Verhaegen, Ames Research Center, Moffett Field, California 
March 1987 
National Aeronautics and 
Space Administration 
Ames Research Center 
Moffett Field, California 94035 
ON THE RELIABLE AND FLEXIBLE SOLUTION OF PRACTICAL 
SUBSET REGRESSION PROBLEMS 
n. Verhaegen' 
Aerospace Engineer 
Ames Research Center, bffett Field, California 
ABSTRACT 
A new algorithm for solving subset regression 
problems is described. The algorithm performs a QR 
decomposition with a new column-pivoting strategy, 
which permits subset selection directly from the 
originally defined regression parameters. This, in 
combination with a number of extensions of the new 
technique, makes the method a very flexible tool 
for analyzing subset regression problems in which 
the parameters have a physical meaning. 
INTRODUCTION 
The subset regression problem analyzed in this 
paper is formulated for the following least-squares 
problem: 
The system matrix 
contain a number of columns that are "nearly" lin- 
early dependent on an independent set of columns 
of A .  The task now is to find the minimal number 
of columns comprising that independent set. Fur- 
thermore, the adjective "nearly" is used to indi- 
cate that small perturbations of A ,  that is, of 
the order of the inaccuracies ( O ( E ) )  on the entries 
of A, establish that linear dependency. 
A E Rmxn (m ? n) is assumed to 
A numerically reliable technique used so far 
to detect near dependencies is that of singular- 
value decomposition (SVD) [ l ] .  This technique 
determines a minimal set of columns. Say that the 
dimension of this set is k. But, each of these 
selected columns now is a linear combination of the 
original columns of A. Hence, that minimal set 
w i l l  generally contain more than k individual 
columns of A. This is a maJor drawback for prac- 
tical subset regression problems, such as aerody- 
namic model identification [21,[31. 
This paper describes a new technique that 
allows the direct selection of k individual col- 
umns of A.  The technique performs a QR decomposi- 
tion with a new column-pivoting strategy. Further- 
more, it is not necessary to compute the full SVD 
of A in order to determine k. This is because 
.accurate estimates of the singular values, which 
allow such a determination, are also obtained. 
This new technique is described first, and 
then a number of extensions and modifications of it 
are presented that allow its application to practi- 
cal problems in subset regression analysis. These 
problems are in the order to be discussed: 
between the defined independent and dependent col- 
umns of A by the new technique. 
2. The inclusion of a priori information 
about the structure of the regression model, as 
well as about estimates of the individual compo- 
nents of 
uncertainties. 
of A that are contaminated by noise and others 
that are noise-free without using scaling 
techniques. 
4. 
related total least-squares problem, described in 
1 .  The determination of the interrelationship 
in ( 1 )  and their corresponding 
3.  The capability to Jointly process columns 
The efficient solution of the closely 
r41. 
NEW COLUHN-PIVOTING STRATEGY IN QR DECOMPOSITION 
The QR decomposition of the matrix A in ( 1 )  
is denoted as 
wbere R c Rnxn 
Q Qo = I . This transformation differs from the 
S8D in &at no right orthogonal transformation on 
A is applied. It is precisely this right- . 
transformation that obscures the selection of k 
individual columns of A. To avoid this problem, 
the only right-transformation on A allowed is a 
column permutation. 
is upper triangular and 
*National Research Council Research Associate 
The idea of colmn pivoting is to exploit the 
Freedom introduced by this column permutation to 
solve the rank-deficient least-squares problem 
(1). Based on the new pivoting strategy, it will 
be shown that this can be done as reliably aa with 
the SVD method without relying on complete SVD. 
Let us clarify this for the rank-one deficient 
case, that is, For k = n - 1, implying that 
an = O ( c )  and > >  on. (For the sake of 
brevity, we restrict without loss of generality to 
this case throughout the whole paper.) 
singular values oi of A are ordered in 
decreasing magnitude. 
Here the 
From [ 5 ]  (for example) it is known that the 
magnitude of the last diagonal element of R in 
(2) is an upper bound for on. Therefore, the 
rank-one deficiency can be revealed with column 
pivoting if we can Find a column permtation 
matrix rrn such that 
c 
where Ir:2a is "small" of O ( E ) .  The existence 
of such a p rmutation is indicated by exercise 
P-6-4-4 of [6]. 
Furthermore this exercise gives a constructive 
way to find this permutation. The key information 
is the computation of the right singular vector 
corresponding to the smallest singular value 
This can be done, for example, by the inverse 
iteration method 161. 
an, 
The decomposition (3) combined with (2) 
produces : 
r 'I 
and 
T Therefore, if we explicitly accumulate Qn (which 
can be done efficiently because we are working only 
in the n-dimensional parameter pace as opposed to 
the sample space) and perform Qo with 
Householder transformations such as described in 
[7], the basic solution x+, the residual 2 ,  and 
the estimate of ( 1 )  can be computed in a manner 
that is completely similar to that in [7]. These 
quantities can now be given as follows: 
4 
(8 )  
where (hn)k denotes the k selected columns 
OF A. 
DETERMINING INTERRELATIONSHIPS IN THE 
REGRESSION MODEL 
The above algorithm rearranges the columns 
Here the of A as [ac1 ... a Ck a Ck+l . . . acn 1. 
vectors [a,, ... ack] correspond to the so-called 
identifiable or independent components of E ,  
denoted previously as 
corresponds to the dependent ones. 
xb, and [ack+l ... a,] 
For a number of applications it is necessary 
In to specify these dependencies more precisely. 
this case, it is necessary to know on which minimal 
set of vectors from [acl ... ack] each vector from 
[ack+l ... a, ] depends. n 
The decomposition obtained in (4) reveals the 
First, focus on the following least-squares 
information necessary to answer this question. 
problem: 
Elements of the solution 
dinates of the projection of acn in the space 
spanned by [acl ... a 1. Of course, "large" 
components in g suggest a strong correlation 
between the corresponding columns of 
[ac1 e.. ac -1 ] and acn. However, it is more 
convenient €0 normalize these coordinates. 
therefore define the Following cosines: 
= R;ir12 are the coor- 
. en- 1 
Let us 
with pi denoting the ith column of R l l  in 
(9). This set OF numbers ranges from -1 to + 1  and, 
therefore, can be interpreted completely anal- 
ogously to the comnonly used cross-correlation 
coefficients [ 5 ] .  Normally, the latter are 
I 
retr ieved from the covariance m a t r i x  o f  the least- 
squares estimates o f  5 i n  ( 1 ) .  Since A here is 
assumed t o  be rank-deficient, t h i s  covariance 
matrix cannot be computed unless use is made o f  
pseudoinverses. 
Based on the fol lowing theorem we can however 
gain an understanding of t h i s  covariance matrix 
because of the inclusion o f  a dependent column. 
Theorem 1: 
IF an upper t r iangular matrix R c R"'" is 
given as 
w i t h  R l l  o f  rank n-1 and rZ2 arb i t rary ,  then 
1 2 
22' wi th I Rilr12 and a = l/r 
- Proof: The proof o f  Theorem 1 could be obtained 
d i rec t l y  from the d e f i n i t i o n  o f  the inverse o f  a 
m a t r i x .  
For the decomposition obtained i n  ( 6 ) .  
Theorem 1 allows us t o  exactly calculate the in f lu-  
ence o f  any r22 on the  estimated va iances o f  x+ 
(given by the diagonal elements o f  (Rl1Rll)-'). f 
When we make the p rac t i ca l  observation that "a 
component o f  E is cal led not i den t i f i ab le  i f  i t s  
corresponding estimated variance becomes unaccepta- 
b ly  large," then another procedure fo r  determining 
the dependencies more precisely corresponds t o  
determining whether an unacceptable increase i n  the 
variances o f  the components o f  xb occurs. 
Both o f  the procedures described above supply 
parameters and insights that have c o m n l y  been 
used to analyze i d e n t i f i a b i l i t y  problems i n  pract i -  
ca l  regression problems. 
INCLUSION OF A PRIORI INFORHATION 
For the physical applications addressed i n  
th i s  paper, information i s  of ten avai lable about 
(1) the structure o f  the regression model, that  is, 
which o f  the components o f  E i n  ( 1 )  have t o  be 
i n  xb, and (2) a p r i o r i  estimates o f  indiv idual  
components o f  x w i th  corresponding variance. For 
the sake o f  brev i ty ,  only the capabi l i t ies  of the 
new procedure i n  dealing wi th  the f i r s t  a p r i o r i  
information source are presented here. 
We can include that  information a f t e r  the 
decomposition given i n  (4) has been obtained. Here 
ue would then use the cosines defined i n  (10) or 
the revealed influences of the dependent terms on 
the estimated variances o f  the independent terms 
o f  E .  This information would allow us t o  i n te r -  
change columns o f  [acl ... ackJ wi th  columns o f  
p e r k l a t i o n  ["C 'n can again be done e f f i c i e n t l y ,  since we 
can remain working on the compressed m a t r i x  struc- 
ture defined i n  (4). Furthermore, information 
could be revealed about the e f fec t  o f  t h i s  i n te r -  
change by computing the change i n  residual  v, 
e f f i c i e n t l y  computed as given i n  ( 7 ) .  
... a 1 .  Algorithmically, the necessary 
JOINTLY PROCESSING NOISE-FREE AND NOISY 
COLUMNS I N  A 
I n  many regression models there is a so-called 
o f f se t  term. This resul ts  i n  a column o f  A of 
a l l  ones. The other columns might resu l t  from 
observations and are, therefore, contaminated by 
errors, of ten of d i f f e ren t  magnitudes fo r  each 
column. A rough pa r t i t i on ing  o f  the A-matrix i n  
( 1 )  is therefore 
A = [  A2 1 (12) 
\ 
Noise-free Contaminated by errors 
t J  'I 
J 
The described algorithm requires, as the SVD, 
a threshold a t o  determine the rank o f  A i n  ( 1 )  
[51.  The use o f  such a s ingle o requires column 
scal ing [ 5 ] ,  so that the magnitude o f  the errors on 
the columns o f  t h i s  scaled m a t r i x  i s  o f  the same 
order. The mixing as given i n  (12.) resul ts  i n  a 
very bad conditioning o f  t h i s  scaled matrix and, 
therefore, should be avoided. The new algorithm 
might handle the s i t ua t i on  i n  (12 )  without scaling 
- 
when modified i n  the fo l lowing stages. 
Stage 1: 
Process Al  by the algorithm given 
second section to  produce 
This investigates the dependencies among 
free columns and, therefore, a threshold 
i n  the 
(13) 
the noise- 
- 
o 1  = O(machine precision) should be used. 
continuing to  the second stage, check whether 
columns o f  A2 are dependent on the columns 
represented i n  R l l  i n  (13). This corresponds t o  
checking whether 
Before 
as outlined in [SI. 
ith column of R in (13) and ai  the corre- 
sponding standard2aeviation of the errors on that 
column. -This procedure reduces R 13, R23 in (13) 
to R13,R23, respectively. 
Stane 2: 
Here t 3 ( i )  denoted the 
Applx the algorithm given in the second sec- 
tion to R23: 
using a threshold G2 = 0 (magnitude of errors on 
A$. 
Remark: 
to the following three-stage procedure, where A 
is partioned as 
This two-stage procedure might be extended 
A = [  A1 I A2 I A3 , 1
u-'7 
Noise-free Mild errors Large errors 
EFFICIENT SOLUTION OF TOTAL LEAST-SQUARES PROBLEM 
This problem is not addressed here; the reader 
is referred to [8] for that discussion. 
CONCLUDING REMARKS 
A new scheme was presented for efficiently 
solving subset regression problems. The scheme is 
as reliable as the singular value decomposition 
technique, but does the subset selection directly 
from the defined set of regression parameters. 
Furthermore, a number of interesting extensions and 
modifications of this new procedure have been pre- 
sented thrt allow the flsxlble solutton of subset 
regression problem where the regression parameters 
are of physical Interest. 
REFERENCES 
( 1 1  G.H. Golub, V. Klma, and C.W. Stewart, "Rank 
Degeneracy and Least Squares Problems," 
TR-456, Department of Computer Science, U.  of 
Maryland, College Park, Md., 1976. 
[2] M.H. Verhaegen, "A New Class of Algorithms in 
Linear System Theory: 
Real-Time Aircraft Model Identification," 
Ph.D. dissertation, Catholic University, 
Leuven, Belgium, Nov. 1985. 
With Applications to 
[3] V .  Klein, J.G. Batterson, and P.C. Murphy, 
"Determination of Airplane Structure from 
Flight Data by Using Modified Stepwise Regres- 
sion." NASA TP-1916, 1981. 
[4] C.H. Colub and C.F. Van Loan, "An Analysis of 
the Total Least Squares Problem," SIAM Journal 
on Numerical Analysis, Vol. 1 1 ,  1980, 
pp. 472-479. 
[5] C.L. Lawson and R.J. Hanson, Solving Least 
Squares Problems, Prentice-Hall, Englewood 
Cliffs, N.J., 1974. 
[6] C.H. Colub and C.F. Van Loan, Matrix Computa- - tions, The Johns Hopkins University Press, 
Baltimore, 1983. 
[7] J.J. Dongerra, C. Moler, J. Buch, and 
C. Stewart, LINPACK User's Guide, SIAM, Phila- 
delphia, 1979. 
[8] M.H. Verhaegen, ."The Minimal Residual OR 
Decomposition for Reliably Solving Rank Defi- 
cient Least Squares Problems," submitted to 
International Symposium on Mathematical Theory 
of Networks and Systems, Phoenix, June 1987. 
Report Documentation Page 
I 15. Supplementary Notes 
Poin t  of con tac t :  M. Verhaegen, Ames Research Center ,  MS 210-9, Moffet t  F i e l d ,  
C a l i f o r n i a  94035 (415) 694-5983 o r  FTS 446-5983 
1. Report No: 
iUASA TM 89442 
17. Key Words (Suggested by Author(s)) 
2. Government Accession No. 
18. Distribution Statement 
I 
4. Title and Subtitle 
On t h e  Re l i ab le  and F lex ib l e  Solut ion of 
P r a c t i c a l  Subset Regression Problems 
19. Security Classif. (of this report) 
Unclas s i f i ed  
7. Author(sl 
M. Verhaegen 
3. Security Classif. (of this page) 
Unclassif ied 
9. Performing Organization Name and Address 
Ames Research Center 
Moffet t  F i e ld ,  CA 94035 
12. Sponsoring Agency Name and Address 
National  Aeronautics and Space Administration 
Washington, DC, 20546 
3. Recipient's Catalog No. 
5. Report Date 
March 1987 
6 .  Performing Organization Code 
8. Performing Organization Report No. 
A87158 
10. Work Unit No. 
505-66-11 
11. Contract or Grant No. 
13. Type of Report and Period Covered 
lkchni  cal  Memo r and um 
14. Sponsoring Agency Code 
Subset r eg res s ion  
QR f a c t o r i z a t i o n  
Column p ivo t ing  
Unlimited - Unclas s i f i ed  
Subjec t  Category: 66 
21. No. of pages 22. Price 2- 
I I I I 
NASA FORM 1626 OCT 06 


On the reliable and flexible solution of practical subset regression problems

Abstract

Similar works

Full text

Available Versions

NASA Technical Reports Server