IDENTIFICATION OF MULTI-DIMENSIONAL SYSTEM BASED ON

A NOVEL CRITERION

Yue Zhao, Kueiming Lo

School of Software, Tsinghua University, Key Lab for ISS, MOE China, Tsinghua University, Beijing 100084, P. R. China

Wook-Hyun Kwon

School of Electrical Engineering, Seoul National University, Seoul 151-742, Korea

Keywords:

Multi-dimensional system, identiﬁcation criterion, ARMAX model, recursive algorithm.

Abstract:

Most system recursive identiﬁcation algorithms are based on the prediction error (PE) criterion. Such a re-

cursive algorithm only considers the present estimation residual error instead of all estimation residuals. It

would result in large estimation error when the signal noise disturbs strongly. In this paper, a new iden-

tiﬁcation criterion is proposed. It considers both the errors between the actual outputs and the estimation

result and the difference of each estimation error. Under this criterion, a new recursive algorithm MSDCN

(Multi-dimensional System Disturbed by Color Noise) is proposed. For multi-dimensional systems, weight-

ing different values on the estimation errors and the difference of each error, MSDCN could both decrease the

estimation errors and got smooth prediction curves. Several simulation examples are given to illustrate the

method’s anti-disturbance performance.

1 INTRODUCTION

There already have many recursive algorithms for

modeling systems disturbed by white noises (Ander-

sson and Broman 1998, Ljung 1999, Grifﬁth Jr 1999,

Mershed 2000, etc.). The most characteristic fea-

ture of disturbance, however, is that its value is not

known beforehand. A physical system is affected by

many factors such as color noise disturbance and an

un-modeled structure. These types of algorithms en-

counter difﬁculties when the measurements are dis-

turbed (Kuo 2000, Trump 2001). For multi-discrete

systems disturbed by color noise, current identiﬁca-

tion algorithms cannot give precise estimate results

either (Schoukens 1991, Ljung 1985).

Efﬁciency modeling depends on the choice

of identiﬁcation criterion. Combining frequency-

domain method and time-domain method the GPE

criterion was proposed for interference systems (Lo

and Kwon 2002, 2003). Some extended recursive al-

gorithms (ERA) were put forward based on the GPE

criterion. (Lo and Kimura 2003, Lo et al. 2006, and

Lo and Huang 2006) Usually, a GPE criterion con-

tains a weighting matrix, in which there are many pa-

rameters for free choice. When the weighting matrix

is the identity, the ERA becomes the recursive least

squares algorithm (Lo and Kimura 2003). In addition,

the standard Yule-Walker method uses identity ma-

trix as weighting matrix. It may have poor accuracy,

and increasing the dimension of the weighting ma-

trix may well degrade the accuracy (Stoica, Friedlan-

der and Soderstrom 1987). Stoica and Jansson (2001)

proposed another method which derived the optimal

weight in a simple way and guaranteed the optimal

weighting matrix to be consistent and non-negative

deﬁnite, while still the choice of weighing matrix is

hard. Furthermore, for all the above raised implemen-

tations, a positive deﬁnite weighting matrix must be

weighted out in order to get a reliable estimate. The

optimal weight in general depends on unknown quan-

tities and hence must be itself estimated before its use

become possible.

Considering both the prediction errors and the dif-

ference of each error, performancecriterion in this pa-

per applies different weights. One part is the errors

between the actual outputs and the other is the esti-

mation result and the difference of each estimation er-

ror. Further, this paper develops a recursive weighting

matrix for fast calculation in recursive algorithms. In

this matrix, the values of each element depends on the

143

Zhao Y., Lo K. and Kwon W. (2008).

IDENTIFICATION OF MULTI-DIMENSIONAL SYSTEM BASED ON A NOVEL CRITERION.

In Proceedings of the Fifth International Conference on Informatics in Control, Automation and Robotics - SPSMC, pages 143-148

DOI: 10.5220/0001492001430148

Copyright

c

SciTePress

weights of the performance function. Based on this

new performance function, a new recursive algorithm

for multi-system identiﬁcation, MSDCN, is proposed.

MSDCN algorithm is new not only because it is

based on a new performance function, but also be-

cause it is based on a new concept of estimating the

noise part separately. Its two-step estimation feature

