Open Collections

MAXIMUM LIKELIHOOD IDENTIFICATION OF LINEAR DISCRETE-TIME SYSTEMS by mchel^De Glas Tngenieur, Ecole Superieure d" Informatique Electronique et Automatique, France 1973 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENT FOR THE DEGREE OF MASTER OF APPLIED SCIENCE in the Department of E l e c t r i c a l Engineering We accept this thesis as conforming to the required standard THE UNIVERSITY OF BRITISH COLUMBIA July, 1976 0 Michel De Glas, 1976 In p r e s e n t i n g t h i s t h e s i s i n p a r t i a l f u l f i l m e n t o f t h e r e q u i r e m e n t s f o r an a d v a n c e d d e g r e e a t t h e U n i v e r s i t y o f B r i t i s h C o l u m b i a , I a g r e e t h a t t h e L i b r a r y s h a l l make i t f r e e l y a v a i l a b l e f o r r e f e r e n c e and s t u d y . I f u r t h e r a g r e e t h a t p e r m i s s i o n f o r e x t e n s i v e c o p y i n g o f t h i s t h e s i s f o r s c h o l a r l y p u r p o s e s may be g r a n t e d by t h e H e a d o f my D e p a r t m e n t o r by h i s r e p r e s e n t a t i v e s . I t i s u n d e r s t o o d t h a t c o p y i n g o r p u b l i c a t i o n o f t h i s t h e s i s f o r f i n a n c i a l g a i n s h a l l n o t be a l l o w e d w i t h o u t my w r i t t e n p e r m i s s i o n . D e p a r t m e n t o f E I €. chr * <L &\ £ t*- J? / **- ^ ~g- v The U n i v e r s i t y o f B r i t i s h C o l u m b i a 2075 Wesbrook P l a c e Vancouver, Canada V6T 1W5 D a t e T c J v Z l ^ j >Q7f> ABSTRACT / RESUME Abs t r a c t The t h e o r e t i c a l p r o p e r t i e s of the Maximum L i k e l i h o o d e s t i m a t o r , f o r both s i n g l e i n p u t - s i n g l e output and m u l t i v a r i a b l e systems, are considered. New r e s u l t s r e l a t i v e to convergence p r o p e r t i e s of some i d e n t i f i c a t i o n methods of s i n g l e i n p u t - s i n g l e output systems are obtained. A u n i f i e d approach to the Maximum L i k e l i h o o d i d e n t i f i c a t i o n method of m u l t i v a r i a b l e systems i s proposed. Numerical t e s t s on a computer are performed. Resumg Nous considerons l e s p r o p r i e t e s theoriques de l ' e s t i — mateur du Maximum de Vraisemblance dans l e cas de systemes mono-v a r i a b l e s et de systemes m u l t i v a r i a b l e s . Nous e'tudions des methodes d ' i d e n t i f i c a t i o n de systemes monovariables, et de nouveaux r e s u l t a t s r e l a t i f s a l a convergence de ces methodes sont obtenus. Nous proposons une approche g l o b a l e de 1 ' i d e n t i f i c a t i o n des systemes m u l t i v a r i a b l e s au sens du Maximum de Vraisemblance. Cette procedure est i l l u s t r e e par des exemples numeriques. i i i TABLE OF CONTENTS Page INTRODUCTION ' 1 Chapter 1 SINGLE INPUT - SINGLE OUTPUT SYSTEMS . . . . 6 1 Introduction 6 2 Preliminaries 1 3 Moving Average noise model 12 3.1 Introduction 12 3.2 Global minimum points 14 3.3 Local extremum points 16 4 Autoregressive noise model 20 4.1 Introduction . . . 20 4.2 Global minimum points ......< 21 4.3 Local extremum points 25 4.4 Comments 30 5 Conclusion 31 Chapter 2 MULTIVARIABLE SYSTEMS . . . . . . . 32 1 Introduction 32 2 Canonical structures . 34 2.1 Introduction 34 2.2 Preliminaries 35 2.3 Canonical forms 36 3 Canonical input - output equations 42 3.1 Introduction 42 3.2 Noise-free version 43 3.3 Noisy version 55 3.4 Comments • • • • $9 Page 4 Maximum Likelihood estimator 61 4.1 Introduction . 61 4.2 Equation of error 61 4.3 St a t i s t i c a l properties of the estimates . . . 65 5 Minimizing the loss functions 71 5.1 Introduction 71 5.2 Minimization procedures using derivatives . . 72 5.3 Minimization procedures without using derivatives 73 6 Conclusion 75 Chapter 3 EXPERIMENTS . . . 78 1 Introduction 78 2 Experiment 1 78 3 Experiment 2 79 CONCLUSIONS 84 REFERENCES 85 V LIST OF TABLES Table P a g e 3.1 I d e n t i f i c a t i o n r e s u l t s . T h i r d order system (S/N = ») 80 3.2 I d e n t i f i c a t i o n r e s u l t s . T h i r d order system (S/N = 10) . . 81 3.3 I d e n t i f i c a t i o n r e s u l t s . Fourth order system (S/N = 10) 82 3.4 I d e n t i f i c a t i o n r e s u l t s . Fourth order system (S/N =3) • 83 v i ACKNOWLEDGEMENT I would like to express my gratitude to Dr. E.V. Bohn, my supervisor, for his guidance during the research and writing of my thesis. I wish to thank Dr. A.C. Soudack for his helpful suggestions. I am thankful to the Canada Council for the financial support afforded me in the form of a Graduate Scholarship from 1974 to 1976. INTRODUCTION P r e l i m i n a r i e s During the past decade, i n c r e a s i n g a t t e n t i o n has been devoted to the d i f f e r e n t aspects of system i d e n t i f i c a t i o n . However, although these t o p i c s have been discussed i n a m u l t i t u d e of papers and a r a t h e r l a r g e number of survey papers [3,8,9,25] have been p u b l i s h e d , the f i e l d of i d e n t i f i c a t i o n does not appear as a u n i f i e d s u b j e c t . I t i s then of importance to have i n mind the b a s i c con-cepts c h a r a c t e r i z i n g an i d e n t i f i c a t i o n problem. One of the b a s i c i n g r e d i e n t s i n the f o r m u l a t i o n of the i d e n t i f i c a t i o n problem i s the choice of the c l a s s of models [21] (parametric, non parametric) and of the model s t r u c t u r e ( l i n e a r -non l i n e a r , c o n t i n u o u s - d i s c r e t e t i m e , . . . ) . Such choices cannot be done s y s t e m a t i c a l l y and depend on the a p r i o r i knowledge of the case t r e a t e d and on the purpose of the i d e n t i f i c a t i o n . Due to the f a c t t h a t , i n most r e a l i s t i c s i t u a t i o n s , the measurements made on the system considered are corrupted by random di s t u r b a n c e s , one may i n v e s t i g a t e the i d e n t i f i c a t i o n problem through s t a t i s t i c a l methods, l e a d i n g to an e s t i m a t i o n problem. The estima-t i o n problem can be formulated as the choice of an estimator (Least-Squares, Markov, Bayes, Maximum L i k e l i h o o d ) . Such a f o r -mulation makes i t p o s s i b l e to d e r i v e mathematical p r o p e r t i e s of the estimates and provides a r i g o r o u s framework to the f i e l d of i d e n t i -f i c a t i o n . I t must be n o t i c e d that the p r o b a b i l i s t i c i n t e r p r e t a t i o n can be eluded by t a k i n g i n t o account the e f f e c t s of the dist u r b a n c e s from an e m p i r i c a l p o i n t of view. In any s i t u a t i o n ( p r o b a b i l i s t i c or d e t e r m i n i s t i c ) , the i d e n t i f i c a t i o n problem can co n c e p t u a l l y be s t a t e d as the f i n d i n g of the m o d e l — s u b j e c t e d to the same input s i g n a l s as the s y s t e m — t h a t i s o p t i m a l i n some sense. The o p t i m a l i t y has to be defined u s i n g a c r i t e r i o n w i t h respect to the output s i g n a l s , which are r e l a t e d to the parameters values through a f u n c t i o n a l r e l a t i o n s h i p . The c r i t e r i o n i s g e n e r a l l y defined by means of a l o s s f u n c t i o n , l e a d i n g to an o p t i m i z a t i o n problem. The choice of an o p t i m i z a t i o n method ( d e r i v a t i v e - t y p e methods, search methods) i s , once ag a i n , s t r o n g l y r e l a t e d to the problem under study. A l l these aspects are e q u a l l y important steps of the i d e n -t i f i c a t i o n problem which then cannot be considered as a simple e s t i m a t i o n of a set of parameters or as a simple o p t i m i z a t i o n problem. Aim of Thesis The present work i s concerned w i t h i d e n t i f i c a t i o n of l i n e a r d i s c r e t e - t i m e systems. As o u t l i n e d above, the choice of an estimator i s a c r u c i a l step of the i d e n t i f i c a t i o n problem. A very popular choice i s t h a t of the Least-Squares (LS) esti m a t o r . However, i t i s w e l l known th a t the LS estimator i s g e n e r a l l y b i a s e d . In order to overcome the problem of c o r r e l a t e d n o i s e , one may p o s t u l a t e a system w i t h c o r r e l a t e d n o i s e , by i n t r o d u c i n g a noise model. Depending on the a p r i o r i knowledge of the system and of the n o i s e , the f o l l o w i n g e s t i m a t o r s may be chosen: - the Markov estimator f o r which the knowledge of the covariance matrix of the n o i s e i s needed - the Bayes' estimator f o r which the knowledge of the p r o b a b i l i t y d e n s i t y f u n c t i o n of the no i s e and of the parameters values i s needed - the Maximum L i k e l i h o o d estimator f o r which the proba-b i l i t y d e n s i t y f u n c t i o n of the n o i s e i s needed. The assumption of known covariance m a t r i x of the noise or of known p r o b a b i l i t y d e n s i t y f u n c t i o n of the parameters values s e v e r e l y l i m i t s the p r a c t i c a l a p p l i c a b i l i t y of the two f i r s t e s t i m a t o r s . The aim of Chapter 1 i s to analyze, f o r s i n g l e i n p u t - s i n g l e output systems, the mathematical p r o p e r t i e s of the Maximum L i k e l i h o o d estimator. Two new r e s u l t s w i l l be obtained: - i n case of a Moving Average n o i s e model (Astrom [2 ] ) , i t w i l l be shown that the Maximum L i k e l i h o o d estimates may con-verge to wrong va l u e s . Counterexamples to general convergence w i l l be given ( S e c t i o n 3 ) . - i n case of an Autor e g r e s s i v e n o i s e model, i t w i l l be e s t a b l i s h e d that - under s u i t a b l e assumptions -the Maximum L i k e l i h o o d estimates converge to the true values of the paraneters (Section 4 ) . This r e s u l t provides a convergence proof of the w e l l -know Generalized-Least-Squares (GLS) method, pro-posed by Clarke In Chapter 2, m u l t i v a r i a b l e systems w i l l be c o n s i -dered. The i d e n t i f i c a t i o n of m u l t i v a r i a b l e systems cannot be viewed as a simple g e n e r a l i z a t i o n of the s i n g l e i n p u t - s i n g l e output case and, before formulating the i d e n t i f i c a t i o n problem the following considerations must be taken i n t o account : the choice g e nerally made f o r the c l a s s of the models (state-space representation) implies the determination of a s u i t a b l e canbni c a l form and the d e r i v a t i o n of a canonical set of input-output r e l a t i o n s . The problem of f i n d i n g state-space canonical forms "have been i n v e s t i g a t e d by Luenberger[l5 J , Gopinath [ 12 J and K a y n e [ l 6 J . I t can be shown that the use of a s e l e c t o r matrix or of a set of indexes allows one to describe the st r u c t u r e of a canonical model. In Section 2, the r e s u l t s r e l a t i v e to thJS s problem w i l l be reported. Although the problem of d e r i v i n g a set of canonical input-output r e l a t i o n s has been studied by many authors, no general approach i s yet a v a i l a b l e i n l i t e r a t u r e . The methods proposed by Gopinath[15 ], Ackernan [ 1 ], Zuercher L28J or Guidorzi[l3 ] are subject to strong r e s t r i c t i o n s . In Section 3, the s o l u t i o n to t h i s problem in the most general case i s propo-sed and comparison with the previous methods i s e s t a b l i s h e d . The basic d i f f e r e n c e s between the approach presented i n t h i s Section and the previous ones l i e i n the f a c t that i t allows a d i s c r i m i n a t i o n of the input-output data and that i t i s ex-tended to the noisy case. In Section k, using the r e s u l t s of Chapter 1, a Ma-ximum L i k e l i h o o d estimator i s defined and the consistency of the estimates, i n case of c o r r e l a t e d noise, i s proved. In. S e c t i o n 5, the o p t i m i z a t i o n problem i s discussed. The various- types of minimization, algorithms are described and compared. F i n a l l y , i n Chapter 3, experimental r e s u l t s are pre-sented. The t h e o r e t i c a l r e s u l t s of Chapter 1 and Chapter 2 are tested on numerical examples. 6 Chapter 1: SINGLE INPUT - SINGLE OUTPUT SYSTEMS 1 - INTRODUCTION As pointed out i n the introductory Chapter, a type of i d e n t i f i c a t i o n problem may be obtained by embedding i t i n a p r o b a b i l i s t i c framework, leading to an estimation problem. Such a point of view allows a rigorous approach to the f i e l d of i d e n t i f i c a t i o n and c o n s t i t u t e s one of i t s main aspects. The aim of t h i s Chapter i s to i n v e s t i g a t e the mathe-matical properties of the Maximum L i k e l i h o o d estimation method fo r various types of noise model s t r u c t u r e s . ( i ) In Section 3, the Moving Average noise model w i l l be considered, ( i i ) In Section 4, the Autoregressive noise mo-"el w i l l be studied. In both cases, two cl a s s e s of Maximum L i k e l i h o o d e s t i -mators w i l l be i n v e s t i g a t e d and new r e s u l t s r e l a t i v e to the convergence properties of some e x i s t i n g methods w i l l be e s t a -b l i s h e d . 2^  P R E L I M I N A R I E S The c l a s s o f l i n e a r d i s c r e t e - t i m e s i n g l e i n p u t — single o u t p u t systems c o n s i d e r e d i s r e p r e s e n t e d by t h e d i f -f e r e n c e e q u a t i o n ACz" 1 ) y(k) = B ( z _ 1 ) u(k) +• n(k) ( 1 . 1 ) where |y(k), : k=l , . . . . , N | i s the system output Ju(k), k=l,....,NJ i s the system input. |n(k), k = l , . . . , N j i s an additive zero-mean noise and n a ACz"1) = 1 + a.z"1 (1.2) i = l BU"1) = Yl b.sT1 (1.3) i=0 1 It. i a assumed that, the noise process can be expressed as a process driven by a white noise e(k) . n(k) = H(z~1) e(k) (1.4) where HU-lj ~Sl£bi (1.5) D(z~i) CCz" 1 ) = 1 + c A z r ± < 1 . 6 ) i = l DCz" 1) = 1 d - z _ i ( 1 . 7 ) 1 = 1 The Likelihood, function L can then be defined as the probabi-l i t y density function of £(k), where the numbers £(k) are the s o - c a l l e d r e s i d u a l s defined by H ( z _ 1 ) £(k) = A(z~ 1) y(k) - B(z- 1) u(k) ( 1 . 8 ) 2 Assuming that e i s (0,; G~e) gaussian, the r e s i d u a l £ ( k ) i s 2 a sequence of independent and gaussian ( 0 , °~ ) v a r i a b l e s and the log-Likelihood function A takes the form H C (k) - N log tr - N l o g 2 TT ( 1 . 9 ) 2 <rd k=l 6 e The Maximum Likeli h o o d estimation procedure can then be i n -terpreted as the f i n d i n g of a model A ( z _ 1 ) y(k) = BU" 1) u(k) + H ( z - 1 ) G(k) ( 1 . 1 0 ) i n such a way that the lo g - L i k e l i h o o d function A i s maximized. From ( 1 . 9 ) , . i t i s c l e a r that maximizing A Is equivalent to 9 minimizing the l o s s function N" V N(a,b,c,d) = ^ rr !C £ 2(k) " k=l (1.11) where £(k) i s obtained through the equation of error C ( k ) = H-^z" 1) A ( z - 1 ) y(k) - i T ^ z " 1 ) k z " 1 ) u(k) (1.12) In the fo i l o w i n g , there i s an advantage i n put t i n g T £ = T , T T , a b e d T ". where and — I an a ~ . ..a 1 2 n c T = (1.13) 1 2 n d, d0»...d yd. n The Maximum Li k e l i h o o d estimator can thus be defined as follows the Maximum L i k e l i h o o d estimate(s) of p_, say p_, i s ( a r e ) the ab-solute minimum point(s) o f V J J ( E ) E ] L s £ £ S x = • j i J J : V N(pJJ) = min V N ( £ ) P (1.14) However, from equation (1,12), i t i s c l e a r that £ i s non l i n e a r * * *M i n c and d„ Conseauently, the f i n d i n g of p has to be done i t e r a t i v e l y using a search r o u t i n e . I t i s then of importance to know i f V,T(T>) has l o c a l extremum points This leads us to redefine the Maximum L i k e l i h o o d estimator i n the fo l l o w i n g way 1 0 p £ Sp = 'm d V ( £ ) = 0 ( 1 . 1 5 ) m The estimator E± (resp. E £ ) w i l l be co n s i s t e n t i f and only i f S 1 ( r e s p . S 2) = j C l e a r l y S x Q S £ and E £ consistent, i m p l i e s E± consistent. ( 1 . 1 6 ) From a p r a c t i c a l point of view,, i t i s hig h l y d e s i r a b l e that E 2 i s c o n s i s t e n t , E 2 i n c o n s i s t e n t means that, even i f the gl o b a l minimum point of the l o s s function c o i n c i d e s with the true value of the parameter vector £ (E^ c o n s i s t e n t ) , the Maximum L i k e l i -hood estimation procedure may not converge i n t o and may then give a wrong estimate. The equation (1,-k) can be r e w r i t t e n as DU" 1) n(k) = C U x) e(k) - 1 , ( 1 . 1 7 ) which im p l i e s that, the noise process i s modeled as an autore-gressive moving average (ARMA) process. In the next Sections, the f o l l o w i n g cases w i l l be i n v e s t i g a t e d 11 a) MA noise model (D=l), This case has been treated by Astrom {_2. ], Astrom and Soderstrom [ ^  ] and Soderstrom b) AR noise model (0=1).. This choice of noise s t r u c -ture i s the basic ingredient i n the formulation of the Generalized Least-Squares algorithm proposed by Clarke [ 6 K In both cases, the s t a t i s t i c a l properties of the estimates w i l l be considered. I t w i l l be shown that — for a MA' model (Section 3) , the e s t i r e l a t o r i s consistent while i s i n c o n s i s t e n t . for a A 3 noise, model (Section ^ ) , the estimators E^ and E 2 are c o n s i s t e n t . The concept of p e r s i s t e n t l y e x c i t i n g s i g n a l s w i l l be used i n the f o l l o w i n g , A s i g n a l u(k) i s said to be p e r s i s t e n t l y e x c i t i n g of order n [ 2 ] i f N l i m i Yl u(k) u(k+3) = R ( o ) (1.18a) N — • 0 0 k=l u e x i s t s and A n = [a U-j) 1 i=l,....,n (1.18b) Is p o s i t i v e definite.. I f u(k) i s p e r s i s t e n t l y e x c i t i n g of o r d e r n n T (1.19a) n X I h. u ( k - i ) = 0 1=0 i m p l i e s h± = 0 i=0,..».,n (1.19b) 3. MOVING AVERAGE NOISE MODEL 3«~1» I n t r o d u c t i o n If a Movi n g A v e r a g e (MA) n o i s e model i s u s e d , t h e o u t p u t d a t a a r e g o v e r n e d by t h e d i f f e r e n c e e q u a t i o n A(z"" 1) y ( k ) = B ( z " 1 ) u ( k ) + C ( z _ 1 ) e ( k ) (1.20) The model s t r u c t u r e i s t h e n g i v e n by M z ' 1 ) y ( k ) = B ( z 1 ) u ( k ) + C ( z _ 1 ) € ( k ) (1.21) To a p p l y Maximum L i k e l i h o o d (ML) method, i t i s assumed t h a t - A l . A l l t h e p r o c e s s e s a r e e r g o d i c - A2> The p o l y n o m i a l s A ( z ~ 1 ) and C(z~^~) have a l l t h e i r z e r o s i n s i d e t h e u n i t c i r c l e £this c o n d i -t i o n i m p l i e s t h a t t h e s y s t e m (1*20) i s s t a b l e ] " A3 . The p o l y n o m i a l s A ( z _ 1 ) , B ( z - 1 ) and C(z~^~) a r e r e l a t i v e l y prime [ t h i s c o n d i t i o n i m p l i e s t h a t t h e system (1.20) i s c o n t r o l l a b l e either, from u or from e ] . [ T T T"lT a _b c l , say p =|^ a b c J , i s obtained by minimization of the l o s s function VN(£) = ^  YZ £ 2 ( k ) (1.22) k=l where the r e s i d u a l £(k) depends on p v i a (1.21) * -1 -1 £(k) = ^ - r 1 y(k> - ? ( 2 > u(k) (1.23) C ( z _ ± ) C ( z - ± ) A s a f i r s t r e s u l t , the gl o b a l minimum points of the l o s s func-t i o n w i l l be considered and co n d i t i o n s w i l l be given for the estimator to be consistent.. This w i l l be the object of Section 3.2. A s mentioned above, i t i s of importance to know i f there i s a unique l o c a l minimum point of V^(£).. In p r a c t i c e , there are cases f o r which the l o s s f u n c t i o n can have more than one l o c a l minimum. The reasons f o r such d i f f i c u l t i e s are the following — the model structure may not be appropriate — the number of data, N, may be too small — the model structure may have the inherent pro-perty that V^ (_I>) bas several l o c a l minima. In order to avoid the f i r s t two p i t f a l l s , the following addi-tional assumptions w i l l be made - M,It i s assumed that ( 1 . 2 0 ) r e a l l y holds - A5.The asymptotic loss function V(£).= lim V„(£) = E [ £ 2 ] ( 1 . 2 4 ) w i l l be considered instead of V^(£). The case covered by the third explanation w i l l be the main object of Section 3.3., where i t w i l l be shown that the model structure considered implies that,- in general, a unique loc a l minimum cannot be obtained and consequently that i s inconsis-tent. 3 - 2 . Global minimum points Let us define - £ K the global minimum point(s) of V^(p) : V%> = m ^ V£) <l-25) P — £ the global minimum point(s) of V( £) V(p M) = min V(jg) (1.26) Astrora has shown [ 2 ] t h a t £ M = £ (1-27) l i m jpjj = 2n with p r o b a b i l i t y 1 (1.28) N-* oo provided u i s p e r s i s t e n t l y e x c i t i n g o f order n + n,. These two * a D r e l a t i o n s establish, that the g l o b a l minimum p o i n t ( s ) of the l o s s function V N ( j 3 . ) converged) as. N — t o the unique g l o b a l minimum point of the asymptotic l o s s f u n c t i o n V ( ] D ) , , which c o i n c i d e s with the true value of the parameters, ja. This means that, under the assumptions A l - A5, = | £j and that Theorem 1,1 The ML estimator E.^  i s asymptotically c o n s i s t e n t . The r e s u l t s of t h i s Section do not give any informa-t i o n about l o c a l extremum points and thus do not prove that JD converges to j). The convergence p r o p e r t i e s of the ML method, c r u c i a l l y depends on the existence of multiple l o c a l extrema of the l o s s f u n c t i o n . I t must be noticed that, although the ML method has been extensively a p p l i e d to systems whose model i s given- by (1.20), no general r e s u l t r e l a t i v e to l o c a l extremum points i s av a i l a b l e i n l i t e r a t u r e - except i n the case of a pure ARMA pro-cess ( B = 0) [ if j . In the next Section' a c l a s s of counter exam-ples to a unique l o c a l minimum i s given : t h i s i m p l i e s the i n -\ 16 consistency of E^ and the non-convergence of the method d e s c r i -bed above, which i s one the most commonly used i n system iden-t i f i c a t i o n . 3-3. L o c a l extremum points Combining (1.20) and (1.21) and dropping the argu-ments z-"*" and k, for convenience, y i e l d s £ = AB - A3 u + AC. e Cl.29) AC AC Since u and e are assumed to be uncorrelated, the asymptotic l o s s function takes then the form 2V E AB - AB AC 2 AC u + E — e AC which can be rewri t t e n , using Parseval s theorem, as 2 2V = °e 2 7Ti (1-30) AC AC ^ 7 ' 1 ' AC AC where the i n t e g r a t i o n path i s the u n i t c i r c l e and 3* u(z) I s the di s c r e t e power spectrum of u. We can now consider the st a t i o n a r y points of V, i . e . the points which are solutions of ^ = 0 (1.32) L e t u s f i r s t c o n s i d e r t h e s o l u t i o n s o f dV/dc . = 0 , j = l l . . | n . A f t e r s i m p l e c a l c u l a t i o n s , i t c a n be w r i t t e n Dc . 2 i r t / ACC AC . 1A&45 ( i ) A B ^ C . - I , 4> u(z) z ^ - 1 d » 3=1,..., 2 i r i J ACC AC (1.33a) o r aS. J A*C* U J - l , . . . , n c ( 1 . 3 3 b ) where i 2 AC "* * A3-; A'B * * * * * * f<z> = - °t - A C ~ ( A B -A B ) * (1.3 '+a) ACC ACC a * n * A ( z ) = z a A ( z ~ L ) (1.3*»b) B * ( z ) = z b B ( z ~ 1 ) (1.3^+c) C ( z ) = z C C ( z _ 1 ) (1.3ifd) S i n c e t h e i n t e g r a n d o f 3 V / 3 c . has n +n p o l e s ( f ( z ) i s a n a l y -3 a c t i c ) , the e q u a t i o n dV/c£ = 0 has s e v e r a l s o l u t i o n s ' i n £_. 18 L e t u s now c o n s i d e r dV/da and DV/Db. U s i n g ( 1 . 3 0 ) , one may w r i t e B AC z" Ju AB-AB-AC + E A AC A C AC j = l , . . . , n. 3V E 1 - j — z J u u AC j= 0 , . . . , n . ( 1 . 3 5 ) or av .1! i = l a. E 1 n T '3 - i -— z u + E \ z- Je" "C - i " —T z e AC AC AC AC E 1=0 B - j — z °u AC 1 - i — z u ( 1 . 3 6 a ) 2>b n v4 E i = l i z - j u B - i — z u AC i=0 i z ^ u 1 - i — z u c)b . ( 1 . 3 6 b ) The e q u a t i o n *a -21 = 0 SV 3b = 0 ( 1 . 3 7 a ) can t h e n be r e w r i t t e n as M - N ' -If T ( 1 . 3 7 b ) 19 where M,- IT,. F,. £ and r are defined by M. . N . . 1 J = E TB - i — z u [AC E AC 5- z-Ju A AC i z~ Ju C - i — z e A AC P. . = E 1 3 1 - i — z u . C C A AC (1 . 3 9 b ) (1 . 3 9 c ) ( 1 . 3 9 a ) da (1.40a) r = - M ( 0 ) ob (1.40b) For any C, the matrices M and P are p o s i t i v e d e f i n i t e . The equa-t i o n ( 1 . 3 8 ) has then, at l e a s t , one solution i n a and b, for any £. Since ( 1 . 3 3 ) ' has several solutions i n £, the gradients of V may vanish for several values of The points JD s a t i s f y i n g ( 1 . 3 3 ) and ( 1 . 3 8 ) belong to S 2. This means that i s inconsistent.. Remark The above analysis i s v a l i d for any (A,B,C). The non-uniqueness of the estimates i s thus an inherent property of the model considered. 20 We h a v e t h u s shown t h e e x i s t e n c e o f c o u n t e r - e x a m -p l e s t o a u n i q u e l o c a l minimum a n d c o n s e q u e n t l y t o g e n e r a l c o n -v e r g e n c e o f t h e ML m e t h o d i n t h e c a s e o f a MA n o i s e m o d e l . T h i s d r a w b a c k i s a n i n h e r e n t p r o p e r t y o f t h e m o d e l s t r u c t u r e c o n s i -d e r e d . . I n t h e n e x t s e c t i o n , , a s i m i l a r a n a l y s i s f o r t h e A u t o r e g r e s s l v e n o i s e m o d e l c a s e i s p r o p o s e d . 4 . AUT05EGRSSSIVE, NOISE MODEL 4.1. I n t r o d u c t i o n I n t h e c a s e o f a A u t o r e g r e s s l v e ( A R ) n o i s e m o d e l , t h e e q u a t i o n s ( 1 . 8 ) a n d ( 1 . 1 0 ) t a k e t h e f o r m A(z.~ 1 ) y(k) = B ( z -1) u(k) + 1 — e ( k ) 1 \ ( 1 . 4 D cC^ _ x) Mz" 1) y(k) = rkz" 1) u(k) + - 1 e<k) (1.42) S i m i l a r l y w i t h t h e p r e v i o u s c a s e , t h e £ ,. i s o b t a i n e d by m i n i m i z i n g t h e l o f u n c t i o n V N ( P ) = 2 I 51 £ 2 ( k ) (1.43) w h e r e t h e r e s i d u a l £ ( k ) i s g i v e n b y C ( k ) = C ( z - 1 ) A C z " 1 ) y ( k ) - C ( z _ 1 ) B C z " 1 ) u ( k ) (1.44) A s s u m i n g t h e c o n d i t i o n s A 1 - A 5 o f S e c t i o n 3 s t i l l h o l d , we s h a l l a n a l y z e f i r s t , t h e g l o b a l minimum p o i n t s ( S e c t i o n 4-2) a n d s e c o n d , t h e l o c a l e x t r e m u m p o i n t s ( S e c t i o n A 4-3) o f t h e a s y m p t o t i c l o s s f u n c t i o n V ( p ) . T h e a s y m p t o t i c c o n -s i s t e n c y o f t h e e s t i m a t o r s a n d E^ w i l l be p r o v e d . 4.2- - G l o b a l minii ' iUUii p o i n t s T h e p u r p o s e o f t h i s s e c t i o n i s t o p r o v e t h e a s y m p t o t i c c o n s i s t e n c y o f t h e ML e s t i m a t o r E ^ . To do s o , we s h a l l p r o c e e d a s f o l l o w s : — F i r s t , we s h a l l p r o v e t h a t t h e g l o b a l minimum p o i n t j2 o f t h e a s y m p t o t i c l o s s f u n c t i o n V(p_) c o i n c i d e s w i t h t h e t r u e v a l u e o f t h e p a r a m e t e r s v e c t o r , p_ - S e c o n d , we s h a l l show t h a t t h e g l o b a l minimum p o i n t ( o f V , T ( p ) , p „ c o n v e r g e ( s ) t o v a s N—» co . 11 ii C o m b i n i n g (1.41) a n d (1.42) a l l o w s u s t o r e w r i t e t h e e q u a t i o n (I . 4 4 ) a n d bo o b t a i n t h e e q u a t i o n o f e r r o r A A ' A A c ~ AB - A S \ AC , , . C K 22 I t then follows 2V = E A A AM AB - AB 1 , J A C j L (1.46) and 2 V . > E AC AC 2 111 T A C v Z j A C U > z " V (1.47) where the i n t e g r a t i o n path i s the u n i t c i r c l e |z| = 1. We can now postulate the Lemma 1.1 2V > o-: (1.48) The e q u a l i t y holds i f and only i f p_ = p_ (A=A, B-B and C=C), Proof r Let us consider the i n t e g r a l 2 I = 2 i r i AC AC ( z ) - 1 A C U ) dz z > 0 (1.49) Then I = V - 2 2 TTi A C l Z ; Z 2 7 T i dz z (1.50) C l e a r l y , since A ( z _ 1 ) and C ( z 1 ) have a l l .their zeros inside the unit circle,(assumption A 2 ) and A ( 0 ) = C ( 0 ) = 1 , i t follows that 2 3 27Ti AC, v dz A C ( Z ) T A(0) 0 ( 0 ) A ( 0 ) C ( 0 ) =• 1 1 I dz 2 l T i f z = 1 (1.51) and I = V - Cr" > 0 e (1.52) Hence 2 Y > 0 - : (.1.53) The e q u a l i t y h o l d s i f and o n l y i f ( i ) 1 = 0 w h i c h i m p l i e s AC = AC ( i i ) E (1.5if) = 0, w h i c h i s e q u i v a l e n t , p r o v i d e d u i s p e r s i s t e n t l y e x c i t i n g o f o r d e r n b + n c , t o A B = AB (1.55) C o m b i n i n g (1,5*4-) and (1,55) and u s i n g A s s u m p t i o n A3 y i e l d s A- = A, B = B , C = C (1.%) Q , E , .D„ F o r a l l N, V"N(_£) and V(jo) a r e p o l y n o m i a l s . Then, the a s s e r t i o n o f co n v e r g e n c e o f V"N t o V as N-»<*> I m p l i e s t h a t V" N(p) c o n v e r g e s u n i f o r m l y t o V(r>) on e v e r y compact s e t . Hence, c o n s i -d e r i n g t h e e s t i m a t e s i n a compact s e t G a l l o w s u s t o e s t a b l i s h t h e Lemma 1 .2 l i m TL! = v'1 w i t h p r o b a b i l i t y 1 TIT ( 1 . 5 7 ) P r o o f L e t « b e an a r b i t r a r y s m a l l p o s i t i v e r e a l number and be the s e t £ - £ M ( 1 . 5 8 ) S i n c e V i s c o n t i n u o u s , t h e r e e x i s t s fi > 0 s u c h t h a t min V ( J ) > V ( £ M ) + /? G-G ( 1 . 5 9 ) S i n c e V" c o n v e r g e s t o V u n i f o r m l y , t h e r e e x i s t s an i n t e g e r M If s u c h t h a t i f IT>M t h e n (1.60) f o r a l l p 6 G, Thus B i n V M ( £ ) > min v ( £ ) - / 3 > V ( P M ) + 2 / 3 G-G T h i s means t h a t £ N i s i n G^. Q,E..D. (1.61) F i n a l l y , the lemmas 1 . 1 and 1 . 2 y i e l d s t h e T h e o r e m 1.2 r T h e ML e s t i m a t o r i s a s y m p t o t i c a l l y c o n s i s t e n t . 4-3 — L o c a l e x t r e m u m p o i n t s T h e g o a l o f t h i s S e c t i o n i s t o show t h a t t h e e s t i m a -t o r i s a s y m p t o t i c a l l y c o n s i s t e n t , i . e . . t h a t - j T ° ^ ° so, we s h a l l p r o v e t h a t t h e g r a d i e n t s o f V(jo) v / i t h r e s p e c t t o V v a n i s h o n l y a t _p_ =• T ) -I n a r e c e n t p a p e r [23], S o d e r s t r o m h a s shov/n t h e e x i s t e n c e o f c o u n t e r - e x a m p l e s t o a u n i q u e l o c a l minimum o f V, f o r " s m a l l " s i g n a l - t o - n o i s e r a t i o s : t h i s i s a m a j o r d r a w b a c k to t h e i m p l e m e n t a t i o n o f t h i s method.. H o w e v e r , i t w i l l be shown t h a t , t h e i n t r o d u c t i o n o f one a d d i t i o n a l a s s u m p t i o n a l l o w s u s to o v e r c o m e t h i s p r o b l e m * V/e s h a l l t h u s i n v e s t i g a t e t h e c o n d i t i o n s f o r t h e • g r a d i e n t s o f V t o v a n i s h i . e . . ^ 7 = 2 (1.62) D e f i n i n g t h e p o l y n o m i a l s n E = I S = 1 + £ f z ~ k (1.63) k = l tn ^ = ^ = | ^ k z ~ k (1 - 6 4 ) where n = n + n, a b m = n, + n D C confers to V" the form 2 V " = E A -, 2 r F + E A G 6 F B - H A A (1 . 6 5 a ) (1 . 6 5 b ) (1.66) and allows us to rev/rite (1.62) as n da. k=l "~J df 3V I . ^ ~ = o k - J m Sb. k=0 0 c i ^ — r ~ = 0 j=l,-..,n (1.67a) 0=0,....n, (1.67b) n m = ^ > a, . ~ ~ = / \ b 5V Be. k=l ° k - ; 5 df. k=0 5 V 0 5 = 1 , . . .,n ( 1 . 6 7 c ) In the following, i t w i l l be assumed that A6 The polynomials A and C, and B and C are r e l a t i v e l y prime. With t h i s assumption, (1.67) im p l i e s j=l,...,n (1.68a) Df. av _ = 0 j=0,... ,m ( 1 . 6 8 b ) 27 Using ( 1 . 6 6 ) , ( 1 . 6 8 ) can be rew r i t t e n as E i c e + E FB-HA A F B - H A , A -u o r n YL h E - 1 i=l i A C z~^u , - 3 = 0 0 (1 . 6 9 a ) (1 . 6 9 b ) A C n E i = l B - i A 2 U f z ^ u ro E i = 0 - l z u + ^ ( 0 ) = 0 (1.70a) n YL E m i=l K - i A Z U r ^ u l - Z h i = 0 E — - X z u z 3 U + dv dh. Let us define the matrices (1 . 7 0 b ) M l = E - X z \z~i 1 A C 6 _ A C 6 i = l , 3=1, . ,n (1.7D H 2 = B - i jz u £ z- Ju A 1=1, 3 = 1, . ,n .. ,n (1.72) N =- E r - i z u 7 A i=l,....,n 3=1,. . . , m (1.73) E ^ z - 1 u z'^uj i=1, • . . , EI , i = l , . . , , ra (1.7'+) and the vectors f = |VV"?n] U- ? 5 ) „T = 4-(o> ( 1 - 7 7 ) The equation (1.70) takes then the form Defining y i e l d s H 2 (1.78) 0 ^ + M 2) f_ - N h + £ = 0 (1.79a) N T f - P h + r = 0 (1.79b) or eq u i v a l e n t l y ^ M l + M 2 ~ N P _ 1 n T ) £ = N p _ 1 H ~ £ (1.80a) h = P" 1 N T f + P - 1 r (1.80b) v = | u (1.81) [ R y ( i - j ) ] (1.82a) N = | R u v ( i - j ) j (1.82b) N T = [ R v u ( i - j ) ] (1.82c) P = | ^ R u ( i - j ) ] (1.82d) -1 T C l e a r l y , the matrix M N P N i s p o s i t i v e d e f i n i t e . Since Is p o s i t i v e d e f i n i t e , the matrix + M 2 - N N^, as a sum of two symmetric p o s i t i v e d e f i n i t e m a t r i c e s , i s p o s i t i v e d e f i n i t e and t h e n n o n - s i n g u l a r . The e q u a t i o n (1.80) has t h u s a u n i q u e s o l u t i o n , g i v e n by 1 = ( M 1 + M 2 - I f P ~ 1 N T ) ~ 1 ( N P ~ 1 r - q ) ( 1 . 8 3 a ) h = P ~ 1 N - T ( M 1 + H 2 - N P " 1 N T ) - 1 ( N p - l p _ a ) + p - l j . ( 1 . 8 3 b ) T h i s p r o v e s t h a t the g r a d i e n t s o f V v a n i s h f o r one and o n l y one v a l u e o f _p_. From Theorem 1.2, t h i s v a l u e c o i n c i d e s w i t h J D . We have t h u s e s t a b l i s h e d t h a t under a s s u m p t i o n A1-A6, £>2 = j JE] a n c * t h a t Theorem 1.3 The M.L. e s t i m a t o r E^ i s a s y m p t o t i c a l l y c o n s i s t e n t . A d i r e c t consequence o f the above s t u d y i s t h a t -p r o v i d e d the a s s u m p t i o n s A1-A6 a r e v e r i f i e d - the i i L method, i n c a s e o f an AR n o i s e model, i s c o n v e r g e n t . I t must be n o t i c e d t h a t t h e G e n e r a l i z e d . L e a s t - S q u a r e s (GLS) i d e n t i f i c a t i o n a l g o r i t h m , p r o p o s e d by C l a r k e [ 6 ] i s no-t h i n g e l s e t h a n a ML method f o r w h i c h the l o s s f u n c t i o n i s m i n i m i z e d v i a a g r a d i e n t a l g o r i t h m . The theorems 1.2 and 1.3 t h u s e s t a b l i s h the convergence p r o o f o f the GLS method. B e f o r e c o n c l u d i n g we s h a l l p r o c e e d t o make some comments about the v a r i o u s a s s u m p t i o n s w h i c h have been made-/{..Iv, Comments The r e a s o n why a f i n i t e N c a n c a u s e d i f f i c u l t i e s i s t h a t t h e v a l u e s o f VN(_p_)» f o r f i x e d J D, a r e s t o c h a s t i c v a r i a b l e s w h i l e V(_p_), f o r f i x e d £, i s a d e t e r m i n i s t i c f u n c t i o n . The a s s u m p t i o n A6 i s a r e s t r i c t i v e assumption s i n c e no m i n i m i z a t i o n ' a l g o r i t h m c a n g u a r a n t e e t h a t t h i s c o n d i t i o n i s f u l f i l l e d a t e a c h s t a g e o f t h e p r o c e d u r e . I t i s however e a s y , when a s t a t i o n a r y p o i n t i s r e a c h e d , t o t e s t i f i t i s so or n o t . In t h e v/hole s e c t i o n , i t was assumed t h a t t h e n o i s e p r o c e s s may be m o d e l l e d as an p r o c e s s CCz" 1 ) n ( k ) = e ( k ) (1.8*0, T h i s means t h a t t h e r e e x i s t s a n o i s e w h i t e n i n g f i l t e r w h i c h can be r e p r e s e n t e d by I t s t r u n c a t e d i m p u l s e r e s p o n s e | , i = 0 , . . . , n c T h i s a s s u m p t i o n has t h e drawback t h a t t h e r e a r e no s y s t e m a t i c r u l e s f o r c h o o s i n g t h e o r d e r n c o f t h e a u t o r e g r e s s i o n s u c h •that t h e c. s, f o r i > n , c a n be n e g l e c t e d . However, v/e can JL O r e a s o n a b l y a r g u e t he f o l l o w i n g p o i n t s ; . - an a n a l y s i s of the n o i s e can f a c i l i t a t e the choice of the order of the a u t o r e g r e s s i o n - n c can be determined i t e r a t i v e l y by i n c r e a s i n g i t u n t i l the r e d u c t i o n of the l o s s f u n c t i o n i s no longer s i g n i f i c a n t . 5 - CONCLUSION In t h i s Chapter, the ML e s t i m a t o r , f o r l i n e a r s i n g l e i n p u t - s i n g l e output systems, has been i n v e s t i g a t e d . The s t a t i s t i c a l p r o p e r t i e s of the ML estimator have been e s t a b l i s h e d : i t has been shown to depend on the noise model s t r u c t u r e considered. For t h i s purpose, two new r e s u l t s have been obtained - f o r a MA noise model, the ML method may not converge. Counterexamples to general convergence can be found. - f o r an AR noise model, the ML method i s convergent. In the next Chapter, the above study w i l l be extended to the m u l t i v a r i a b l e case where i t w i l l be shown t h a t , once a c a n o n i c a l s t r u c t u r e i s d e f i n e d , the extension i s s t r a i g h t f o r w a r d . Chapter 2: MULTIVARIABLE SYSTEMS i - INTRODUCTION A s o u t l i n e d i n t h e p r e v i o u s c h a p t e r , t h e g e n e r a l d e f i n i t i o n ' o f s y s t e m s i d e n t i f i c a t i o n a l l o w s many d e g r e e s o f f r e e d o m i n t h e f o r m u l a t i o n o f t h e i d e n t i f i c a t i o n p r o b l e m , l e a d i n g t o t a k e i n t o a c c o u n t t h e f o l l o w i n g p o i n t s : - T h e c h o i c e o f a c l a s s o f m o d e l s : i m p u l s e r e s p o n s e - t r a n s f e r f u n c t i o n - s t a t e - s p a c e r e p r e s e n t a t i o n . . - T h e d e t e r m i n a t i o n o f a c l a s s o f i n p u t s i g n a l s - T h e c h o i c e o f t h e e s t i m a t o r - The c h o i c e o f t h e o p t i m i z a t i o n a l g o r i t h m . When f o r m u l a t i n g a n d s o l v i n g s u c h p r o b l e m s , i t i s i m p o r t a n t t o h a v e t h e p a r t i c u l a r t y p e o f t h e s y s t e m t r e a t e d a n d t h e f i n a l g o a l o f t h e i d e n t i f i c a t i o n p r o c e d u r e i n m i n d . T h i s c a n n o t t h e n be done f r o m a p u r e l y m a t h e m a t i c a l p o i n t o f v i e w . I t i s , h o w e v e r , h i g h l y d e s i r a b l e t o a c h i e v e a u n i f i e d a p p r o a c h o f t h e s e p r o b l e m s a n d , by e m b e d d i n g them i n a n a b s t r a c t f r a m e w o r k , o b t a i n g e n e r a l r e s u l t s . ( i ) T h e c h o i c e made f o r t h e c l a s s o f t h e m o d e l s o f m u l t i v a r i a b l e s y s t e m s i s g e n e r a l l y t h a t o f s t a t e - s p s c e r e p r e s e n t a t i o n ^ S u c h a c h o i c e i m p l i e s - t h e d e t e r m i n a t i o n o f t h e s y s t e m o r d e r a n d o f a s u i t a b l e c a n o n i c a l f o r m - the d e r i v a t i o n of a set of input-output r e l a t i o n s . The use of a s e l e c t o r m a t r i x S, as suggested by Gopinath [12], or of a set N of indexes, as shown by Mayne [16] a l l o w s one to d e s c r i b e the s t r u c t u r e of a c a n o n i c a l model. These approaches have the drawback that there i s no method f o r the choice of S or N. However, s i n c e such a method should be a p p l i e d to input-output sequences, t h i s problem cannot be performed i n r e a l i s t i c cases - i . e . when no p r e c i s e a p r i o r i knowledge of the n o i s e c h a r a c t e r i s t i c s i s a v a i l a b l e . The v a r i o u s aspects of the determination of a c a n o n i c a l s t r u c t u r e w i l l be t r e a t e d i n S e c t i o n 2. Although the problem of d e r i v i n g a set of input-output r e l a t i o n s from the c a n o n i c a l state-space r e p r e s e n t a t i o n has been i n v e s t i g a t e d by many authors, no u n i f i e d approach can be found i n l i t e r a t u r e . The methods proposed by Gopinath [12], Ackerman [1], Zuercher [28] or G u i d o r z i [13] are subject to strong r e s t r i c t i o n s . In S e c t i o n 3, i t i s shown how to work out t h i s problem i n the most general case, p e r m i t t i n g the decomposition of the system i n t o subsystems. Moreover t h i s approach i s extended to the n o i s y case i n which a s u i t a b l e noise model i s introduced i n such a f a s h i o n as to preserve the advantages of the decomposition -even i n the case of c o r r e l a t e d n o i s e . ( i i ) In Section k, a Maximun Likelihood estimator i s described and the s t a t i s t i c a l properties of the Maximum Like-lihood estimates are investigated.. ( i i i ) F i n a l l y , i n Section 5, the computational aspects are discussed. Two classes of optimization algorithms are des-cribed and compared. 2 - CANONICAL STRUCTURES 2,1.- Introduction The c l a s s of l i n e a r discrete-time raultivariable systems considered i s represented by the state equations ' z(k+l) = Fz(k) + Gu(k) (2.1a) £(k) = Hz(k) (2.1b) where z ( k ) f R n (state vector), u_(k)CRIn (input v e c t o r ) , ^(k)£R r (out-put vector), F£R n x. R n (state matrix), G £ R n x R m ( c o n t r o l matrix) a n d H£R X-R (observation matrix). I t i s w e l l known (see e,g, [ 15 ] ) that. (F^, G-^ , H^) and (F^, G-,,, H-,) are s i m i l a r - i . e , they lead to the sane transfer function K(z) =. H(zl - H) -"^ — i f there e x i s t s a non-singular p n matrix TER" x R such that F 2 " " l ' " 1 <2.2a> ° 2 = T G 1 (2.2b) "2 - V 1 (2.2c) The non-uniqueness of (F, G, H) r e s u l t s i n a degree of a r b i t r a r i n e s s i n the r e a l i z a t i o n algorithm K — G , H ) , which cannot be removed by simple n o r m a l i z a t i o n . The next sections show that each representation ( F r G, II) i s characterized by a set of i n t e g e r s (n^, n 2 >.-..,n s) and that, i f (n, ,n,,, ,.n ) i s knov/n, the representation (F,. G T H ) of (2 -1 ) i s unique.. 2..2.. P r e l i m i n a r i e s Let r = [ H T (HF) T....- ( H F N - 1 ) T ] T A =.[G FG...- F n _ 1 c - ] (2.3) be, r e s p e c t i v e l y , the o b s e r v a b i l i t y and c o n t r o l l a b i l i t y matrices. (F, G) and (F, H ) are assumed to be, r e s p e c t i v e l y , c o n t r o l l a b l e and observable, which i s equivalent to rank(T) = rank (A ) = n ( 2 . 4 ) ; T Let h^ denote the i - t h rpw of H and l e t I and J denote the sets I = { 1,2 r j - (2 .5 ) J = { J 2 , - . - . J s ) ( 2 , 6 ) where l ^ j ^ ^ r for a l l k = 1,2,... .s, s^Tr. J i s a permutation of any subspace of I . Let T be the non-singular (n x n) matrix defined by T = T . T. ^2 T . 0. h . T —1 h . T F — l h-> X F 1 — x i € J ( 2 . 7 ) S i n c e t h e p a i r ( F , H ) i s o b s e r v a b l e , T n o n - s i n g u l a r i m p l i e s ( n ± + l ) = s ( 2 . 8 ) i € J T = P where P i s a p e r m u t a t i o n m a t r i x ( 2 . 9 ) From ( 2 . 7 ) and ( 2 . 9 ) i t i s . c l e a r t h a t t h e t r a n s f o r m a t i o n m a t r i x • i x T i s c o m p l e t e l y d e f i n e d by e i t h e r t h e s e t j n . ,n. , . . . , n . \ ' 3 1 J 2 ° s J or t h e m a t r i x . P, T h i s m a t r i x i s n o t h i n g e l s e t h a n t he s e l e c t o r m a t r i x , w h i c h i s t h e b a s i c c o n c e p t o f t h e a p p r o a c h s u g g e s t e d by G o p i n a t h [ l 2 ] , A l t h o u g h t h e s e two d e f i n i t i o n s o f T a r e e q u i v a l e n t , t h e s e t o f i n d e x e s f n. ,n. , . . . , n . 1,rather t h a n P, w i l l be ' J l J 2 J s > u s e d i n t h e f o l l o w i n g . 2 . 3 - C a n o n i c a l forms D e f i n i n g t h e new s t a t e v e c t o r x ( k ) - T z ( k ) y i e l d s t h e s t a t e - s p a c e e q u a t i o n s x ( k + l ) = A x ( k ) + B u ( k ) ( 2 a 0 a ) y_(k) = C x ( k ) ( 2 . 1 0 b ) where (A, B, C) i s o b t a i n e d from ( F , G, H) t h r o u g h ( 2 . 2 ) . 37 We can now e s t a b l i s h the Lemma 2 . 1 : Let T - a. — l be the i - th rov/ of A . - b . T — i be the i - th row of B - c . T — l be the i - th row of C T — e. — l be the i - th u n i t vec tor i ' i = 1 , 2 , . s, the indexes i - 1 ^ = 1 v. = n . + i 1 1 k = l J k ' Then (i> % T = <L±+l i £ { l , 2 , . " . . , n J - j % 1 1 , v 2 , . . . ^ s J ( 2 . 1 1 a ) a . r = a . . T i £ j W " - ' 9 s ! ( 2 . 1 1 b ) —x —x T T ( i i ) c . = e . . i = 1 , 2 , . . . , s ( 2 . 1 2 a ) - i . " -v. J i i T T c . = c . — X — X i € Jl,2,...,rJ - j j ± . , . . . J s ) ( 2 . 1 2 b ) Before i n v e s t i g a t i n g these r e s u l t s f u r t h e r , i t i s necessary to prove the Le mraa 2 . 2 : The t r i p l e t ( A , B, C ) , d e f i n e d by ( 2 . 1 1 ) - ( 2 . 1 2 ) i s u n i q u e . Proof : Let us assume there e x i s t two represe 'ntat ions ( A ^ , B ^ , C and ( A 2 > B 2 , C 2 ) of ( 2 . 1 0 ) 38 r - F T 1 2 1 Let and r*2 be the corresponding o b s e r v a b i l i t y matrices. Since the representations are s i m i l a r , (2.13) and there e x i s t s a non singular permutation matrix P such that r I. P r i = n M , P T " 2 -•where I i s the n x n u n i t matrix, n (2 . -14) Combining (2,13) and (2.14) y i e l d s T PF± = Pr £ T = M 2T and T I . Q E D n * ^ (2.15) The lemma 2;.l and 2.2 e s t a b l i s h that the knowledge of a set .Cvf indexes j n^,n2,.. .n } implies the uniqueness of the t r i p l e t (A, B, C) c h a r a c t e r i z i n g the state representation of a multivariable system. In the following, we can assume, without.loss of g e n e r a l i t y , that 0\ = i i = 1 s (2.16) i . e . J = |1,2,...sj l < s < r (2.17) 39 Hence, the canonical form of the t r i p l e t (A, 3, C) , given i n ( 2 . 1 1 ) and ( 2 . 1 2 ) can be rewritten.. The matrices A, B and C exhi-b i t the following structure ( i ) A matrix n^+1 A = n^+l XL (n.+l) = n 1 n +1 s ( 2 . 1 8 ) • n where A ^ , i = l , 2 , . r s , i s the £(n^+l) x n J matrix n 0 0 . . . . 0 1 0 0 0 0 0 0 1 0 a. —x T n ± + l ( 2 . 1 9 ) ( i i ) B matrix B = Bn B s «-m,-n^+1 n^+l n n +1 s 4-( 2 . 2 0 ) v/here , i = 1 , 2 , ,s, i s the (n i+l) x m matrix B_. = T x T ^.+1 x T b x — -o. +n. l l n.^+1 ( 2 . 2 1 ) ( i i i ) C matrix T T c T — r - s r - s ( 2 , 2 2 ) n The state .equations take then the form 41 (k+1) 1 x v 1 + n i ( k + 1 ) x v (k+1) 2 x , A (k+1) V 2+n 2 •x (k+1) s x v +n ( k + 1 ) S 8 A x v (k) x.v (k) 1 n l *v (k) ( k ) 2 2 *S> Ck) s X V + n < k ) s s — ^ „ 1 ~ V n i — 2 ^ T + b ~ V n 2 • b T — s b T — V +n s s u(k) ( 2 . 2 3 a ) Za 0 0 Z 2 (k) T •v.. l J i = 1 s x (k) ( 2 . 2 3 b ) where i ^ k ) = [y±(k) y 2 ( k ) . . . y s ( k ) ] T ^ 2 < k ) = t y s + i ( k ) y s + 2 ( k ) . . . y r ( k ) ] T (2.24a) (2.24b) Having obtained a c a n o n i c a l s t r u c t u r e f o r the matrices A, B and C, a set of input-output equations can now be deduced from (A, B, C). Remarkable f e a t u r e s of t h i s f o r m u l a t i o n are that i t e x h i b i t s the s t r u c t u r e of the system and that i t shows the l i n k between the i n t e r n a l (state-space) r e p r e s e n t a t i o n and the e x t e r n a l (input-output) equations. T h i s extends the r e s u l t s of Gopinath [12], Zuercher [28] and G u i d o r z i [ 1 3 ] . 3 - CANONICAL INPUT-OUTPUT EQUATIONS 3.1 I n t r o d u c t i o n A d i r e c t consequence of the s t r u c t u r e of the matrices A and C, d e r i v e d i n the above S e c t i o n , i s that the o r i g i n a l system has been decomposed i n t o s int e r c o n n e c t e d subsystems. The main goal of t h i s S e c t i o n i s to d e r i v e an input-output r e -p r e s e n t a t i o n of the system such that the subsystems can be se p a r a t e l y i d e n t i f i e d . S e c t i o n 3.2 w i l l be concerned w i t h the n o i s e - f r e e c a s e w h i l e t h e n o i s y c a s e w i l l be t r e a t e d i n S e c t i o n 3-3. 3.2.. N o i s e - f r e e v e r s i o n A c c o r d i n g t o (2.23) t h e o u t p u t v e c t o r £ ( k ) c a n be d i v i d e d i n t o two c o m p o n e n t s ^ ( k ) a n d ^ ( k ) , £ ^ 0 0 £ r e s p . y_p(k)] b e i n g r e l a t e d t o x(k) t h r o u g h a r e l a t i o n whose p a r a -m e t e r s a r e known [ r e s p . u n k n o w n ] . T h e a i m o f t h i s s e c t i o n i s t h e d e r i v a t i o n o f an i n p u t - o u t p u t e q u a t i o n u i n v o l v i n g t h e p a r a m e t e r s o f t h e m a t r i c e s A a n d B a n d t h e n a n i n p u t - o u t p u t e q u a t i o n u_ »-2L2 C O n^ a^ n^- nS ^ n e u n k n o w n p a r a m e t e r s o f C. 3.2.1. p a r a m e t e r s o f A a n d 3 F r o m (2.19) a n d (2.21), i t f o l l o w s t h a t t h e i - t h s u b s y s t e m , i = 1,2,... , s , i s d e s c r i b e d by t h e s e t o f e q u a t i o n . 1 xv C k ) = x v + 1 ( k ) + b / u ( k ) i 1 i *v + 1 ( k ) = x y + 2 ( k ) + b / + u ( k ) i i i V + n . - l < k > = X v . + n ( k ) + * V n , - l * ( k ) X X X X X X x v + n ( k ) = a . T x ( k ) + b „ T + n u ( k ) i i i i (2 C o m b i n i n g ( 2 . 2 2 ) a n d ( 2 . 2 5 ) l e a d s t o x:v (k) = y (k) i x (k) = zy (k) - b T u(k) i i x v _ + 2 ( k ) = z 2 y ; L ( k ) - D ^ T + 1 u . ( k ) . - b v T z u ( k ) x v + n ( k ) = ^ V 1 0 ~ b7+n -1 u ( k ) J ^ 1 1 1 " 1 H.(k) i i i i i ( 2 . 2 6 ) Remark I t must be emphasized that i s the time delay operator, defined by z~ J f(k) = f ( k - j ) Then, since x(k) = ^ • + 1 1 ( k ) ' Vn 1 I ' x.. (k) x , A . (k) v +n s s the whole state vector can be expressed as ( 2 . 2 7 ) x(k) =. H(z) Z l ( k ) - N(z) u(k) ( 2 . 2 8 ) 45 w h e r e M(z) = 1 1 1 | 1 z 1 1 1 1 1 0 1 . . . | 0 Z | 1 1 1 [ 1 1 z • ° i • • 1 i | L_ i 1 l 1 | i l j z o | i 0 j *-1 l 1 1 ! 2 n= (n x s) ( 2 . 2 9 ) and 0 1 T n-i-1 m b V 1 + n 1 - l + Z b V n i - 2 + —• + Z b V x 0 b T "*2 N(z) = T Z ^ 2 + n 2 - 2 + * n 2 - l T + 2 V 2 • * 0 b T • * v + n -1 + s s Z^V +-n -2 + ' ' s s -• • ."'"V s (2.30) The s u b s t i t u t i o n of (2.28) i n t o (2.10a) leads to the input-output r e l a t i o n [ ( z l - A) H-J^Oc) = [(zl -A) IT + B ] u(k.) (2.3D In (2.3D, howovjr, only the ( v 2 - l ) t h , (N), - l ) t h , .... (\^ s - l ) t h , n-th equations are s i g n i f i c a n t , the remaining ones 47 b e i n g i d e n t i t i e s . By r e m o v i n g t h e s e i d e n t i t i e s , ( 2 . 3 D t a k e s t h e form P ( z ) v x ( k ) = Q ( Z ) U ( k ) where P and Q a r e p o l y n o m i a l m a t r i c e s i n z Q(z) = P ( z ) = ^ l jCz ) j i 1, . • > s1, .. , s 1 , . . , s 1, . . ,10 ( 2 . 3 2 ) ( 2 . 3 3 a ) ( 2 . 3 3 b ) The p o l y n o m i a l s P . . ( z ) and 0, .(z) a r e o b t a i n e d by s i m p l e i n s -p e c t i o n o f ( 2 . 3 D - I t f o l l o w s i j (z) = </, , zn j + 1 n. a. , z. - a. _ z - a. r , V l x,». i = 1,.. . , s j = 1,. .. , s (2 . 3 4a ) Q* . ( z ) =/3 . z 1 " 1 + i = l , . . . , s . j = 1,...,m (2.3*+b) where t h e / ^ ' s a r e r e l a t e d t o the a*s and b's t h r o u g h t h e e q u a t i o n 1 p 21 n l fii. ftim = L b l l b 1 2 '21 J l m J2m nl nm (2.35) 48 w i t h L = h i ) \ : \ ' -(2.36) L i i a i , v . + l " ^ ,^+2 ' ' - a i,V ± +2 x ,V.+n. ' x x - a i , V n i o (2.37a) a i , - * . + l ~ ai,S» +2 J J - a i'Y2 0 - a 0 0 (2.37b) the m a t r i c e s b e i n g [ ( n i + 1) x (n^. + 1) J . Remark The m a t r i x L i s a l w a y s non s i n g u l a r : d e t ( L ) = 1 , s i n c e d e t ( L . . ) = cf. .. x j x j 49 We then obtain s subsystems whose input-output d e s c r i p t i o n i s p i i ( z ) y i ( k ) + + p i s ( z ) y s ( k ) = Q 1^(z) U ; L ( k ) + ... + Q*m(z)• u f f l(k) i = l , . . , s (2.38) Defining the integers n\ = max ( n 1 , n £, ... ,. n i _ 1 , nj,+l, n±+1, . , n Q) i = l , . . . , s (2.39) ,and the parameters Pjk = " ai,-> .+n.-k+l 3=1,...s k=l,..,n J + l q J k s s A 1 + i r 1 - k + I , J - 3=1,..,S k=l,..,n.+l i = l , . . , s (2.40) allows us to define the r e c i p r o c a l polynomials , • n.-n.+l . n.-n. . -n. " if + P j l Z + + P j , n + 1 Z n..+l 3 n .-n. -k+1 l . 0 l 13" ' P 3 k Z n.-n.+l V 7 = cf. ,z J 1 + Z_ i - l , _ , , S j = l , . . , s (2.41) n.-n. . n.-n.+l . - n . X „ X X , X X X , x l =• q j l Z + q j 2 Z + — • + q o , n . + l Z n .+1 a q i k z k=l J k n . - r u - k + l 1 x i  i = l , . . . , s j = l , . . . , m (2.42) and t o r e v / r i t e (2.38) i n t h e form P i l ( z ~ 1 ) y x ( k ) + ... + P ^ U " 1 ) y s ( k ) = Q i l ( z ~ 1 ) u x ( k ) + ... + Q i f f i ( z ~ 1 ) u m ( k ) i = l , . . , s (2.43) I n summary, a s e t o f s i n p u t - o u t p u t e q u a t i o n s o f t h e form s m Yl P^A^1) y M = X ^ 0 . . ( z _ 1 ) u ( k ) (2.44a) j = l 0 J 0=1 J J n.+l J , n.-n +1 \ . n.-n.-k+l P i 0 ( z ~ } = ^ i j 2 + ^ p j k Z i = l , . . , s (2.44b) n. +1 x 1 / • " ^ - m - k + l Q i 0 ( Z " ) = ^ 1 q ^ k Z J=r,-.-,s (2.44c) 1 = 1, • j s has been o b t a i n e d . . R e c a l l i n g ( 2 . 4 0 ) , ( 2 . 4 4 ) c a n be u s e d i n o r d e r t o e s t i m a t e t h e p a r a m e t e r s o f A a n d B. 3 . 2 . 2 , P a r a m e t e r s o f C I n e q u a t i o n s - ( 2 . 4 4 ) , which, a l l o w u s t o i d e n t i f y t h e s y s t e m d y n a m i c s , o n l y t h e i n f o r m a t i o n c o n t a i n e d i n u_ a n d i s used.. I n o r d e r t o i d e n t i f y t h e m a t r i x C, we s h a l l now d e r i v e a s e t o f e q u a t i o n s m a k i n g u s e o f t h e i n f o r m a t i o n p r o -v i d e d by y_ 2. A c c o r d i n g t o ( 2 . 2 3 ) , i t f o l l o w s i 2 ( k ) = C x ( k ) ' ( 2 . 4 5 ) w h e r e ( 2 . 4 6 ) T h e n , c o m b i n i n g ( 2 . 4 6 ) a n d ( 2 . 2 8 ) y i e l d s 2 2 ( k ) = C M ( z ) v ^ k ) - C N ( z ) u ( k ) o r , e q u i v a l e n t l y P*(z) v x ( k ) - £ 2 ( k > = Q.«(z) u ( k ) C = T T ^2 c T — r - s ( 2 . 4 7 ) ( 2 . 4 8 ) 52 where P ' ( z ) and Q'(z) a r e t h e m a t r i c e s P ( z ) •= CM(z) = F i J ( z ) i = l , . j = l , . . , r - s ( 2 . 4 9 a ) Q ( z ) = CN(z) = Ci ( z ) i = l , . J = l , . , r - s ( 2 . 4 9 b ) whose c o n s t i t u t i v e element P _ and a r e o b t a i n e d from s i m -p l e c a l c u l a t i o n s • * n. n - l P . . ( - z ) = c . , z J + c. - , z J x 0 x , V n . i . V V 1 " + . . . + C . . n z + c . X, \> . + 1 X , V . J J i = l , . . . , c - s j = l , . . . , s ( 2 , 5 0 a ) t * 0 - . ( z ) = . z 10 j n , x n - l v n - 2 j n - l , i .+ . + t . ^ . Z• + ^. . ^ , A/ , . ( j - l ) n , x ( j - l ) n + l , x i = l , . . . , r - s 'j = l,...,m ( 2 . 5 0 b ) i n w h i c h n = max(n.) x ( 2 . 5 D n - k + 1 s u * 7 7 ( 3 - l ) n + k , i = ^ £ a ° i . - V + v + k - 1 b v u + v ~ 1 ^ ( 2 . 5 2 ) V/e t h e n o b t a i n r - s s u b s y s t e m s whose i n p u t - o u t p u t d e s c r i p t i o n i s 5 3 p i i ( a ) . y ' i ( k ) + — + p ± s y s ( k ) - W k ) = ; Q i l ( z ) u : ( k ) + ... + ^ ( z ) u m ( k ) i = l , . . . , r - s ( 2 . 5 3 ) R e c a l l i n g ( 2 . 5 D and i n t r o d u c i n g t h e p a r a m e t e r s P j k = C i , v . + n . - k + l J=l>-->s k=l,..,n.. + l 3 3 q3:, = #\ , . j = l , . . . , m k = l , , . , n H j k j n - k + l , i o » > , 5 i = l , . . . , r - s (2.34) a l l o w s u s t o d e f i n e t h e r e c i p r o c a l p o l y n o m i a l s P. ( z L ) = z n p ( z ) J 3 . n .—n . n .-n+1 «v = p.,z J + p ^ z J + ... + p. . , , z " n - j l 32 - j . n ^ + l y ± n .-n-k+1 = p * z J i = l , . . , r - s j = l , . . , s (2.55) k = l J k t _1 _ « »* Q. .(z ) = z Q. .(z) v x j v ^13 i -1 i -2 , J i -n q .,z + q . _ z +.... + q . z j l " ^32" " * \ j , n n q ^ z " " ^ i = l , . . , r - s 3=1,..,in (2.56) •k=l J K and t o rev/ri-.e ( 2 , 5 3 ) i n t h e form P^C*-'1) y i ( k ) + ... + F^Cz" 1) y s ( k ) - y e + j L ( k ) = Q i l ( z " 1 } u l ( k ) + + Q i n ( z " 1 ) u m ( k ) i = l , . . r » r - s ( 2 . 5 7 ) Remark U s i n g ( 2 . 5 4 ) and ( 2 . 5 2 ) , t h e q 1, 's car. be e x p r e s s e d i n t e r m o f t h e p3:. 's J k ~ n -n+k s u q* = p 1 « , . b , • . ( 2 . 5 8 ) J k u=l v=l ^u,n u-v-n+k+l ^.+v-.l,j so t h a t t h e number o f unknown p a r a m e t e r s i n ( 2 . 5 3 ) i s n. I n summary, a set. o f r - s i n p u t - o u t p u t e q u a t i o n s o f th e form s . m ^ P . ' ^ - 1 ) y ( k ) - ^ + - ; ( k ) = YL Q T ( z " 1 ) u . ( k ) ( 2 . 5 9 a ) n a + i * r .1 i n.-n-k+1 P i j ( z ? = ^ p j k z J = l , ^ - , s (2.59b) C^Cz" 1) = q^z k j=l,,....,m (2,59c) i = l , ; , , . , r - s "has been o b t a i n e d . R e c a l l i n g (2.54), (2,59) a l l o w s u s t o e s t i m a t e t h e p a r a m e t e r s o f t h e m a t r i x CT, E q u a t i o n s (2.44) and ( 2 . 5 9 ) c o n s t i t u t e t h e c a n o n i c a l i n p u t - o u t p u t d e s c r i p t i o n o f ( 2 . 1 0 ) . The e q u i v a l e n c e between t h e c l a s s o f c a n o n i c a l s t a t e - s p a c e r e p r e s e n t a t i o n s ( 2 . 1 0 ) and the i n p u t - o u t p u t d e s c r i p t i o n (2,44)-(2..59) has t h u s been o b t a i n e d . I n t h e n e x t S e c t i o n , t h e r e s u l t s o b t a i n e d above a r e e x t e n d e d t o t h e c a s e o f d a t a c o r r u p t e d by n o i s e . 3.3.- N o i s y . v e r s i o n I f we assume t h a t t h e i n p u t and o u t p u t s e q u e n c e s a r e c o r r u p t e d by an a d d i t i v e n o i s e y ± ( k ) =3^ 00 + v ± ( k ) i = l , . . . , r (2.60) u ± ( k ) = u ± ( k ) + w ± ( k ) i = l , . . . , m ( 2 . 6 1 ) t h e e q u a t i o n s (2.44) and (2,59) t a k e t h e f o r m p I i = Q ». + p I x + Q I ( 2 , 6 2 ) p'zi ~ Z2 = Q'ii + -v 2-+ q'w (2.63) where v 1 ( k ) = [ v ^ k ) v 2 ( k ) v 2 ( k ) = [ v s + l ( k > v s + 2 ( k > v ( k ) = [ v ^ ( k ) ' v£(k)" T w(k) = [w-^k) w 2 ( k ) ... T (2.64) (2.65) J (2.66) \W]T (2.67) The i - t h component n ^ ( k ) o f t n e n o i s e v e c t o r n ( k ) = P v + Q w ( 2.63) d e f i n e d by s m n ( k ) = P . ^ z - 1 ) v. ( k ) + O.Az'1) w, ( k ) ( 2.69) 0 = 1 J j = l J has a p o w e r - d e n s i t y s p e c t r u m g i v e n by s m Q, S(S) 0 . .(r 1) <* (<?) ( 2 . 7 0 ) 0 = 1 ^ a j W i where and <|> a r e t h e p o w e r - d e n s i t y s p e c t r a o f v. and w. , v^ ± 1 1 S i n c e $ n (%) can be e x p r e s s e d as t h e p r o d u c t 4> i%) = * n V%) r (r1) (2 .71) n i n. " " n. 1 n. i i l t h e n , - 1 ^ H.(rX) = i — r - ( 2 . 7 3 ) T n <r ) c o n s t a n t i s t h e t r a n s f e r f u n c t i o n o f a r e a l i z a b l e f i l t e r . S i n c e , i n a d d i t i o n , ( 2 . 7 3 ) i s e q u i v a l e n t t o , - 1 , * N ( ? ) H . ( ^ ) \(%~X) = cl (2 .74) 57 H.(^ -^) i s the transfer function of a 'whitening f i l t e r . I t follows immediately that n ±(k) = H ±(z - 1 ) x e,(k) (2.75) where e^(k) i s a sequence of independant random v a r i a b l e s E [ e j = 0 E[e.(j) e.(k)] = <r\ cf (2.76a) (2.76b) I f the whitening f i l t e r i s represented by i t s impulse response, the process (2.?5) takes the form of an auto-regressive pro-cess where n (k) = e (k) D . U " 1 ) i V D. ( z " 1 ) = 1 + d^z"' j=l 3 ( 2 . 7 7 ) (2.78) form The equation ( 2 . 6 2 ) can, then, be re w r i t t e n i n the P y 1 = Q u + D e (2.79) where e(k) = [e-^k) e 2(k) ... e s ( k ) J T ( 2 . 8 0 ) ( 2 , 8 1 ) DizT1) = diag ^(z - 1 ) V 2 ' ^ 58 The e q u a t i o n ( 2 . 6 3 ) , ' i n t h e sane manner, becomes where u + D e i r • ( k ) e 2 ( k ) . ( 2 . 8 2 ) ' e r - s ( k ) ] T ( 2 * 8 3 ) D ' ( Z _ 1 ) = d i a g 1 D 1 ( Z _ 1 ) D ( z X ) r - s ( 2 . 8 4 ) V/e t h e n o b t a i n the n o i s y v e r s i o n o f t h e e q u a t i o n s ( 2 . 4 4 ) and (2.59) m Yl P, A z - 1 ) y n o = Yl Q. . 5=1 10 .. ( z _ 1 ) u . ( k ) + ^ - e ( k ) 0=1 1 J 3 D , ( z - X ) 1 ( 2 . 8 5 a ) m Yl P.'.CZ " 1 ) 5=1 J 3(k) - y s + i ( k ) = Yl oj - - - - 1 1— s ( 2 . 8 5 b ) These r e s u l t s g e n e r a l i z e t h o s e o b t a i n e d i n S e c t i o n 3.2- The e q u i v a l e n c e between t h e s t a t e - s p a c e r e p r e s e n t a t i o n and t h e i n p u t - o u t p u t d e s c r i p t i o n o f a m u l t i v a r i a b l e s y s t e m has t h u s been e s t a b l i s h e d i n t h e c a s e o f d a t a c o r r u p t e d by n o i s e . I t must be p o i n t e d o u t t h a t , a s r e v e a l e d by ( 2 . 8 5 ) , t h e a d v a n t a g e s o f t h e d e c o m p o s i t i o n of t h e s y s t e m i n t o s u b s y s t e m s has been p r e s e r v e d — I n t h e sense t h a t t h e s e s u b s y s t e m s c a n be s e p a r a t e l y i d e n t i f i e d - although, no r e s t r i c t i v e a s s u m p t i o n on n o i s e c h a r a c t e r i s t i c s has been made.. I t i s now a s i m p l e m a t t e r t o f o r m u l a t e t h e i d e n t i f i -c a t i o n p r o b l e m . However, b e f o r e s t a t i n g and s o l v i n g t h i s p r o -blem,, we s h a l l p r o c e e d t o some comments. 3.4- Comments The s t r u c t u r e o f the m a t r i c e s A, 3 and C, g i v e n i n S e c t i o n 2 , and, c o n s e q u e n t l y , t h e s t r u c t u r e o f the i n p u t - o u t p u t e q u a t i o n s , d e r i v e d i n S e c t i o n 3 , r e s u l t f rom an a r b i t r a r y c h o i c e o f t h e i n d e x e s n^. I t would be o f i m p o r t a n c e t o propose an e f -f i c i e n t method f o r the f i n d i n g o f t h e s e t j n^, n^, .... , n 0 j » o p t i m a l w i t h r e s p e c t t o a g i v e n c r i t e r i o n . However, v/hatever c r i t e r i o n we c h o o s e , such a method s h o u l d be based on the measured i n p u t - o u t p u t sequences and would l e a d t o the f o l l o w i n g p r o b l e m : how can one t a k e i n t o a c c o u n t t h e n o i s e whose s t a t i s t i c s a r e unknown ? F o r example, t h e m i n i m a l r e p r e s e n t a t i o n case l e a d s us t o c o n s i d e r i n p u t -o u t p u t d a t a m a t r i c e s : how can we a s s e r t t h a t t h e s e m a t r i c e s a r e s i n g u l a r or n o t ? ( I s a m a t r i x whose d e t e r m i n a n t i s 10~ 1 0 s i n g u l a r ?"). I t t h e n seems t h a t t h i s p r o b l e m cannot be s o l v e d i n r e a l i s t i c c a s e s . N e v e r t h e l e s s , t h e c h o i c e f o r t h e n.'s can be g u i d e d by some p h y s i c a l c o n s i d e r a t i o n s : — T h e s t r u c t u r e may be a s s i g n e d by t h e p h y s i c a l m e a n i n g o f t h e s t a t e v a r i a b l e s . — I f t h e j - t h o u t p u t o f t h e s y s t e m i s m a r k e d l y more f r e e f r o m n o i s e t h a n t h e r e m a i n i n g o n e s , a m a x i m i -z a t i o n o f t h e c o r r e s p o n d i n g n . c o u l d be a d v a n t a -g e o u s . — T h e number o f unknown p a r a n e t e r s i s i n (2.85a) : n + m ( n i + l ) + n ^ i = l , . . . , s i n (2.85b) : n + n d i = l , . . . , r - s I f n i s a m u l t i p l e o f s, i t c a n be j u d i c i o u s , i n o r d e r t o m i n i m i z e t h e c o m p u t i n g r e q u i r e n e n t s , t o c h o o s e n^= n / s , s i n c e t h e s u b s y s t e m s a r e s e p a r a t e l y i d e n t i f i e d . A n o t h e r p r o b l e m a r i s e s w i t h r e g a r d t o t h e n o i s e m o d e l : i n t r o d u c i n g a n o i s e f i l t e r i n c r e a s e s t h e number o f unknown p a r a m e t e r s which a f f e c t s t h e I d e n t i f i c a t i o n p r o c e -d u r e , s i n c e t h e e f f i c i e n c y o f a n y o p t i m i z a t i o n a l g o r i t h m d e -c r e a s e s when t h e number o f p a r a m e t e r s i n c r e a s e s . H o w e v e r , t h e n o i s e m o d e l a l l o w s u s — to r e d u c e t h e e f f e c t o f t h e c h o i c e o f t h e n ^ s on t h e i d e n t i f i c a t i o n r e s u l t s s i n c e t h i s e f f e c t i s d i r e c t l y r e l a t e d t o t h e amount o f n o i s e c o r r u p -t i n g t h e d i f f e r e n t o u t p u t s . - to i m p r o v e t h e e f f i c i e n c y o f some s i m p l e m i n i m i z a -t i o n a l g o r i t h m s , w h i c h a r e s e n s i t i v e t o n o i s e . M o r e o v e r • s u c h an a p p r o a c h p r o v i d e s a p r o b a b i l i s t i c framework w h i c h p e r m i t s us t o e s t a b l i s h t h e m a t h e m a t i c a l p r o p e r t i e s o f th e e s t i m a t e s . 4 - MAXIMUM LIKELIHOOD ESTIMATOR i + . l . I n t r o d u c t i o n The o b j e c t o f t h i s S e c t i o n i s t o p r e s e n t a Maximum L i k e l i h o o d e s t i m a t o r . The e q u a t i o n o f e r r o r and t h e l i k e l i h o o d f u n c t i o n w i l l be d e r i v e d ( S e c t i o n 4 . 2 ) and t h e c o n s i s t e n c y o f t h e Maximum L i k e l i h o o d e s t i m a t e s w i l l be p r o v e d ( S e c t i o n 4 . 3 ) . 4 . 2 , E q u a t i o n o f e r r o r R e c a l l i n g ( 2 . 8 5 ) l e a d s u s t o choose t h e model s m 1=1 , , • • , s ( 2 , 8 6 a ) s m i i — • » > 3?" — S ( 2 , 8 6 b ) 62 where n .+1 , n,-n,+l n , - n , - k + l (2.87) n.+l n . - n . - k + l Q^OT1) =- YL q ^ z " 1 ^ ( 2 . 8 8 ) n .+1 ,. , n .-n-k+1 P ( z — ) ^ Z _ p ^ z J (2,89) 1 3 k = l J n i n d D.Cz"' 1) = 1 + YL d^z" 3' (2.9D 1 i = i * i n d D.(z X ) = 1 +: ZL, (T- z J (2.92) 1 j = l 3 The r e s i d u a l s £^(t) [ ^ ( t ) ] 3 X 6 i n d e p e n d a n t . and G a u s s i a n ( 0 , ; r±) [ ( 0 , ^ ) ] . R e w r i t i n g (2,86) i n t h e fo r m 63 e±(t) = D. S \—' ^ P^ , y , ( t ) 0=1 3=1 0. . u . ( t ) "13 3 1—1 y »#-»•) S L 3=1 ( 2 . 9 3 a ) ^—• Q_. • u 3=1 1 — 1 , , , » j T - S y i e l d s t h e s o - c a l l e d e q u a t i o n s o f e r r o r ( 2 . 9 3 b ) £ ( t ) = Z _ a l h=0 n.,+1 y f~* - i <^ —1 p., y .(t+n .-n. -k-h+1) 3=1 K = l 3 3 0 + y i ( t + n i _ n j L - h + l ) -n. +1 m x 0=1 k = l q 3 k u o ( t + n . - n . - k - h + l ) i — 1 j » « • j s ( 2 . 9 4 a ) n € , ( t ) = E?: h=0 s n .+1 3=1 k = l m n P,, - y , ( t + n -n-k-h+1) 3 - ^ 3 3 E C 1 " I -s + 1 - ' ^ ^ i k " i 8 1 0=1 k = l J k J q , v u . ( t - k - h ) 3.—1 j » « • i X*~"S ( 2 . 9 4 b ) 64 The G a u s s i a n ess o f c ± and €^ a l l o w s , us t o d e f i n e t h e l o g - L i k e l i h o o d f u n c t i o n s N L i £ 2 ( t ) + Nlogcr + NTLo S2Tr i = l , . . , s ( 2.95a) 2 * 2 t = l i 1 L. = -l 2 cr ^ x~x i » G i 2 ( t ) + Nlog<r! + Nlog27r i = l , . . . , r - s ( 2 . 9 5 b ) where N i s t h e number of a v a i l a b l e d a t a . The Maximum L i k e l i h o o d e s t i m a t e s o f Pj_-pQjLji a n d Qj_j a r e o b t a i n e d by m a x i m i z a t i o n o f L^ and L^ or e q u i v a l e n t l y by m i -n i m i z a t i o n o f the l o s s f u n c t i o n s N W. = ~ Z^cf(t) i=l,....,s (2.96a) 1 ^ t=l N w| = ~~ H £-2(t) l=l,...,r-s (2.96b) 1 ^ t=l H a v i n g d e f i n e d t h e Maximum L i k e l i h o o d e s t i m a t i o n i p r o b l e m as t h e f i n d i n g o f t h e g l o b a l m i n i m a - o f V/. and V/., x x we c a n d e t e r m i n e t h e s t a t i s t i c a l p r o p e r t i e s o f t h e e s t i m a t e s by e x a m i n i n g t h e s e minima.. 4-3.- S t a t i s t i c a l p r o p e r t i e s o f t h e e s t i m a t e s I n t h e whole S e c t i o n , t h e f o l l o w i n g a s s u m p t i o n s a n made : . - - -- A l . To ensu r e t h a t t h e model s t r u c t u r e i s a p p r o p r i a t e , i t i s assumed t h a t ( 2 . 8 5 ) r e a l l y h o l d s . —1 ' —1 _A2 v F o r a l l ( i , j ) , t h e p o l y n o m i a l s P. . ( z ) ? . . ( z ) , — i • — i D. ( z ) and D. ( z ) a r e assumed t o have a l l t h e i r x x z e r o s i n s i d e t h e u n i t c i r c l e . T h i s i m p l i e s t h a t , t h e e q u a t i o n s (2.85) a r e s t a b l e . _ A3 . F o r a l l i = l , . . . . ,s and j = l , . . . , min (m, s) £i=l,.. . , r -and j = l , . . . . ,min(m,r-s) 1 t h e p o l y n o m i a l s P. .(z~^~) and Q ij(z""' 1') ^ ^ ( z " 1 ) and Q i ^ . ( z ~ " 1 " ) J a r e r e l a t i v e l y p r i m e . These c o n d i t i o n s a r e f u l f i l l e d p r o v i d e d t h e sy s t e m i s c o n t r o l l a b l e from u. _ M .- A l l s t o c h a s t i c p r o c e s s e s a r e e r g o d i c . M o r e o v e r — S i n c e , f o r a f i n i t e number o f d a t a , N, the p r o p e r t i e s t o f and may change d r a s t i c a l l y f rom one e x p e r i m e n t t o a n o t h e r , t h e a s y m p t o t i c l o s s f u n c t i o n s V . = l i m W i = l , . . . . , s ( 2 . 9 7 a ) v ! = l i m W\ i = l , . . . , r - s ( 2 . 9 7 b ) w i l l be c o n s i d e r e d i n t h e f o l l o w i n g . - For convenience, E f ( t ) v / i l l denote N J 77 ^ » lim N _ » o o " t=l f ( t ) (2.98) provided the l i m i t e x i s t s . I f f ( t ) i s an ergodic s t o c h a s t i c process, E ^ f ( t ) j i s the expected value of f ( t ) . The output data are governed by P = Q u + D e i i P ii - 12' = Q H. + D 1 (2.99a) (2.99b) while the model i s defined by P zx = Q u + D e (2,100a) (2.100b) Let u s , • f i r s t , consider (2.99a) and (2.100a). Combining these two equations y i e l d s £ = D • P P D e + D P P -1 * - l ~ Q - D Q (2.101) Then, assuming that u and £ are uncorrelated and s e t t i n g (2.102) ~ - l " -1 F = D P P D A -J A 1 A "I A Ki = D P P Q — B Q (2.103) allows us to write E [ i £ T ] = E [F. e e T F T J + E [ G u u T G T J ( 2 . 1 0 4 ) w h i c h c l e a r l y i m p l i e s t h a t E ['£ £ T ] > E [ F £ £ T p T ] = V ( 2 . 1 0 5 ) where A > B means t h a t A-B i s p o s i t i v e d e f i n i t e . A p p l y i n g P a r s e v a l ' s theorem, V becomes ( A s s u m p t i o n A4) ^ = M f TO) R p T ^ " 1 ) i f ' ( 2 . 1 0 6 ) ' 1*1=1 where R = E e T J ( 2 . 1 0 7 ) L e t us c o n s i d e r t h e i n t e g r a l J = ~ r i [?{%) - i j R[F T ( r 1 ) - i ] -d# > 0 2 l f i ( 2 . 1 0 8 ) T h e n J = V f J x - 2 J 2 ( 2 . 1 0 9 ) where J , =<fc» R -Tr = R (2.110) 7 m=i 5 J - =cb F(S) R ~ = R (2.111) 2 J |%| =1 * The l a t t e r e q u a l i t y i s due t o th e f a c t t h a t t h e e l e m e n t s o f F ( ^ ) have t h e i r zeros o u t s i d e the u n i t c i r c l e (Assumption A2) and that F ( 0 ) = D~ 1(0) P(0) P~ 1(0) D(0) = I s i n c e D ±(0) = D ±(0) = 1 and P (0) = P (0) = J H e n c e J = V - R > 0 (2.112) and, r e c a l l i n g (2.105) and (2.107) E [ i £ T ] - E [ e £ T ] > 0 (2.113) A p p l y i n g S y l v e s t e r ' s theorem i m p l i e s t h a t t h e d i a g o n a l e l e m e n t s must be p o s i t i v e V i = E [ € i ] ^ E [ e i ] = °i '1=1,....8 (2.114) T h e e q u a l i t y , i n (2.114) o c c u r s i f and o n l y i f (1) J = V - R = 0 i . e . F = D"1 P P - 1 D = I (2.115) 69 ( i i ) t h e second term o f the r i g h t - h a n d s i d e o f ( 2 . 1 0 4 ) v a n i s h e s i . e . - p r o v i d e d u i s p e r s i s -t e n t l y e x c i t i n g o f o r d e r 2n + n ^ + 1 ( s e e chap -t e r 1 ) - i f G .= D"1 P P" 1 Q - D" 1 Q = 0 ( 2 . 1 1 6 ) C o m b i n i n g (2.115) and (2.116) y i e l d s ' D P = D P ( 2 . 1 1 7 ) - D Q = D Q ( 2 . 1 1 8 ) w h i c h i m p l i e s ( A s s u m p t i o n A3) A P = P Q = Q A D = D ( 2 . 1 1 9 ) ( 2 . 1 2 0 ) ( 2 . 1 2 1 ) We have t h e n e s t a b l i s h e d t h e Lemma 2 . 3 . V ± ^ > crjr i = l , . . . , 8 ( 2 . 1 2 2 ) A A The e q u a l i t i e s a r e o b t a i n e d i f and o n l y i f P = P, Q = Q and D = D. Tha lemma 2 . 3 s t a t e s t h a t t h e e s t i m a t e s o f t h e p a r a -m e t e r s o f t h e m a t r i c e s A and B c o n v e r g e t o t h e t r u e v a l u e s as N »-=»:> . 70 Let us now consider ( 2 . 9 9 b ) and ( 2 . 1 0 0 b ) .Remarking, from ( 2 . 1 0 0 a ) , t h a t v., i s defined by ^ = P" 1 Q u + P - 1 D e ( 2 . 1 2 3 ) the equations ( 2 . 9 9 b ) and ( 2 . 1 0 0 b ) can be re w r i t t e n as X 2 = ( p' P" 1 Q - Q f) u +• P* P" 1 D e - D* e' ( 2 . 1 1 4 a ) £ 2 = (F* P" 1 Q - Q ) u + P f P" 1 D e - D' ( 2 . 1 1 4 b ) fr-am which we can deduce € = [D 1 (P - P ) P" 1 Q - D - 1 ( Q " - Qr)j u 4W.,Pr)P-1D]e+[;r-1D,]er ( 2 . 1 1 5 ) Proceeding as above lead us to the Lemma 2 . 4 . V i > * l ' 2 i = l , . . . , r - s (2.116) The e q u a l i t i e s are obtained i f and only i f p - p Q = Q and D" = D . The lemma 2 . V p r o v e s t h a t t h e e s t i m a t e s o f t h e p a r a m e t e r s o f t h e m a t r i x C c o n v e r g e t o t h e t r u e v a l u e s as N—>~ <>=> . We t h e n o b t a i n t h e Theorem 2.1 The Maximum L i k e l i h o o d e s t i m a t e o f t h e t r i p l e t (A, B, C) i s a s y m p t o t i c a l l y c o n s i s t e n t . The r e m a i n i n g p r o b l e m i s the m i n i m i z a t i o n o f the l o s s f u n c t i o n s . . I n s p e c t i n g (2.96) and (2.94) shows t h a t m i n i m i z i n g i W^  and i s a n o n - l i n e a r o p t i m i z a t i o n p r o b l e m . The o b j e c t o f t h e n e x t s e c t i o n i s co p r e s e n t t h e v a r i o u s c l a s s e s o f n o n - l i n e a r o p t i m i z a t i o n a l g o r i t h m s . 5 - MINIMIZING-THE LOSS FUNCTIONS 5.1. I n t r o d u c t i o n The n o n l i n e a r m i n i m i z a t i o n p r o b l e m can be f o r m a l l y ' s t a t e d as M i n i m i z e V(p_) p € E n (2.12?) 7 2 A large number of methods have been proposed to solve the general nonlinear minimization problem [ 2 0 , 1 ^ ] . These methods can be divided i n two classes : — The methods that use derivatives - The methods that do not use derivatives In the following, the most commonly used methods are br i e f l y described 5.2... Minimization procedures usinr; derivatives We f i r s t consider how to solve problem ( 2 . 1 2 7 ) by algorithms that make use of f i r s t and, possibly, second d e r i -vatives of V(p_). 5.2.1. Gradient methods i l k ] At the k-th stage, the transition from a point - (k) . . " . . (k+1) . JD to another point, i s given by ^(k+D s £ ( k ) _ x v ( p ( k ) } ( 2 a 2 8 ) The negative of the gradient gives the direction for optimiza-tion but not the magnitude of the step, so that various steepest doscent algorithms are possible, depending on the choice of A . Many methods of selecting A are available. It can be shown that, under suitable conditions, the method converges as k • 0 0 . However this, theoretical.result i s of l i t t l e i n -terest i n practice because the rate of convergence can be intolerably slow. 5 - 2 . 2 - Newton s m e t h o d s [ 1 9 ] 1 The s e c o n d - d e r i v a t e m e t h o d s , among w h i c h Newton s m e t h o d i s t h e b e s t known, o r i g i n a t e f r o m t h e q u a d r a t i c a p p r o -x i m a t i o n o f V(p_).. The t r a n s i t i o n f r o m t o _p_^ k + 1^ i s V ( p ( k > ) nn • ( 2 . 1 2 9 ) T h e c o n v e r g e n c e i s g u a r a n t e e d i f t h e i n v e r s e o f t h e H e s s i a n m a t r i x o f V i s p o s i t i v e d e f i n i t e . T h i s i s a m a j o r d r a w b a c k t o t h i s m e t h o d s i n c e f o r f u n c t i o n s v / h i c h a r e n o t s t r i c t l y c o n v e x , t N ewton s m e t h o d c a n d i v e r g e . 5 . 2 . 3 - Q u a s i - N e w t o n m e t h o d s [ l 4 j Q u a s i - N e w t o n m e t h o d s a p p r o x i m a t e t h e H e s s i a n m a t r i x o r i t s i n v e r s e b u t -use i n f o r m a t i o n f r o m o n l y f i r s t - o r d e r d e r i -v a t i v e s . A t s t a g e ( k+ 1 ) , p ^ k + 1 ^ i s c o m p u t e d f r o m J 3 ^ K ^ t h r o u g h £ ( k + l ) = £ ( 1 ° - * k A ( p ( k ) ) V £(p< k>) ( 2 . 1 3 0 ) w h e r e A i s a n a p p r o x i m a t i o n o f t h e i n v e r s e o f t h e H e s s i a n . The p r o b l e m o f a p p r o x i m a t i n g t h e H e s s i a n h a s b e e n i n v e s t i g a t e d by P e a r s o n f l ? ] r , D a v i d o n - F l e t c h e r - P o w e l l [ " 1 4 ] a n d F l e t c h e r [ 1 0 ] -5.3«- M i n i m i z a t i o n p r o c e d u r e s w i t h o u t u s i n g d e r i v a t i v e s T h i s S e c t i o n i s c o n c e r n e d w i t h d e r i v a t i v e - f r e e t y p e o f m e t h o d s ( s e a r c h m e t h o d ) . As a g e n e r a l r u l e , f i r s t - d e r i v a t i v e a n d s e c o n d - d e r i v a t i v e m e t h o d s c o n v e r g e f a s t e r t h a n s e a r c h m e t h o d s . H o w e v e r , i n p r a c t i c e , t h e d e r i v a t i v e - t y p e m e t h o d s h a v e two m a i n d r a w b a c k s t o t h e i r i m p l e m e n t a t i o n . F i r s t , i t c a n be l a b o r i o u s o r i m p o s s i b l e t o p r o v i d e a n a l y t i c a l f u n c t i o n s f o r t h e d e r i v a t i v e s . S e c o n d , t h e d e r i v a t i v e - t y p e m e t h o d s r e q u i r e a l a r g e amount o f p r o b l e m p r e p a r a t i o n a s c o m p a r e d w i t h s e a r c h m e t h o d s . Two o f t h e many e x i s t i n g s e a r c h a l g o r i t h m s a r e b r i e f l y d e s c r i b e d , A c o m p l e t e d e s c r i p t i o n o f t h e s e two m e t h o d s i s g i v e n i n [14]. 5 . 3 . I . R o s e n b r o c k s m e t h o d [ 2 2 ] 1 ( k + 1 ) R o s e n b r o c k s m e t h o d l o c a t e s _t> ' by s u c c e s s i v e u n i d i m e n s i o n a l s e a r c h e s f r o m a n i n i t i a l p o i n t a l o n g a s e t o f o r t h o n o r m a l d i r e c t i o n s v T . ^ ^ ' * * ' » V n ^ ^ g e n e r a t e d by Gram -S c h m i d t p r o c e d u r e . 5 , 3.2. P o w e l l s m e t h o d 1 P o w e l l s m e t h o d l o c a t e s t h e minimum o f V b y s u c c e s -s i v e u n i d i m e n s i o n a l s e a r c h e s f r o m an i n i t i a l p o i n t a l o n g a s e t of c o n j u g a t e d i r e c t i o n s . C h o o s i n g an a l g o r i t h m i s n o t a s i m p l e m a t t e r s i n c e no one m e t h o d a p p e a r s t o be f a r s u p e r i o r t o a l l t h e o t h e r s . T h e c h o i c e cf a p a r t i c u l a r a l g o r i t h m r e s t s on t h e f o r m u l a t i o n of t h e p r o b l e m a n d t h e e x p e r i e n c e o f t h e p r a c t i t i o n e r s . For the i d e n t i f i c a t i o n of m u l t i v a r i a b l e systems, search methods seem more app r o p r i a t e than d e r i v a t i v e - t y p e methods because of the g e n e r a l l y l a r g e number of parameters. However, i f the com-p u t a t i o n a l (storage) requirements are considered, a quasi-Newton method can be chosen. In Chapter 3, the v a r i o u s aspects of these problems w i l l be t r e a t e d on some examples. 6 - CONCLUSION In t h i s chapter a u n i f i e d approach to the i d e n t i f i c a t i o n of m u l t i v a r i a b l e systems has been presented. Although the problem of d e r i v i n g the input-output de-s c r i p t i o n of a m u l t i v a r i a b l e system from i t s state-space represen-t a t i o n has been deeply i n v e s t i g a t e d by many authors, no general r e s u l t i s a v a i l a b l e i n l i t e r a t u r e . The r e s u l t s of S e c t i o n 3 b r i d g e t h i s gap. Compared w i t h the previous methods, the present one has the f o l l o w i n g advantages. - This i s the most general approach, s i n c e a l l the pre-v i o u s methods are s p e c i a l cases of the above method. a) I f the number of subsystems i s equal to the number of outputs (s = r ) , one obtains the input-output r e p r e s e n t a t i o n de-r i v e d by Gopinath [12] and Zuercher [28]. b) I f the n^'s are chosen as l a r g e as p o s s i b l e , one obtains the s o - c a l l e d m i n i m i z a t i o n r e a l i z a t i o n d e r i v e d by Acker-man [ 1 ]. c) I f the two above c o n d i t i o n s can be s a t i s f i e d simultaneously, the input-output r e p r e s e n t a t i o n , suggested by G u i d o r z i [13], i s obtained. - The i d e n t i f i c a t i o n of the system dynamics (matrices A and B) i s based on the i n f o r m a t i o n provided by u. and y_ This allows a d i s c r i m i n a t i o n of the input-output data which i s of i n t e r e s t i f some outputs of the system are markedly more f r e e from . n o i s e than the other ones. - Equations (2.28) and (2.25) provide a very simple s t a t e estimator. Having estimated the parameters of P and Q ( i . e . A and B), equation (2.99a) gives an estimate of the output v e c t o r y_^, say y_, from which the noise has been removed Py^ = Qu (2.131) Then, provided i i i s n o i s e - f r e e , the s t a t e estimator x(k) = Mx 1(k) - Nu(k) (2.132) gives an unbiased estimate of x ( k ) . This f o r m u l a t i o n leads to an a l t e r n a t i v e approach to the i d e n t i f i c a t i o n of the parameters of the o b s e r v a t i o n matrix: the parameters of C can be estimated through equation (2.45). The Maximum L i k e l i h o o d e s t i m a t o r , introduced i n Se c t i o n 4, all o w s us to e s t a b l i s h the con s i s t e n c y of the estimates i n case of c o r r e l a t e d n o i s e . A l l these r e s u l t s c o n s t i t u t e a g l o b a l approach to the i d e n t i f i c a t i o n of l i n e a r m u l t i v a r i a b l e systems. The whole procedure described i n the previous Sections can be extended, i n an e n t i r e l y obvious way, to systems whose input-output s t r u c t u r e i s known a p r i o r i ( t r a n s f e r f u n c t i o n s , impulse response). Due to the r e c u r s i v e s t r u c t u r e of the computation scheme, the procedure can be extended, w i t h only minor m o d i f i c a t i o n s , to an o n - l i n e i d e n t i f i c a t i o n procedure. 7 8 C h a p t e r 3 : EXPERIMENTS 1 - INTRODUCTION The aim o f t h i s C h a p t e r i s t o p r e s e n t some examples o f Maximum L i k e l i h o o d i n d e n t i f i c a t i o n o f m u l t i v a r i a b l e s y s t e m s . F o r t h i s p u r p o s e , t e s t s w i l l be made on s i m u l a t e d d a t a — f o r d i f f e r e n t m i n i m i z a t i o n s a l g o r i t h m s - f o r d i f f e r e n t n o i s e powers 2 - EXPERIMENT 1 We c o n s i d e r t h e f o l l o w i n g t h i r d - o r d e r - syscem 0 1 0 x ( k + l ) =• a l l a l 2 a l 3 a 2 1 a 2 2 a 2 3 1 0 0 i ( k ) = 0 0 1 x ( k ) + orde b b 11 12 '°21 b 2 2 b 3 i b,, u ( k ) ( 3 . 1 ) x ( k ) ( 3 . 2 ) w i t h 2 i n p u t s and 2 o u t p u t s . T h i s s y s t e m has been s i m u l a t e d on an IBM 370 computer. The in p u t s and the o u t p u t s vere sequences o f l e n g t h N = 2 5 0 . The l o s s f u n c t i o n was m i n i m i z e d v i a the F l e t c h e r a l -g o r i t h m £ lo ] ..... T a b l e s 3-1- and 3 - 2 . g i v e t h e r e s u l t s f o r 79 d i f f e r e n t v a l u e s o f the s i g n a l - t o - n o i s e r a t i o , N N " s i = l If k = l 3. EXPERIMENT 2 The f o l l o w i n g system has been s i m u l a t e d on an IBM 370 computer: x ( k + l ) = zoo = 0 1 0 0 b l l b 1 2 a l l a 1 2 a13 a 1 4 x ( k ) + b 2 1 b 2 2 u ( k ) ( 3 . 3 ) 0 0 0 1 b 3 l b 3 2 a 2 1 a 2 2 a 2 3 a 2 4 Al V 1 0 0 0 0 0 1 0 x ( k ) ( 3 . 4 ) c l l C 1 2 C 1 3 c 1 4 The i d e n t i f i c a t i o n o f t h i s s y s t e m has been a c h i e v e d by u s i n g t h e R o s e n b r o c k m i n i m i z a t i o n a l g o r i t h m f o r a s e t N = 250 i n p u t - o u t p u t data.. The r e s u l t s a r e r e p o r t e d i n T a b l e s 3 . 3 and 3 . 4 . The r e s u l t s , g i v e n i n T a b l e s 3 . 1 — 3 . 4 show the u s e -f u l n e s s , , i n c a s e o f s m a l l s i g n a l - t o - n o i s e - r a t i o s , of the n o i s e model w h i c h a l l o w s one t o e l i m i n a t e the b i a s on t h e p a r a m e t e r s . T h i s c o n f i r m s t h e t h e o r e t i c a l r e s u l t s o f t h e p r e v i o u s c h a p t e r s . p a r a m e t e r s t r u e v a l u e s e s t i m a t e d v a l u e s a l l - 0.5 - 0 . 5 0 0 0 a 1 2 - 0 . 2 5 - 0 . 2 4 9 9 a 1 3 0.25 0 . 2 5 0 0 a 2 1 1.375 1.3749 a 2 2 i-5 1 . 5 0 0 0 a 2 3 - 0.75 - 0 . 7 5 0 0 b l l 0 . 2 0.1999 b 1 2 0.3 0.2999 b 2 1 0 . 1 0.0999 b 2 2 - 0 . 2 7 5 - 0 . 2 7 5 0 b 3 l 0.85 0.8499 b 3 2 0 . 9 5 0 . 9 5 0 0 T a b l e 3.1 : T; = e o 81 parameters true va3.ues estimated values n d = 0 n a = 3 - 0.5 0.083 - 0.489 a i l - 0.25 - 0.236 - 0.252 a12 a 13 0.25 0.018 0.249 1.328 a ~ 1.375 1.005 d21 a ~ 1.5 1.435 1-506 . 22 a — - 0.75 - 0.595 - 0.751 23 b 0.2 0.181 0.210 D l l b 0.3 0.329 0.296 12 b 0.1 0.084 0.097 D21 b - 0.275 - 0.237 - 0.259 22 b , 0.85 0.426 0.879 D 3 l 0.898 b 0.95 0.874 32 Table 3-2 ': | = 10 (n, : noise model order) estimated values parameters true values n , = 0 d n d = 3 a l l - 1. -1 .274 - 1.061 a l 2 0.5 0.067 0.506 a13 2 2.186 2.003 a U 1.5 1.324 1.491 a21 1. 0.645 0.979 a22 2.5 2.484 2.485 a23 0.23 - 0.016 0.261 a24 1 0.715 1.083 b l l 0.5 0.674 0.501 b12 0.5 0.323 0.518 b21 2 2.711 2.215 b22 - 1 - 1.617 • - 1.184 V 4 3.P77 3.899 • b32 3.2 1.006 2.994 V 2.7 - 0.199 2.645 V - - 0.2 —1.157 - - 0.247 c l l 1 - 0.874 1.012 C12 - 1 - 1.015 - 1.024 C13 - 1 - 1.637 - 1.127 c 14 1 0.421 0.919 T a b l e 3.3 : 77 = 1 0 (rt, : noise model order) p a r a m e t e r s t r u e v a l u e s e s t i m a t e d n . = 0 d v a l u e s n , = 3 a l l - 1 1 . 4 2 9 - 1 . 0 6 4 a l 2 ' 0.5 - 0 . 1 2 6 0 . 512 a 13 2 2 . 1 1 4 2 . 0 0 4 a14 1.5 1.677 1 . 4 9 4 a 2 1 •1 0 . 4 2 9 0 . 8 7 5 a 2 2 2 . 5 2 . 7 3 4 2 . 5 1 5 a 23 a 2 4 b l l 0 . 2 5 1 0 . 5 - 0 . 1 2 6 0 . 8 3 4 0 . 7 2 2 0 . 2 9 9 1 . 0 7 7 0 . 4 9 4 b12 0 . 5 0 . 2 9 9 0 . 5 7 7 b21 . 2 2 . 6 1 2 2 . 1 2 9 b22 - 1 - 1 . 7 3 9 - 1 . 2 2 7 4 2 . 9 8 4 3 . 7 2 5 b32 3.2 1 . 0 0 8 2 . 9 9 9 2.7 - 0.237 2 . 8 2 5 - 0 . 2 - 0 . 9 7 5 - 0 . 2 5 1 - c l l 1 0 . 8 2 7 1 . 1 0 8 • C12 - 1 - 1 . 1 1 5 - 1 . 0 0 9 C13 - 1 - 1 . 7 2 5 - 1 . 1 1 8 1 i c. . ! 14 1 0 . 6 2 2 0 . 9 0 5 T a b l e 3 . 4 . * N ( n , : n o i s e model o r d e r ) CONCLUSIONS The purpose of t h i s t h e s i s i s to i n v e s t i g a t e the Maximum L i k e l i h o o d i d e n t i f i c a t i o n of l i n e a r d i s c r e t e - t i m e systems. The s t a t i s t i c a l p r o p e r t i e s of the Maximum L i k e l i h o o d e s t i m a t o r , f o r s i n g l e i n p u t - s i n g l e output systems, are analyzed. Two d i f f e r e n t types of n o i s e model s t r u c t u r e are considered and two new r e s u l t s r e l a t i v e to the convergence p r o p e r t i e s of i d e n -t i f i c a t i o n methods are obtained. A general method f o r d e r i v i n g the input-output d e s c r i p t i o n o f a m u l t i v a r i a b l e system from i t s state-space r e p r e s e n t a t i o n i s proposed. I t i s shown th a t t h i s approach i s s u p e r i o r , from d i f -f e r e n t p o i n t s of view, t o those proposed i n l i t e r a t u r e . Based on the r e s u l t s obtained i n the s i n g l e i n p u t - s i n g l e output case, a Maximum L i k e l i h o o d e s t i m a t i o n procedure f o r m u l t i -v a r i a b l e systems i s des c r i b e d and convergence p r o p e r t i e s are d e r i v e d . I t i s shown th a t the r e s u l t s can be g e n e r a l i z e d i n v a r i o u s d i r e c t i o n s . F i n a l l y , numerical r e s u l t s of Maximum L i k e l i h o o d i d e n -t i f i c a t i o n a re obtained. 85 REFERENCES [1] Ackermann, J.E. and Bucy, R.S., "Canonical minimal r e a l i z a t i o n of a m a t r i x of impulse response sequences", Inform. Contr. 19, pp. 224-231, 1971. O It [2] Astrom, K.J, and B o h l i n , T., "Numerical i d e n t i f i c a t i o n of l i n e a r dynamic systems from normal o p e r a t i n g r e c o r d s " , IFAC Symp. on theory of s e l f - a d a p t i v e systems, Teddington, England, 1965. In theory of s e l f - a d a p t i v e c o n t r o l systems (ed. by P.H. Hammond), Plenum Press, New York, 1966. O II [3] Astrom, K.J., and Eykhoff, P., "System i d e n t i f i c a t i o n - A survey", Automatica 7, pp. 123-162, 1971. o I t I t t l [4] Astrom, K.J., and Soderstrom, T. , "Uniqueness of the maximum l i k e l i h o o d estimates of the parameters of an ARMA model", IEEE Trans. Aut. C o n t r o l , S p e c i a l i s s u e on i d e n t i f i c a t i o n and time s e r i e s a n a l y s i s , AC-19, pp. 769-773, 1974. [5] B o h l i n , T., "On the problem of a m b i g u i t i e s i n maximum l i k e l i h o o d i d e n t i f i c a t i o n " , Automatica 7, pp. 199-210, 1971. [6] C l a r k e , D.W., "Generalized Least-Squares e s t i m a t i o n of the para-meters of a dynamic model", IFAC Symp. on i d e n t i f i c a t i o n i n autom. . c o n t r o l systems, Prague, paper 3.17, 1967. [7] Cramer, H., "Mathematical methods of s t a t i s t i c s " , P r i n c e t o n U n i v e r s i t y P ress, 1946. [8] Cuenod, M. and Sage, A., "Comparison of some methods used f o r process i d e n t i f i c a t i o n " , Automatica 4, pp. 235-269, 1968. [9] Eykhoff, P., "System i d e n t i f i c a t i o n , parameter and s t a t e e s t i m a t i o n " , Wiley, 1974. [10] F l e t c h e r , R., "A new approach to v a r i a b l e m e t r i c methods", Computer J . 13, pp. 317-322, 1970. [11] F l e t c h e r , R. and Powell, M.J.P., "A r a p i d l y convergent descent method f o r m i n i m i z a t i o n " , Computer J . 6, pp. 163-168, 1963. [12] Gopinath, B., "On the i d e n t i f i c a t i o n of l i n e a r time i n v a r i a n t systems from input-output data", B e l l . Syst. Tech. J . 48, pp. 1101-1113, 1969. [13] G u i d o r z i , R., "Canonical s t r u c t u r e s i n the i d e n t i f i c a t i o n of m u l t i v a r i a b l e systems", Automatica 11, pp. 361-374, 1975. [14] Himmelblau, P.M., "Applied n o n l i n e a r programming", McGraw H i l l 1972. 86 [15] Luenberger, D.G., "Canonical forms f o r l i n e a r m u l t i v a r i a b l e systems", IEEE Trans. Aut. C o n t r o l , AC-12, pp. 290-293, 1967. [16] Mayne, D.Q., " P a r a m e t r i z a t i o n and i d e n t i f i c a t i o n of l i n e a r m u l t i v a r i a b l e systems", Lecture notes i n mathematics, Springer V e r l a g , 1972. [17] Pearson, J.D., "Unconstrained m i n i m i z a t i o n of a f u n c t i o n of s e v e r a l v a r i a b l e s " , Computer J . 13, pp. 171-178, 1969. [18] P o w e l l , M.J.D., "An e f f i c i e n t method of f i n d i n g the minimum of a f u n c t i o n of s e v e r a l v a r i a b l e s without c a l c u l a t i n g d e r i v a t i v e s " , Computer J . 7, pp. 155-162, 1965. [19] P o w e l l , M.J.D., A survey of numerical methods f o r unconstrained o p t i m i z a t i o n , SIAM Rev. 12, pp. 79-97, 1970. [20] R i c h a l e t , J . , R a u l t , A. and Pouliquen, R. , "Process i d e n t i f i -c a t i o n by the model method" ( i n French), Gordon and Breach, 1971. [21] Rosenblueth, A. and Wiener, N., "The r o l e of models i n s c i e n c e " , Philosophy of Science 12, pp. 316-321, 1945. [22] Rosenbrock, H.H., "An automatic method f o r f i n d i n g the g r e a t e s t or l e a s t value of a f u n c t i o n " , Computer J . 3, pp. 175-184, 1970. [23] SciderstrHm, T., "Convergence p r o p e r t i e s of the Generalized Least-Squares i d e n t i f i c a t i o n method", Automatica 10, 1974. [24] SBderstrcltn, T., "On the uniqueness of Maximum L i k e l i h o o d i d e n t i f i c a t i o n " , Automatica 11, pp. 193-197, 1975. [25] Soudack, A.C, Suryanarayanan, K.L. and Rao, S.G. , "A u n i f i e d approach to d i s c r e t e - t i m e systems i d e n t i f i c a t i o n " , I n t . J . C o n t r o l , V o l . 14, No. 6, pp. 1009-1029, 1971. [26] Wilde, D.J., "Optimum seeking methods", P r e n t i c e H a l l , 1964. [27] Woo, K.T., "Maximum L i k e l i h o o d i d e n t i f i c a t i o n of n o i s y systems." Second IFAC Symp. i d e n t i f i c a t i o n and process parameter estima-t i o n , Prague, paper 3.1, 1970. [28] Zuercher, H., "Linear systems i d e n t i f i c a t i o n a lgorithms f o r a minicomputer", M.A.Sc. T h e s i s , Dept. of E l e c t . Engg.,^University of B r i t i s h Columbia, 1974.