Ultra wideband OFDM channel estimation ... - Wiley Online Library
EUROPEAN TRANSACTIONS ON TELECOMMUNICATIONS Eur. Trans. Telecomms. 2008; 19:761–771 Published online 11 September 2008 in Wiley InterScience (www.interscience.wiley.com) DOI: 10.1002/ett.1324
Ultra wideband OFDM channel estimation through a wavelet based EM-MAP algorithm† S. M. S. Sadough1,2∗ , M. Ichir2 , P. Duhamel2 and E. Jaffrot3 1 UEI, ENSTA, 32 boulevard Victor, 75015 Paris, France des Signaux et Syst`emes, CNRS-Sup´elec, Plateau de Moulon, 91190 Gif–sur–Yvette, France 3 Universidad Nacional de San Martin, B1650ANQ San Martin, Pcia. de Buenos Aires, Argentina
2 Laboratoire
SUMMARY Ultra wideband (UWB) communications involve very sparse channels, since the bandwidth increase results in a better time resolution. This property is used here to propose an efficient algorithm jointly estimating the channel and the transmitted symbols. More precisely, this paper introduces an expectation-maximisation (EM) algorithm within a wavelet domain Bayesian framework for semi-blind channel estimation of multiband orthogonal frequency-division multiplexing (MB-OFDM) based UWB communications. A prior distribution is chosen for the wavelet coefficients of the unknown channel impulse response (CIR) in order to model a sparseness property of the wavelet representation. This prior yields, in maximum a posteriori (MAP) estimation, a thresholding rule within the EM algorithm. We particularly focus on reducing the number of estimated parameters by iteratively discarding ‘insignificant’ wavelet coefficients from the estimation process. Simulation results using UWB channels issued from both models and measurements show that under sparsity conditions, the proposed algorithm outperforms pilot based channel estimation in terms of mean square error (MSE) and bit error rate (BER). Moreover, the estimation accuracy is improved, while the computational complexity is reduced, when compared to traditional semi-blind methods. Copyright © 2008 John Wiley & Sons, Ltd.
1. INTRODUCTION An ultra wideband (UWB) radio signal is defined as any signal whose bandwidth is larger than 20% of its centre frequency or greater than 500 MHz [1]. In recent years, UWB system design has experienced a shift from the traditional ‘single-band’ radio that occupies the whole 7.5 GHz allocated spectrum to a ‘multiband’ design approach [2]. In these multiband UWB systems, the available spectrum is divided into several subbands, each one occupying approximately 500 MHz. Multiband orthogonal frequency-division multiplexing (MB-OFDM) [3] is an appropriate candidate for multiband UWB which enables high data rate UWB transmission to inherit all the strength of OFDM that is widely used in numerous communication systems. This approach uses a
conventional coded OFDM system [4] together with bit interleaved coded modulation (BICM) [5] and frequency hopping over different subbands to improve the diversity and to enable multiple access. It is well known that an efficient detection of the transmitted symbols requires a good knowledge of the channel realisations at the receiver. The required channel estimation is often based on the transmission of known training symbols (also called pilots) [6]. However, the length of the required training sequence may become large when the channel impulse response (CIR) is highly resolved, such as in UWB OFDM systems. Besides, obtaining an accurate channel estimate in highly mobile environments only through the use of pilots, would require inserting multiple training symbols per frame, which can result in a waste of bandwidth and power.
* Correspondence to: Seyed Mohammad Sajad Sadough, UEI, ENSTA, 32 boulevard Victor, 75015 Paris, France. E-mail: [email protected]. † A previous version of this paper was presented in the 13th European Wireless Conference (EW 2007), Paris, France.
Copyright © 2008 John Wiley & Sons, Ltd.
Accepted 16 June 2008
762
S. M. S. SADOUGH ET AL.
In order to reduce the number of pilot symbols, the information carried by the observations corresponding to the data symbols can also be exploited for channel estimation. An attractive technique in this area is that based on the so-called expectation-maximisation (EM) algorithm [7] which provides a framework to iteratively estimate the unknown channel thanks to the a posteriori probabilities on the unknown data. Moreover, by using the EM algorithm, one can naturally combine the iterative process of channel estimation with the decoding operation of encoded data. Several papers have recently reported promising results on the combination of EM-based channel estimation and data decoding process. For instance, iterative or ‘turbo’ processing that includes the channel estimation into the iterative process of decoding turbo like codes is addressed in Reference [8–10]. In particular, Reference [9] addresses a turbo channel estimation based on a Karhunen–Loeve expansion of the unknown channel for reducing the number of estimated parameters. However, the proposed scheme requires the channel covariance matrix as an a priori information which can not be available in realistic situations. In Reference [11], the authors present several EM based algorithms for time or frequency domain channel estimation of an OFDM system that is subject to slow timevarying frequency-selective fading. In Reference [12], Lu et al. proposed an iterative receiver for space-time blockcoded OFDM systems based on a reduced complexity EM algorithm. Though iterative joint channel estimation and data detection scheme outperforms receivers using a pilot-only based channel estimation, it has a higher complexity that may be a critical concern for its practical implementations. This complexity is mainly driven by the number of estimated parameters for channel updating and the decoding algorithm within each iteration. In this work, we consider a semi-blind joint channel estimation and data detection scheme for MB-OFDM systems based on the EM algorithm, with the objective of minimising the number of estimated parameters and enhancing the estimation accuracy. This is achieved by expressing the unknown UWB CIR in terms of its discrete wavelet expansion, which has been shown to provide a parsimonious representation [13, 14]. Our wavelet domain framework enables us to choose a particular a priori distribution for the channel wavelet coefficients that renders the maximum a posteriori (MAP) channel estimation equivalent to a hard thresholding rule at each iteration of the EM algorithm. The latter is then exploited to reduce the estimator’s computational complexity by discarding ‘insignificant’ wavelet coefficients from the Copyright © 2008 John Wiley & Sons, Ltd.
estimation process. Moreover, the adopted wavelet domain prior model adapts itself to the actual channel scenario since its parameters are learnt from the observed data during the channel estimation process. More precisely, the prior model is not restricted to a specific propagation environment. We show that when the channel has a sparse wavelet expansion, the proposed algorithm outperforms the classical EM estimators and provides a significant reduction in the number of estimated parameters. Hence, the proposed receiver has not to modify the prior information for different channel scenarios. This paper is organised as follows: Section 2 introduces the MB-OFDM system and its equivalent model for wavelet domain channel estimation. Section 3, first presents the derivation of the MAP channel estimator based on the EM algorithm and then shows how the number of estimated parameters can be reduced through the EM iterations. Section 4 presents the combination of the channel estimation with the decoding operation, as well as some implementation issues. Section 5 illustrates, via simulations, the performance of the proposed receiver in different realistic UWB channel environments and Section 6 concludes the paper. Notational conventions are as follows: Dx is a diagonal matrix with diagonal elements x = [x1 , . . . , xN ]T , Ex [·] refers to expectation with respect to x, IN denotes an (N × N) identity matrix; CN(m, ) denotes complex Gaussian distribution with mean m and covariance matrix ; Card{·} denotes set cardinality; (·)T and (·)H denote matrix or vector transpose and Hermitian transpose, respectively.
2. SYSTEM MODEL AND WAVELET DOMAIN PROBLEM FORMULATION 2.1. MB-OFDM transmission We consider the MB-OFDM transmission proposed in Reference [3]. This scheme divides the spectrum between 3.1 and 10.6 GHz into several non-overlapping subbands each one occupying 528 MHz of bandwidth [3]. Information is transmitted using OFDM modulation over one of the subbands in a particular time slot. As shown in Figure 1, the transmitter architecture for the MB-OFDM system is very similar to that of a conventional wireless OFDM system. The main difference is that MB-OFDM uses a timefrequency code (TFC) to select the centre frequency of different subbands. This code is used to provide frequency diversity as well as to distinguish between multiple users (see Figure 2). Here, we consider MB-OFDM in its basic Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
763
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
frequency domain as† Ym = DSm Hm + Zm
m = 1, . . . , Msym
(2)
where Ym = [y1,n , y2,n+1 , y3,n+2 ]T , Sm = [s1,n , s2,n+1 , ¯ 1,n , H ¯ 2,n+1 , H ¯ 3,n+2 ]T and Zm = [z1,n , s3,n+2 ]T , Hm = [H T z2,n+1 , z3,n+2 ] are (M × 1) vectors, with M = 3N and Msym = Nsym /3. In the remainder, unless otherwise mentioned, we will not write the time index m for notational convenience. Figure 1. TX architecture of the multiband OFDM system.
2.2. A proper model for wavelet domain channel estimation
Figure 2. Example of time-frequency coding for the multiband OFDM system: TFC = {1, 3, 2, 1, 3, 2, . . .}.
mode, i.e. employing the first three subbands with N data subcarriers over each subband. Let us consider a single-user MB-OFDM transmission. At the receiver, assuming a cyclic prefix (CP) longer than the channel maximum delay spread and perfect carrier synchronisation, OFDM converts a frequency-selective channel into N parallel flat fading subchannels [4]. Under these conditions, the transmission of the nth OFDM symbol (inside a frame of size Nsym ) over the ith subband can be written as [4] ¯ i,n + zi,n yi,n = Dsi,n H
(1)
for i ∈ {1, 2, 3} and n = 1, . . . , Nsym ; (1 × N) vectors ¯ i,n denote respectively the received and yi,n , si,n and H transmitted symbols, and the channel frequency response; the noise vector zi,n is assumed to be a zero-mean circularly symmetric complex Gaussian random vector with distribution CN(0, σ 2 IN ). In what follows, we group the data and the observation corresponding to three subbands of the MB-OFDM system in a single vector. This operation can be written in the Copyright © 2008 John Wiley & Sons, Ltd.
In order to take advantage of the wavelet based estimation, the CIR is expressed in terms of its orthogonal discrete wavelet coefficients. Let FM,L be the truncated fast Fourier transform (FFT) matrix constructed from the (M × M) FFT matrix by keeping the first L columns where L is the length of the CIR over a group of three subbands. We define W as the (L × L) orthogonal discrete wavelet transform (ODWT) matrix. The unknown channel frequency response can be expressed as H = FM,L WH g, where g is the (L × 1) vector of the CIR wavelet coefficients. The observation model (2) is rewritten as Y = DS Tg + Z
(3)
where T = FM,L WH . In this model, although the channel is practically used (by the transmitter) by slices of 528 MHz bandwidth (corresponding to a single subband), on the receiver side, three received OFDM symbols are grouped for estimating the wavelet coefficients of the CIR, taken over all three subbands (1.584 GHz bandwidth). This is motivated by the fact that estimating the channel over a wider bandwidth leads to a sparser representation in the wavelet domain. Most of wavelet based estimation algorithms, rely on an observation model in which the unknown wavelet coefficients are corrupted by an AWGN [13]. Since the model (3) does not provide such a framework, our first step consists in enforcing this property. In order to do so, the AWGN in Equation (3) is split into two independent Gaussian terms as suggested in Reference [15] Z = DS Z1 + Z2
(4)
† For
the sake of notational brevity in Equation (2), we have assumed that the TFC is equal to {1, 2, 3, . . .}.
Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
764
S. M. S. SADOUGH ET AL.
where Z1 and Z2 are (M × 1) independent Gaussian noise vectors such that p(Z1 ) = CN(0, α2 IM ) and p(Z2 ) = † CN(0, σ 2 IM − α2 DS DS ). Since we are using power † normalised QPSK symbols, DS DS = IM and the covariance matrix 2 of Z2 reduces to 2 = (σ 2 − α2 )IM . We define the positive design parameter ρ α2 /σ 2 , (0 ρ < 1) as the proportion of noise that is assigned to Z2 . Obviously, the value assigned to parameter ρ is important and will have to be tuned. Note that setting ρ = 0 leads to Z1 = 0 and is equivalent to working with the initial model (3). However, for 0 < ρ < 1, the above noise decomposition allows the ˜ and a two-stage introduction of a hidden channel vector H observation model defined as
˜ = Tg + Z1 H ˜ + Z2 . Y = DS H
3. THE EM-MAP CHANNEL ESTIMATOR Starting from Equation (5), one has to estimate g in the ˜ are unknown. Obviously, this MAP sense while S and H MAP estimation has not a closed-form solution. In such situations, the EM algorithm can be used to iteratively solve this problem. Here, the EM algorithm is used for MAP estimation of g (we call it the EM-MAP algorithm) by considering the observation model (5). ˜ be the complete-data set in the EM Let X = {Y, S, H} algorithm terminology. Note that the observation set Y determines only a subset of the space ᐄ of which X is an outcome. The EM-MAP algorithm searches for g that maximises log p(g|X). After initialisation by a short pilot sequence at the beginning of the frame, the EM algorithm alternates between the following two steps (until some stopping criterion is fulfilled) to produce a sequence of estimates {g(t) , t = 0, 1, . . . , tmax }. Expectation Step (E-step): The conditional expectation of the complete log-likelihood given the observed vector and the current estimate g(t) is computed. This quantity is called to a left-multiplication by the matrix T† .
Copyright © 2008 John Wiley & Sons, Ltd.
˜ Y, g(t) Q g, g(t) = ES,H˜ log p(Y, S, H|g)
(6)
Maximisation Step (M-step): The vector of estimated parameters is updated according to g(t+1) = arg max Q g, g(t) + log π(g)
(7)
g
where π(g) a prior distribution for the wavelet coefficients which ensures a certain percentage of coefficients to be set to zero (see Subsection 3.2). When applied to Equation (5), each step can be calculated as follows:
(5)
˜ which provides This model introduces a hidden vector H ‡ us with a direct relation between the true and the estimated wavelet coefficients corrupted by an AWGN, even if the two-stage observation model (5) is equivalent to Equation (3).
‡ Up
the auxiliary or the Q-function and defined as
3.1. E-step: computation of the Q-function The complete likelihood is ˜ ˜ g)p(S|H, ˜ g)p(H|g). ˜ p(Y, S, H|g) = p(Y|S, H, ˜ Y and g are According to Equation (5), given H, independent. Furthermore, S which results from coding and ˜ and g. Since interleaving of bit sequence is independent of H Z1 is a complex white Gaussian noise, the complete loglikelihood can be simplified to
˜ ˜ ˜ log p(Y, S, H|g) = log p(Y|S, H)p(S)p( H|g) ˜ = log p(H|g) + cst.1 ˜ g† T† Tg − 2Re g† T† H + cst.2 =− α2 where cst.1 and cst.2 are constant terms that do not depend on g. According to Equation (6), we have
Q g, g
(t)
˜ g† T† Tg − 2Re g† T† H = ES,H˜ − α2 (t) + cst.2Y, g =−
˜ (t) − Tg2 H + cst.3 α2
(8)
where the latter equation is obtained by adding and ˜ (t) E ˜ [H|Y, ˜ ˜ (t) 2 , with H subtracting H g(t) ] and S,H cst.3 representing a constant term. Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
765
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
From Equation (8), it is obvious that the E-step involves ˜ (t) , as follows: the computation of H ˜ (t) = H
S∈Ꮿ
(t)
µH˜ (S)p(S|Y, g(t) )
(9)
˜ EH˜ H|Y, g(t) =
˜ Ᏼ H∈
˜ H|Y, ˜ ˜ Hp( g(t) , S)dH.(10)
† (t) µH˜ = Tg(t) + ρDS Y − DS Tg(t) .
(11)
By introducing Equation (11) in Equation (9) we obtain †
˜ (t) = (1 − ρ)Tg(t) + ρDS Y H
(12)
where DS = s∈Ꮿ DS p(S|Y, g(t) ). Here, we assume that a part of the receiver called soft-input soft-output (SISO) decoder, is able to provide the vector of probabilities p(S|Y, g(t) ) required in Equation (12) (see Section 4 for more details). The E-step is then completed by inserting ˜ (t) into Q(g, g(t) ) of Equation (8). H 3.2. M-step: wavelet based MAP estimation In this step the estimate of the parameter g is updated as given in Equation (7) where Q(g, g(t) ) is given by Equation (8). We have ˜ (t) − Tg2 H = arg max − + log π(g) . (13) α2 g
g
It can be easily shown that the expression (13) is equivalent to
g(t+1)
(15)
Equation (9) results from the independence between S ˜ belonging respectively to the sets Ꮿ and Ᏼ. and H ˜ (t) in Equation (9), we first have to In order to evaluate H (t) ˜ in Equation (10). evaluate the conditional mean µH˜ (S) of H (t) ˜ ˜ S)p(H|g ˜ (t) ). To this end we write p(H|Y, g , S) ∝ p(Y|H, (t) ˜ ˜ Since both p(Y|H, S) and p(H|g ) are Gaussian densities, it is well known that their product remains Gaussian. We have to evaluate the mean of this Gaussian density. After some algebra we get [16]
(t+1)
˜ (t) = (1 − ρ)g(t) + ρ(DS T)† Y. g˜ (t) = T† H
In fact, g(t+1) in Equation (14) is not more than the MAP estimate of g from the observation model:
where (t) µH˜ (S)
where
2
˜g(t) − g = arg max − + log π(g) α2 g
Copyright © 2008 John Wiley & Sons, Ltd.
(14)
g˜ (t) = g + Z1
(16)
where Z1 ∼ CN(0, α2 IL ). Note that Equation (16) is very important since it shows that the initial observation model is equivalent to an observation model which involves a direct relation between the unknown wavelet coefficient g and g˜ (t) , and this direct relation is corrupted by an AWGN. This is the reason of using the two level observation model in Equation (5). As clear from Equation (16), our channel estimation problem can now be viewed as a standard wavelet domain denoising problem. In what follows, we derive the update formula of our wavelet domain channel estimator. From the Bayes theorem, the posterior distribution of g is given by p gg˜ (t) ∝ p g˜ (t) g π (g)
(17)
where from Equation (16), p(˜g(t) |g) is the Gaussian likelihood, g˜ (t) ∼ CN(g, α2 IL ). In this approach, π (g) is a prior distribution, chosen for the wavelet coefficients g of the unknown CIR as follows. Choice of the Bernoulli–Gaussian Prior Model: Here, the parsimonious characteristic of wavelet basis is imposed through the following prior model: each wavelet coefficient is assumed to have a probability λ to be zero and a probability 1 − λ to be Gaussian distributed as CNgj (0, τ 2 ). This corresponds to an independent and identically distributed (i.i.d.) Bernoulli–Gaussian [13] prior model for the probability density of the jth wavelet coefficient gj as π(gj ) = λδ(gj ) + (1 − λ)CNgj (0, τ 2 )
(18)
for j = 1, . . . , L. The parameters λ and τ (hyperparameters in the Bayesian wording) are estimated from the observed data. In order to deal with that particular model, we introduce an additional state variable (or indicator) βj ∈ {0, 1} such Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
766
S. M. S. SADOUGH ET AL.
that we can express the above prior model conditionally as
(gj |βj = 0) ∼ δ(gj ) (gj |βj = 1) ∼ CNgj
0, τ 2
with probability λ, withprobability 1 − λ.
(t) p βj = 0|˜gj = λ N 0, α2 /c (t) p βj = 1|˜gj = (1 − λ) N 0, α2 + τ 2 /c
(20)
where the parameter c = λN(0, α2 ) + (1 − λ)N(0, α2 + τ 2 ). From this set of equations, we notice that the indicator variable βj allows us to discriminate between the noise coefficients (for βj = 0) and the effective channel wavelet coefficients (for βj = 1), eventually corrupted by noise. The indicator variables βj are estimated in the MAP sense by (t+1)
βj
=
0, 1,
elsewhere.
=
α2
(t+1)
=0
(t+1)
= 1.
if βj τ2 (t) g˜ , + τ2 j
if βj
(21)
(22)
Updating the Prior Parameters τ and λ: The prior parameters τ and λ stand respectively for the significantwavelet coefficients energy and insignificant coefficient probability. Since τ is a scale parameter, a Chi-square prior is chosen for its inverse while a Dirichlet prior is chosen for the probability λ. We point out that these two priors are widely used in Bayesian estimation theory where scale and discrete probabilities are in concern. The main advantage for such priors is that they generally result in relatively simple update algorithms. Here, the update rules for these two parameters are MAP based rules derived from assigning conjugate priors to these parameters [17]. After Copyright © 2008 John Wiley & Sons, Ltd.
˜ = Card{j | βj = 0} and η = where L
(23)
βj =1
(t+1) 2 g . j
Reduction of the Number of Estimated Parameters: The thresholding procedure derived in this section, provides an automatic framework for reducing the number of estimated coefficients. This can be achieved by discarding at each iteration, the elements of g(t+1) that are replaced by zero in Equation (22). The underlying assumption is as follows: whenever the estimator assimilates an unknown wavelet coefficient to noise (i.e. replaces it by zero), this coefficient will always be considered as noise, hence it will not be estimated in future iterations. We verified by simulations that incorporating this scheme into the EM algorithm reduces the number of estimated parameters without any significant performance degradation. This operation is shown in Figure 3 and can be modelled as (t+1)
Therefore, the MAP estimates of the channel wavelet coefficients are obtained by a simple denoising/thresholding rule as
(t+1) gj
˜ τˆ 2 = η/(L − L)
gtr
(t) if p βj = 0g˜ j 0.5
0,
˜ − 1/2)/(L − 1), λˆ = (L
(19)
This prior model, conditionally on the state variable βj , leads to a Gaussian posterior for gj which makes the (t) estimation explicit; from the direct observation model g˜ j = (Equation (16)), we can express these posterior gj + Z1,j probabilities of βj as [16]
some algebra we obtain [16]
= g(t+1) ,
Ttr = T
(24) (t+1)
where the truncation operator (·) gathers in gtr the components of g(t+1) that must be kept and the operator (·) constructs Ttr from T by keeping the columns corresponding to the kept indexes. During the first iteration (t = 0), the algorithm does not perform any truncation and the EM algorithm estimates all the coefficients. However, after each M-step, the number of unknown parameters to be estimated in the next iteration is reduced according (t+1) and Ttr in the updating to Equation (24) by using gtr formula of the E-step (12).
Figure 3. EM-MAP channel estimation combined with the decoding process. Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
767
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
4. DECODING METHOD AND IMPLEMENTATION ISSUES The block diagram of the receiver is shown in Figure 3. Besides the channel estimation part, the rest of the receiver principally consists of the combination of two sub-blocks that exchange soft informations with each other. The first sub-block, refered to as soft demapper (also called detector), produces bit metrics (probabilities) from the input symbols and the second one is a SISO decoder. Here, SISO decoding is performed using the well known forward-backward algorithm [18]. Let ck,i (i = 1, . . . , B) be the ith§ coded and interleaved bit corresponding to the kth constellation symbol Sk and let Yk be the corresponding received symbol (k = 1, . . . , MMsym ). Actually, according to Equation (9), we have to use the information on transmitted symbols, obtained from the SISO decoder through the probability ˆ (t) ), to update the channel estimate at each P(Sk |Yk , H k iteration. Furthermore, the soft demapper requires an estimate of the channel in order to provide the probability of encoded bits (see Figure 3). Hence, the proposed semi-blind channel estimation algorithm is naturally combined with the ˆ (t) ) of process of data decoding. The probability P(Sk |Yk , H k the unknown symbols Sk is calculated using the a posteriori probabilities provided by the SISO decoder at the end of the tth iteration as B ˆ (t) = P Sk |Yk , H Pdec (ck,j ) k
(25)
j=1 j =i
where Pdec (ck,j ) are prior probabilities coming from the SISO decoder. At the first iteration, where no a priori information is available on bits ck,j , the probabilities Pdec (ck,j ) are set to 0.5. The main steps of the EM-MAP estimator are summarised as follows:
r Initialisation (t = 0)
r § In
– Set all probabilities of coded bits Pdec (ck,i ) to 0.5 and ˆ (0) ) for all k according to Equation derive P(Sk |Yk , H k (25) and then derive p(S|Y, H(0) ) = p(S|Y, g(0) ). – Initialise the unknown vector g by g(0) obtained from pilot symbols. for {t = 1, . . . , tmax }
this paper, B = 2 since we are using a QPSK modulation.
Copyright © 2008 John Wiley & Sons, Ltd.
r
– Use the previous estimate g(t−1) and p(S|Y, g(t−1) ) to calculate g˜ (t−1) according to Equation (15). – Use g˜ (t−1) to obtain the updated channel parameters g(t) by using Equation (22). – Discard the wavelet coefficients that are replaced (t) by zero in g(t) by evaluating gtr and Ttr from Equation (24). (t) if {t = tmax }: Use the current estimate gtr to update the probability of encoded bits Pdec (ck,i ) and derive ˆ (t) ) from Equation (25) and then derive P(Sk |Yk , H k (t) p(S|Y, H ) = p(S|Y, g(t) ). else: Decode the information data by thresholding the uncoded bit probabilities. end
5. SIMULATION RESULTS In this section we present a comparative performance study of the proposed EM-MAP algorithm. The performance comparison is made in terms of mean square error (MSE) for channel estimation and bit error rate (BER) for the combination of channel estimation with the decoding process. The binary information data are encoded by a nonrecursive non-systematic convolutional encoder with rate R = 1/2 and constraint length 3, defined in octal form by (5, 7)8 . Throughout the simulations, each frame is composed of Nsym = 9 OFDM symbols with N = 128 subcarriers each. Channel coefficients are kept constant during each fading block and changed to new independent realisations (measures) from one frame to the next. We devote one OFDM pilot symbol for each subband in order to initialise the EM algorithm. Data and pilot symbols belong to the QPSK constellation with gray labelling. The interleaver is pseudo-random, operating over the entire frame of size NI = Nsym NB bits (excluding pilots, obviously). Among different wavelet families, ‘symmetric’ wavelet basis functions [19] providing a sparser representation [14] have been considered. Unless otherwise mentioned, the BER and MSE curves correspond to the fourth iteration of the algorithm. The parameter ρ affecting the variance of the noise in Equation (16) is set to 0.4 which minimises the MSE between the true and estimated wavelet coefficients. Different propagation environments are considered for performance comparison. In all environments, the Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
768
S. M. S. SADOUGH ET AL.
bandwidth of the channel frequency response is 1584 MHz with 384 coefficients and the CIR has a total number of 96 taps. For each transmitted frame, a different measurement/realisation of the channel has been used. First, a sparse channel model where only 20 wavelet coefficients out of total 96 have non-zero values, is considered. This propagation environment provides the best adequacy between the sparseness prior assumption and the actual propagation environment. The second channel is a line of sight (LOS) scenario issued from realistic UWB indoor channel measurements [20] where the receive and transmit antennas are located in a corridor separated by 9 meters. We also consider the non-LOS UWB channel model CM2 specified by the IEEE 802.15.3a channel modelling subcommittee report [21]. The two latter channels characterise a rich scattering environment where the sparseness assumption in the wavelet domain is not necessarily satisfied. For the sake of comparison, we consider two pilot-only based approach using ML and minimum mean square error (MMSE) channel estimation, refered to as pilot-ML and pilot-MMSE. We also compare the proposed algorithm with two semi-blind channel estimation based on the EM algorithm. The first approach called EM-Freq, consists of estimating the 384 frequency coefficients of the channel over all three subbands, using the model (3), similar to the approach adopted in Reference [11]. The second semi-blind method, called EM-Wav, is a wavelet domain EM based channel estimator where the prior model is set to have a uniform distribution. First, consider the case of the sparse channel. Figure 4 depicts the MSE as a function of Eb /N0 . It can be observed that, although the pilot-MMSE approach improves the estimation accuracy for low SNR values, the performance of pilot based channel estimation methods are very far from the family of semi-blind methods. Comparing the wavelet domain semi-blind approach (EM-Wav) and the frequency domain approach (EM-freq), shows that significant gain is achieved by the former method. As shown, the best performance is achieved by the EM-MAP method. We see that by using EM-MAP, a gain of almost 4 dB in SNR is achieved at a MSE of 2 × 10−3 , as compared to the EM-Wav method. This clearly shows the adequacy of the EM-MAP method for the case where the unknown channel has few non-zero wavelet coefficients, which is in perfect agreement with the prior model. Figure 5 shows the BER results along with the BER for the case of perfect channel state information (CSI). It can be seen that at a BER of 10−3 , the pilot-ML and the EM-Freq approaches are respectively 3.9 and 2 dB of SNR far from Copyright © 2008 John Wiley & Sons, Ltd.
Figure 4. Mean square error between the true and estimated coefficients for the sparse channel model.
Figure 5. BER performance of different channel estimation methods over the sparse channel model.
the BER obtained with the perfect channel. Furthermore, the performance of the Pilot-MMSE approach is not shown since it was very close to that of Pilot-ML. Also, we observe that wavelet based semi-blind methods perform closely to the perfect CSI case. For example, at a BER of 10−4 , the EM-MAP and EM-Wav method have respectively about Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
0.2 and 0.5 dB of SNR degradation from the performance obtained with perfect CSI. We now evaluate the performance of the EM-MAP algorithm for the case of the Corridor and CM2 channels. Figures 6 and 7 show that wavelet based methods again outperform pilot based and EM-Freq methods in terms
769
Figure 6. Mean square error between the true and estimated coefficients over the Corridor channel.
of MSE and BER. Furthermore, we observe that the performance of the EM-MAP method is now comparable to that of the EM-Wav method. This can be explained by noting that when the channel is not very sparse (i.e. for low λ values in (18)), the algorithm assigns small values to λ. This leads to a Gaussian prior model with a large variance, which can be approximated to a uniform prior. As a result, the prior becomes ‘less informative’ and the EM-MAP performs close to the EM-Wav, as shown in Figures 6 and 7. Thus, the EM-MAP algorithm proposed here is globally able to adapt its prior model parameters to any propagation environment. We now compare the above semi-blind algorithms with respect to the average number of estimated parameters at each iteration of the EM algorithm. This is shown in Figure 8 for different channel scenarios. First, recall that the EMFreq and EM-Wav methods have to continuously estimate 384 and 96 coefficients at each iteration, respectively. As explained in Subsection 3.2, by discarding the coefficients that are replaced by zero, the EM-MAP approach is able to reduce the number of estimated parameters, specially for sparse channels. This can be seen in Figure 8 for the sparse channel where the number of estimated parameters is reduced from 96 down to 20 parameters after the second iteration. Furthermore, for the non-sparse Corridor and CM2 channels, we observe that the EM-MAP method is to be preferred to the EM-Wav, due to its lower computational load, since it estimates about 40 coefficients after the second iteration. Although we observed that in this case
Figure 7. BER performance of different channel estimation methods over the CM2 channel.
Figure 8. Reduction of the number of estimated parameters through iterations, Eb /N0 = 8 dB.
Copyright © 2008 John Wiley & Sons, Ltd.
Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
770
S. M. S. SADOUGH ET AL.
a larger bandwidth is likely to separate the paths and to provide a sparser wavelet expansion. It was observed that when the channel has a very sparse wavelet expansion, the prior model parameters which are estimated from the observed data, carry this sparseness information to the EM-MAP algorithm. Moreover, we showed that in this case, the EM-MAP method provides significant reduction in the number of estimated parameters and outperforms all considered pilot based and semi-blind methods. Under non-sparse channels, although both EM-MAP and EM-Wav methods perform closely, EM-MAP takes the advantage over the EM-Freq and EM-Wav schemes due to its lower computational complexity.
REFERENCES
Figure 9. Convergence of the MSE with respect to the number of iterations over the CM2 channel.
these two methods exhibit close performance, the EM-MAP algorithm brings a reduction of about 60% on the number of estimated parameters when compared to the EM-Wav approach. Finally, we analyse the number of iterations that the EMMAP algorithm requires for convergence. In Figure 9, the MSE performance of the EM-MAP algorithm is presented as a function of the number of iterations. It is obvious from these curves that the MSE performance of the proposed algorithm converges within two to four iterations, depending on the average SNR. This justifies our choice of four iterations in the results presented above.
6. CONCLUSION We proposed an EM-based channel estimation and data detection algorithm that integrates the advantages of wavelet based parameter estimation. By expressing the unknown UWB channel in terms of its discrete wavelet coefficients and choosing a proper prior distribution for them, we were able to capture the sparseness property of UWB channels in the wavelet domain. This led to a MAP estimator equivalent to a hard thresholding procedure at each iteration of the EM algorithm, which is used to reduce the number of estimated coefficients. Note that the sparseness condition is almost built in by considering all subbands of the multiband transmission at the same time: Copyright © 2008 John Wiley & Sons, Ltd.
1. FCC. First report and order, revision of part 15 of the commission’s rules regarding ultra-wideband transmission systems. Technical report, February 2004. 2. Roy S, Foerster JR, Somayazulu VS, Leeper DG. Ultrawideband radio design: the promise of high-speed, short range wireless connectivity. Proceedings of the IEEE 2004; 92:295–311. 3. Batra A, Balakrishnan J, Aiello GR, Foerster JR, Dabak A. Design of multiband OFDM system for realistic UWB channel environments. IEEE Transactions Microwave Theory and Techniques 2004; 52:2123–2138. 4. Prasad R. OFDM for Wireless Communications Systems. Artech House Publishers, Norwood, MA, USA, 2004. 5. Caire G, Taricco G, Biglieri E. Bit-interleaved coded modulation. IEEE Transactions and Information Theory 1998; 44:927–945. 6. van de Beek JJ, Edfors O, Sandell M. On channel estimation in OFDM systems. In Proceedings of Vehicular Technology Conference, September 1995; 815–819. 7. McLachlan GJ, Krishnan T. The EM Algorithm and Extensions. Wiley Series in Probability and Statistics, Wiley, Chichester, UK, 1997. 8. Zhuang X, Vook FW. Iterative channel estimation and decoding for a turbo-coded OFDM system via the EM algorithm. In Proceedings of IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), May 2002. 9. Jaffrot E, Siala M. Turbo channel estimation for OFDM systems on highly time and frequency selective channels. In Proceedings of IEEE International Conference on Acoustics, Speech, Signal Processing (ICASSP), June 2000. 10. Lu B, Wang X. Bayesian blind turbo receiver for coded OFDM systems with frequency offset and frequency-selective fading. IEEE Journal on Selected Areas in Communications. 2001; 19:2516–2527. 11. Ma X, Kobayashi H, Schwartz SC. EM-based channel estimation algorithms for OFDM. EURASIP Journal on Applied Signal Processing 2004; 10:1460–1477. 12. Lu B, Wang X, Li YG. Iterative receivers for space-time block coded OFDM systems in dispersive fading channels. IEEE Transactions on Wireless Communications 2002; 1:213–225. 13. Johnstone IM, Silverman BW. Neddles and straw in haystacks: empirical bayes estimates of possibly sparse sequences. Annals of Statistics 2004; 32:1594–1649. 14. Sadough S, Jaffrot E. A wavelet packet based model for an ultrawideband propagation channel. In Proceedings of ECPS, 15–18 March 2005.
Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
15. Nowak RD, Figueiredo MAT. Fast wavelet-based image deconvolution using the EM algorithm. Asilomar Conference on Signals, Systems and Computers, 2001; 371–375. 16. Sadough S. Ultra Wideband OFDM Systems: Channel Estimation and Improved Detection Accounting for Estimation Innacuracies. PhD thesis, Universit´e Paris-Sud 11, January 2008. 17. Bernardo JN, Smith AFM. Bayesian Theory. Wiley and Sons: Chichester, England, 1994.
Copyright © 2008 John Wiley & Sons, Ltd.
771
18. Bahl L, Cocke J, Jelinek F, Raviv J. Optimal decoding of linear codes for minimizing symbol error rate. IEEE Transactions on Information Theory 1974; 20(2):284–287. 19. Mallat S. A Wavelet Tour of Signal Processing. Academic Press, 1999. 20. Bories S, Sibille A, Roblin C. UWB indoor channel measurement. In Proceedings of IWAT, March 2005; 466–469. 21. Foerster J. Channel modeling sub-committee report final. Technical report, IEEE802.15-02/490, 2003.
Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
Ultra wideband OFDM channel estimation through a wavelet based EM-MAP algorithm† S. M. S. Sadough1,2∗ , M. Ichir2 , P. Duhamel2 and E. Jaffrot3 1 UEI, ENSTA, 32 boulevard Victor, 75015 Paris, France des Signaux et Syst`emes, CNRS-Sup´elec, Plateau de Moulon, 91190 Gif–sur–Yvette, France 3 Universidad Nacional de San Martin, B1650ANQ San Martin, Pcia. de Buenos Aires, Argentina
2 Laboratoire
SUMMARY Ultra wideband (UWB) communications involve very sparse channels, since the bandwidth increase results in a better time resolution. This property is used here to propose an efficient algorithm jointly estimating the channel and the transmitted symbols. More precisely, this paper introduces an expectation-maximisation (EM) algorithm within a wavelet domain Bayesian framework for semi-blind channel estimation of multiband orthogonal frequency-division multiplexing (MB-OFDM) based UWB communications. A prior distribution is chosen for the wavelet coefficients of the unknown channel impulse response (CIR) in order to model a sparseness property of the wavelet representation. This prior yields, in maximum a posteriori (MAP) estimation, a thresholding rule within the EM algorithm. We particularly focus on reducing the number of estimated parameters by iteratively discarding ‘insignificant’ wavelet coefficients from the estimation process. Simulation results using UWB channels issued from both models and measurements show that under sparsity conditions, the proposed algorithm outperforms pilot based channel estimation in terms of mean square error (MSE) and bit error rate (BER). Moreover, the estimation accuracy is improved, while the computational complexity is reduced, when compared to traditional semi-blind methods. Copyright © 2008 John Wiley & Sons, Ltd.
1. INTRODUCTION An ultra wideband (UWB) radio signal is defined as any signal whose bandwidth is larger than 20% of its centre frequency or greater than 500 MHz [1]. In recent years, UWB system design has experienced a shift from the traditional ‘single-band’ radio that occupies the whole 7.5 GHz allocated spectrum to a ‘multiband’ design approach [2]. In these multiband UWB systems, the available spectrum is divided into several subbands, each one occupying approximately 500 MHz. Multiband orthogonal frequency-division multiplexing (MB-OFDM) [3] is an appropriate candidate for multiband UWB which enables high data rate UWB transmission to inherit all the strength of OFDM that is widely used in numerous communication systems. This approach uses a
conventional coded OFDM system [4] together with bit interleaved coded modulation (BICM) [5] and frequency hopping over different subbands to improve the diversity and to enable multiple access. It is well known that an efficient detection of the transmitted symbols requires a good knowledge of the channel realisations at the receiver. The required channel estimation is often based on the transmission of known training symbols (also called pilots) [6]. However, the length of the required training sequence may become large when the channel impulse response (CIR) is highly resolved, such as in UWB OFDM systems. Besides, obtaining an accurate channel estimate in highly mobile environments only through the use of pilots, would require inserting multiple training symbols per frame, which can result in a waste of bandwidth and power.
* Correspondence to: Seyed Mohammad Sajad Sadough, UEI, ENSTA, 32 boulevard Victor, 75015 Paris, France. E-mail: [email protected]. † A previous version of this paper was presented in the 13th European Wireless Conference (EW 2007), Paris, France.
Copyright © 2008 John Wiley & Sons, Ltd.
Accepted 16 June 2008
762
S. M. S. SADOUGH ET AL.
In order to reduce the number of pilot symbols, the information carried by the observations corresponding to the data symbols can also be exploited for channel estimation. An attractive technique in this area is that based on the so-called expectation-maximisation (EM) algorithm [7] which provides a framework to iteratively estimate the unknown channel thanks to the a posteriori probabilities on the unknown data. Moreover, by using the EM algorithm, one can naturally combine the iterative process of channel estimation with the decoding operation of encoded data. Several papers have recently reported promising results on the combination of EM-based channel estimation and data decoding process. For instance, iterative or ‘turbo’ processing that includes the channel estimation into the iterative process of decoding turbo like codes is addressed in Reference [8–10]. In particular, Reference [9] addresses a turbo channel estimation based on a Karhunen–Loeve expansion of the unknown channel for reducing the number of estimated parameters. However, the proposed scheme requires the channel covariance matrix as an a priori information which can not be available in realistic situations. In Reference [11], the authors present several EM based algorithms for time or frequency domain channel estimation of an OFDM system that is subject to slow timevarying frequency-selective fading. In Reference [12], Lu et al. proposed an iterative receiver for space-time blockcoded OFDM systems based on a reduced complexity EM algorithm. Though iterative joint channel estimation and data detection scheme outperforms receivers using a pilot-only based channel estimation, it has a higher complexity that may be a critical concern for its practical implementations. This complexity is mainly driven by the number of estimated parameters for channel updating and the decoding algorithm within each iteration. In this work, we consider a semi-blind joint channel estimation and data detection scheme for MB-OFDM systems based on the EM algorithm, with the objective of minimising the number of estimated parameters and enhancing the estimation accuracy. This is achieved by expressing the unknown UWB CIR in terms of its discrete wavelet expansion, which has been shown to provide a parsimonious representation [13, 14]. Our wavelet domain framework enables us to choose a particular a priori distribution for the channel wavelet coefficients that renders the maximum a posteriori (MAP) channel estimation equivalent to a hard thresholding rule at each iteration of the EM algorithm. The latter is then exploited to reduce the estimator’s computational complexity by discarding ‘insignificant’ wavelet coefficients from the Copyright © 2008 John Wiley & Sons, Ltd.
estimation process. Moreover, the adopted wavelet domain prior model adapts itself to the actual channel scenario since its parameters are learnt from the observed data during the channel estimation process. More precisely, the prior model is not restricted to a specific propagation environment. We show that when the channel has a sparse wavelet expansion, the proposed algorithm outperforms the classical EM estimators and provides a significant reduction in the number of estimated parameters. Hence, the proposed receiver has not to modify the prior information for different channel scenarios. This paper is organised as follows: Section 2 introduces the MB-OFDM system and its equivalent model for wavelet domain channel estimation. Section 3, first presents the derivation of the MAP channel estimator based on the EM algorithm and then shows how the number of estimated parameters can be reduced through the EM iterations. Section 4 presents the combination of the channel estimation with the decoding operation, as well as some implementation issues. Section 5 illustrates, via simulations, the performance of the proposed receiver in different realistic UWB channel environments and Section 6 concludes the paper. Notational conventions are as follows: Dx is a diagonal matrix with diagonal elements x = [x1 , . . . , xN ]T , Ex [·] refers to expectation with respect to x, IN denotes an (N × N) identity matrix; CN(m, ) denotes complex Gaussian distribution with mean m and covariance matrix ; Card{·} denotes set cardinality; (·)T and (·)H denote matrix or vector transpose and Hermitian transpose, respectively.
2. SYSTEM MODEL AND WAVELET DOMAIN PROBLEM FORMULATION 2.1. MB-OFDM transmission We consider the MB-OFDM transmission proposed in Reference [3]. This scheme divides the spectrum between 3.1 and 10.6 GHz into several non-overlapping subbands each one occupying 528 MHz of bandwidth [3]. Information is transmitted using OFDM modulation over one of the subbands in a particular time slot. As shown in Figure 1, the transmitter architecture for the MB-OFDM system is very similar to that of a conventional wireless OFDM system. The main difference is that MB-OFDM uses a timefrequency code (TFC) to select the centre frequency of different subbands. This code is used to provide frequency diversity as well as to distinguish between multiple users (see Figure 2). Here, we consider MB-OFDM in its basic Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
763
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
frequency domain as† Ym = DSm Hm + Zm
m = 1, . . . , Msym
(2)
where Ym = [y1,n , y2,n+1 , y3,n+2 ]T , Sm = [s1,n , s2,n+1 , ¯ 1,n , H ¯ 2,n+1 , H ¯ 3,n+2 ]T and Zm = [z1,n , s3,n+2 ]T , Hm = [H T z2,n+1 , z3,n+2 ] are (M × 1) vectors, with M = 3N and Msym = Nsym /3. In the remainder, unless otherwise mentioned, we will not write the time index m for notational convenience. Figure 1. TX architecture of the multiband OFDM system.
2.2. A proper model for wavelet domain channel estimation
Figure 2. Example of time-frequency coding for the multiband OFDM system: TFC = {1, 3, 2, 1, 3, 2, . . .}.
mode, i.e. employing the first three subbands with N data subcarriers over each subband. Let us consider a single-user MB-OFDM transmission. At the receiver, assuming a cyclic prefix (CP) longer than the channel maximum delay spread and perfect carrier synchronisation, OFDM converts a frequency-selective channel into N parallel flat fading subchannels [4]. Under these conditions, the transmission of the nth OFDM symbol (inside a frame of size Nsym ) over the ith subband can be written as [4] ¯ i,n + zi,n yi,n = Dsi,n H
(1)
for i ∈ {1, 2, 3} and n = 1, . . . , Nsym ; (1 × N) vectors ¯ i,n denote respectively the received and yi,n , si,n and H transmitted symbols, and the channel frequency response; the noise vector zi,n is assumed to be a zero-mean circularly symmetric complex Gaussian random vector with distribution CN(0, σ 2 IN ). In what follows, we group the data and the observation corresponding to three subbands of the MB-OFDM system in a single vector. This operation can be written in the Copyright © 2008 John Wiley & Sons, Ltd.
In order to take advantage of the wavelet based estimation, the CIR is expressed in terms of its orthogonal discrete wavelet coefficients. Let FM,L be the truncated fast Fourier transform (FFT) matrix constructed from the (M × M) FFT matrix by keeping the first L columns where L is the length of the CIR over a group of three subbands. We define W as the (L × L) orthogonal discrete wavelet transform (ODWT) matrix. The unknown channel frequency response can be expressed as H = FM,L WH g, where g is the (L × 1) vector of the CIR wavelet coefficients. The observation model (2) is rewritten as Y = DS Tg + Z
(3)
where T = FM,L WH . In this model, although the channel is practically used (by the transmitter) by slices of 528 MHz bandwidth (corresponding to a single subband), on the receiver side, three received OFDM symbols are grouped for estimating the wavelet coefficients of the CIR, taken over all three subbands (1.584 GHz bandwidth). This is motivated by the fact that estimating the channel over a wider bandwidth leads to a sparser representation in the wavelet domain. Most of wavelet based estimation algorithms, rely on an observation model in which the unknown wavelet coefficients are corrupted by an AWGN [13]. Since the model (3) does not provide such a framework, our first step consists in enforcing this property. In order to do so, the AWGN in Equation (3) is split into two independent Gaussian terms as suggested in Reference [15] Z = DS Z1 + Z2
(4)
† For
the sake of notational brevity in Equation (2), we have assumed that the TFC is equal to {1, 2, 3, . . .}.
Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
764
S. M. S. SADOUGH ET AL.
where Z1 and Z2 are (M × 1) independent Gaussian noise vectors such that p(Z1 ) = CN(0, α2 IM ) and p(Z2 ) = † CN(0, σ 2 IM − α2 DS DS ). Since we are using power † normalised QPSK symbols, DS DS = IM and the covariance matrix 2 of Z2 reduces to 2 = (σ 2 − α2 )IM . We define the positive design parameter ρ α2 /σ 2 , (0 ρ < 1) as the proportion of noise that is assigned to Z2 . Obviously, the value assigned to parameter ρ is important and will have to be tuned. Note that setting ρ = 0 leads to Z1 = 0 and is equivalent to working with the initial model (3). However, for 0 < ρ < 1, the above noise decomposition allows the ˜ and a two-stage introduction of a hidden channel vector H observation model defined as
˜ = Tg + Z1 H ˜ + Z2 . Y = DS H
3. THE EM-MAP CHANNEL ESTIMATOR Starting from Equation (5), one has to estimate g in the ˜ are unknown. Obviously, this MAP sense while S and H MAP estimation has not a closed-form solution. In such situations, the EM algorithm can be used to iteratively solve this problem. Here, the EM algorithm is used for MAP estimation of g (we call it the EM-MAP algorithm) by considering the observation model (5). ˜ be the complete-data set in the EM Let X = {Y, S, H} algorithm terminology. Note that the observation set Y determines only a subset of the space ᐄ of which X is an outcome. The EM-MAP algorithm searches for g that maximises log p(g|X). After initialisation by a short pilot sequence at the beginning of the frame, the EM algorithm alternates between the following two steps (until some stopping criterion is fulfilled) to produce a sequence of estimates {g(t) , t = 0, 1, . . . , tmax }. Expectation Step (E-step): The conditional expectation of the complete log-likelihood given the observed vector and the current estimate g(t) is computed. This quantity is called to a left-multiplication by the matrix T† .
Copyright © 2008 John Wiley & Sons, Ltd.
˜ Y, g(t) Q g, g(t) = ES,H˜ log p(Y, S, H|g)
(6)
Maximisation Step (M-step): The vector of estimated parameters is updated according to g(t+1) = arg max Q g, g(t) + log π(g)
(7)
g
where π(g) a prior distribution for the wavelet coefficients which ensures a certain percentage of coefficients to be set to zero (see Subsection 3.2). When applied to Equation (5), each step can be calculated as follows:
(5)
˜ which provides This model introduces a hidden vector H ‡ us with a direct relation between the true and the estimated wavelet coefficients corrupted by an AWGN, even if the two-stage observation model (5) is equivalent to Equation (3).
‡ Up
the auxiliary or the Q-function and defined as
3.1. E-step: computation of the Q-function The complete likelihood is ˜ ˜ g)p(S|H, ˜ g)p(H|g). ˜ p(Y, S, H|g) = p(Y|S, H, ˜ Y and g are According to Equation (5), given H, independent. Furthermore, S which results from coding and ˜ and g. Since interleaving of bit sequence is independent of H Z1 is a complex white Gaussian noise, the complete loglikelihood can be simplified to
˜ ˜ ˜ log p(Y, S, H|g) = log p(Y|S, H)p(S)p( H|g) ˜ = log p(H|g) + cst.1 ˜ g† T† Tg − 2Re g† T† H + cst.2 =− α2 where cst.1 and cst.2 are constant terms that do not depend on g. According to Equation (6), we have
Q g, g
(t)
˜ g† T† Tg − 2Re g† T† H = ES,H˜ − α2 (t) + cst.2Y, g =−
˜ (t) − Tg2 H + cst.3 α2
(8)
where the latter equation is obtained by adding and ˜ (t) E ˜ [H|Y, ˜ ˜ (t) 2 , with H subtracting H g(t) ] and S,H cst.3 representing a constant term. Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
765
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
From Equation (8), it is obvious that the E-step involves ˜ (t) , as follows: the computation of H ˜ (t) = H
S∈Ꮿ
(t)
µH˜ (S)p(S|Y, g(t) )
(9)
˜ EH˜ H|Y, g(t) =
˜ Ᏼ H∈
˜ H|Y, ˜ ˜ Hp( g(t) , S)dH.(10)
† (t) µH˜ = Tg(t) + ρDS Y − DS Tg(t) .
(11)
By introducing Equation (11) in Equation (9) we obtain †
˜ (t) = (1 − ρ)Tg(t) + ρDS Y H
(12)
where DS = s∈Ꮿ DS p(S|Y, g(t) ). Here, we assume that a part of the receiver called soft-input soft-output (SISO) decoder, is able to provide the vector of probabilities p(S|Y, g(t) ) required in Equation (12) (see Section 4 for more details). The E-step is then completed by inserting ˜ (t) into Q(g, g(t) ) of Equation (8). H 3.2. M-step: wavelet based MAP estimation In this step the estimate of the parameter g is updated as given in Equation (7) where Q(g, g(t) ) is given by Equation (8). We have ˜ (t) − Tg2 H = arg max − + log π(g) . (13) α2 g
g
It can be easily shown that the expression (13) is equivalent to
g(t+1)
(15)
Equation (9) results from the independence between S ˜ belonging respectively to the sets Ꮿ and Ᏼ. and H ˜ (t) in Equation (9), we first have to In order to evaluate H (t) ˜ in Equation (10). evaluate the conditional mean µH˜ (S) of H (t) ˜ ˜ S)p(H|g ˜ (t) ). To this end we write p(H|Y, g , S) ∝ p(Y|H, (t) ˜ ˜ Since both p(Y|H, S) and p(H|g ) are Gaussian densities, it is well known that their product remains Gaussian. We have to evaluate the mean of this Gaussian density. After some algebra we get [16]
(t+1)
˜ (t) = (1 − ρ)g(t) + ρ(DS T)† Y. g˜ (t) = T† H
In fact, g(t+1) in Equation (14) is not more than the MAP estimate of g from the observation model:
where (t) µH˜ (S)
where
2
˜g(t) − g = arg max − + log π(g) α2 g
Copyright © 2008 John Wiley & Sons, Ltd.
(14)
g˜ (t) = g + Z1
(16)
where Z1 ∼ CN(0, α2 IL ). Note that Equation (16) is very important since it shows that the initial observation model is equivalent to an observation model which involves a direct relation between the unknown wavelet coefficient g and g˜ (t) , and this direct relation is corrupted by an AWGN. This is the reason of using the two level observation model in Equation (5). As clear from Equation (16), our channel estimation problem can now be viewed as a standard wavelet domain denoising problem. In what follows, we derive the update formula of our wavelet domain channel estimator. From the Bayes theorem, the posterior distribution of g is given by p gg˜ (t) ∝ p g˜ (t) g π (g)
(17)
where from Equation (16), p(˜g(t) |g) is the Gaussian likelihood, g˜ (t) ∼ CN(g, α2 IL ). In this approach, π (g) is a prior distribution, chosen for the wavelet coefficients g of the unknown CIR as follows. Choice of the Bernoulli–Gaussian Prior Model: Here, the parsimonious characteristic of wavelet basis is imposed through the following prior model: each wavelet coefficient is assumed to have a probability λ to be zero and a probability 1 − λ to be Gaussian distributed as CNgj (0, τ 2 ). This corresponds to an independent and identically distributed (i.i.d.) Bernoulli–Gaussian [13] prior model for the probability density of the jth wavelet coefficient gj as π(gj ) = λδ(gj ) + (1 − λ)CNgj (0, τ 2 )
(18)
for j = 1, . . . , L. The parameters λ and τ (hyperparameters in the Bayesian wording) are estimated from the observed data. In order to deal with that particular model, we introduce an additional state variable (or indicator) βj ∈ {0, 1} such Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
766
S. M. S. SADOUGH ET AL.
that we can express the above prior model conditionally as
(gj |βj = 0) ∼ δ(gj ) (gj |βj = 1) ∼ CNgj
0, τ 2
with probability λ, withprobability 1 − λ.
(t) p βj = 0|˜gj = λ N 0, α2 /c (t) p βj = 1|˜gj = (1 − λ) N 0, α2 + τ 2 /c
(20)
where the parameter c = λN(0, α2 ) + (1 − λ)N(0, α2 + τ 2 ). From this set of equations, we notice that the indicator variable βj allows us to discriminate between the noise coefficients (for βj = 0) and the effective channel wavelet coefficients (for βj = 1), eventually corrupted by noise. The indicator variables βj are estimated in the MAP sense by (t+1)
βj
=
0, 1,
elsewhere.
=
α2
(t+1)
=0
(t+1)
= 1.
if βj τ2 (t) g˜ , + τ2 j
if βj
(21)
(22)
Updating the Prior Parameters τ and λ: The prior parameters τ and λ stand respectively for the significantwavelet coefficients energy and insignificant coefficient probability. Since τ is a scale parameter, a Chi-square prior is chosen for its inverse while a Dirichlet prior is chosen for the probability λ. We point out that these two priors are widely used in Bayesian estimation theory where scale and discrete probabilities are in concern. The main advantage for such priors is that they generally result in relatively simple update algorithms. Here, the update rules for these two parameters are MAP based rules derived from assigning conjugate priors to these parameters [17]. After Copyright © 2008 John Wiley & Sons, Ltd.
˜ = Card{j | βj = 0} and η = where L
(23)
βj =1
(t+1) 2 g . j
Reduction of the Number of Estimated Parameters: The thresholding procedure derived in this section, provides an automatic framework for reducing the number of estimated coefficients. This can be achieved by discarding at each iteration, the elements of g(t+1) that are replaced by zero in Equation (22). The underlying assumption is as follows: whenever the estimator assimilates an unknown wavelet coefficient to noise (i.e. replaces it by zero), this coefficient will always be considered as noise, hence it will not be estimated in future iterations. We verified by simulations that incorporating this scheme into the EM algorithm reduces the number of estimated parameters without any significant performance degradation. This operation is shown in Figure 3 and can be modelled as (t+1)
Therefore, the MAP estimates of the channel wavelet coefficients are obtained by a simple denoising/thresholding rule as
(t+1) gj
˜ τˆ 2 = η/(L − L)
gtr
(t) if p βj = 0g˜ j 0.5
0,
˜ − 1/2)/(L − 1), λˆ = (L
(19)
This prior model, conditionally on the state variable βj , leads to a Gaussian posterior for gj which makes the (t) estimation explicit; from the direct observation model g˜ j = (Equation (16)), we can express these posterior gj + Z1,j probabilities of βj as [16]
some algebra we obtain [16]
= g(t+1) ,
Ttr = T
(24) (t+1)
where the truncation operator (·) gathers in gtr the components of g(t+1) that must be kept and the operator (·) constructs Ttr from T by keeping the columns corresponding to the kept indexes. During the first iteration (t = 0), the algorithm does not perform any truncation and the EM algorithm estimates all the coefficients. However, after each M-step, the number of unknown parameters to be estimated in the next iteration is reduced according (t+1) and Ttr in the updating to Equation (24) by using gtr formula of the E-step (12).
Figure 3. EM-MAP channel estimation combined with the decoding process. Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
767
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
4. DECODING METHOD AND IMPLEMENTATION ISSUES The block diagram of the receiver is shown in Figure 3. Besides the channel estimation part, the rest of the receiver principally consists of the combination of two sub-blocks that exchange soft informations with each other. The first sub-block, refered to as soft demapper (also called detector), produces bit metrics (probabilities) from the input symbols and the second one is a SISO decoder. Here, SISO decoding is performed using the well known forward-backward algorithm [18]. Let ck,i (i = 1, . . . , B) be the ith§ coded and interleaved bit corresponding to the kth constellation symbol Sk and let Yk be the corresponding received symbol (k = 1, . . . , MMsym ). Actually, according to Equation (9), we have to use the information on transmitted symbols, obtained from the SISO decoder through the probability ˆ (t) ), to update the channel estimate at each P(Sk |Yk , H k iteration. Furthermore, the soft demapper requires an estimate of the channel in order to provide the probability of encoded bits (see Figure 3). Hence, the proposed semi-blind channel estimation algorithm is naturally combined with the ˆ (t) ) of process of data decoding. The probability P(Sk |Yk , H k the unknown symbols Sk is calculated using the a posteriori probabilities provided by the SISO decoder at the end of the tth iteration as B ˆ (t) = P Sk |Yk , H Pdec (ck,j ) k
(25)
j=1 j =i
where Pdec (ck,j ) are prior probabilities coming from the SISO decoder. At the first iteration, where no a priori information is available on bits ck,j , the probabilities Pdec (ck,j ) are set to 0.5. The main steps of the EM-MAP estimator are summarised as follows:
r Initialisation (t = 0)
r § In
– Set all probabilities of coded bits Pdec (ck,i ) to 0.5 and ˆ (0) ) for all k according to Equation derive P(Sk |Yk , H k (25) and then derive p(S|Y, H(0) ) = p(S|Y, g(0) ). – Initialise the unknown vector g by g(0) obtained from pilot symbols. for {t = 1, . . . , tmax }
this paper, B = 2 since we are using a QPSK modulation.
Copyright © 2008 John Wiley & Sons, Ltd.
r
– Use the previous estimate g(t−1) and p(S|Y, g(t−1) ) to calculate g˜ (t−1) according to Equation (15). – Use g˜ (t−1) to obtain the updated channel parameters g(t) by using Equation (22). – Discard the wavelet coefficients that are replaced (t) by zero in g(t) by evaluating gtr and Ttr from Equation (24). (t) if {t = tmax }: Use the current estimate gtr to update the probability of encoded bits Pdec (ck,i ) and derive ˆ (t) ) from Equation (25) and then derive P(Sk |Yk , H k (t) p(S|Y, H ) = p(S|Y, g(t) ). else: Decode the information data by thresholding the uncoded bit probabilities. end
5. SIMULATION RESULTS In this section we present a comparative performance study of the proposed EM-MAP algorithm. The performance comparison is made in terms of mean square error (MSE) for channel estimation and bit error rate (BER) for the combination of channel estimation with the decoding process. The binary information data are encoded by a nonrecursive non-systematic convolutional encoder with rate R = 1/2 and constraint length 3, defined in octal form by (5, 7)8 . Throughout the simulations, each frame is composed of Nsym = 9 OFDM symbols with N = 128 subcarriers each. Channel coefficients are kept constant during each fading block and changed to new independent realisations (measures) from one frame to the next. We devote one OFDM pilot symbol for each subband in order to initialise the EM algorithm. Data and pilot symbols belong to the QPSK constellation with gray labelling. The interleaver is pseudo-random, operating over the entire frame of size NI = Nsym NB bits (excluding pilots, obviously). Among different wavelet families, ‘symmetric’ wavelet basis functions [19] providing a sparser representation [14] have been considered. Unless otherwise mentioned, the BER and MSE curves correspond to the fourth iteration of the algorithm. The parameter ρ affecting the variance of the noise in Equation (16) is set to 0.4 which minimises the MSE between the true and estimated wavelet coefficients. Different propagation environments are considered for performance comparison. In all environments, the Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
768
S. M. S. SADOUGH ET AL.
bandwidth of the channel frequency response is 1584 MHz with 384 coefficients and the CIR has a total number of 96 taps. For each transmitted frame, a different measurement/realisation of the channel has been used. First, a sparse channel model where only 20 wavelet coefficients out of total 96 have non-zero values, is considered. This propagation environment provides the best adequacy between the sparseness prior assumption and the actual propagation environment. The second channel is a line of sight (LOS) scenario issued from realistic UWB indoor channel measurements [20] where the receive and transmit antennas are located in a corridor separated by 9 meters. We also consider the non-LOS UWB channel model CM2 specified by the IEEE 802.15.3a channel modelling subcommittee report [21]. The two latter channels characterise a rich scattering environment where the sparseness assumption in the wavelet domain is not necessarily satisfied. For the sake of comparison, we consider two pilot-only based approach using ML and minimum mean square error (MMSE) channel estimation, refered to as pilot-ML and pilot-MMSE. We also compare the proposed algorithm with two semi-blind channel estimation based on the EM algorithm. The first approach called EM-Freq, consists of estimating the 384 frequency coefficients of the channel over all three subbands, using the model (3), similar to the approach adopted in Reference [11]. The second semi-blind method, called EM-Wav, is a wavelet domain EM based channel estimator where the prior model is set to have a uniform distribution. First, consider the case of the sparse channel. Figure 4 depicts the MSE as a function of Eb /N0 . It can be observed that, although the pilot-MMSE approach improves the estimation accuracy for low SNR values, the performance of pilot based channel estimation methods are very far from the family of semi-blind methods. Comparing the wavelet domain semi-blind approach (EM-Wav) and the frequency domain approach (EM-freq), shows that significant gain is achieved by the former method. As shown, the best performance is achieved by the EM-MAP method. We see that by using EM-MAP, a gain of almost 4 dB in SNR is achieved at a MSE of 2 × 10−3 , as compared to the EM-Wav method. This clearly shows the adequacy of the EM-MAP method for the case where the unknown channel has few non-zero wavelet coefficients, which is in perfect agreement with the prior model. Figure 5 shows the BER results along with the BER for the case of perfect channel state information (CSI). It can be seen that at a BER of 10−3 , the pilot-ML and the EM-Freq approaches are respectively 3.9 and 2 dB of SNR far from Copyright © 2008 John Wiley & Sons, Ltd.
Figure 4. Mean square error between the true and estimated coefficients for the sparse channel model.
Figure 5. BER performance of different channel estimation methods over the sparse channel model.
the BER obtained with the perfect channel. Furthermore, the performance of the Pilot-MMSE approach is not shown since it was very close to that of Pilot-ML. Also, we observe that wavelet based semi-blind methods perform closely to the perfect CSI case. For example, at a BER of 10−4 , the EM-MAP and EM-Wav method have respectively about Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
0.2 and 0.5 dB of SNR degradation from the performance obtained with perfect CSI. We now evaluate the performance of the EM-MAP algorithm for the case of the Corridor and CM2 channels. Figures 6 and 7 show that wavelet based methods again outperform pilot based and EM-Freq methods in terms
769
Figure 6. Mean square error between the true and estimated coefficients over the Corridor channel.
of MSE and BER. Furthermore, we observe that the performance of the EM-MAP method is now comparable to that of the EM-Wav method. This can be explained by noting that when the channel is not very sparse (i.e. for low λ values in (18)), the algorithm assigns small values to λ. This leads to a Gaussian prior model with a large variance, which can be approximated to a uniform prior. As a result, the prior becomes ‘less informative’ and the EM-MAP performs close to the EM-Wav, as shown in Figures 6 and 7. Thus, the EM-MAP algorithm proposed here is globally able to adapt its prior model parameters to any propagation environment. We now compare the above semi-blind algorithms with respect to the average number of estimated parameters at each iteration of the EM algorithm. This is shown in Figure 8 for different channel scenarios. First, recall that the EMFreq and EM-Wav methods have to continuously estimate 384 and 96 coefficients at each iteration, respectively. As explained in Subsection 3.2, by discarding the coefficients that are replaced by zero, the EM-MAP approach is able to reduce the number of estimated parameters, specially for sparse channels. This can be seen in Figure 8 for the sparse channel where the number of estimated parameters is reduced from 96 down to 20 parameters after the second iteration. Furthermore, for the non-sparse Corridor and CM2 channels, we observe that the EM-MAP method is to be preferred to the EM-Wav, due to its lower computational load, since it estimates about 40 coefficients after the second iteration. Although we observed that in this case
Figure 7. BER performance of different channel estimation methods over the CM2 channel.
Figure 8. Reduction of the number of estimated parameters through iterations, Eb /N0 = 8 dB.
Copyright © 2008 John Wiley & Sons, Ltd.
Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
770
S. M. S. SADOUGH ET AL.
a larger bandwidth is likely to separate the paths and to provide a sparser wavelet expansion. It was observed that when the channel has a very sparse wavelet expansion, the prior model parameters which are estimated from the observed data, carry this sparseness information to the EM-MAP algorithm. Moreover, we showed that in this case, the EM-MAP method provides significant reduction in the number of estimated parameters and outperforms all considered pilot based and semi-blind methods. Under non-sparse channels, although both EM-MAP and EM-Wav methods perform closely, EM-MAP takes the advantage over the EM-Freq and EM-Wav schemes due to its lower computational complexity.
REFERENCES
Figure 9. Convergence of the MSE with respect to the number of iterations over the CM2 channel.
these two methods exhibit close performance, the EM-MAP algorithm brings a reduction of about 60% on the number of estimated parameters when compared to the EM-Wav approach. Finally, we analyse the number of iterations that the EMMAP algorithm requires for convergence. In Figure 9, the MSE performance of the EM-MAP algorithm is presented as a function of the number of iterations. It is obvious from these curves that the MSE performance of the proposed algorithm converges within two to four iterations, depending on the average SNR. This justifies our choice of four iterations in the results presented above.
6. CONCLUSION We proposed an EM-based channel estimation and data detection algorithm that integrates the advantages of wavelet based parameter estimation. By expressing the unknown UWB channel in terms of its discrete wavelet coefficients and choosing a proper prior distribution for them, we were able to capture the sparseness property of UWB channels in the wavelet domain. This led to a MAP estimator equivalent to a hard thresholding procedure at each iteration of the EM algorithm, which is used to reduce the number of estimated coefficients. Note that the sparseness condition is almost built in by considering all subbands of the multiband transmission at the same time: Copyright © 2008 John Wiley & Sons, Ltd.
1. FCC. First report and order, revision of part 15 of the commission’s rules regarding ultra-wideband transmission systems. Technical report, February 2004. 2. Roy S, Foerster JR, Somayazulu VS, Leeper DG. Ultrawideband radio design: the promise of high-speed, short range wireless connectivity. Proceedings of the IEEE 2004; 92:295–311. 3. Batra A, Balakrishnan J, Aiello GR, Foerster JR, Dabak A. Design of multiband OFDM system for realistic UWB channel environments. IEEE Transactions Microwave Theory and Techniques 2004; 52:2123–2138. 4. Prasad R. OFDM for Wireless Communications Systems. Artech House Publishers, Norwood, MA, USA, 2004. 5. Caire G, Taricco G, Biglieri E. Bit-interleaved coded modulation. IEEE Transactions and Information Theory 1998; 44:927–945. 6. van de Beek JJ, Edfors O, Sandell M. On channel estimation in OFDM systems. In Proceedings of Vehicular Technology Conference, September 1995; 815–819. 7. McLachlan GJ, Krishnan T. The EM Algorithm and Extensions. Wiley Series in Probability and Statistics, Wiley, Chichester, UK, 1997. 8. Zhuang X, Vook FW. Iterative channel estimation and decoding for a turbo-coded OFDM system via the EM algorithm. In Proceedings of IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), May 2002. 9. Jaffrot E, Siala M. Turbo channel estimation for OFDM systems on highly time and frequency selective channels. In Proceedings of IEEE International Conference on Acoustics, Speech, Signal Processing (ICASSP), June 2000. 10. Lu B, Wang X. Bayesian blind turbo receiver for coded OFDM systems with frequency offset and frequency-selective fading. IEEE Journal on Selected Areas in Communications. 2001; 19:2516–2527. 11. Ma X, Kobayashi H, Schwartz SC. EM-based channel estimation algorithms for OFDM. EURASIP Journal on Applied Signal Processing 2004; 10:1460–1477. 12. Lu B, Wang X, Li YG. Iterative receivers for space-time block coded OFDM systems in dispersive fading channels. IEEE Transactions on Wireless Communications 2002; 1:213–225. 13. Johnstone IM, Silverman BW. Neddles and straw in haystacks: empirical bayes estimates of possibly sparse sequences. Annals of Statistics 2004; 32:1594–1649. 14. Sadough S, Jaffrot E. A wavelet packet based model for an ultrawideband propagation channel. In Proceedings of ECPS, 15–18 March 2005.
Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
ULTRA WIDEBAND OFDM CHANNEL ESTIMATION
15. Nowak RD, Figueiredo MAT. Fast wavelet-based image deconvolution using the EM algorithm. Asilomar Conference on Signals, Systems and Computers, 2001; 371–375. 16. Sadough S. Ultra Wideband OFDM Systems: Channel Estimation and Improved Detection Accounting for Estimation Innacuracies. PhD thesis, Universit´e Paris-Sud 11, January 2008. 17. Bernardo JN, Smith AFM. Bayesian Theory. Wiley and Sons: Chichester, England, 1994.
Copyright © 2008 John Wiley & Sons, Ltd.
771
18. Bahl L, Cocke J, Jelinek F, Raviv J. Optimal decoding of linear codes for minimizing symbol error rate. IEEE Transactions on Information Theory 1974; 20(2):284–287. 19. Mallat S. A Wavelet Tour of Signal Processing. Academic Press, 1999. 20. Bories S, Sibille A, Roblin C. UWB indoor channel measurement. In Proceedings of IWAT, March 2005; 466–469. 21. Foerster J. Channel modeling sub-committee report final. Technical report, IEEE802.15-02/490, 2003.
Eur. Trans. Telecomms. 2008; 19:761–771 DOI: 10.1002/ett