also make it a novel method. The MSDCN algorithm

has good anti-disturbance. For multi-dimensional

system, using the novel criterion to do the complex

estimation in each dimension of parameters, MSDCN

can give more precise results than current algorithms

do.

In the second part of this paper, the extended re-

cursivealgorithm for multi-dimensionalsystems is in-

troduced; The third part proposes a novel identiﬁca-

tion criterion; In the fourth part a recursive algorithm

for multi-dimensional system, MSDCN, is proposed;

An analysis of the performance of new criterion and

new algorithm is given in the ﬁfth section, and then

its simulation results are compared with other algo-

rithms.

2 EXTENDED RECURSIVE

ALGORITHM FOR

MULTI-DIMENSIONAL

SYSTEM

Consider system:

A(q)y(t) = B(q)u(t) + w(t) (1)

where, y(t) and w(t) are p-dimensional vectors, u(t)

is m-dimensional vector. y(t), u(t), and w(t) are sys-

tem output, input, and noise respectively. A(q), B(q)

are the backward operator q

−1

polynomial expression

A(q) = I +

n

a

∑

k=1

A

k

q

−k

, B(q) =

n

b

∑

k=1

B

k

q

−k

.

Denote:

θ = (A

1

, A

2

, ·· · , A

n

a

, B

1

, B

2

, ·· · , B

n

b

)

T

ϕ

t

= (−y

T

(t − 1), ··· , −y

T

(t − n

a

),

u

T

(t − 1), ·· · , u

T

(t − n

b

))

T

where θ ∈ R

(n

a

·p+n

b

·m)×p

is a matrix formed by the

system part parameters. ϕ

t

∈ R

n

a

·p+n

b

·m

is a regres-

sion vector. Then (1) is rewritten as

y(t) = θ

T

ϕ

t

+ w(t) (2)

Denote ε

t

be the prediction error of system output

y(t):

ε

t

= y(t) − θ

T

ϕ

t

, (3)

Then the identiﬁcation criterion can be expressed

as:

J(N) = tr([ε

1

, ε

2

, ··· , ε

N

]

T

Q(N)[ε

1

, ε

2

, ··· , ε

N

]) (4)

where trA represents the trace of matrix A. The

weighting matrix Q(N) ∈ R

N×N

is a symmetrical pos-

itive deﬁnite matrix and expressed as:

Q(N) =

Q(N − 1) α(N)

α(N)

T

q

N

, t = 1, 2, ··· , N

where α(N) ∈ R

N−1

. Denote:

Φ(N) = (ϕ

1

, ϕ

2

, ··· , ϕ

N

)

T

Y(N) = (y

1

, y

2

, ··· , y

N

)

T

Then the identiﬁcation criterion (4) can be expressed

as:

