Nonlinearity Recovering in Wiener System Driven with Correlated Signal

printed in IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. 49, NO. 10, PP. 1805-1810, OCTOBER 2004

Nonlinearity Recovering in Wiener System Driven with Correlated Signal W odzimierz Greblicki Abstract— The characteristic of the nonlinear part of the Wiener system is estimated. The system is driven by a Gaussian random process which may not be white. Three algorithms are presented, two semirecursive and one of the off-line type. Their pointwise convergence in probability is shown and results of numerical simulation are given. Index Terms— System identi cation, nonparametric identi cation, nonparametric regression, Wiener system.

I. I NTRODUCTION In the block–oriented approach, we identify systems consisting of linear dynamic and nonlinear memoryless subsystems. In particular, in the Wiener system, a linear dynamic part is followed by a nonlinear memoryless one. The objective is to recover descriptions of the subsystems from input-output observations of the whole system. In this note we focus our attention on the nonlinear subsystem. In the literature, the a priori information about both subsystems is usually parametric, see, e.g., [1], [3], [7], [8], [10]. As a result, the nonlinear characteristic is a function known up to some coef cients, [10], or just a polynomial, [3], [7], [8]. The random input signal is Gaussian, [3], [8], but not necessarily, [10], white, [3], or correlated, [8], [10]. Due to serious theoretical dif culty, properties of some algorithms are only suggested by simulation examples, [1], [7]. The characteristic can be recovered also when the a priori information is smaller than parametric, i.e., nonparametric. Nonparametric algorithms work properly even if we know nothing about the nonlinear characteristic which can be any Borel measurable function, [5]. Needless to say that parametric algorithms are then completely useless. Thus, for a user, nonparametric methods have a great advantage, since they can be applied when the knowledge about the system is so poor that its description can't be parameterized. Off–line, [4], [5], and semirecursive, [6], nonparametric kernel algorithms have been examined. In all those papers, the input signal is a Gaussian white random process. The novelty of this note is the fact that we admit also a correlated input signal. To estimate the nonlinear characteristic, we propose three kernel algorithms, of which two are semirecursive and one is off-line. We show that all are pointwise consistent and present results of numerical simulation. W. Greblicki is with the Institute of Engineering Cybernetics, Wroc aw University of Technology, 50-370 Wroc aw, Poland. Email: [email protected]

1

II. S YSTEM AND A LGORITHMS The Wiener system with input Un and output Yn , shown in Fig. 1 with the thick line, consists of two subsystems. First is linear, dynamic while the other nonlinear, memoryless. The dynamic subsystem has an impulse response f n ; n = 0; 1; : : :g which means that n X Vn = (1) n i Ui ; i= 1

where fUn ; n = : : : ; 1; 0; 1; : : :g. The impulse response f n g is unknown but we assume that j

n 1 "1

nj

(2)

with some unknown 1 and "1 , where 0 "1 < 1. Observe that the restriction is satis ed by every stable ARMA system. We assume that the input signal is correlated and is of the following form: Un =

n X

(3)

n i i;

i= 1

where f n ; n = : : : ; 1; 0; 1; : : :g is a stationary white random Gaussian process with zero mean and unknown variance 2 . The signal f n g is not measured. The weighting function f n g is unknown but satis es the following inequality: j

nj

n 2 "2

(4)

with some unknown 2 and "2 , where 0 "2 < 1. Thus, fUn g is output of some system with unknown impulse response f n g driven by f n g, see the thin line in Fig. II. Output of the subsystem is disturbed by additive stationary white Gaussian noise fZn ; n = : : : ; 1; 0; 1; : : :g independent of both fUn g and fVn g and having zero mean and unknown nite variance. Owing to all that, both fVn ; n = : : : ; 1; 0; 1; : : :g as well fWn ; n = : : : ; 1; 0; 1; : : :g are also stationary and Gaussian. The unknown characteristic of the nonlinear part, denoted by , has a continuous derivative such that, for all w 2 R, 0 < j 0 (w)j

(5)

c

with some unknown c . Observe that, owing to that, 1 is strictly monotonous and exists in R (R), where (R) is the image of R under the mapping .

Zn κn ξn

of

λn Un

µ Vn Wn

Yn

Our goal is to estimate the characteristic the nonlinear subsystem from observations

2

printed in IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. 49, NO. 10, PP. 1805-1810, OCTOBER 2004

(U0 ; Y0 ) ; (U1 ; Y1 ) ; : : : taken at input and output of the Wiener system. The basis for our algorithms is the following lemma owing to which, to recover the inverse of (y), we estimate the regression EfU0 jYk = yg. In the lemma and further parts of the paper, 2 k

=

1 X

P1 with i = j=0 j i+j . Lemma 1: Let fUn g be a stationary Gaussian random process with zero mean. For k = 0; 1; : : :, EfU0 jYk = yg = k 1 (y): (6) Proof: Using (1) and (3), one can verify that cov [U0 ; Wk ] = k 2W . Since the pair (U0 ; Wk ) has a normal distribution with zero marginal means, an application of Lemma 2 in Appendix A, leads to EfU0 jWk = wg = k w. Since EfU0 j (Wk ) = 1 yg = EfU0 jWk = (y)g, the proof has been completed. To recover k 1 (y), we apply the following offline algorithm: n X

