Estimation and Identification Using Interval Arithmetic

Copyri ght © IFAC Identificati on and System Parameter Estimation , Budapest, Hungary 199 1

ESTIMATION AND IDENTIFICATION USING INTERV AL ARITHMETIC S. M. Markov and E. D. Popova Division 0/ Mathematical Modeling in Biology, Institute 0/ Biophysics, Bulgarian Academy o/ Sciences, Acad. C . Bonchev str. , Bldg 21 , I I 13 Sofia, Bulgaria

Abstract. We formulate estimation and identification problems in the context of curve fitting and interpolation under measurement uncertainties. We consider the computation of the bounding functions for the solution sets and give two examples of the use of an interval arithmetic approach. Keywords. Estimation; identification; interpolation; interval arithmetic; curve fitting; polynomial approximation; uncertain measurements; unknown but bounded errors.

INTRODUCTION

variable can be neglected, then the traditional strategy is to determine (to identify) the parameter A in such a way that h(z,) = y" i = 1, · · · , n , (3)

We consider the problems of curve fitting and interpolation in the presence of unknown but bounded errors. Let fA : D -> R, D eR!", be a function depending on a real argument E D eR!", and let the model function h depends on a parameter A E A C RI<. In what follows we shall 888ume that f(A; e) = fA(e) is continuous both in A and and that A is a connected compact set.

if such A does exist at all.

e