J(N) = tr([Y(N) − Φ(N)θ]

T

Q(N)[Y(N) − Φ(N)θ]

(5)

Since

J(N) = tr(Y

T

(N)Q(N)Y(N)) −tr(θ

T

Φ

T

(N)Q(N)Y(N))

−tr(Y

T

(N)Q(N)Φ(N)θ) + tr(θ

T

Φ

T

(N)Q(N)Φ(N)θ)

The gradient of J(N) is:

∂J(N)

∂θ

=

∂

∂θ

[tr(Y

T

(N)Q(N)Y(N))

−tr(θ

T

Φ

T

(N)Q(N)Y(N))

−tr(Y

T

(N)Q(N)Φ(N)θ)

+tr(θ

T

Φ

T

(N)Q(N)Φ(N)θ(N))]

=

∂

∂θ

[tr(Y

T

(N)Q(N)Y(N))]−

∂

∂θ

[tr(θ

T

Φ

T

(N)Q(N)Y(N))]

−

∂

∂θ

[tr(Y

T

(N)Q(N)Φ(N)θ)]

+

∂

∂θ

[tr(θ

T

Φ

T

(N)Q(N)Φ(N)θ(N))]

= −(Y

T

(N)Q(N)Φ(N))

T

− Φ

T

(N)Q(N)Y(N)+

(Φ

T

(N)Q(N)Φ(N) + Φ

T

(N)Q(N)

T

Φ(N))θ(N)

and

∂J

2

(N)

∂θ

2

= Φ

T

(N)Q(N)

T

Φ(N)

where Q(N) = Q

T

(N) and

∂J

2

(N)

∂θ

2

is positive and def-

inite. Let

∂J(N)

∂θ

= 0

which minimizes J(N) and yields:

θ(N) = [Φ(N)

T

Q(N)Φ(N)]

−1

Φ(N)

T

Q(N)Y(N)

(6)

ICINCO 2008 - International Conference on Informatics in Control, Automation and Robotics

144

t = 1, 2, · ·· Then at time t, denote:

P

t

= Φ

T

t

Q

t

Φ

t

a

t

= 1+ ϕ

T

t

P

−1

t−1

Φ

T

t−1

α

t

σ

t

= q

t

− α

T

t

Φ

t−1

P

−1

t−1

Φ

T

t−1

α

t

b

t

= a

t

+ a

−1

t

σ

t

ϕ

T

t

P

−1

t−1

ϕ

t

Then we can get:

Theorem 1

P

−1

t

= P

−1

t−1

− b

t

P

−1

t−1

(ϕ

t

β

T

t

Φ

t−1

+ Φ

T

t−1

β

t

ϕ

T

t

)P

−1

t−1

+

1

a

t

b

t

P

−1

t−1

(ϕ

T

t

P

−1

t−1

ϕ

t

Φ

T

t−1

β

t

β

T

t

Φ

t−1

− σ

t

ϕ

t

ϕ

T

t

)P

−1

t−1

θ

t

= θ

t−1

+

1

a

t

b

t

P

−1

t

(β

t

Φ

T

t

β

t

+ σ

t

ϕ

t

)(y

t

− θ

T

t−1

ϕ

t

)

+

1

a

t

b

t

P

−1

t−1

(β

t

ϕ

t

− ϕ

T

t

P

−1

t−1

ϕ

T

t−1

β

t

)β

T

t

(Y

t−1

− Φ

t−1

θ

t−1

)

As Extended Recursive Algorithm for multi-

dimensional system.

3 IDENTIFICATION CRITERION

New identiﬁcation criterion is:

J(N) = λ

N

∑

t=1

k ε

t

k

2

+µ

N−1

∑

t=1

k ε

t+1

− ε

t

k

2

(7)

where ε

t

represents the estimation errors at time t, in

multi-dimensional system, ε

t

is p-dimensional vector,

which dimension is the same as the output y(t).

λ represents the weight on the estimation errors

µ represents the weight on the difference between

each estimation error.

Transforming equation (7), we can get

J(N) =

N

∑

t=1

(λε

T

t

ε

t

) +

N−1

∑

t=1

µ(ε

t+1

− ε

t

)

T

(ε

t+1

− ε

t

))

=

N

∑

t=1

tr(λε

t

ε

T

t

) +

N−1

∑

t=1

µtr(ε

t+1

− ε

t

)(ε

t+1

− ε

t

)

T

= tr

(λ+ 2µ)

N

∑

t=1

ε

t

ε

T

t

− µ

N−1

∑

t=1

ε

t+1

ε

T

t

− µ

N−1

∑

t=1

ε

t

ε

T

t+1

)

Assume E(N) = (ε

1

, ε

2

, ··· , ε

N

)

T

∈ N × p

and

Q(N) =

q

11

q

12

··· q

1N

q

21

q

22

··· q

2N

.

.

.

.

.

.

.

.

.

.

.

.

q

t1

q

t2

··· q

NN

∈ N × N,

in which,

q

i,i

= λ+ 2µ

q

i+1,i

= −µ i = 1, 2, ··· , N.

q

i,i+1

= −µ

Then equivalent with equation (3), we can get

J(N) = trE(N)

T

Q(N)E(N) (8)

as the new expression for the performance criterion.

Thus at time t, the weighting matrix Q

t

can be ex-

pressed in this recursive form:

Q

t

=

Q

t−1

β

t

β

T

t

p

t

, t = 1, 2, · ·· , N

in which, β

t

= [0, 0, · ·· , −µ] ∈ R

t−1

and p

t

= λ+ 2µ.

4 RECURSIVE ALGORITHM

FOR MULTI-DIMENSIONAL

SYSTEM

Consider the ARMAX model