i=1 n X

m ^ k (y) =

y Yi+k hn

Ui K

(7) K

y Yi+k hn

i=1

and two semirecursive ones, namely:

m ~ k (y) =

n X

Ui h1i K

i=1 n X

y Yi+k hi

(8) 1 hi K

y Yi+k hi

Ui K

y Yi+k hi

Ui K

y Yi+k hi

i=1

and

mk (y) =

n X

:

(9)

i=1

In all above algorithms, K and fhn g are a suitably selected kernel function and a positive number sequence, respectively. On the Borel measurable kernel K, we impose the following restrictions: sup jK(y)j < 1;

(10)

Z

(11)

y2R

K(y)dy < 1;

yK(y) ! 0 as jyj ! 1; jK(x)

K(y)j

cK jx

yj

(12) (13)

with some cK , for all x; y 2 R. We denote dK = supy jK(y)j. Depending on algorithm, the positive

(15)

nh2n ! 0 as n ! 1;

hn

n X i=1

n X i=1

(16)

1 ! 0 as n ! 1; h2i

(17)

hi ! 1 as n ! 1:

(18)

Algorithms (7)–(9) just estimate the regression in (6). Estimates (8) and (9) are semirecursive since their numerators and denominators can be calculated recursively. For the sake of the simplicity of the notation, W and Y stand for Wk and Yk , respectively. III. C ONVERGENCE OF ALGORITHMS Theorems 1 and 2 establish convergence of algorithms (8) and (9). Having their proofs, the reader can easily verify convergence of (7), i.e., prove Theorem 3. Theorem 2: Let satisfy (5). Let K satisfy (10)– (13). Let the number sequence satisfy (15) and (17). Then, at every y 2 (R),

m ~ k (y) ! k 1 (y) as n ! 1 in probability. Proof: First of all, observe that from Lemma 4 in Appendix B, it follows that the density f of Y is continuous and 0 < f (y) < 1 at every y 2 (R). Let y be any point in (R). We denote n 1 y Yk+i 1X Ui K ; g~(y) = n i=1 hi hi n

i=1

n X

hn ! 0 as n ! 1;

1 n2

i k i

2 W i=0

number sequence fhn g satis es some of the following conditions: fhn g is monotonous (14)

1X 1 f~(y) = K n i=1 hi

y

Yk+i hn

and observe m ~ k (y) =P g~(y)=f~(y). Applying Lemma1, n we get E~ g (y) = n 1 i=1 m(y; hi ) with m(y; h) = =

1 y Yk E U0 K h h 1 y Y 1 (Y )K k E h h

(19) :

Using (15) and Lemma 3 in Appendix A, we nd Z m(y; h) ! k 1 (y)f (y) K(x)dx as h ! 0;

(20) since Ej 1 (Y )j = EjW j < 1. Thus, E~ g (y) converges to the same limit as n ! 1. Passing to variance, we nd var[~ g (y)] = An (y) + Bn (y) with An (y) =

n 1 X 1 var U0 K n2 i=1 h2i

y

Yk hi

;

printed in IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. 49, NO. 10, PP. 1805-1810, OCTOBER 2004 n n 1 X X n2 i=1

Bn (y) =

j=1;j6=i

cov Ui K and d2K

y

Yk+i hi

y

; Uj K

observe Pn

2 2 Un

that 0 2 h = O(n i=1 i

2

y

Yk+i hi

Yk+j hj

An (y) 2 h i=1 i ) and

y

; Uj K

Yk+j hj

:

By virtue of Lemma 7 in Appendix B, the double sum in the above expression isPbounded in Pn i 1 absolute value by " k (y) i=1 hi 2 j=1 "i j n X P1 n " k (y) hi 2 with = n=0 " . Hence i=1 Pn 2 2 Bn (y) = g (y)] = PnO(n 2 i=1 hi ). Therefore, var]~ 2 O(n h ). In this way, we have nally shown i=1 i R that E~ g (y) ! k 1 (y)f (y) K(x)dx as n ! 1 in probability. Using similar R arguments, we can easily verify that f~(y) ! f (y) K(x)dx as n ! 1 in probability. Since, as we have already noticed, f (y) 6= 0, the proof has been completed. Theorem 3: Let satisfy (5). Let K satisfy (10)– (13). Let the number sequence satisfy (14), (15) and (18). Then, at every y 2 (R), mk (y) ! k 1 (y) as n ! 1 in probability. Proof: Let y 2 (R). Clearly, mk (y) = g(y)=f (y) with 1 g(y) = Pn

i=1

and

n X

hi

1 f (y) = Pn

Yk+i hi

i=1

hi

i=1

y

Ui K

n X

K

y

i=1