In the literature the identification problem is often stated in terms of a model operator F = F(M) : A x R!"x" -> R!', which maps the data z = (Z1>' . • , z,,) E R!"x" into the nvector y(z) = (h(Zl)" '" fA(z,,» by means of fA EM. A commonly used class M is the class £ of all linear on A functions of the form h(e) = E~=l A,'Me). For this cla.ss the operator F F(£) satisfies F(£; A, z) A(z»., where the n x k-matrix A(z) is given by

e

Let us a.ssume that experimental data are provided in the form (Z"II,), z, = (Zl" Z2,,"', zm,)T E Rm, 11, E R, i = 1, ... , n. The problem of curve fitting can be roughly formulated as the problem of finding a parameter A, such that

h(z,)

~ 11"

i = 1, .. ·, n,

=

=

(1)

that is, the computed values h(z,) by means of the model function h should be close to the measured quantities y, . When deciding what to mean by "close" , the following arguments can be considered: 1. the measured variable y may be of stochastic origin; 2. the model function may have to be simple, with a small number of parameters in comparison to the number of data, i. e. k <: n; 3. the measurements may contain additional errors of non-1ltochastic origin (e.g. systematic errors due to the imprecision of the measurement instruments) . When arguments 1. and (or) 2. are prevailing, an usual A E A} of strategy is to look for a class M = {h smoothing model functions, depending on some parameter A E A (which type may be predetermined by the particular modelling situation), and then to determine (to estimate) the parameter A by means of an operator

" " z,,), y = (Y1" '" V,,) :

Both the identification operator F and the estimator e2, . .. , eft) and assume that it is bounded. Then assuming that y, are the true values of the measured variable (which can be in particular of stochastic nature) we may write in (1)

h(z,)

(2)

where 1ej 1< ei,

Alternatively, when the stochastic nature of the measured 221

£

~

=

y; + e;

, i

(£1, £ 0, .. •

=1, · .. , n ,

,en) >

(4)

O. Or, if we denote

the intervals [-~" +~,) = E, we may write

where the set A is given by (6) or (8) depending on the problem considered. Our results are formulated in terms of interval arithmetic. The basic concepts of interval arithmetic used by us will be given below.

The intervals Y; can be considered as measurement intervals, containing with guarantee the true values of the measured quantities. It seems that experimental scientists can provide such intervals in most situations (what is the use of the experiment, if they can not!) .

Interval Arithmetic: Basic Concepts Under the assumption that h(e) = !(A ; ~) is continuous on A, equality (10) defines an interval-valued function L : R!" -+ I R, where by I R we denote the set of all intervals

We now have two alternative ways to interprete relation (5) : JA(z,) ::::: y; , which says that the values JA(z,) should be close to the measurement intervals Y; . Either, "close" does not necessarily mean ''is contained in", or alternatively, "close" means "E" . Following the first interpretation, we can compute the set of all A, produced by the estimator 4>M(Z, V) in (2), whenever the numeric value V varies in the interval measurement Y (YI ," ' , Y.. ), that is the set

Y of the form Y

lL

~

y.

Throughout the paper we shall use the following two simple interval arithmetic operations (Alefeld and Herzberger, 1974; Moore, 1966, 1968) for addition of two intervals X, Y E I R and for scalar multiplication by a real number O'ER:

=

X+Y a~O,

aX

(6) The set A~ is called the estimate uncertainty set (Milanese, 1989; Milanese and Belforte, 1982).

(11)

0'<0.

We shall also make use of the following simple applic&tion of interval arithmetic. Given a real valued vector a = (0'1>"', a .. ) and an interval valued vector Y = (YI, "' , y .. )T we can write

In the alternative situation we do not make use of an estimator and we look for parameters A, such that for a fixed class M , the function JA E M satisfies

JA(z,) E y;, i = 1, .. . ,n.

= [r, y),

{ay : V E Y} {aIVI

(7)

+ a2V2 + ... + a .. v.. :VI

E YI,' "

,v.. E Y.. } (12)

The set of all such parameters Using (11) we obtain for the endpoints of (12) ; is called the feasible parameter set (Milanese, 1989; Milanese and Belforte, 1982), or likelihood set (Walter and Piet-Lahanier, 1988).

wherein V"II'''' = {1!, a < 0 ; y, a ~ O}, v-- = v+ = = v- = lL' Similarly, if 21 is a real valued matrix

y,

V-+

The set of all corresponding solutions 2111,

.1'" = {!~ EM: A E A}

(9)

2112,

(13)

21= (

will be called the feasible solution set or the estimate solution set, according to whether A = AY, or A = A~ .

Zhb Zh2,

then using (12) we may write {ZV ; VEY}

The determination of the sets (6) and (8) is often of considerable interest. However, definition (6) is in no way constructive, since an algorithm computing A = 4>(21, V) can be of limited use for the determination (or description) of the set A~. Here we have to look for a specific new algorithm which can produce the set A~ for interval input data. This is still more valid for the computation of the set AY defined by (8). Although we may have an algorithm for finding A such that JA (21,) = v" i = 1,· . . , n, this algorithm is usually of limited value for the determination of the set (8).

{(ZIV,Z2y, "· ,ZhV) : VEY}

=

(zIY,"', z.Y)

= zY

(14)

We thus see that (12) is a special case of (14). LINEAR ESTIMATION UNDER INTERVAL MEASUREMENTS Consider the case when the estimator 4>dz, y) is linear of the form 4>dz, y) HV, where H E Rh .. , and H H(z) may depend on 21 but not on V. Then formula (6) can be

=

=

writen as

In the next sections we consider the sets (6) and (8) for a linear estimation and for a linear identification problem respectively. We then study the determination for a fixed ~ E D of the set

A~

= {A E Rh , A = 4>c(z, V) ; V E Y} {Hy : yEY} = HY,

=

where by the last equality relation (14) has been used.

(10)

222

(15)

which delivetll bounds for the feasible solution set (9) under the 888umption that the latter is not empty.

Let us further assume that £ is a cl888 of linear on A functiODl! of the form J(Aj e) = L~;~ A;tMe) = t/I(e)A, where t/I(e) = (t/Io(e), ···, t/lk-l(e)), A = (Ao, ···, Ak_l)T . Formula (10) can be written in the form

L(e)

= = = =

U(A je) : A E

INTERPOLATION UNDER INTERVAL MEASUREMENTS

An

{t/I(ep, A=Hy : yEY}

Let us 888ume as before that the set M of approximating (interpolating) functions h is fixed. Denote the interval function (10) "comprising" all h E £ such that h(Zi) E Y; by L(z,Y je) and denote Y; = L(Z ,Y jZi). It is a trivial observation that Y; C Y; and that the interval vectors Y = (Y1 , Y" ... ,Yn) and Y = (Yi , l';, .. . ' Yn) generate same feasible solution sets. However, the intervals Y; have the property of poIIIIe8IIing no "excess points" , that is such points, through which no single-valued solutions h p888 . We give the following definition with respect to this property.

{t/I(e)H y : yE Y}

(t/I(e)H)Y

= r(e)Y .

(16)

Note that the interval- valued function (16) gives the exact bounds for the solution set (9). We now consider the case of a least~uare estimator. In the case of multiple linear regression, let us denote = (1,6 ,···, e.-d and 888ume t/li(e) = ei, i = O,··· , k -1 , so that Jpj e) t/I(e)A Ao + AI6 + ... + Ak-le.-l eA. Denoting

e

=

X=

=

=

=

Dd. For a fixed Z (ZI ' ·· · ' zn) and a cl888 M of approximating functions a system Y = (Y1 , ··· , Yn) of intervals is called M~ompatible , if for any point (z, y) from these intervals there is an element of M, whose graph P888es through (z, y) and interpolates all Y;.

1 ZI1

. . . Zlk-l) .. . , ( 1 Znl . .. Zn.-I

we obtain the matrix H in the form H = (XT X)-l XT . Substituting in (15) and (16) we obtain

Notationally, Y is M~ompatible if for any i = 1,· ·· , n and Yi E Y;, there exist a hEM, such that !l(Zi) = Yi and h(zj) E Yj for j :F i, j = 1, ·· ·, n . Two systems Y, Y, generating same feasible solution sets are called equivalent. Practically, the problem of finding the feasible solution set, can be divided into two steps: 1) to find an equivalent Mcompatible system Y of (measurement) intervals, and 2) to find the feasible solution set generated by Y.

A~ = HY=(XTX)-l)XTy, (17) T L(e) = r(e)Y = (eH)Y (e(X X)-l XT) Y, (18)

=

where r(e)

= e(X TX)-l X T = (-yl(e), · · ·, 'Yn(e).

In the case k = 2, the approximating function obtaiDl! the form J(A j e) = Ao + Ale· For the components 'Yi(e) of the n--dimensional vector r(e) we obtain 'Yi(e) =

(e(X TX)-l XT)i

:;;; (:e, - ¥)(e - :f')/S... + l/n , i = 1,···, n

=

=

=

In the situation when M = £ is a set of linear on A functions we can formulate an algorithm, using interval arithmetic, which can produce the bounds of the interval function at any e. For simplicity we shall cOllllider the polynomial case £ = {L~;~ Aiek- 1 - i }; the generalisation for an arbitrary linear family of the form L~;~ Ait/li(e) is evident .

,

where Z Li.l Zi/n, Su Li. l z?-nz Li.. I(Zi- Z )'. The boundary functions of the interval function L(e) r( e)Y are lines in each interval with end-points two neighbouring ei, i = 1, · · · , n, where ei are defined by 'Yi(ei) = 0, that is:

=

Let us fitllt consider the case k = n . In this case the computation of the interval function, producing the bounds of the feasible solution set is straightforward (Crane, 1975j Garloff, 1979j Herzberger, 1978j Krawczyk, 1969j Markov, 1986j Markev, 1991). We have

ei = z+S."./(n(z-zi)) . For more details and the geometric interpretation see Markev (1990).

k

The polynomial and multinomial cases produce similar results under the corresponding choice of the matrix X .

LW

=Lk-l(Yje) = Llk-l,;(e)y;,

(19)

1=1

where the polynomials Ik_1 •i of (k - l)st degree are the coefficients in the Lagrange interpolation formula. Using interval arithmetic (19) can be effectively computed at any E R. For a generalisation of the case k = n towards consideration of interval values for the argument Z the reader may consult Rockne (1972) .

Let us note that formulae (17) and (18) are constructive in the sense that they can be used for effective computation in any programming environment which supports interval arithmetic . Also we should like to stress on the fact that the interval function (18) produces the exact bounds for the estimate solution set.

e

Consider now the case k < n. If we need the value of the interval function L at we can compute at the interval functions for all subsets of k intervals out of the given n

e

In the next section we demonstrate the interval arithmetic approach for the computation of the interval function (10)

223

e

(by meaIlll of (19» and to take the intersection of these interval values. Such an algorithm can be used to produce, in particular, the values of L at the mesh points z,' that is the C-;:ompatible intervals Y;.

Crane, M.A. (1975). A bounding technique for polynomial functions. SIAM J. Appl. Math., Vol. 29, NoA. Garloff, J. (1979). Optimale Schranken bei Intervall interpolation mit Polynomen und mit Funlr.tionen az b • Z. Angew. Math. Mech., Vo!. 59, T59--T60.

It is an open problem to find effective algorithlDll for the computation of the boundary functions of L, which make use 1) of the C-;:ompatible intervals y;, and 2) of simple properties of the boundary functioDB . In this connection we should like to formulate the following hypothesis:

Herzberger, J. (1978). Note on a bounding technique for polynomial functions. SIAM J. Appl. Math., Vo!. 34, NoA.

Hypothesis. The boundary functions of L are elements of C in every interval between two neighbouring mesh points.

Krawczyk, R. (1969). Approximation durch Intervallfunktionen. Interner Bericht de8 In8t. f. Informatile, 69/7, Univ. Karlsruhe.

The above hypothesis is true for le = 2 and le = 3. For le = 2 we can write down the interval function (19) in the following simple form

Markov, S.M. (1986). On the usage of interval and computer arithmetic in numerical analysis. In BI. Sendov (Ed.), Optimal algorithm8, Publ. House of the Bulg. Acad. of SeL, Sofia, pp. 183-186.

if

z, ~ ~ ~

Markov, S.M. (1990). Least-square aproximatioDB under interval input data. In Ch. Ullrich (Ed.), Contribution8 to computer arithmetic and 6el/-validating numerical method8, J.C. Baltzer AG, Sci. Pub!. Co., Basel, Switzerland, pp. 133-147.

Zi+l!

i = 1,···,n - 1j

if

~ ~ Z1, ~ ;:::

z",

wherein the mesh points z, are assumed to be in increasing order: Z1 < Z2 < ... < z".

Markov, S.M. (1991). Polynomial interpolation of vertical segments in the plane. In E. Kaucher, S. Markov, G. Mayer (Eds.), Computer arithmetic, validated computation and mathematical modelling, to appear.

CONCLUSION AB we already mentioned, the interval arithmetic expressioDB for the computation of the interval-valued functioDB can be effectively computed in a software environment which suppotrs interval arithmetic. Such environments provide computer operatioDB with directed roundings, so that the computed interval bounds are automatically rounded tc:r wards outside and contain with guarantee the true results. Thus, the computed bounds comprise all possible kinds of input and computational errors. This fact opeDB a new way to the practical implementation and interpretation of the computed results. For example, assume that the measurement intervals contain with guarantee the true values of the measured quantities and that we have computed the interval function L{Cj z, Y) from these measurements under the assumption that the model function fA belongs to C. Assume that the experimental scientist provides us a new measurement (z, Y) such that L{Cj z, Yj ~ Y . The only conclusion then can be that the class C of model functioDB is unadequate for the description of the data.

Milanese, M. (1989). Estimation and prediction in the presence of unknown but bounded uncertainty: a survey. In M. Milanese, R. Tempo and A. Vicino (Eds.), Robwtne88 in identification and control, Plenum Press, New York, pp. 3-24. Milanese, M., G. Belforte (1982). Estimation theory and uncertainty intervals evaluation in presence of unknown but bounded errors. Linear families of models and estimators. IEEE Tran8. Autom. Control, Vo!. AC-27, pp. 4OH14. Moore, R.E. (1966). Interval anal1l8i8. Prentice Hall, Englewood Cliffs, N. J. Moore, R.E. (1968). Practical aspects of interval computation. Aplileace mathematilell, pp. 53-92.

e)

Rokne J. (1972). Explicit calculation of Lagrangian interval interpolating polynomial. Computing, Vo!. 9, pp. 149-157.

ACKNOWLEGEMENTS. The authors wish to thank Prof. M. Milanese for his attention, useful advises and support during the preparation of this note. This work was partially supported by the Ministry of Science and Higher Education under contract 755/89.

Schweppe, F.C. (1968) . Recursive state estimation: unknown but bounded errors and system inputs. IEEE Tran8. Autom. Control, Vol. AC-13, pp. 22-28. Schweppe, F.C. (1973) . Uncertain DlInamic SlI8tem8. Prentice Halla, Englewood Cliffs N.J.

REFERENCES

Walter, E., H.Piet-Lahanier. (1988). Estimation of parameter bounds from bounded-error data: a survey. Proc. l£th IMACS World Congreu, Paris, pp.467472.

Alefeld, G., J. Herzberger (1974). Einfuhrung in die Intervallrechung. Bibliographisches IDBtitut, Mannheim /Wien/ Zurich, pp. 2~30 .

224