A(q)y(t) = B(q)u(t) +C(q)w(t) (9)

where, y(t) and w(t) are p-dimensional vectors, u(t)

is m-dimensional vector. y(t), u(t), and w(t) are sys-

tem output, input, and color noise respectively. A(q),

B(q), C(q)is the backward operator q

−1

polynomial

expression

A(q) = I +

n

a

∑

k=1

A

k

q

−k

, B(q) =

n

b

∑

k=1

B

k

q

−k

,

C(q) = I +

n

c

∑

k=1

C

k

q

−k

.

where A

k

, C

k

∈ p × p, and B

k

∈ p× m. Denote:

θ = (A

1

, A

2

, ··· , A

n

a

, B

1

, B

2

, ··· , B

n

b

)

T

,

ϕ

t

= (−y

T

(t − 1), · ·· , −y

T

(t − n

a

), u

T

(t − 1), · ·· , u

T

(t − n

b

))

T

,

ρ = (C

1

, C

2

, ··· , C

n

c

)

T

,

φ

t

= (w(t − 1), w(t − 2), ··· , w(t − n

c

))

T

.

in which θ ∈ R

[n

a

·p+n

b

·m]×p

is the matrix formed by

the system part parameters; and ρ ∈ R

[n

c

·p]×p

is the

matrix formed by the system’s noise part parameters.

Here we introduce a two-step algorithm to do the

estimation for system (9).

Step 1: Transform (9) into:

y(t) − ρ

T

t−1

φ

t

= θ

T

t

ϕ

t

+ w(t) (10)

IDENTIFICATION OF MULTI-DIMENSIONAL SYSTEM BASED ON A NOVEL CRITERION

145

Estimate the system parameters θfor system (10). For

the identiﬁcation of the parameters of system θ, use

the following measurements:

J(t) = tr(λ

N

∑

t=1

k ε

t

k

2

+µ

N−1

∑

t=1

k ε

t+1

− ε

t

k

2

). (11)

Step 2: Transform (10) into:

y(t) − θ

T

t

ϕ

t

= ρ

T

t

φ

t

+ w(t). (12)

Estimate the parameter ρ of ﬁlter C(q) of system (9).

Denote

Φ

t

= (ϕ

1

, ϕ

2

, ··· , ϕ

t

)

T

, Y

t

= (y

1

, y

2

, ··· , y

t

)

T

,

g

t

= y

t

− ρ

T

t−1

φ

t

, G

t

= (g

1

, g

2

, ··· , g

t

)

T

And according to equation (8), the novel criterion J(t)

can be expressed as:

J(t) = trE(t)

T

Q

t

E(t)

= tr[ε

1

, ε

2

, ··· , ε

t

]

T

Q

t

[ε

1

, ε

2

, ··· , ε

t

]

= tr[G

t

− Φ

t

θ]

T

Q

t

[G

t

− Φ

t

θ]

If Φ

T

t

Q

t

Φ

t

is not singularthen minimize J

t

to get the

system’s parameter vector θ

t

’s optimal solution:

θ

t

= [Φ

T

t

Q

t

Φ

t

]

−1

Φ

T

t

Q

t

G

t

, t = 1, 2, ··· , N. (13)

Similar with Extended Recursive Algorithm, we

can get MSDCN algorithm for ARMAX model (9):

θ

M

(t) = θ

M

(t − 1)

+

1

a

t

b

t

P

−1

t

(β

t

Φ

T

t

β

t

+ σ

t

ϕ

t

)(y

t

− ρ

M

(t − 1)

T

φ

t

− θ

M

(t − 1)

T

ϕ

t

)

+

1

a

t

b

t

P

−1

t−1

(β

t

ϕ

t

− ϕ

T

t

P

−1

t−1

ϕ

T

t−1

β

t

)β

T

t

(G

t−1

− Φ

t−1

θ

M

(t − 1))

P

−1

t

= P

−1

t−1

− b

t

P

−1

t−1

(ϕ

t

β

T

t

Φ

t−1

+ Φ

T

t−1

β

t

ϕ

T

t

)P

−1

t−1

+

1

a

t

b

t

P

−1

t−1

