OPTIMAL ESTIMATION IN UNCERTAINSTOCHASTIC DISCRETE-TIME SYSTEMS ANDREI V. BORISOV Department of Industrial & Operations Engineering University of Michigan Ann Arbor, MI 48109-2117 ALEXEI R. PANKOV Department of Applice Mathematics Moscow Aviation Institute Moscow, 127080, USSR Technical Report 91-27 October 1991

OPTIMAL ESTIMATION IN UNCERTAIN-STOCHASTIC DISCRETE-TIME SYSTEMS Andrei V. Borisov Department of Industrial & Operations Engineering University of Michigan Ann Arbor, MI 48109 and Alexei R. Pankov Department of Applied Mathematics Moscow Aviation Institute Moscow, 127080, USSR ( Revised on September 30, 1991 ) Abstract In this paper we consider estimation problems in uncertainstochastic discrete-time dynamic systems. To solve these problems, we prove necessary and sufficient conditions for the identifiability and optimality of parameter estimates in uncertain-stochastic linear regression. Using these results we derive optimal filtering and smoothing algorithms. We also present a suboptimal filtering algorithm for nonlinear uncertain-stochastic systems. Keywords: uncertain-stochastic system, filtering, two-filter smoothing. Section 1. Introduction. The Kalman filter is known as the best mean square state estimator if true values of the input signal characteristics are known. These estimates are used widely in the such applied areas as a navigation, radar tracking, data processing in complicated electronic systems, numerous control problems etc. [ 1 ], [ 2 ], 1

[ 3 ],[ 4 ],[ 5 ]. In many practical situations a priori information about the input signal characteristics is incomplete. In order to overcome this serious obstacle one may use adaptive or robust filters [ 6 ], [ 7 ], [ 8 ]. All of these filters use a limited amount of available a priori information. So the corresponding estimates may be very inaccurate if the lack of a priori information is serious. It is however well-known that very often it is possible to observe some components of the input vector. Our aim is to derive the recursive filtering and smoothing algorithms which allow us to obtain acceptable results even in the case when we know nothing about the characteristics of the input signal. That algorithm was presented in [ 9 ]. The optimal filtering algorithm described in this paper is a generalization of the previous one in the case of complex observations containing the information about both the state vector and the input signals. The results below make the construction of an optimal twofilter smoothing algorithm for uncertain system processes possible. Similar results for stochastic systems are given in [ 10 ]. The suboptimal recursive filter for the uncertain-stochastic nonlinear system is also considered The results of the numerical experiments are given. Using input signal observations in the presented algorithm allows us to obtain accuracy of the new filtering estimates under a priori input uncertainty close to ones of the Kalman filter with known input signal characteristics. It also turns out that the smoothing estimates are significantly more accurate than the filtering ones. Thus, smoothing is preferable than filtering in the case of postexperimental data processing. Section 2. Problem formulation. Consider the dynamic system given by the following difference equation Xt = atxt-h + btut + tt v t = h,2h,....; x = u ( 1 ) where xt ~ Rn is a state vector; vu Rn is an uncertain initial 2

condition; ut E Rm is an uncertain input vector, t R" is a random input vector with known characteristics E ( t ) = 0; cov( St,as ) = 8s C; h > 0 is a time increment. Here and further we suppose v and ut to be unrestricted and nonrandom. Let us consider the following observation model Yo = pou + o t ( 2) l yt = t"u +, txt + qtet + ct t t t = h,2h,..., where,t E R" is a vector of random disturbances with known characteristics E (wt ) = 0; cov( tcos ) = 8 Qt Qt > 0 for all t - 0. Matrices a t, bt, ot, t ) of appropriate dimensions are known. We suppose ( t ), ( wt ) to be independent. Our aim is to construct a linear estimate xt given observations ( ya, 0 s: - t ), which minimizes the following criterion t = s u p E (( x - xt )T ( x - xt ), ( uT )O, v where 2 > 0 is known symmetric weight matrix. (3) Section 3. Parameter identifiability of the uncertainstochastic linear regression model. In this section we state some preliminary results concerning the general type of the linear regression model. Consider the multidimension linear regression Y= U + H +, (4) where U E Rq is an uncertain ( nonrandom ) vector; s. ~! Q E~ R1 are independent random vectors with known characteristics E( Z ) = 0; Ef 0 } = 0; cov( Z, s ) = C; cov(, Q ) = Q > 0. Our aim is to find the best linear minimax estimate ( BLME ) x = aY of x = A U + D i, which minimizes the following criterion Ax) J(x) = s u p E ( ( x - x ) x - x ), (5) U Let us call x identifiable if there exists at least one estimate 3