Yk+i hn

:

Applying LemmaP 1, we obtain Eg(y) = Pn n h m(y; h )= h with m(y; h) as in i i=1 i i=1 i P1 (19). Since i=1 hi = 1, see R (18), using (20), we nd Eg(y) ! k 1 (y)f (y) K(x)dx as n ! 1. To examine variance, observe that var[g(y)] = An (y) + Bn (y) with n X 1 An (y) = Pn var U0 K 2 ( i=1 hi ) i=1

Bn (y) =

(

i=1

Yk+i hi

2

hi )

n n X X

Yk hi

;

i=1 j=1;j6=i

Yk+j : hj Pn 2 2 2 Now, 0 An (y) dP K U n=( i=1 hi ) P n n 2 2 dK U (nhn = i=1 hi )(1=hn i=1 hi ). Since fhn g is cov Ui K

y

Pn

1

y

; Uj K

y

1 and An (y) = n

i 1

XX 2 Bn (y) = Pn 2 ( i=1 hi ) i=1 j=1 y Yk+i y Yk+j cov Ui K ; Uj K hi hj

Pn

n i 1 2 XX 1 Bn (y) = 2 n i=1 j=1 hi hj

cov Ui K

Pn monotonous, Pn nhn = i=1 hi O(1=hn i=1 hi ). In turn,

1 hi hj

3

:

By virtue of Lemma 7 in Appendix B, the double sum in the expression bounded in abPn is P i 1 solute value by (y)" k i=1 j=1 (hj =hi )"i j . Since fhn g is monotonous, P the obtained double n Pi 1 sum is not greater than hn 1 i=1 j=1 hj "i j = Pn Pn Pn i j hn 1 P = hn 1 j=1P hj with j=1 hj i=j+1 " 1 n n = " . Hence, B (y) = O(1=h n n n=0 i=1 hi ) Pn and, thus, var[g(y)] = O(1=hn i=1 hRi ). We have thus shown that g(y) ! k 1 (y)f (y) K(x)dx as n ! 1 in probability. Using similar R arguments, we can easily verify that f (y) ! f (y) K(x)dx as n ! 1 in probability. Since f (y) 6= 0, the proof has been completed. Theorem 4: Let satisfy (5). Let K satisfy (10)– (13). Let the number sequence satisfy (15) and (16). Then, at every y 2 (R), m ^ k (y) !

k

1

(y) as n ! 1 in probability.

IV. S IMULATION E XAMPLE In the simulation example, the dynamic part of the system has a transfer function z= (z 0:5). The nonlinearity is of the following form: (w) = sign(w)[(jwj+1)2 1]. The input signal is produced by a dynamic system with a transfer function z=(z 0:2) driven by a white Gaussian process with zero mean and variance 2 = 1. Disturbance has variance 2Z = 0:1. In the example, k = 0 and 0 0:64. The kernel K(y) = exp( y 2 ) has been applied. The global R 20 E( m ~ error is measured with MISE = 0 (y) 20 1 2 (y)) dy. 0 For n = 500, observations and the resulting off–line estimate obtained for h500 minimizing the MISE is shown in Fig. 1. In Fig. 2, for each n; n 2 [10; 1280], the hn minimizing the error has been selected. For n = 1280, MISE 0:6 which means that the average squared pointwise error measured on the interval [ 20; 20] is not greater than 0:015. The error gets small rapidly as the number of observations increases, especially for n counted in dozens, and becomes quite small for n reaching a few hundreds. V. F INAL R EMARKS Within the nonparametric approach, new in this note is that the input signal is correlated. For input being a white process, an algorithm similar to (7) has been already examined in [4] and [5] while (8) and (9) in [6]. Applying arguments applied in [5], one can extend our results to any Borel function .

4

printed in IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. 49, NO. 10, PP. 1805-1810, OCTOBER 2004

3

Lemma 6: Let X be a random variable with a density function f and let ' be a Borel function such that Ej'(X)j < 1. If a Borel measurable kernel K satis es (10)–(12), then Z x X 1 E '(X)K ! '(x)f (x) K(y)dy h h

n =500

"0 µ !1(y) estimate

2 1

as h ! 0 and

u 0

1 sup E '(X)K h h>0

-1

x

X

< 1 dK h0 (y; h0 ); > :

with some

2 dK

1,

'(y; h) =

+

2,

1 E h

1 3 hn cK c 3,

(

h0 "n k '(y; h0 ); for 0 k < n;

and 1

(Y )) K

y

Y h

n k

2 dK +

'(y);

1 3 hn cK c

for 0 < n

k;

h0 "n k '(y); for 0 k < n:

R EFERENCES

Using (4) again, we nally obtain jSn j

1 d K h0 "

Since supn jhn j < 1, the proof has been completed.

Pk which is bounded by c' "n j= 1 "k j Ej j j c' 2 "n . Hence, using Lemma 5, we get jBni j c' 2 "n k Ef (W ) (W )g. This and (26) give jAni c' "n

8 >