(ϕ

T

t

P

−1

t−1

ϕ

t

Φ

T

t−1

β

t

β

T

t

Φ

t−1

− σ

t

ϕ

t

ϕ

T

t

)P

−1

t−1

ρ

M

(t) = ρ

M

(t − 1) +

R

−1

t−1

φ

t

1+ φ

T

t

R

t−1

φ

t

(y

t

− θ

M

(t)

T

ϕ

t

− ρ

M

(t − 1)

T

φ

t

)

R

t

= R

t−1

+ φ

t

φ

T

t

bw(t) = y(t) − ρ

M

(t)

T

φ

t

− θ

M

(t)

T

ϕ

t

(14)

in which,

(

p

t

= λ+ 2µ

β

t

= [0, 0, · ·· , −µ] ∈ R

t−1

Initially, θ and ρ can be zero matrix; R

0

, P

0

to-

gether constitute the identity matrix.

5 SIMULATIONS

All the simulations were conducted in the same

computation environment. The main criterions of the

computer were: CPU 1.66GHz, RAM 1G bytes and

with Windows XP OS.

Experiment. A Black-box model is as follows:

y(t) =

B

1

q

−1

I + A

1

q

−1

u(t) +

I +C

1

q

−1

I + A

1

q

−1

w(t)

The real parameters of the system were:

A

1

=

0.325 −1

0.5 −1.1

,

B

1

=

1.9 3.7

−0.4 −0.8

,

C

1

=

−0.7 −0.1

0.8 0.9

.

The color noise sequence {w(t)}

N

1

was a two di-

mensional vector, of which each dimension was com-

posed of different sawtooth wave and random num-

ber generator with variance 2 and the input signal

u(t) = [1 2], each dimension of which was gener-

ated by a square generator. The experiment was con-

ducted using the LS method and the algorithm MS-

DCN method with the sample number N = 1000. The

results are shown in Figures 1,2,3,4,5,6 and the statis-

tics are as follows:

(1)Least-Squares method

The system statistics results are:

θ =

0.2889 0.0361

−0.9823 −0.0177

0.3855 0.1145

−1.2366 0.1366

1.7990 0.1010

3.5981 0.1019

−0.1692 −0.2308

−0.3384 −0.4616

−0.7499 0.0499

0.0907 −0.1907

0.9270 −0.1270

0.8506 0.0494

.

The simulation results are shown in Figure 1, 2

and 3.

(2)MSDCN method

if we choose λ = 0.8andµ = 0.2,

then according to (7), β

t

= [0, ··· , −0.2, −0.2] and

p

t

= 1.2.

Then the system statistics results are:

ICINCO 2008 - International Conference on Informatics in Control, Automation and Robotics

146

0 100 200 300 400 500 600 700 800 900 1000

−2

−1.5

−1

−0.5

0

0.5

1

1.5

2

0.325

−1.1

Figure 1: Estimation results of parameter A under ELS.

0 100 200 300 400 500 600 700 800 900 1000

−4

−3

−2

−1

0

1

2

3

4

5

6

3.7

1.9

−0.4

−0.8

Figure 2: Estimation results of parameter B under ELS.

0 100 200 300 400 500 600 700 800 900 1000

−2

−1.5

−1

−0.5

0

0.5

1

1.5

2

−0.1

−0.7

0.8

0.9

Figure 3: Estimation results of parameter C under ELS.

θ =

0.2620 0.0630

−0.9602 −0.0398

0.4448 0.0552

−1.2425 0.1425

1.7595 0.1405

3.5189 0.1811

−0.4267 0.0267

−0.8534 0.0534

−0.6499 −0.0501

0.0409 −0.1409

0.8232 −0.0232

0.9362 −0.0362

The simulation results are shown in Figure 4, 5

and 6.

0 100 200 300 400 500 600 700 800 900 1000

−2

−1.5

−1

−0.5

0

0.5

1

1.5

2

0.325

−1.1

Figure 4: Estimation results of parameter A under MSDCN.

0 100 200 300 400 500 600 700 800 900 1000

−4

−3

−2

−1

0

1

2

3

4

5

6

3.7

1.9

−0.4

−0.8

Figure 5: Estimation results of parameter B under MSDCN.

0 100 200 300 400 500 600 700 800 900 1000

−2

−1.5

−1

−0.5

0

0.5

1

1.5

2

0.9

0.8

−0.1

−0.7

Figure 6: Estimation results of parameter C under MSDCN.

Through a comparison of the estimation results

of ELS and MSDCN, we can see that the new

performance function has proved to be efﬁcient for

the multi-dimensional systems. For the example

given above, we weight 0.8 on the estimation errors

and 0.2 on the difference between each estimation

error, which means we want the estimation curves

close to the actual values more than make the curves

smooth. The weighting matrix is deﬁnite after β

t

and

p

t

are ﬁxed by λ and µ. Weighting different on λ and

µ will result in different estimation results. Those

ﬁgures show that MSDCN has both decreased the

estimation errors and got smooth prediction curves.

IDENTIFICATION OF MULTI-DIMENSIONAL SYSTEM BASED ON A NOVEL CRITERION

147

6 CONCLUSIONS

This paper proposes a new identiﬁcation criterion for

multi-dimensional system disturbed by color noise,

and further develops a recursive algorithm, MSDCN,

based on it. Weighting both on the prediction er-

rors and the difference of each prediction error, the

identiﬁcation criterion makes the weighting matrix

deﬁnite in calculation. Based on this performance

criterion, MSDCN is developed using a two-step

method to estimate both the system parameters and

the noise part. The MSDCN algorithm has high anti-

disturbance performance in the prediction of multi-

dimensional systems disturbed by color noise. It both

decreased the estimation errors and got smooth pre-

diction curves. The performance of the MSDCN al-

gorithm was demonstrated by simulations.

ACKNOWLEDGEMENTS

This work is supported by the Funds NSFC60672110,

NSFC60474026,and the JSPS Foundation.

REFERENCES

Andersson, A. and Broman, H. (1998), A second-order re-

cursive algorithm with applications to adaptive ﬁlter-

ing and subspace tracking. IEEE Trans. Signal Pro-

cessing, vol. 46, pp. 1720C1725, June.

Grifﬁth Jr, D. W. and Arce, G. R. (1999), A partially decou-

pled RLS algorithm for volterra ﬁlters. IEEE Trans.

Signal Processing, vol. 47, pp. 579C582, Feb.

Kuo, C. J. , Jr, R. D. , Lin, C. Y. and Tsai, Y. C. (2000), Set-

theoretic estimation based on a priori knowledge of

the noise distribution. IEEE Trans. Signal Processing,

vol. 48, pp. 2150C2156, July.

Lo, K. M. , et al. (2006), Empirical frequency-domain

optimal parameter estimate for Black-box processes.

IEEE Transactions on Circuits and Systems-I: Regu-

lar Papers, vol.53, No.2, 419-430.

Lo, K. M. , Kwon, W. H. (2003), New identiﬁcation ap-

proaches for disturbed models. Automatica.vol.39,

No.9, 1627-1634.

Lo, K. M. , Kimura, H. (2003), Recursive Estimation Meth-

ods for Discrete Systems. IEEE Trans. On Circruits

and Systems.vol.49, NO.6, No.6, 439-446.

Ljung, L. (1985), On the estimation of transfer function.

Automaticavol. 21, 677-696.

Ljung, L. (1999), System Identiﬁcation: Theory for the

User. Upper Saddle River,NJ: Prentice-Hall, 1999.

Mershed, R. and Sayed, A. J. (2000), Order-recursive RLS

Laguerre adaptive ﬁltering. IEEE Trans. Signal Pro-

cessing, vol. 48, pp. 3000C3010, July.

Stoica, P., Friedlander, B.and Soderstrom, T.(1987) Op-

timal Intrumental Variable Multistep Algorithms for

Estimation of the AR Parameters of an ARMA Pro-

cess, Int. J. Control 45(1987), 2083-2107.

Stoica, P. and Jansson, M.(2001) Estimating Optimal

Weights for Instrumental Variable Methods. Digital

Signal Processing - A Review Journal vol. 11, no. 3,

pp. 252-268, Jul. 2001.

Trump, T. (2001), Maximum likelihood trend estimation

in exponential noise. IEEE Trans. Signal Processing,

vol. 49, pp. 2087C2095, Sept.

ICINCO 2008 - International Conference on Informatics in Control, Automation and Robotics

148