x = p(Y) ( (.) is some measurable function of Y ) for which ^# AV J( x) < Proposition 1. i) If x is identifiable then E A = ZA (6) where + is the matrix 0 pseudoinverse; ii) if x is identifiable then x is BLME if and only if E cov( x - x, Y- U (7) where U is an arbitrary estimate of U, which satisfies EE(x-AU)=0; iii) BLME x is also the best estimate within the class of Ah J -optimal unbiased estimates, if Z and C2 are Gaussian random vectors. The proof is given in the Appendix Section 4. Optimal filtering algorithm. First we find the conditionally-optimal solution within the class of all linear unbiased recursive filters xt = atxt-h + gtyt, E( xt - t ) = 0. Denote,Lt = tbt + ~t, et = at + )t, and let At xt - Xt be the estimate error, where t= xt-X = (at-at-ttat) Xt-h +t +att-h. ( bt-bLtgt) Ut -9t(ctet+wt) The following conditions are necessary and sufficient for xt to be unbiased at = ( I - Stat )at t = btt + Zt( I - / ). 8 ) btgt+ = bt where zt is an arbitrary matrix of appropriate dimension. Then At= -t zt( I - ity )vt, where v t = +tt +t *tatAt-h, tt = atat-h +t -btgvt, cov(vtt)=At 4

