Signal Processing 80 (2000) 311}331
Fixed-order H and H optimal deconvolution "lter designs 2 = Bor-Sen Chen*, Jui-Chung Hung Department of Electrical Engineering, National Tsing-Hua University, Hsin-Chu 300, Taiwan, ROC Received 11 February 1999; received in revised form 7 September 1999
Abstract For the simplicity of implementation and saving of operation time, the "xed-order optimal deconvolution "lter design is appealing for engineers in signal processing from practical design perspective. In this study, a design method based on genetic algorithms is proposed to simultaneously treat with H and H optimal signal reconstruction design problem 2 = with prescribed "lter order. Genetic algorithms are optimization and machine learning algorithms, initially inspired from the processes of natural selection and evolutionary genetics. They tend to "nd the global optimum solution without becoming trapped at local minima. Both IIR and FIR cases are discussed in this study. The convergence property of our design algorithm is also discussed by Markov chain method. Two design examples of H and H optimal deconvolution 2 = with "xed-order "lters are given to illustrate the design procedure and the performance of the proposed design methods, respectively. ( 2000 Elsevier Science B.V. All rights reserved. Zusammenfassung Vom Standpunkt des praktischen Entwurfs aus betrachtet ist ein optimales Entfaltungs"lter fester Ordnung wegen seiner einfachen Implementierung und E$zienz attraktiv fuK r Anwendungen in der Signalverarbeitung. In dieser Studie wird eine auf genetischen Algorithmen beruhende Entwurfsmethode vorgeschlagen. Damit lassen sich Entwurfsprobleme fuK r H - und H - optimale Signalrekonstruktionen mit vorgeschriebener Filterordnung gleichzeitig behandeln. Gen2 = etische Algorithmen sind Algorithmen zur Optimierung und zum maschinellen Lernen, welche anfangs von Prozessen der natuK rlichen Selektion und evoleuionaK ren Genetik inspiriert wurden. Sie tendieren dazu, global optimale LoK sungen zu "nden, ohne durch lokale Minima gefangen zu werden. Sowohl der IIR- als auch der FIR-Fall werden in dieser Studie diskutiert. Weiters wird das Konvergenzverhalten unseres Entwurfsalgorithmus mittels einer auf Marko!-Ketten beruhenden Methode untersucht. Die Entwurfsprozedur und die LeistungsfaK higkeit der vorgeschlagenen Entwurfsmethoden werden durch zwei Entwurfsbeispiele zur H - und H - optimalen Entfaltung mit Filtern fester Ordnung 2 = illustriert. ( 2000 Elsevier Science B.V. All rights reserved. Re2 sume2 Pour la simpliciteH d'impleH mentation et l'eH conomie de temps d'opeH rations, la conception de "ltres de deH convolution optimaux d'ordre "xe est attrayante pour les ingeH nieurs en traitement de signaux dans une perspective de conceptions pratiques. Dans cette eH tude, nous preH sentons une meH thode de conception reposant sur des algorithmes geH neH tiques pour traiter simultaneH ment le proble`me de la conception de reconstructions optimales H et H de signaux avec un ordre de 2 =
* Corresponding author. Tel.: 886-35-731155; fax: 886-35-715971. E-mail address:
[email protected] (B.-S. Chen) 0165-1684/00/$ - see front matter ( 2000 Elsevier Science B.V. All rights reserved. PII: S 0 1 6 5 - 1 6 8 4 ( 9 9 ) 0 0 1 3 0 - 9
312
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
"ltre prescrit. Les algorithmes geH neH tiques sont des algorithmes d'optimisation et d'apprentissage, inspireH s initialement des processus de la seH lection naturelle et de la geH neH tique de l'eH volution. Its tentent de trouver l'optimum global sans tomber dans des minima locaux. Nous discutons dans cette eH tude a` la fois des "ltres IIR et FIR. Nous discutons aussi des proprieH teH s de convergence de notre algorithme de conception au moyen d'une meH thode de cham( nes de Markov. Deux exemples de conception pour la deH convolution optimale H et H avec des "ltres d'ordres "xes sont donneH s pour 2 = illustrer respectivement la proceH dure de conception et les performances de la meH thode de conception proposeH e. ( 2000 Elsevier Science B.V. All rights reserved. Keywords: Genetic algorithm; Deconvolution "lter; H criterion; H criterion 2 =
1. Introduction The deconvolution problem is to reconstruct a signal embedded in noise with the addition of removing the e!ect of any distortion in the signal transmission channel. The deconvolution problem is widely explored in the engineering literature, especially in many signal processing applications [18,6,7,4,22,2,12,17,5,13,11,20], including seismology, equalization, image restoration, etc. Deconvolution problems are generally solved via H optimal method by Wiener "ltering technique 2 from frequency domain perspective [22,2,12,17,5, 13,11,20], or by Kalman "ltering technique from time domain perspective [4,6,18]. Recently, minimax deconvolution "lter design based on H the= ory has received great attention for its robustness property [8,11,20]. Furthermore, if the orders of signal models and the channel dynamics are higher, the di$culty and the complexity on the design of the H and H optimal deconvolution "lters are 2 = increased so that the cost of realization increases and the operation time of the optimal deconvolution "lters become longer. In practical application, the "xed-order IIR or FIR H and H optimal 2 = deconvolution "lter designs are appealing for the simplicity of implementation and saving of operation time. However, to design "xed-order optimal deconvolution "lter is will new up to the knowledge of the authors. Since the "xed-order deconvolution "lter is used, conventional H and H optimization techniques 2 = cannot be employed to obtain a closed form of H and H optimal deconvolution "lters. How to 2 = specify the coe$cients of "xed-order deconvolution "lter to achieve the H or H minimization of 2 =
estimation error is a highly nonlinear minimization problem, in which many local minima may exist, and a local minimum may be easily reached via conventional algorithms. In this study, the design procedure of "xed-order IIR H and H deconvolution "lter is divided into 2 = two steps. In the "rst step, based on Jury's stability criterion, the stability domain of the coe$cients of denominator in a "xed-order deconvolution "lter is speci"ed. In the second step, the minima in the stability domain of the coe$cients (parameter space) will be searched via genetic algorithm. So the stability domain speci"cation not only can guarantee the stability of the deconvolution "lter but also can improve the convergent rate of the genetic algorithm. In FIR case, without the stability problem, the "rst step can be neglected and the coe$cients of "xed-order FIR deconvolution "lter are speci"ed via genetic searching in the parameter space. Since the stability of "xed-order FIR optimal deconvolution "lter is always guaranteed under channel parameter variation and noise uncertainty, it is more appealing from practical design perspective. However, the FIR H and H optimal "lter designs are not easy to 2 = obtain and still need further studies. Genetic algorithms are of parallel, global optimal search technique that emulates natural genetic operations because it simultaneously evaluates many points in the parameters space, it is more likely to converge toward the global solution [9,14]. It does not need to assume that the search space is di!erentiable or continuous, and can also iterate several times on each data received. The genetic algorithms apply operations inspired by the mechanics of natural selection to a populations of binary strings encoding the parameter space. At
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
each generation, it explores di!erent areas of the parameter space, and then directs the search to regions where they are of high probability to improve the performance. By working with a population of solutions, the algorithms can in e!ect search for many local minima and, thereby increase the likelihood of "nding the global minimum. Global optimization can be achieved via a number of genetic operators, e.g., reproduction, mutation, crossover and extinction. [9,3]. Based on the genetic algorithms, the coe$cients (parameters) of "xed-order deconvolution "lter will be tuned in the stability domain to achieve the H or H optimal design. Genetic algorithms are 2 = more suitable to solve the iterative "xed-order H and H optimization problem than other 2 = major searching method such as gradient-based algorithm and random searching algorithm due to the following reasons. First, the search space is large. Second, the performance surface does not need to be di!erentiable with respect to the change of parameters. Since, the gradient-based searching algorithm that depends on the existence of derivatives is ine$cient. Third, the likely "t terms are less likely to be destroyed under genetic operator, thereby often leading to faster convergence.
2. Problem description of deconvolution systems Consider a deconvolution system as shown in Fig. 1, the signal u(t) is generated by ¸(q~1) driven
313
by white noise w(t). The signal u(t) is transmitted through a channel system whose transform function is H(q~1) and is corrupted by additive colored noise r(t), which is generated by D(q~1) driven by white noise n(t). The received signal y(t) is given by y(t)"x(t)#r(t)"H(q~1)u(t)#D(q~1)n(t),
(1)
where q~1 is a backward shift operator, w(t) and n(t) are assumed to be zero-mean white noises with covariance E[w2(t)]"Q, E[n2(t)]"R, E[w(t) ) n(t)]"0,
(2)
respectively, where Q and R are assumed to be positive scales, and the system is assumed to have reached a statistical steady state. Let uL (t) denote the estimate of u(t) in the presence of both signal distortion through channel H(q~1) and the corrupted noise r(t). The objective of this paper is to "nd a "xed-order IIR deconvolution "lter h #h q~1#2#h q~i#2#h q~l 1 i l F(q~1)" 0 1#f q~1#2#f q~i#2#f q~l 1 i l (3) so that the H and H minimization of the estima2 = tion error e(t) is achieved with prescribed "lter order l by designer. The stability of the "xed-order deconvolution "lter F(q~1) in Eq. (3) is the most important thing. By Jury's stability test [1], the stability domain S of the coe$cients ( f ,2 f ) can 1 l be easily speci"ed. In order to guarantee the
Fig. 1. The block diagram of a linear discrete-time deconvolution system.
314
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
stability of F(q~1), the following constraint is imposed: ( f ,2 f )3S. 1 l
(4)
Remark 1. In FIR case, set f "f "2"f "0 1 2 l and neglect the stability condition in Eq. (4).
achieve the following minimization of the peak of power spectrum (i.e., minimax) [8,11,13,20]. min EU(z~1)E ": min sup DU(e~+w)D = F(z~1) F(z~1) w|*0,p+ ": min sup U(e~+w). F(z~1) w|*0,p+
From Fig. 1, the estimation error is given by
Remark 3. Since U(e~+w) is positive and symmetric within [!p, 0] and [0, p], we only consider the peak within [0, p] in Eq. (8).
e(t)"u(t)!u( (t), "[1!F(q~1)H(q~1)]¸(q~1)w(t) #F(q~1)D(q~1)n(t).
(5)
Remark 2. If the transmitted signal u(t) and channel noise r(t) are all white, then ¸(q~1)" D(q~1)"1 so that Eq. (4) is reduced to e(t)"[1!F(q~1)H(q~1)]u(t)!F(q~1)n(t). The power spectral density of e(t), U(z~1) is given by U(z~1)"(1!F(z~1)H(z~1))¸(z~1) ]Q¸H(z~1(1!F(z~1)H(z~1))H #F(z~1)D(z~1)NDH(z~1)FH(z~1),
(6)
where AH(z~1) denotes the complex conjugate of A(z~1), i.e., AH(z~1)"A(z). In the "xed-order H optimal deconvolution "l2 ter design case, as table "lter F(z~1) in (3) will be speci"ed to minimize the following least-mean square error [19].
Q
(8)
1 dz min Ee2(t)" min U(z~1) , (7) 2pi z @z@/1 F(z~1) F(z~1) where U(z~1) is de"ned in Eq. (6). In general, the calculus variation technique with the aid of spectral factorization and Cauchy theorem can be employed to solve the above H op2 timization design problem to obtain a full-order H optimal deconvolution "lter [5,12,17]. How2 ever, in "xed-order H optimal deconvolution "lter 2 design case, a closed-form solution is impossible and iterative searching methods are necessary to treat this problem. In the "xed-order H optimal = deconvolution "lter design case, the coe$cients of the "lter F(z~1) in Eq. (3) will be speci"ed to
Some techniques, such as inner}outer factorization, Nehari's theorem, and Riccati equation, etc. have been utilized to solve the above H decon= volution problem to obtain a full-order H opti= mal deconvolution after [11,13,20]. However, in "xed-order H deconvolution "lter case, the above = techniques cannot be applied to solve this problem. In this paper, genetic algorithms are also employed to treat this design problem. Genetic algorithms are optimization and machine learning algorithms, initially inspired from the processes of natural selection and evolution genetics. Therefore, in this study, genetic algorithms will be employed to specify the coe$cients of "xed-order "lter F(z~1) in Eq. (3) to solve the H optimal deconvolution "lter design problem in 2 Eq. (7) and the H optimal H optimal deconvolu= = tion "lter design problem in Eq. (8), respectively.
3. Fixed-order H and H optimal deconvolution 2 = 5lter designs In this section "xed-order H and H optimal 2 = deconvolution "lter designs will be introduced in the following two subsections, respectively. 3.1. Fixed-order H optimal deconvolution xlter 2 Consider the "rst-order H optimal deconvolu2 tion in Eq. (7). Let 1 I " 2 2pi
Q
@z@/1
U(z~1)
dz . z
(9)
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
And let us perform the following decomposition:
A
B
b b l c 1 i # i U(z~1)" 0 # + , (10) z#a z z~1#a z i i i/1 where a , b and c are all the function of the coe$i i i cients h for i"0, 1,2, l and f for i"1, 2,2, l. i i Suppose the poles z"!a are inside DzD"1 and i the poles z"!1/a are outside DzD"1 for all i i"1, 2,2, l. Substituting Eq. (10) into Eq. (9), we get 1 I " 2 2pi 1 " 2pi
Q Q A @z@/1
B(z)B(z~1) dz A(z)A(z~1) z
A
b b l c 0#+ i # i z z#a z~1#a i @z@/1 i i/1
BB
dz
l "+ b , (11) j j/0 where b , j"0, 1,2, l are the function of the coe$j cients. To evaluate the integral in Eq. (7), we employ the Cauchy's residue theorem
Q
1 B(z)B(z~1) dz min I "min 2 2pi A(z)A(z~1) z @z@ fi ,hi fi ,hi "min + Res(a ) k fi ,hi k l "min + b (h ,2, h , f ,2 f ), (12) j 0 l 1 l fi ,hi j/0 where the points a are the poles of U(z~1) inside k the unit circle. In general, I is a very highly nonlinear function 2 of the coe$cients h for i"0, 1,2, l and f for i i i"1, 2,2, l. There may exist many local minima. It is very di$cult to "nd the global minimum of I in Eq. (12), by the conventional searching 2 methods. From the above analysis, "lter is divided into the following these steps. 3.1.1. Design procedure I Step 0. Given the order l of the deconvolution "lter F(z~1), the channel model H(z~1), signal model ¸(z~1), noise model D(z~1), and noise covariances Q and R.
315
Step 1. Use Jury's stability test [1] to specify the parameter domain S of the coe$cients ( f 2 f ) of 1 l denominator of F(z~1) to guarantee the stability of deconvolution "lter. Step 2. Perform the factor decomposition
A
B
1 b b l c i # i U(z~1)" 0 # + . (13) z z#a z z~1#a i i i/1 Step 3. Find I ( f 2 f , h 2h ) from Eq. (11). 2 1 l 1 l Step 4. Solve the following minimization via genetic algorithm min I ( f 2 f , h 2h ). 2 1 l 0 l fi ,hi
(14)
3.2. Fixed-order H optimal deconvolution xlter = design Consider the "xed-order H optimal deconvolu= tion design problem in Eq. (8), let I ( f 2 f , h 2h ) = 1 l 1 l "EU(z~1)E " sup U(e~+w). (15) = w|*0,p+ Since we only need to "nd the peak of U(z~1) within the internal [0, p] of the frequency axis in Eq. (15), in this paper EU(z~1)E will be obtained = via searching the maximum value of U(e~+w) the within w3[0, p]. To solve the H minimization = problem min EU(z~1)E is to specify the coe$cient, = h for i"0, 1,2 l, and f for i"1, 2,2 l so that i i the maximum value of U(e~+w) is as small as possible. Therefore, the "xed-order H optimal "lter = design problem is equivalent to solving the following minimization min I ( f 2 f , h 2h ). = 1 l 0 l fi ,hi In general, I ( f 2 f , h 2h ) is also a highly = 1 l 0 l nonlinear function of the parameters ( f 2 f , 1 l h 2h ). At present, there is no good method to 0 l solve this problem. Therefore, genetic searching algorithm will be employed to solve the H minim= ization problem in Eq. (16). Based on the above analysis, the design procedure of "xed-order H optimal "lter design is divided into the follow= ing steps.
316
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
3.2.1. Design procedure II Step 0. Given the order l of the deconvolution "lter F(z~1), the channel model H(z~1) signal model ¸(z~1) and noise model D(z~1), and noise covariances Q and R. Step 1. Use Jury's stability test [1] to specify the parameter domain S of the coe$cients ( f 2 f ) of 1 l denominator of F(z~1) to guarantee the stability of deconvolution "lter F(z~1). Step 2. Find I ( f 2 f , h 2h ) from Eq. (15) = 1 l 0 l Step 3. Solve the following H minimization via = genetic algorithm min I ( f 2 f , h 2h ). = 1 l 1 l fi ,hi
(16)
From the analysis in the above two subsection, how to search the coe$cients of F(z~1) via genetic algorithm to solve the H minimization in Eq. (14) 2 or the H minimization in Eq. (17) is the key point = in our design.
4. Genetic algorithm in 5xed-order H and H 2 = optimal deconvolution 5lter designs 4.1. A simple description of genetic algorithm Genetic algorithms are optimization and machine learning algorithms. Their initial inspiration comes from the processes of natural selection (Darwinism) and evolutionary genetics. The underlying principles of genetic algorithm were "rst published by Holland in 1962 [14]. The mathematical framework was developed in the late 1960s, and has been presented in Holland's pioneering books, adaptive in Natural and Arti"cial Systems [14]. Genetic algorithms have been proven to be e$cient in many areas [3,9,10,14,16]. Recently, genetic algorithm are used to learn membership functions of fuzzy parameters and to select IF}THEN rules in fuzzy systems [16]. The more details about genetic algorithm could be found in Goldberg [3,9,21]. Genetic algorithms are powerful search algorithms based on the mechanics of natural selection and natural genetics. Some di!erences arise between genetic algorithm and conventional searching algorithms.
They are summarized as follows: (1) The algorithms work with a population of strings, searching many peaks in parallel, not a single point. (2) Genetic algorithms work directly with strings of characters representing the parameter set, not the parameters themselves. (3) Genetic algorithms use probabilistic transition rules to guide their search, not deterministic rules. (4) Genetic algorithms use objective function information instead of derivatives of other auxiliary knowledges. Genetic algorithms work with a population of binary strings, not the parameters themselves. The coding that has shown to be the optimal one is binary coding. So the parameter coe$cient set of the "xed-order deconvolution "lter, f ,2, f , h 2, h , in this paper could be coded as 1 l 0 l binary strings of 0's and 1's. For example, f "1, f "0, h "2, could be coded 0001, 0000, 1 2 0 0010. Note that these parameters are required to be tuned. The choice of bit number for each parameter is concerned with a `resolutiona we want in the search space. The binary strings } called chromosomes } then explore the search space, and each chromosome represents one possible solution to the problem. Genetic algorithms require only information concerning the quality of the solution produced by each parameter set (cost function values or evaluation function values). This di!ers from many optimization methods which require derivative information or, worse yet, complete knowledge of problem structure and parameters. Since genetic algorithm do not require such problem-speci"c information, they are more #exible than most searching algorithm. 4.2. A simple genetic algorithm for xxed-order H and H optimal deconvolution xlter designs 2 = A simple genetic algorithm is composed of four operators: (1) reproduction, (2) crossover, (3) mutation, and (4) extinction [14,3,10,16]. These operators are implemented by performing the basic tasks of coping strings, exchanging portion of strings, and changing the state of bit from 1's to 0's.
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
These operators ensure that the best members of the population will survive, and their information contents are preserved and combined to generate better o!springs. That is to improve the performance of next generation. It is shown in Schema theorem [9] that the genetic search algorithm will converge exponentially from the view point of schema. We describe the genetic algorithm as follows. 4.2.1. Coding and decoding Genetic algorithms work with a population of strings, not parameter themselves. For simplicity and convenience, binary coding is used in this article. Based on binary coding method, every element I i of the parameter I is coded as string of length ¸, k k which consists of `0'sa and `1'sa, for the desired resolution R . In general, we have i ; !¸ i, R" i i 2L!1 where ; and ¸ are the upper and lower bounds of i i the parameter I i . For example, suppose the parak meter domain of ( f , h ) is 0 0 0.05)f )0.9, 0 0)h )0.63 0 and the desired resolutions are R "0.06 and 1 R "0.01. Then the bit numbers need for f and 2 0 h are four and six, respectively. Therefore, we can 0 code f and h according to the following 0 0 f code h code 0 0 0 0000 0 000000 0.06 0001 0.01 000001 F
F
F
F
0.9
1111 0.63 111111
minimize in Eq. (14) or Eq. (17) and the cost function is only de"ned for those parameters ( f 2 f , h 2 h )3S. 1 l 0 l Our objective is to search ( f 2 f , h 2 h ) with 1 l 0 l ( f 2 f , h 2 h )3S to achieve the minimization 1 l 0 l of Eq. (18). Then the cost function E( f 2 f , 1 l h 2 h ) takes a chromosome (a possible 0 l ( f 2 f , h 2 h )) and returns a value. The value 1 l 0 l of the cost is then mapped into a "tness value F( f 2 f , h 2 h ) so as to "t into the genetic 1 l 0 l algorithm. The "tness value is a reward based on the performance of the possible solution represented by the string, or it can be though of as how well a "xed-order deconvolution "lter can be tuned according to the string to actually minimize the cost function. The better the solution encoded by a string (chromosome), the higher the "tness. To minimize E( f 2 f , h 2 h ) is equivalent to 1 l 0 l getting a maximum "tness value in genetic searching process, a chromosome that has a lower E( f 2 f , h 2 h ) should be assigned a large "t1 l 0 l ness value. Then genetic algorithm tries to generate better o!springs to improve the "tness. Therefore, a better deconvolution "lter could be obtained by better "tness in genetic algorithms. So are let the "tness value as 1 F( f 2 f , h 2 h )J . (18) 1 l 0 l E( f 2 f , h 2 h ) 1 l 0 l There are a number of methods to perform this mapping known as "tness techniques. In this paper, we use so-called windowing [16] as described in Fig. 2. In Fig. 2, E is the largest cost value in the 80345 generation being evaluated, and E is the smallest "%45 cost value. F and F are the corresponding "t8 " ness values. Furthermore, the relation between F( f 2 f , h 2 h ) and E( f 2 f , h 2 h ) could 1 l 0 l 1 l 0 l be expressed as a linear equation F( f 2 f , h 2 h )"mE( f 2 f , h 2 h )#b, 1 l 0 l 1 l 0 l (19)
4.2.2. Fitness and cost function In this paper, the cost function is de"ned as E( f 2 f , h 2 h )"I ( f 2 f , h 2 h ), 1 l 0 l k 1 l 0 l ( f 2 f , h 2 h ),3S, k"2, R, 1 l 0 l
317
(17)
where k"2 for H optimal deconvolution case 2 and k"R for H optimal deconvolution case. = The right-hand side of Eq. (18) is what we want to
where the constant m and b are computed by F ,F ,E and E in each generation as " 8 "%45 80345 F !F F !F " 8 , b"F !E " 8 . m" " "%45 E !E E !E "%45 80345 "%45 80345 (20)
318
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
Fig. 2. The relation between F( f ,2 f , h ,2, h ) and F( f ,2 f , h ,2, h ). 1 l 0 l 1 l 0 l
4.2.3. Reproduction Reproduction is based on the principle of survival of the "tness. The "tness of the ith string, F is i assigned to each individual string in the population where higher F means better "tness. These strings i with large "tness would have large number of copies in the new generation. For example, in Roulette wheel selection (see Fig. 3), the ith string with high "tness values F is given a proportionately i high probability of reproduction, P , according to i this distribution F P" i . i +F i
(21)
Once the strings are reproduced or copied for possible use in the next generation, they are reproduced in a mating pool where they await the action of the other three operators, crossover, mutation and extinction (see Fig. 3). 4.2.4. Crossover If the chromosomes are operated only by reproduction, they search toward the best existing individes but does not create any new individuals. By
the second operator (i.e, crossover) strings exchange information via probabilistic decision. Crossover provides a mechanism for strings to mix and match their desirable qualities through a random process. After reproduction simple crossover proceeds in three steps. First, two newly reproduced string are selected from the mating pool produced by reproduction. Second, a position along the two strings is selected uniformly at random. This is illustrated below where two binary coded strings, ( f 2 f , h 2 h ) and 1 l 0 lA ( f 2 f , h 2 h ) are shown aligned for cross1 l 0 lB over: ( f 2 f , h 2 h ) "11000010 . . . . . . 1111, 1 l 0 lA ( f 2 f , h 2 h ) "11000001 . . . . . . 1101, 1 l 0 lB Csplice point,
the third step is to exchange all characters following the crossing sit. For example, the two strings ( f 2 f , h 2 h ) and ( f 2 f , h 2 h ) with 1 l 0 lA 1 l 0 lB a crossover sit at 4 become: ( f 2 f , h 2 h ) "11000011 . . . . . . 1101, 1 l 0 l A{ ( f 2 f , h 2 h ) "11000000 . . . . . . 1111. 1 l 0 l B{
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
319
Fig. 3. Flow chart of design procedure.
Although crossover uses random choice, it should not be thought of as a random walk through the search space. When combined with reproduction, it is an e!ective means of exchanging information and combing portions of high-quality solutions. 4.2.5. Mutation Reproduction and crossover give genetic algorithms most of their search power. The third operator, mutation, enhances an ability of genetic algorithms to "nd near optimal solution. Mutation is the occasional alternation of a value at a particular string position. An insurance policy against the
permanent loss of any simple bit and it is applied with a low probability that is chosen so that on the average one string in the population is mutated: Bmutation ( f 2 f , h 2 h )"1100001001 . . . . . . 1111, 1 l 1 l ( f 2 f , h 2 h )"1100001101 . . . . . . 1111. 1 l 1 l In the case of binary coding, the mutation operator simply #ips the state of a bit from 0 to 1 or vice versa. Mutation should be used sparingly because it is a random search operator, and with high mutation rates, the algorithm should become little more than a random search.
320
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
4.2.6. Extinction Extinction eliminates all of the chromosomes in the current generation except the chromosome corresponding to the minimum estimation error. N!1 (N is number of chromosomes in each generation) chromosomes are randomly generated to "ll out the population. M!1 (M is number of chromosomes chosen as parents for mating) chromosomes associated with the least estimation errors among these immigrants are then selected as the parents. Together with the surviving chromosome, these are allowed to mate as usual to form the next generation. For convenience, we say that another era begins. Extinction is analogous to a particular time varying mutation rate in probability is close to 1 at the beginning of each era and then small for the remaining generations within the era.
4.3. Convergence analysis A genetic algorithm with elitist strategy is chosen in this paper to guarantee the convergence of the parameter estimation algorithm. The elitist strategy always survives the best chromosome intact into next generation. The convergence of genetic algorithm with elitist strategy has been proved in [21] from viewpoint of Markov-chain. It has also been proved that the genetic algorithm converges to the global optimum as the generation approaches to in"nity. In this paper, the convergence analysis of the proposed parameter estimation algorithm can be discussed by following the path of proof of genetic algorithm in [21]. Let all of the chromosomes in a generation be arranged into a binary string with length N¸ (N is t the population size, and ¸ is the length of each t chromosome). The size of the state space is denoted by S"2NLt . We denote the Xt"(Xt , Xt ,2, Xt ), N 1 2 in which Xt denotes the kth chromosome in k a population at generation t. A simple genetic algorithm uses three operations, i.e., reproduction, crossover, and mutation, to change state, from the state of the operations of the population at generation t, say state j, to another state at generation t#1, say state k. The evolution of a simple genetic algorithm is regarded as a "nite ergodic Markov-
chain [21]. The probabilistic changes of the genes within the population caused by genetic operations are captured by the transition matrix ¹ of a Markov-chain with dimension N¸ ]N¸ , in which the t t matrix ¹"R ) C ) M (R , C and M denote the e r u e r u intermediate transitions caused by reproduction, crossover, and mutation, respectively). The "tness F(Xt ) of the kth chromosome at genk eration t is a positive function. We wish to "nd the global maximum FH of F. Let us denote the best chromosome at generation t as XK t"MXt : 1)k)NDF(Xt )" max F(Xt )N. h k k 1xhxN Next, we will show the transition matrix ¹ of a simple genetic algorithm is positive (i.e., every element ¹ '0 for i, j"1, 2,2, N¸ ). i,j t The intermediate transition matrices of ¹ (i.e., R , C , and M ) are stochastic, because each of e r u them is a state transition matrix, and each state of S is translated probabilistically to another state of S. The transition state j changes to state k by mutation, whose probability can be written as m "pHj,k (1!p )NLt ~Hj,k '0 for all j, k3S, where j,k m m p is the mutation probability, and H is the m j,k Hamming distance [21]. Thus, M is positive (if u p O0). Next, let ;"C M , ¹"R ;. Since C is m r u e r a stochastic transition matrix, there exists at least one positive element in each row of C , then r ; "+NLt C i,k M k,j '0 for all i, j31, 2,2, N¸ . i,j k/1 r u t Similarly, R is stochastic, so ¹ is positive. e After the positiveness of the transition matrix ¹ is shown, the convergence of estimation algorithm is discussed as following. In an ergodic Markov-chain the expected transition time between initial state j and any other state k is "nite [15]. Therefore, a simple genetic algorithm can reach the best "tness (i.e., the global maximum "tness value) in "nite generations. However, the best "tness may not survive in the next generation because in a simple genetic algorithm the best chromosome is not reserved in each generation. In this paper, in order to avoid the drawback of the simple genetic algorithm, the best chromosome will remain intact in the generation, i.e., the genetic algorithm contains the elitist strategy. In this situation, if we assume that FH is found at generation t, the
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
F(Xt)"FH for i*t. Therefore, we can prove lim F(Xt)"FH. t?= Actually, the property that the proposed algorithm can converge to the global optimum as time approaches in"nity is not of practical interest because the algorithm is generally performed over "nite generations. By using extinction in the genetic algorithm, we can increase the convergence rate to the global maximum of the likelihood function.
321
mutation probability, and ¸ is the length of every parameter. Example 1. Consider a non-minimum phase deconvolution system in Fig. 1 with the received signal corrupted by colored noise such that y(t)"H(z~1)u(t)#D(z~1)n(t), u(t)"¸(z~1)w(t), where the signal, channel, and noise models are described as follows: 1#0.8z~1 ¸(z~1)" , H(z~1) 1#0.5z~1 (1#0.9z~1)(1#0.8z~1) " , (1#0.1z~1)(1!0.1z~1)
5. Design examples In this section, two "xed-order deconvolution "lter design examples based on H and H opti2 = mal criteria are given to illustrate both the design procedure and the performance of the proposed method. In the following example, the parameters of genetic algorithm are set as follows: N"320, M"80, P "0.01, ¸"17 bit, . where the N is the number of population, M is the number of chosen parents for mating, P is the .
(1#0.9z~1)(1#0.8z~1) . D(z~1)" (1#0.1z~1)(1!0.1z~1) The driving signal Mw(t)N and disturbance noise Mn(t)N are also assumed to be independent, stationary and white with zero mean and variances as given by Q"1, R"0.1, respectively. By the proposed design algorithms, the "rst order H optimal 2 convolution "lters F (z~1), the second-order 1 H and H optimal deconvolution "lters F (z~1) 2 = 2
Fig. 4. The undisturbed signal x(t) and the colored noised r(t) in example 1.
322
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
Fig. 5. The cost function I ( f , h , h ) of the "rst-order H optimal deconvolution "lter F (z~1) via genetic algorithm in example 1. 2 1 0 1 2 1
Fig. 6. u( (t) by the "rst-order H optimal deconvolution "lter F (z~1) in example 1. 2 1
and F (z~1), and the "fth-order FIR H optimal 3 = deconvolution "lters F (z~1) are given in order 4
1!0.2802z~1#0.1642z~2 F (z~1)"0.635] , 3 1#1.0109z~1#0.1562z~2
1!0.4601z~1 F (z~1)"0.6831] , 1 1#0.8765z~1
F (z~1)"1#0.8333z~1!5.6667z~2!0.6667z~3 4 #2.8334z~4#z~5.
1#0.0745z~1!0.0317z~2 , F (z~1)"0.9374] 2 1#1.6248z~1#0.6249z~2
The poles and zeros by the proposed method are shown in Table 2. The impulse response of the
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
323
Fig. 7. The cost function I ( f , f , h , h , h ) of the second-order H optimal deconvolution "lter F (z~1) via genetic algorithm in 2 1 2 0 1 2 2 2 example 1.
Fig. 8. The cost function I ( f , f , h , h , h ) of the second-order H optimal deconvolution "lter F (z~1) via genetic algorithm in = 1 2 0 1 2 = 3 example 1.
channel followed by the deconvolution "lters are shown in Figs. 13}16. The undisturbed signal x(t) and colored noise r(t) of the convolution system are shown in Fig. 4. The cost function I ( f , h , h ) of 2 1 0 1 the "rst-order H optimal "lter via genetic algo2
rithm is given in Fig. 5. The cost function converges exponentially. The input signal u(t) and the reconstructed signal u( (t) of the "rst-order H 2 deconvolution "lter F (z~1) are shown in Fig. 6, 1 the estimation error I K0.3. The cost function 2
324
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
Fig. 9. u( (t) by the "rst-order H optimal deconvolution "lter F (z~1) in example 1. 2 2
Fig. 10. u( (t) by the "rst-order H optimal deconvolution "lter F (z~1) in example 1. = 3
I ( f , f , h , h , h ) and I ( f , f , h , h , h ) of the 2 1 2 0 1 2 = 1 2 0 1 2 second-order H and H optimal deconvolution 2 = "lter F (z~1) and F (z~1) in example 1 are shown 2 3 in Figs. 7 and 8, respectively. It can be seen that it converges exponentially. The input signal u(t) and the reconstructed signal u( (t) by F (z~1) and 2
F (z~1) are shown in Figs. 9 and 10, respectively. 3 The results indicate that we can obtain a good reconstruction performance by the second-order H and H , the estimation error I K0.1 and 2 = 2 I K0.27. The cost function I (h ,2, h ) of the = = 0 5 "fth-order FIR H optimal deconvolution "lter =
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
325
Fig. 11. The cost function I (h , h , h , h , h h ) of the "fth-order H optimal deconvolution "lter F (z~1) via genetic algorithm in = 0 1 2 3 4 5 = 4 example 1.
Fig. 12. u( (t) by the "fth-order FIR H
=
optimal deconvolution "lter F (z~1) in example 1. 4
F (z~1) is shown in Fig. 11, and the reconstructed 4 signal u( (t) by F (z~1) is given in Fig. 12. The simu4 lation results indicate that we can obtain good reconstruction performance by low order FIR
H deconvolution "lter. Since FIR deconvolution = "lter has not stability problem under parameter variation and exogenous noise, it is more appealing for practical design applications (Figs. 13}16).
326
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331 Table 1 The performance of the various case as some parameters are changed of example 1 Error variance Filter Channel (1#0.9z~1)(1#0.8z~1) H(z~1)" (1#0.1z~1)(1!0.1z~1)
F 2 (Second-order H 2 deconvolution "lter)
F 3 (Second-order H deconvolution "lter) =
0.12
0.27
0.165
0.283
0.14
0.268
(Original channel) (1#0.95z~1)(1#0.85z~1) H(z~1)" (1#0.1z~1)(1!0.1z~1) (Change channel) (1#0.9z~1)(1#0.7z~1) H(z~1)" (1#0.1z~1)(1!0.1z~1) (Change channel)
Fig. 13. Impulse response of the channel of example 1 followed by the order-"rst H deconvolution "lter. 2
Remark 4. In addition, for illustrating the robustness of the H criterion, we assume that the para= meters of channel H(z~1) has been changed due to variation in the parameters. The various cases are shown in Table 1. From Table 1, it is found, while
Fig. 14. Impulse response of the channel of example 1 followed by the order-two H deconvolution "lter. 2
channel's parameters have some variation, the deconvolution "lter based on the H criterion is = more robustness than the H criterion. 2 Example 2. In the deconvolution system of Example 1 if the transmission channel H(z~1) is
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
Fig. 15. Impulse response of the channel of example 1 followed by the order-two H deconvolution "lter. =
Fig. 16. Impulse response of the channel of example 1 followed by the order-"fth H deconvolution "lter. =
changed as (0.8#z~1)(1#0.5z~1) H(z~1)" (1#0.3z~1)(1!0.3z~1)(1#0.1z~1) and the other system parameters are same as Example 1. The undisturbed signal x(t) and colored noise w(t) are shown in Fig. 19.
327
Fig. 17. Impulse response of the channel of example 2 followed by the order-"rst H deconvolution "lter. 2
Fig. 18. Impulse response of the channel of example 2 followed by the order-two H deconvolution "lter. =
By the proposed design method based on genetic algorithm, the "rst-order H and the second-order 2 H optimal deconvolution "lters are obtained in = order 1!0.3839z~1 F (z~1)"0.9353] , 5 1#0.9243z~1 1!0.083z~1#0.0143z~2 . F (z~1)"0.5468] 6 1#1.0771z~1#0.2147z~2
328
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
Fig. 19. The undisturbed signal x(t) and the colored noise r(t) in example 2.
Fig. 20. The cost function I ( f , h , h ) of the "rst-order H optimal deconvolution "lter F (z~1) via genetic algorithm in example 2. 2 1 0 l 2 5
The poles and zeros by the proposed method are shown in the Tables 2 and 3. The impulse response of the channel followed by the deconvolution "lters are shown in Figs. 17 and 18. The cost function
I ( f , h , h ) of the "rst-order H optimal decon2 1 0 1 2 volution "lter F (z~1) in the genetic search pro5 cessing is shown in Fig. 20, it again converges exponentially. The reconstructed signal u( (t) by
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
F (z~1) is shown in Fig. 21. The estimation error 5 I K0.17. The convergence of cost function 2 I ( f , f , h , h , h ) of the second-order H opti= 1 2 0 1 2 = mal deconvolution "lter F (z~1) is shown in 5
329
Fig. 22, and the reconstructed signal u( (t) is given in Fig. 23. Since the reconstruction performances I and 2 I of several orders of optimal deconvolution can =
Fig. 21. u( (t) by the "rst-order H optimal deconvolution "lter F (z~1) in example 2. 2 5
Fig. 22. The cost function I ( f , f , h , h ) of the "rst-order H optimal deconvolution "lter F (z~1) via genetic algorithm in example 2. = 1 2 0 l = 6
330
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
Table 2 Example 1
H(z~1) (channel)
F (z~1) 1 ("rst-order H 2 deconvolution "lter)
F (z~1) 2 (second-order H 2 deconvolution "lter)
F (z~1) 3 (second-order H = deconvolution "lter) !0.8205 !0.1904
Poles
0.1 !0.1
!0.8765
!0.6251 !0.9997
Zeroes
!0.8
!0.4601
0.1446
0.1401#0.3802i
!0.2191
0.1401!0.3802i
!0.9
F (z~1) 4 ("fth-order FIR H deconvolution = "lter)
!2.7098 1.8714 !0.4487#0.1311i !0.4487!0.1311i
Fig. 23. u( (t) by the second-order H optimal deconvolution "lter F (z~1) in example 2. = 6
Table 3 Example 2
Poles
Zeroes
F (z~1) 5 H(z~1) ("rst-order H 2 (channel) deconvolution "lter)
F (z~1) 6 (second-order H = deconvolution "lter)
0.3 !0.3 !0.1
!0.8130
!0.5 !1.25
!0.9243 !0.2641 0.3839
0.0415#0.1122i 0.0415!0.1122i
be also obtained in the design procedure, these also gives an important reference for designer to accept this order of optimal deconvolution "lter or not. In general, if the statistical properties of signals and the channel H(z~1) are exactly known, H 2 optimal deconvolution "lters are better than H optimal deconvolution "lters of the same = order. However, as the statistical properties of signals are uncertain and channel H(z~1) has parameter variations, the H deconvolution "lters are = most robust than H deconvolution "lters 2 [11,13,20].
B.-S. Chen, J.-C. Hung / Signal Processing 80 (2000) 311}331
6. Conclusion In this paper, design methods of "xed-order H and H optimal deconvolution "lters have 2 = been introduced via genetic algorithm. The near H and H optimizations are obtained by genetic 2 = operations such as reproduction, mutation, crossover and extinction. The simulations results indicate that the genetic-based design algorithms converge exponentially and the reconstruction performance is acceptable even if the order of the deconvolution "lter is lower. There is a tradeo! between the order of deconvolution "lter and the reconstruction performance. The proposed design methods are suitable for lower order IIR and FIR optimal deconvolution "lter design with simplicity of implementation and saving of operation time and are useful for practical application in the signal reconstruction problems with high order channel and signal model.
References [1] K.J. Astrom, B. Wittenmark, Computer Controlled Systems: Theory and Design, Prentice-Hall, Englewoods Cli!s, NJ, 1984. [2] A.B. Berkhout, Least-square inverse "lter and wavelet inverse, Geophysics 42 (1987) 1369}1383. [3] B.P. Buckles, P.E. Petry (Eds.), Genetic Algorithms, IEEE Computer Society Press., 1992. [4] Y.L. Chen, B.S. Chen, Mini-max robust deconvolution "lter under stochastic parametric and noise uncertainties, IEEE Trans. Signal Process. 42 (1994) 32}45. [5] B.S. Chen, S.C. Peng, Optimal deconvolution "lter design based on orthogonal principle, Signal Processing 25 (3) (1991) 361}372. [6] C.Y. Chi, J.M. Mendel, Performance of minima-variance deconvolution "lter, IEEE Trans. Acoust. Speech, Signal Process. ASSP-35 (2) (1989) 217}226.
331
[7] K. De Jong, Learning with the genetic algorithm: an overview, Machine Learning 3 (1988) 121}137. [8] B.A. Francis, A Course in H Control Theory, Springer, = Berlin, 1987. [9] D.E. Goldberg, Genetic Algorithms in Search, Optimization, and Machine Learning, Reading, Addision-Wesley, Reading, MA, 1989. [10] J.J. Grefensteete, Optimization of control parameters for genetic algorithms, IEEE Trans. Syst., Man, Cybern. 6 (10) (1986) 122}128. [11] M.J. Grimble, H "xed-lay smoothing "lter for scale sys= tems, IEEE Trans. Signal Process. 39 (9) (1991) 1955}1963. [12] M.J. Grimble, H inferential "ltering prediction and 2 smoothing with application to rolling mill gauge estimation, IEEE Trans. Signal Process. 42 (8) (1994) 2078}2093. [13] M.J. Grimble, A.E. Sayed, Solution of the H optimal = "ltering problem for discrete time systems, IEEE Trans. Acoust., Speech, Signal Process. ASSP-38 (7) (1990) 1092}1104. [14] J.H. Holland, Adaptive in Natural and Arti"cial Systems, University of Michigan Press, Ann Arbor, MI, 1975. [15] M. Iosifescu, Finite Markov Processes and Their Applications, Wiley, Chichester, 1980. [16] C.L. Kary, E.J. Gentry, Fuzzy control of pH using genetic algorithms, IEEE Trans. Fuzzy System (1) (1993). [17] D.G. Lainiotis, S.K. Katsikas, S.D. Likothanasis, Adaptive deconvolution of seismic signals-performance, computational analysis, parallelism, IEEE Trans. Acoust. Speech, Signal Process. ASSP-36 (11) (1988) 1715}1734. [18] J.M. Mendel, Optimal Seismic Deconvolution: An Estimation-based Approach, Academic Press, New York, 1983. [19] A. Papoulis, Probability, Random Variable and Stochastic Processes, 2nd Edition, McGraw-Hill, New York, 1984. [20] S.C. Peng, B.S. Chen, A deconvolution "lter for multichannel non-minimum phase systems via mini}max approach, Signal Processing 36 (1) (1994) 71}91. [21] G. Rudolph, Convergence analysis of canonical genetic algorithms, IEEE Trans. Neural Net. 5 (1994) 96}101. [22] M.T. Silvia, E.A. Robinson, Deconvolution Underline of Geophysical Time Series in the Exploration for Oil and Natural Gas, Elsevier, New York, 1979.