Mathematics and Computers North-Holland
PARAMETER
in Simulation
BOUNDING
527
32 (1990) 527-534
FOR TIME-VARYING
SYSTEMS
J.P. NORTON School of Electronic and Electrical Engineering, Birmingham El5 2TT, UK
University of Birmingham, P. 0. Box 363,
S.H. MO Jaguar Cars Ltd., Abbey Road, Whitley, Coventry CV3 4LF, UK
The computation of parameter bounds for models of linear, time-varying systems is considered. Three methods are discussed, all with counterparts in conventional least-squares or statistical parameter estimation: fixed-memory bounding, scalar bound inflation and bound incrementing. Their computational needs are examined, and some results from artificial and real records are presented.
1. Introduction Parameter-bounding identification [l-3] employs a deterministic description of uncertainty, computing bounds on the model parameters instead of point parameter estimates and their error covariance. The parameter bounds are derived from input-output observations by specifying bounds on the model-output errors, reflecting either prior knowledge of the largest possible output error (from an understanding of the sources of error) or requirements on model performance. The parameter bounds define the set of all values consistent with the observations, the model structure and the output-error bounds: the feasible parameter set. The feasible set may turn out to be empty, indicating that no model with the specified structure can represent the output to within the specified maximum error. No attempt is made to derive point parameter estimates; all feasible values are equally valid. This relatively new approach has particular appeal in early stages of identification, when model-output error may be due largely to a tentative and less-than-ideal model structure, and may thus be heavily structured and difficult to characterize statistically. In the same circumstances parameters may vary with time. An example is when high-order and/or non-linear, time-invariant dynamics are represented by a lower-order and/or linear model whose parameters must vary to accommodate the observed behaviour. The efficacy of time-varying models in suggesting model refinements has been pointed out in conventional parameter estimation [3-S]. Several ways of calculating bounds on time-varying parameters are therefore examined in this paper. Parallels can be drawn with parameter estimation, but the computational problems are quite different. An attractive feature of parameter bounding is that the computation often has a simple geometrical interpretation. 037%4754/90/$03.50
0 1990 - Elsevier Science Publishers
B.V. (North-Holland)
528
J. P. Norton, S. H. MO / Parameter bounding for time-varying
systems
2. Methods of tracking time-varying parameters The basic ideas will be discussed in this section, leaving details for later. Three methods are considered: fixed-memory bounding, scalar bound inflation and bound incrementing, corresponding respectively to fixed-memory filtering [6] in state estimation, exponential weighting-out of past errors in parameter estimation [7] and random-walk modelling of parameter variation [5]. The model considered is &=C$;e,+e,,
t=l,
2 )...)
N,
(1)
where y is the scalar observed output, 8 is the p-vector of parameters, C#Iis the known explanatory-variable vector and e is the output error, all at sample instant t. Bounds ]e,]
t=l,2
,...,
N,
Each observation
(2) then provides
Yt - I; G G? =GY, + r,
bounds (3)
on 9,. Fixed-memory bounding computes parameter bounds from L successive observations y,_ L+ 1 to yt, treating the model as constant-parameter over the L-sample window. The window advances by one sample as each new observation arrives. Successive feasible parameter sets, locally valid over L samples, are thereby found. The window width L is the only item to be tuned to suit the variation of the system. This idea originated to statistical state estimation (optimal filtering), where recursive fixed-memory filtering is little more complicated than ordinary expanding-memory filtering; application to parameter estimation merely requires that 0 be treated as a state vector, with dynamics /3,= O,_ 1 and with (1) as the observation equation. An alternative way to drop old information (in parameter estimation as distinct from bounding) is to decrease the weights attached to older errors in weighted least-squares estimation. Recursive weighted least squares with exponentially decreasing weights [3] may be written as (4) (5) (6) where P,’ is w,P,, P,Yl is w,P,_,, P, is the error covariance of 4, and the forgetting factor p, is the ratio of weights wI_r and wr attached to the squared output errors for times t - 1 and t. The only difference from recursive ordinary least squares is the resealing of P to P’ and P”; algorithmically, only a single resealing, (4), is new. The algorithm is popular in self-tuning control schemes, because it is simple and only one item, II. (normally constant), has to be tuned. Simplicity is achieved, as in fixed-memory filtering, only by failing to admit that parameters may change at different rates and require different memory lengths. The bounding counterpart of recursive ordinary least squares is the ellipsoidal outer-bounding algorithm of Schweppe [8] and Fogel and Huang [9], which updates the matrix P, and centre 0, of an ellipsoid a,=(eI(e-e))TPr-l(e-Q
(7)
J. P. Norton, S. H. MO / Parameter bounding for time-oaTing
systems
529
containing all parameter values consistent with the observations up to time c. The algorithm differs from recursive least squares by an extra data-dependent scaling of P. If we were to inflate P by a further resealing like (4), with a fixed scaling factor, we should have a bounding counterpart of exponentially weighted least squares. Geometrically, the bounding ellipsoid is expanded about its centre by the same factor in all directions. A third way to treat time-varying parameters is as variables, modelled explicitly by a set of state equations. The state equations can embody any prior information about the average behaviour of the parameters, including correlation between changes in different parameters. The simplest useful such state-space model is a simple random walk [3,5] r=1,2
e,=e,_,+WI_,,
,...)
N.
(8)
In statistical parameter estimation, the increments { W} are assumed white, with specified mean and covariance Q, _ 1. In bounding, w,_ 1 is specified as belonging to a set %‘_ i. In both cases, time updates, using the state equation to account for the evolution of 8, alternate with observation updates, bringing in the information in the new observation. For model (8) the statistical effect is to add Q,_, to the covariance Pt_l of 4-i. For parameter bounds, the time update is vector (Minkowski) summing of the set ~9-i of feasible a,_, and set W’_, of possible increments, forming gt’, say. Next we look in more detail at the bounding version of each method.
3. Fixed-memory bounding If computational expense is unimportant, any parameter-bounding method may be used to process separate batches of L observations, staggered by one. A recursive method is more economical. The results reported later were given by two recursive algorithms, modified ellipsoidal bounding (m.e.b.) [lo] and the exact polytope algorithm (e.p.a.) [l]. Recursive processing requires that any bounds due to y,_ L+ 1 to yt_ 1 made redundant by observation ,v_ L be recorded, to be reimposed when the effects of Y,__~ are dropped. Even in e.p.a., substantial extra computing is needed to determine and reimpose all such bounds. In m.e.b., expansion of the ellipsoid E,_, to encompass any later bounds formerly made redundant by the deleted bound is so complicated that no method more economical than batch processing of L observations at a time has yet been found. Similar remarks apply to other recursive algorithms which approximate 9, by some simpler set. The following method of choosing L has been found satisfactory. For each parameter, the maximum possible variation per sample interval, assessed from background knowledge or guessed, is used to plot maximum variation over L intervals against L (a straight line). The mean half-range of each parameter, obtained from a number of L-observation parameter-bounding runs, is also plotted against L. The value of L where the plots cross is noted for each parameter, and the largest adopted. The reasoning is that the memory should be long enough for the parameter range to be attributable equally to time variation and memory length, on average. Figure 1 shows how 35 is chosen as L for the model (used also in Section 6) yt = -n,y,_,
- a,~,_,
+ big,-i
+ b2q2
+ e, = ~$0 + e,
(9)
530
J.P. Norton, S.H. MO / Parameter bounding for time-varying systems
Fig. 1. Selection of memory length L for fixed-memory bounding: mean half-range of each parameter (over 5 runs in this simulation example) decreases as L increases, and L is taken as max value at which half-range equals max expected variation over L intervals.
with {e} white and bounded by f 1, and actual 8 equal to [1.1 - 0.28 1 0.6(1 + 2 cos Trt/300)lT. The mean-square signal : noise ratio was 10 and the true bounds on { e} were specified.
4. Scalar bound inflation As noted above, inflation of ellipsoidal feasible-parameter set 8 by a factor l/p is algorithmically a close bounding counterpart to exponentially weighted least squares. The results reported below are for scalar expansion of ellipsoidal bounds. Any other approximate parameter-bounding figure with a centre of symmetry, such as the orthotopes (boxes) employed by Milanese and Belforte [ll], could heuristically be inflated in the same way. An alternative procedure is motivated by interpreting weighted least squares as minimum-covariance-of-parameter-error estimation when the weights are proportional to the inverses of the variances of corresponding error samples and the error sequence is white. Exponentially weighted least squares may then be regarded as implying that if 8, were used retrospectively in _Yt_k=+;_&+et_k,
ka.0,
(10)
the standard deviation of et-k would be l/p’/* times that of e,, increasing exponentially with k. The bound counterpart of the standard deviation of statistical noise is the output-error bound r,. The distance between the two bounds (3) resulting from observation y, is 2rJ II& 11,proportional to r,, so exponential increasing of the model-output-error bounds into the past increases the distance between the two parameter bounds exponentially as they get older. Such bound expansion is conceptually simple but computationally a heavy task, since it requires shifting and reimposition of all earlier bounds, not just those currently active. Expansion of current active bounds alone may cause a newly active bound to be ignored, missing a reduction of the feasible parameter set. Dasgupta and Huang [12] discuss ellipsoidal-bound updating for time-varying parameters. They renormalize the updating equation of Fogel and Huang [9] for updating the centre of the ellipsoid and redefine notation so as to express the updating in terms of a forgetting factor. It can, in fact, be shown that the rewritten updating equations are exactly equivalent to the originals, with the updated ellipsoid parameterized by a scalar quantity which may be chosen at will. The original and revised versions differ in how this quantity is chosen. Of two choices
J. P. Norton, S. H. MO / Parameter bounding for time-varying systems
531
discussed in [9], the most popular minimizes the hypervolume of the updated ellipsoid. The choice in [12] is claimed to minimize a bound on the parameter estimation error, but the quantity minimized is not a bound on the parameter error, nor does it bear a simple relation to it. The resulting algorithm does, however, allow a relatively simple analysis of convergence properties. The ability of the ellipsoid centre to track staircase variations in all parameters of a 4th-order autoregressive process is found to exceed that of recursive-least-squares estimates (apparently without a forgetting factor) in a simulation reported in [12].
5. Bound incrementing based on explicit parameter-variation Consider time-updating the ellipsoid ,g;‘-, = contains
of the parameters
(e l(e -
to accord
4_l)T(2(P,-1 + Q))-'(0
model
with model (8). Schweppe
- ‘t-1)G 1)
[8] shows that
(11)
the vector sum of gt_.i (given by (7) with c - 1 for t) and an ellipsoid W= ( wt-i (~;r_iQ-~w~_i < 1)
02)
containing all possible parameter increments from time t - 1 to time t. Thus time updating could merely involve incrementing P,_i to 2( P,_i + Q,_,), almost exactly parallelling the update from P,_ i to P,_ i + Q,_ 1 in statistical estimation. However, (11) generally gives a loose bound on the vector sum, and cumulatively makes the ellipsoid far too large. A tighter bound is (1 + a)P,_ 1 + (1 + l/a) Q,_,, with (Yfound as by Chernous’ko [13]. It is the unique positive root of &&&=
a(1 p+ a)
03)
where diag( d,, d,, . . . , d,) is the diagonal matrix AQ,AT obtained by the unique transformation [14] that also gives AP,_,AT = I. The algorithm may still give a loose outer bound on 8,. For instance, when a linear combination of the parameters varies rapidly, W is a line, long compared with the width of a,__, in the same direction. The vector sum of E,_, and W is then two half-ellipsoids connected by a long ellipsoidal prism, not another ellipsoid. Time-updating of box bounds on the parameters, with W a box with the same axes, merely requires relaxation of each prior parameter bound by the corresponding increment bound. Boxes are also likely to give loose bounds, though, and vector summing is not simple for boxes with non-coincident axes or for general polytopes. Time-updating a polytope g,_i by a polytope K_i gives a much more complicated polytope LYE’in general, since the vector sum of two polytopes is the convex hull of the vertex set formed by all sums of one vertex from each polytope [15]. One possible way to retain the flexibility of polytope bounds but avoid excessive complexity is to move each ( p - l)-face of g,_ i by the extreme point of K_, in the direction normal to that face, then extend the faces as required to close the new polytope. The result is an outer bound to the vector sum, since the extensions include infeasible points. This technique may be viewed as satisfying the relation [15] 04)
532
J. P. Norton, S. H. MO / Parameter
bounding for time-varying
systems
between the support functions of P’,__1, K_ 1 and P:‘, for directions 9 normal to the ( p - l)faces of P,_r. Some results of incrementing ellipsoidal bounds by ellipsoidal increments are given below. 6. Computational
results
Tests were performed with records generated by (9), all details being as in Section 3. Figure 2 shows fixed-memory bounds for the time-varying b, in (9), computed by e.p.a. with memory L fixed at 20. The bounds reveal the variation of b, clearly, but fluctuate strongly and almost clash at times; a model constant over 20 samples is barely able to explain the observation to within the specified error. The bounds do clash for any L above 25. Little spurious long-term variation occurs in the bounds on the other parameters (not shown).
t Fig. 2. Fixed-memory exact-polytope bounding: bounds on varying parameter b, for memory = 20.
Fig. 3. Fixed-memory modified ellipsoidal bounding: bounds on b, for memory = 20 (tighter bounds) and 80.
x lo-’ 3or.
1
25 01
"1 -10.
20 15 10
-15.
5
-20.
0 -5
-25. -3OL 0
100
t
200
300
-10’
100
t
I
200
300
200
I 300
x10-’ 16 16
12
a2
'8 6 L ; 2 r"6 -8 -10 -12 -IL
100
200 t
300
-3’
.
100
t
Fig. 4. Bounds on all parameters of (9) found by scalar inflation of ellipsoid: l/p = 1.016 (tighter bounds) and 1.03.
J. P. Norton, S. H. MO / Parameter bounding for time-varying systems
533
LOO 300 200 b2 100 0 -100 -200
Fig. 5. Bounds
on b, found by modified
100
t
200
ellipsoidal bounding with incrementing + 0.005 (tighter bounds) and + 0.0141.
300
in b, direction:
increment
bounds
Fixed-memory b, bounds yielded by m.e.b. are shown in Fig. 3. For L = 20 they are so loose and erratic that the actual variation is obscured. A longer memory, L = 80, smooths them but imposes noticeable lag. The m.e.b. bounds do not clash until L = 100, evidence of how loose they are. Figure 4 gives bounds on all parameters of (9), found by m.e.b. with scalar ellipsoid inflation. The inflation factor l/p is 1.016 and 1.03, bracketing the value required for long-term constancy of the bound range. These bounds track medium-term variation better than those of Figs. 2 and 3, and overall variation at least as well. Figure 5 shows the bounds on b, obtained by ellipsoidal increment bounding, for increment bounds kO.005 and + 0.0141 on b, and zero on all other parameters. Again spurious overall variation of bounds on the other parameters was small. Incrementing P,_ i suboptimally to 2( P,_ i + Q,_ i), as suggested by Schweppe [8], was found to loosen the bounds too much.
7. Conclusions Fixed-memory bounding is inconvenient and offers no advantage over bound inflation or incrementing. Experimental results show little to choose in performance between the last two methods, but incrementing is more flexible and can exploit prior knowledge. Exact updating of the parameter-bounding polytope does not seem practicable for a time-varying model of realistic size. Reasonably neat and economical incrementing of ellipsoidal outer bounds is possible. Similar conclusions have been reached in extensive trials of these methods on hydrological records with strongly time-varying characteristics [16].
References [l] S.H. MO and J.P. Norton, Fast and robust algorithm to compute exact polytope parameter bounds, Math. Comput. Simulation 32 (5&6) (1990) 481-493 (this issue). [2] J.P. Norton, Identification and application of bounded-parameter models, Automatica 23 (1987) 497-507. [3] J.P. Norton, An Introduction to Identification (Academic Press, New York, 1986). [4] P.C. Young, Recursive Estimation and Time-series Analysis (Springer, Berlin, 1984). [5] J.P. Norton, Optimal smoothing in the identification of linear time-varying systems, Proc. ZEE 122 (1975) 663-668.
534
J.P. Norton, S.H. MO / Parameter bounding for time-varying systems
[6] A.H. Jazwinski, Limited memory optimal filtering, IEEE Trans. Automat. Control AC-13 (1968) 558-563. [7] R.H. Jones, Exponential smoothing from multivariate time-series, J. Roy. Statist. Sot. Ser. B 28 (1966) 241-251. [8] F.C. Schweppe, Recursive state estimation: unknown but bounded errors and system input, IEEE Trans. Automat. Control AC-13 (1968) 22-29. [9] E. Fogel and Y.-F. Huang, On the value of information in system identification-bounded noise case, Automatica 18 (1982) 229-238. [lo] G. Belforte and B. Bona, An improved parameter identification algorithm for signals with unknown but bounded measurement uncertainty, in: Proc. ZFAC/ZFORS Symp. on Identification and System Parameter Estimation, York (1985) 1507-1512. [ll] M. Milanese and G. Belforte, Estimation theory and uncertainty intervals evaluation in presence of unknownbut-bounded errors, linear families of models and estimators, IEEE Trans. Automat. Control AC-27 (1982) 408-414. [12] S. Dasgupta and Y.-F. Huang, Asymptotically convergent modified recursive least-squares with data-dependent updating and forgetting factor for systems with bounded noise, IEEE Trans. Znform. Theory IT-23 (1987) 383-392. [13] F.L. Chemous’ko, Optimal guaranteed estimates of indeterminacies with the aid of ellipsoids I, Engrg. Cybernetics 18 (1980) l-9. [14] G. Hadley, Linear AIgebra (Addison-Wesley, Reading MA, 1961). [15] B. Griinbaum, Convex PoZytopes (Wiley, London, 1967). [16] S.H. MO, Computation of bounded-parameter models of dynamical systems from bounded-noise recorak, Ph.D. Thesis, School of Electronic & Electrical Eng., Univ. of Birmingham, 1989.