Minimizing criterion Jt and taking into account ( 8 ) we obtain Z = (akt haT T +CtT -b At) (_-gtg) ((I-tt) At( I-t kt = cov( At,t ) = coV( Ctct ) - cov( <tVt )( I-ttII ) T (( I-Agt) A t( I-Agt) T)+( I-_tt) COV( vtt ) ( 9) From (9) and the matrix Schwartz inequality [ 7 ] it follows that A+ tLTA-1 +. -1 t= ( t t ) t )+ t1 In the case when o(To > 0, the initial conditions are T -1 -1 T,,.~-I XO = ( Po0 o 0) PoQ Yo ko =( O 0 )-' ( 10) Hence the filtering algorithm is as follows Xt = atxt-h kt = atkt-haT + C t = t b t + st, = p t + t Pt = ktt Ct i +t A= tTatkthaITt tctt + Qt ~ t = ( A ) >TA- ( 11 ) t +[ bt + PtA ) ] C Yt - taxt 1 t k= k +btL At ( btlLf)T -b t p, -( b p )T -pT A-1p Pt = ktEt + t Cr; t 7 t t t.-ttt tt~ t Proposition 2. i) The following conditions are sufficient for xt to be the Jt- optimal estimate I bIf = b t, for all T E ( O,t ], l(P T0 > ~ i ii) the estimate xt given by 10 )( 11 ) is the solution of the filtering problem for ( 1 )-( 3) not only within the recursive estimates, but also within all linear est imates, iii) x is the best unbiased estimate if ( w),( ~t) and wo 5

are Gaussian. The proof immediately follows from the fact, that xt given by (10),(11) satisfies the conditions ( 6 ),( 7 ) of Proposition 1. Section 5. Fixed interval smoothing algorithm. Consider the dynamic system ( 1 ) without the random component ( ~t ) of the input signal Xt = atXt-h + btut, t = h,2h,...; xo = v ( 12 ) In this case ( 12 ) defines a pure uncertain system. Let us also consider an observation model of the special type Yo= POP + lo t yt = put + )t, t = h,2h,...,T-h,T. Zt= Vtxt+ 2, ( 13 ) YT.= fPTXT+ (OT' where VT > 0, t 0,...,T; cov( o1 ) = 6t; Qt > 0 cov( W<,( ) = P P > 0; cov( u2 ) =0. Consider the problem of finding a Js-optimal linear estimate xs of t, t E [ 0,T ] given all observations ( 13 ) where Jt s u p E (x Xt )T- ( xt - x t t t t) T ( u )}oP Using ( 12 ) we may easily derive the reversed-time dynamic system Xrh= atxrt - au, = x, t = T, T-h,...,0 (14) which is pathwise-equivalent to (12 ). By applying ( 10 ),( 11 ) we obtain from ( 12 ),( 13 ) the forward-time estimate T -1 T -Iy0 X o = ( '0 o )Qo Yo ko = ( oQ'Io )-1 (15) 6

Xt = atxt h + bt( t'Qt 1 -1 Ttt =ak T T 1 I T (1) k t = atkt haTt bt ( (tQ ~t b) t xt = x ktt( ktt + Pt)- zt - t ] 17 k = kt -ktT( tktT + Pt) 1t, Analogously from ( 14 ) it follows that the reversed-time estimate ir is given by the equations xT ( TQT T T ) 'TQ t ' r IT -1 -1 kT= ( ^pTQT T )-, ( 18 ) t-ht - at t t b) tQt y, kr =a-lkr ) + ab kt T t Q )b (a ( )T 19) t-h t t (at) t Pt ) ^r _ r + kr T T + - r Xt-h- t-h kt-h<-h( /t-h t-ht -h t- [ Zt-h- *t-hXt h] ' k r = k k r T r T. P _ h k-l r/ 20 t t t -h t-h( tt-hk t-ht h t-h) -h ( 20 ) From ( 15 )-( 20 ) it follows that xr is the BLME of xt given all " past " observations and jr is the BLME of xt given all " future" observations,i.e. all observations obtained after the moment t up to T. The following proposition describes the connection between xt, xt and the best linear smoothing estimate xt Proposition 3. X = k ( kt )-1 Xt ( kt], 21 kt. = [ ( kt )-1 +( ktr-1 ]-1 where kt, kl and k; are error covariance matrices of xt, and x respectively. The proof of Proposition 3 also follows from Proposition 1. Section 5. Nonlinear filtering algorithm. Consider the nonlinear discrete-time uncertain-stochastic 7

system Xt = a Xt_ht) + b(t)Ut + et ' Xo = ( 22 ) where a( x,t ) is a known differentiable function with respect to x; v and ut are the uncertain initial condition and input signal respectively; ( et ) is the same as in ( 1 ). The observation model is given by Yo = PoV + ~o ' t ' yt = ptut + Wt, t = h,2h,... Z = ^xt, t) + t ( 23 ) where Tp1 >0, t = 0,h,2h,...; ^ x,t ) is a known differentiable function with respect to t. We consider the problem of suboptimal estimation of the state vector in ( 22 ) given observations ( 23 ). Using the linearization and the optimal linear uncertain-stochastic filtering equations ( 10 ),( 11 ) we obtain the following suboptimal nonlinear recursive filter T1 T -1 0TIo31 0= ( OQo '0 ) o k PT IgT -1 T-1 Xt a( xt_-ht) + b(t)( vtQ tl1t ) -pTQt yt kt = AtkthAT b(t)( tQt t -'bT( t) + Ct, _x^(T+ - ktk, - <(Pt t t ( 24 ) ( 25 ) (26 ) where A = ( aa/ax )l _ t '"'t-h,t = ( a(/ax )I Equations ( 24 )-( 26 ) may be considered as the analog of the extended Kalman filter for the nonlinear difference dynamic system ( 22 ),( 23 ) with the uncertain-stochastic input signals. 8

Section 7. Numerical example. 1.Let us consider a motion of the aircraft mass center x[ I f x 0 1 v =2 0- + + 2 '2 X _ [ xt] [ ~ i ] [ x - ] [ ] [ XO2 ] V[ " where xi is a distance between the aircraft and the origin of coordinates, x is a speed, ul is an unknown component of the acceleration which may be interpreted as uncertain input signal ( control ). Unobservable random disturbances et we assume to be a zero-mean discrete-time white noise with the covariance matrix C. Information concerning ( xt ),( x2 ) and ( ut ) is given by the observations: yI I 1; I2 + (28) y0=v +. I; = (28) t= t + ot; Y = Xt+ where ( wt ) and ( it ) are zero-mean discrete white noises with the covariances Q and P respectively, coo = (co,~ T is a zero-mean Gaussian vector with the covariance matrix Q = diag ( qO2 ),h is a time increment. wo I ( et )' ( wt and ( t ) are independent. The filtering problem is to calculate recursively the Jt-optimal estimate of the vector (x,x2)T given observations ( 28 ). The following parameter values were selected for the numerical tests: Q = 225; P = 400; q= 18 108; T = 20;t = 0.; X0 = 5000; x- 100; u(t) =20 + 10 sn( O.1lt ); C 33.3 50.0 0 C = L 50.0 100.0 J Selected parameter values Q and P correspond to accuracy of the on board acceleration sensor and radar respectively, noise et simulates the atmospheric disturbances, values xl, x and function u(t) describe one of the possible modern aircraft manoeuvres. We use the algorithm ( 10 ),( 11 ). Table 1 gives the estimation 9

results: ol and a2 are the estimation error standard deviations of the distance and the speed respectively. oa and o2 are the estimation error standard deviations obtained by the Kalman filter in the ideal situation when the full information about ut is a priori available. From the Table 1 it follows that the estimates given by the algorithm ( 10 ),( 11 ) are close enough to those of the Kalman filter and may be used in the case of unknown but observable input signal when the Kalman filter is useless. Time s o1 o2 2 0 104 104 104 104 1 20.0 20.0 56.85 57.08 3 15.12 15.49 14.63 18.82 5 14.39 15.17 14.09 18.70 7 14.39 15.17 14.07 18.68 9 14.38 15.17 14.06 18.68 15 14.38 15.17 14.06 18.68 20 14.38 15.17 14.06 18.68 A, Table 1. Standard deviations of the distance ( and speed ( 02 ) estimation errors for the filter ( 10 ),( Al ) and the corresponding values a1 and r2 for the Kalman filter. 2. Let us now consider some control system given by x - 0.4x + 0.16x = uI + 0.4 uI + u2, x(0) = xo,x(0)= x ( 29) where ( u1 ),( u2 ),xo and xo are unknown input signals and initial conditions. Equation ( 29 ) may be rewritten as the first order differential system 10

X1 =X2 + U1 a.~~~~~ ~~( 30) x2 = -0,16x1 - 0,4x2 + U2 where x(t) = x1(t). The observations are given by Y1 = U + U2+ ()1 Y2 = U1 - UI + 32 t t E [ 0,3 ] ( 31 ) Z= X + 1 T where ( ) } is the standard Wiener process, ( (co1,w ) is the zero-mean Wiener process with the differential covariance cov( O, ) = 1 0,7 L0,7 1 J The known parameter L characterizes the accuracy of the observat ions. The fixed-interval smoothing problem is to calculate Jt-optimal estimate xt for all t E [ 0,3 ]. First we discretize the model ( 30 ),( 31 ) with some small time increment h and then use the algorithm ( 21 ).Table 2 gives the estimation results: caF, 2, 3 are the filter estimation error standard deviations obtained by gL = 0,1; jt = 0,5; tL = 4 and Ioas,3 are the 1 2' 3 smoothing estimation error standard deviations respectively. From the Table 2 it follows that the smoothing estimate is significantly accurate than the corresponding filter estimate. 11

Cr a$ |I 1 | I 1 2 2 |3 0.0 o 1.305 1B 1.426 o 1.736 0.3 3.511 1.047 3.515 1.098 3.532 1.236 0.6 2.474 0.843 2.486 0.897 2.536 1.076 0.9 1.995 0.704 2.017 0.796 2.108 1.041 1.2 1.701 0.637 1.736 0.759 1.874 1.036 1.5 1.496 0.637 1.546 0.760 1.732 1.037 1.8 1.340 0.681 1.409 0.783 1.641 1.038 2.1 1.218 0.747 1.306 0.822 1.582 1.045 2.4 1.118 0.820 1.227 0.884 1.542 1.075 2.7 1.030 0.891 1.165 0.978 1.515 1.182 3.0 0.957 0.957 1.117 1.117 1.497 1.497 Table 2. Standard deviations of the filtering estimate error ( al, 02.03 ) and the corresponding values a,2 oa3 for the smoothing. Appendix. Proof of the Proposition 1. i) Let x be identifiable,i.e. there exists an estimate x = p(Y): J = s u p E( II x - (1) 12 ) = J~ < o, U where II a 12 = aTE a. Then the corresponding error A is A=x- (Y) =AU +D - p(Y). Then J = s u p E( 11 m Il2 + 1 A - m 112 ), where m = E( A ) is the U biase of q(Y). Obviously E( p(Y) ) = 8( V), where V = E( Y) = * U and 9( V) is some nonrandom function. Then m = A U - ( V). Denote R( i ) the linear span of the vector-columns of 0. Suppose V* E R( ) and U = ( U: h U = V* ). Then U = U*+(I-+) Z 12

where U* ~ U and Z is an arbitrary vector from Rq. Let m*= A U*- O(V*) and G = A( I-+0 ), then for every U E U we have m = A(U*+ ( I-0r)Z ) - 8(V*) = m - GZ. E( It x - (Y) 112 ) > 11 m 112 = I m - GZ 112 for any Z E Rq. If EG * 0 then 1(Z) = I* - GZ 12 = ZT( G Z GT)Z - 2(m*)TEGZ + (m*)TE m* is Zo 2 jo unrestricted from above.So there exists Z~: I| m* - G Z012 > J Thus, for U~= U* + ( I-0* )Z~ and Y0= 0 U~ + H " + Q we have E( II x - p(Y) 112 ) > 1 - m* - G Z~l2 > Jo, which is impossible. From this contradiction, it follows that 2 G = E A( I-0+0 ) = 0. ii) Suppose that Z* is the optimal matrix coefficient,i.e. x = A+Y + Z*(I - +)Y is the BLME. For arbitrary Z = Z+ 3Z we have the corresponding estimate x = A+Y + Z( I - 1+ )Y = x + 5Z( I - 01+)Y. J(x) = J(Z*) + J2(Z) - J3(Z"Z), where J1(Z*) = E( 1 x- x 112) = (x), J2(Z) = tr ( cov(3Z( I - *+)Y,YZ( I - *+)Y) } 0 ( Al ) J(Z*,Z) = 2tr( cov(( I - +)Y,x - x )S5Z ). Note that ( I - 4+)Y = Y - +U, where U = ( 0+ + F( I - 4+))Y is the estimate of U and E( x -A U ) = 0, F is an arbitrary matrix. Then J3(Z*wZ) 2tr( cov(Y- 0 U,x - x )5Z ) ( A2 ) If Scov( Y- 0 U, x - x ) O0 then from ( A2 ),( A3 ) it follows that there exists 6Z~: 6s J 2(Z - J3(Z*,IZ~) 0. Hence J(x) = J(x) + J < J(x) which means that x is not optimal. So Ecov( Y- U, x - x ) = 0 is necessary for x to be BLME. 13

Now suppose that for some Z* cov( Y - * U, x - x )E = 0, then J3( Z,Z ) = 0 for all 6Z, and for x = x + 8Z( I - 4)Y J(x) JI(Z*) + J2(SZ) - Jl(Z) J1(X). Hence x is the BLME. Now it is easy to check that ( 7 ) is fulfilled if x = ( A +Y + D C HTP-( I - +) )Y ( A3 ) where P = cov(Y,Y) = H C HT + Q. iii) Let the unbiased estimate x = qt(Y) be given by ~(Y) = x + (Y), where E { p(Y) ) = 0. Denote v = D Z - A + D C HT-1( I - +) ][ H, w = [ A + + D C HTP-( I - ) ] y, then A = x - x = v - (Y). The bias is E(p(Y) )= y( y)[ (27r) det(P) ]1/2exp( -1/2( y-w) TP-( y-4w) )dy = 0 R' for all w E Rq.Then the function qp(Y) has the following property: (y0~) = ( Yo + w ), w E Rq, YO - R1 ( A4 ) Let Y = ( 1-4 O+)Y + 40 I+Y, 0 +Y E R( 0 ). Using ( A4 ) it may be shown that p(Y) = -( ( (I-0 0 +)Y). So A = v - jp( ( I-I 4+)Y ) and the optimal nonlinear x estimation problem is equivalent to the vector v estimation problem given the observations ( I-4 *+)Y, i.e. q(Y) = E ( v ( I-4+)Y ). But cov( v, ( I-00+)Y ) = O, hence y and ( I-0+)Y being Gaussian are independent and E v I ( I-+ )Y } = E ( v ) = 0. So there is no unbiased estimate that is more accurate than the linear one in the case of Gaussian disturbances. 14

References [ 1 ] Sage, A.P. and Melsa J.L " Estimation theory with applications to communications and control," New York: McGraw-Hill, 1971. [ 2 ] Hashemipour, H.R., Roy, S. and Laub A.J., " Decentralized structures for parallel Kalman filtering " IEEE Trans. Automat. Control, vol. 33, pp. 88-94, 1988. [ 3 ] Grewal, M.S., Henderson, R.S. and Miyasako R.S., "Applications of Kalman filtering to the colibration and alignment of inertial navigation systems " IEEE Trans. Automat. Control, vol.36, pp. 4-14, 1991. [ 4 ] Cortina E., Otero D. and D'Attelis C.E. " Maneuvering target tracking using extended Kalman filter " IEEE Trans. Aerospace and Electronic Systems, vol.27, pp. 155-159, 1991. [ 5 ] Gibson, J.D., Koo, B. and Gray, S.D " Filtering of colored noise for speech enhancement and coding " IEEE Trans. Acoust.,Speech,Signal Processing, vol.39, pp.1732-1743, 1991. [ 6 ] Goodwin, G.C. and Sin, K.S., " Adaptive filtering, prediction and control," New Jersey: Prentice-Hall Inc., 1984. [ 7 1 Fomin V.N., " Recursive estimation and adaptive filtering," Moscow, Nauka, 1984 ( in Russian ). [ 8 ] Kassam, S.A. and Poor, H.V.," Robust techniques for signal processing," Proc. IEEE, vol.73, pp 433-481, 1985. [ 9 ] Borisov, A.V., Pankov, A.R., and Sotsky N.M.," Filtering and smoothing in uncertain-stochastic systems with partially observable input signals," Autom. and Remote Control, N 3, 1991 ( in Russian ). [ 10 ]Wall, J.E., Willsky, A.S. and Sandell, N.R., " On the fixed-interval smoothing problem," Stochastics, vol.5, pp 1-42, 1981. 15