Particle Filtering Equalization Method for a Satellite Communication Channel

We propose the use of particle ﬁltering techniques and Monte Carlo methods to tackle the in-line and blind equalization of a satellite communication channel. The main di ﬃ culties encountered are the nonlinear distortions caused by the ampliﬁer stage in the satellite. Several processing methods manage to take into account these nonlinearities but they require the knowledge of a training input sequence for updating the equalizer parameters. Blind equalization methods also exist but they require a Volterra modelization of the system which is not suited for equalization purpose for the present model. The aim of the method proposed in the paper is also to blindly restore the emitted message. To reach this goal, a Bayesian point of view is adopted. Prior knowledge of the emitted symbols and of the nonlinear ampliﬁcation model, as well as the information available from the received signal, is jointly used by considering the posterior distribution of the input sequence. Such a probability distribution is very di ﬃ cult to study and thus motivates the implementation of Monte Carlo simulation methods. The presentation of the equalization method is cut into two parts. The ﬁrst part solves the problem for a simpliﬁed model, focusing on the nonlinearities of the model. The second part deals with the complete model, using sampling approaches previously developed. The algorithms are illustrated and their performance is evaluated using bit error rate versus signal-to-noise ratio curves.


INTRODUCTION
Telecommunication has been taking on increasing importance in the past decades and thus led to the use of satellitebased means for transmitting information. A major implementation task to deal with such an approach is the attenuation of emitted communication signals during their trip through the atmosphere. Indeed, one of the most important roles devoted to telecommunication satellites is to amplify the received signal before sending it back to Earth. Severe technical constraints, due to the lack of space and energy available on board, can be solved thanks to special devices, namely, traveling-wave-tube (TWT) amplifiers [1]. A common model for such a satellite transmission chain is depicted in Figure 1.
Although efficient for amplifying tasks, TWT devices suffer from nonlinear behaviors in their characteristics, thus implying complex modeling and processing methods for equalizing the transmission channel.
The very first approaches for solving the equalization problem of models similar to the one depicted in Figure 1 were developed in the framework of neural networks. These methods are based on a modelization of the nonlinearities using layers of perceptrons [2,3,4,5,6]. Most of these approaches require a learning or training input sequence for adapting the parameters of the equalization algorithm. However, the knowledge or the use of such sequences is sometimes impossible: if the signal is intensely corrupted by noise at the receiver stage or for noncooperative applications, for instance.  Figure 2: Identification of the model depicted in Figure 1 with a Volterra filter.
Blind equalization methods have thus to be considered. These methods often need precise hypothesis with the emitted signals: Gaussianity or circularity properties of the probability density function of the signal, for instance [7]. Recently, some methods make it possible to identify [8] or equalize [9,10,11] blindly nonlinear communication channels under general hypothesis. These blind equalization methods assume that the transfer function of the system can be modeled as a Volterra filter [12,13].
However, for the transmission model considered here, a Volterra modelization happens to be only suitable for the task of identification and not for a direct equalization. For instance, a method based on a Volterra modelization of the TWT amplifier and a Viterbi algorithm at the receiver stage is considered in [14]. Such an identification method can be easily implemented through a recursive adaptation rule of the filter parameters with a least mean squares approach (cf. Figure 2). The mean of the quadratic error (straight line) and its standard deviation (dotted lines) are depicted in Figure 3 for 100 realizations of binary phase shift keying (BPSK) symbol sequences, each composed of 200 samples. Similarly, the equalization problem of the transmission chain 1 can be considered with a Volterra filter scheme, adapted with a recursive least squares algorithm as depicted in Figure 4. However, in  Figure 4: Equalization of the model depicted in Figure 1 with a Volterra filter.
this case, the error function happens not to converge showing that the Volterra filter is unstable and that the system is not invertible with such modelization. It is then necessary to consider a different approach for realizing blindly the equalization of this communication model. The aim of this paper is thus to introduce a blind and sequential equalization method based on particle filtering (sequential Monte Carlo techniques) [15,16]. For realizing the equalization of the communication channel, it seems interesting to fully exploit the analytical properties of the nonlinearities induced by TWT amplifiers through parametric models of these devices [1]. Sequential Monte Carlo methods, originally developed for the recursive estimation of nonlinear and/or non-Gaussian state space models [17,18], are well suited for reaching this goal. The field of communication seems to be particularly propitious for applying particle filtering techniques, as shown in the recent literature of the signal processing community [15,19] and this issue and of the statistics community (see [20,21], [16,Section 4]). Such Monte Carlo approaches were successfully applied to blind deconvolution [22], equalization of flat-fading channels [23], and phase tracking problems [24], for instance.
The paper is organized as follows. Firstly, models for the emitted signal, the TWT amplification stage, and the other parts of the transmission chain are introduced in Section 2. A procedure for estimating the emitted signal is considered in a Bayesian framework. Monte Carlo estimation techniques are then proposed in Section 3 for implementing the computation of the estimated signal under the assumption of a simpler communication model, focusing on the nonlinear part of the channel. This approach uses analytical formulae of the TWT amplifier model described in Section 2 and sampling methods for estimating integral expressions. The method is then generalized in Section 4 for building a blind and recursive equalization scheme of the complete transmission chain. The sequential simulation algorithm proposed is based on particle filtering techniques. This approach makes it possible to process the data in-line and without the help of a learning input sequence. The performance of the algorithm is illustrated by numerical experiments in Section 5. Finally, some conclusions are drawn in Section 6. Details of the Monte Carlo approach are given in Appendix A.

MODELING OF THE TRANSMISSION MODEL
The model of the satellite communication channel depicted in Figure 1 is roughly the same as the one considered for various problems dealing with TWT amplifiers devices (cf., e.g., [2]). The different stages of this communication channel are detailed below.

Emission stage
The information signal to transmit is denoted by e(t). It is usually a digital signal composed of a sequence of N e symbols (e k ) 1≤k≤Ne . The signal is transmitted under the analog form where T denotes the symbol rate and I Ω (·) is the indicator function of set Ω. Symbols e k are generated from classical modulations used in the field of digital telecommunication, like PSK or quadratic amplitude modulation (QAM), for instance. In the following, the case of 4-QAM symbols is considered. Each symbol can be written as where the sequence of samples (φ k ) 1≤k≤Ne is independently and identically distributed from where U Ω denotes the uniform distribution on the set Ω.
The signal is emitted through the atmosphere to the satellite. The emission process is modeled by a Chebyshev filter. This class of filters admits an IIR representation and their parameters, particularly their cutoff frequency, depend on the value of symbol rate T [2]. A detailed introduction to Chebyshev filters is given in [25], for instance. In the present case, the emission filter is assumed to be modeled with a 3 dB bandwidth equal to 1.66/T. The emitted signal is altered during its trip in the atmosphere by disturbance signals. These phenomena are modeled by an additive corrupting noise n e (t), which is assumed to be Gaussianly, independently, and identically distributed: where N CC 0, σ 2 e is a complex circular Gaussian distribution, with zero-mean and variance equal to σ 2 .
Remark 1. The amplitude of signal (1) is adjusted in practice at the emission stage in order to reach a signal-to-noise ratio (SNR) roughly equal to 15 dB during the transmission.

Amplification
After being received by the satellite, the signal is amplified and sent back to Earth. This amplification stage is processed by a TWT device. A simple model for TWT amplifier is an instantaneous nonlinear filter defined by where r denotes the modulus of input signal. Amplitude gain and phase wrapping can be modeled by the following expressions: These formulae have been shown to model various types of TWT amplifier device with accuracy [1]. Figures 5 and 6 represent functions (6) and (7) for two sets of parameters estimated in [1, Table 1] from real data and duplicated in Table 1. Curves with straight lines represent functions obtained with the set of parameters of the first row of Table 1. The ones with dashed lines represent functions obtained with the other set of parameters. A drawback of model (5) is that it is not invertible in a strict theoretical sense, as drawn in Figure 5. However, only the amplificative and invertible part of the system, represented above the dotted line on Figure 5, will be considered.    Signal processing in the satellite also performs the task of multiplexing. The devices used for this purpose are modeled by Chebyshev filters. Tuning of their parameters is given in [2], for instance. In the present case, filters at the input and at the output of the amplifier are assumed to have bandwidths equal, respectively, to 2/T and 3.3/T.

Reception
The transmission of the signal back to Earth is much less powerful than at the emission stage. This is mainly due to severe technical constraints because of the satellite design. The influence of the atmospheric propagation medium is then modeled by a multipath fading channel [26,Section 11], with one reflected path representing an attenuation of 10 dB in this case: z(t) → z(t) + αz(t − ∆). Moreover, the signal is still corrupted by disturbance signals, modeled by an additive noise signal n r (t), Gaussianly, independently, and identically This noise is always much more intense than at the emission stage. This is mainly due to the weak emission power available in the satellite. The received signal, denoted as r(t), is sampled at rate T s .

Equalization
The goal of equalization is to recover emitted sequence (e k ) 1≤k≤Ne from the knowledge of sampled sequence (r( jT s )) 1≤ j≤Nr . The equalization method proposed in this paper consists in estimating symbol sequence (φ k ) 1≤k≤Ne by considering its posterior distribution conditionally to sequence (r( jT s )) 1≤ j≤Nr of samples of the received signal: To reach this goal, a Bayesian estimation procedure is considered with the computation of maximum a posteriori (MAP) estimates [27].
Remark 2. Bayesian approaches have already been successfully applied in digital signal processing in the field of mobile communication. In [28], for instance, autoregressive models and discrete-valued signals are considered.
The computation of the estimates is implemented via Monte Carlo simulation methods [16,29]. As the complete transmission chain is a complex system, a simpler model focusing on the nonlinear part of the channel is considered in the following section, where Monte Carlo estimation approaches are introduced. These estimation techniques will be used in the equalization algorithm for the global transmission chain in Section 4.

MONTE CARLO ESTIMATION METHODS
As a first approximation, to focus on the nonlinearity of the model, only a TWT amplifier is considered in a transmission channel corrupted with noises at its input and output parts as shown in Figure 7.
The received signal r(t) is assumed to be sampled at symbol rate T. The problem is then to estimate a 4-QAM symbol φ a priori distributed from (3) with the knowledge of the model depicted in Figure 7 (cf. relations (4), (6), (7) and (8)), and information of a received sample r. A Bayesian approach is developed [27] by considering the posterior distribution n e (t) Figure 7: Simple communication channel with TWT amplifier. and its classical MAP estimate. The method proposed in the following consists in estimating values of distribution (10) thanks to Monte Carlo simulation schemes [29] using relations (6) and (7) which model nonlinearities of the TWT amplifier in a parametric manner.

Estimation with known parameters
In order to further simplify the study, parameters of the transmission channel depicted in Figure 7 are firstly assumed to be known. In the sequel, the coefficients of expressions (6) and (7) are denoted by the symbol TWT. This information is taken into account in the posterior distribution (10) thus becoming where A denotes the amplitude of the emitted signal. From Bayes' formula, the probability density function of this distribution is proportional to The prior distribution at the right-hand side of the above expression reduced to p(φ), which is given by (3). The problem is then to compute the likelihood p r A, φ, σ e , TWT, σ r .
Indeed, this formula can be viewed (cf. Appendix A) as the following expectation: with respect to the random variable x which is Gaussianly distributed: Considering a sequence of samples (x ) 1≤ ≤N independently and identically distributed from (15), a Monte Carlo approximation of (14) is given by which is accurate for a number N of samples large enough.
References [29,30] provide detailed ideas and references about Monte Carlo methods. To illustrate such an approach, approximation (16) is computed for the emitted symbol φ = π/4 and the values of TWT amplifier parameters given by  (15) of 100 samples is considered. Table 2 gives mean values obtained from (16) and their standard deviations. The error of the estimated values (16) of probabilities (11) might seem quite large as the standard deviations can be reduced providing a larger number of samples. In the sequel, we are only interested in obtaining rough estimates of (11), enabling comparison of mean values for different φ as shown in Table 2. Thus, even with a reduced number of samples (15), it is possible to estimate accurately the MAP estimate of (11).
Performance of the Monte Carlo estimation method is then considered with respect to SNR at the input and output of the amplifier (cf. Figure 7). The bit error rate (BER) is computed by averaging the results obtained with a MAP approach. Statistics of the Monte Carlo estimates (16) of distribution (11) are computed with 100 realizations of symbol sequences composed of 1, 000 samples each. For each estimate, sequences (15) composed of 100 samples are considered. The results of these simulations for SNR e taking values 10, 12, and 15 dB are depicted in Figure 8 and curves from the bottom to the top are associated to decreasing SNR e .
The Bayesian approach and its Monte Carlo implementation make it possible to estimate the emitted signal with accuracy for a wide range of noise variances (cf. Figure 8). However, the estimation method described previously requires the knowledge of the model parameters. For many applications in the field of telecommunication, it is necessary to assume these parameters unknown. It is the case for nonstationary transmission models and for communication in noncooperative contexts like passive listening, for instance. The equalization problem of the simplified model depicted in Figure 7 is now tackled in the case where the parameters (A, σ e , TWT, σ r ) of the transmission channel are assumed to be unknown.

Estimation with unknown parameters
If the parameters are unknown, there are, at least, two Bayesian estimation approaches to be considered with respect to posterior distribution (10). A first method consists This method makes it possible theoretically to jointly estimate the emitted symbols and the parameters of the transmission channel by implementing MAP and/or posterior mean approaches. The probability density function of distribution (17) being generally very complex, Markov chain Monte Carlo (MCMC) simulation methods [29,30] can be used to perform these estimation tasks. Such an approach is developed in [31] particularly for equalizing the complete transmission chain depicted in Figure 1. However, results obtained with this method happen not to give accurate estimates of the model parameters in practice. Indeed, MCMC methods are generally useful for estimating various models in the field of telecommunication [28,32]. Another approach consists in considering a marginalized version of distribution (10) with respect to the parameters of the model: p(φ) = p φ, A, σ e , TWT, σ r r d A, σ e , TWT, σ r . (18) Such a technique, called Rao-Blackwellization in the statistics literature, for example, [33,34], makes it possible to improve the efficiency of sampling schemes (see [20,21], [16,Section 24]). From Bayes' formula, the integrand of expression (18) is proportional to p r φ, A, σ e , TWT, σ r × p φ, A, σ e , TWT, σ r .
Assuming that symbols and the model parameters are independent, expression (18) From the study of the previous case, the likelihood term in the integrand can be computed via a Monte Carlo estimate of expression (14) with a sequence of samples (15). An approach to estimate (18) is then to consider the integral expression in (20) as the expectation E p(A,σe,TWT,σr ) p r φ, A, σ e , TWT, σ r (21) which is estimated via a Monte Carlo approximation of the following form: where (A k , σ e (k), TWT k , σ r (k)) k=1,...,Np is a sequence of samples independently and identically distributed from the prior distribution p A, σ e , TWT, σ r .
Remark 3. The algorithm for sampling distribution (17) introduced in [31] requires also the setting of prior distribution (23).
The model of the parameters includes generally prior information thanks to physical constraints. For instance, the TWT amplifier is assumed to work in the amplificative part of its characteristic (cf. Figure 5). Thus, as a first rough approximation, it can be assumed that A ∼ U [0,1] a priori. This parameter is also tuned such that SNR e equals 15 dB during the emission process (cf. Remark 1) implying the constraint σ e = 0.2A. In a less strict case, it is sufficient to assume that σ e ∼ U [0.01,0.5] . The parameters (α a , β a , α p , β p ) of the TWT amplifier are supposed to be independent of other variables of the system and also to be mutually independent. From the values introduced in Table 1, an adequate prior distribution is [2,10] . (24) The extremal values of the downlink transmission noise variance σ r can be estimated with respect to prior ranges of values defined above. A uniform U [0.1,1.1] prior distribution for σ r is thus chosen. Once all prior distributions have been defined, it is possible to implement a Monte Carlo estimation procedure for (20) with the help of approximations (22) and (16). Such an approach is tested for the computation of values of posterior distribution for an emitted symbol φ = π/4 considering that the values of the TWT amplifier are given by the first row of Table 1, that the amplitude of the emitted signal is given by A = 0.5, and that noise variances are such that are drawn from distribution (23). For each sequence, as in the case where the model parameters are known, approximations (16) are computed from sequences (25) composed of 100 samples each. Table 3 shows the mean values of the estimated (22) and their standard deviations computed from these simulations. As for the previous case, where parameters (A, σ e , TWT, σ r ) are known, we are only interested in obtaining rough mean values of Monte Carlo estimates of the MAP expression (10) and thus do not consider larger sample sizes for reducing the standard deviation of these estimates.
Performance of this Monte Carlo estimation method is now considered with respect to uplink and downlink SNR. As previously, BERs are computed by averaging results obtained with a MAP approach for the Monte Carlo estimate (22) of posterior distribution (10). One hundred realizations of 1 000-symbol sequences are considered for each value of SNR. For every Monte Carlo estimate, sequences (25) and (15) are composed of 100 samples. The results of simulations for SNR e equal to 10, 12, and 15 dB are depicted in Figure 9. Curves from the bottom to the top are associated to a decreasing uplink SNR. As a comparison, the estimated mean values of BER in the case where the parameters of the TWT amplifier are known are represented with dashed lines.
Performance is not much corrupted in the case where model parameters are unknown. Thus, considering the posterior distribution of interest (18), marginalized with respect to these parameters, seems to be a good strategy for tackling the equalization problem. An in-line simulation method based on the Monte Carlo estimation techniques previously developed is proposed hereinafter for realizing the equalization of the complete transmission chain depicted in Figure 1

Transmission model
Equalizing the complete satellite communication channel depicted in Figure 1 is a difficult problem as it requires taking  (1) effects of the filters modeling, emission, and multiplexing stages; (2) attenuation of the received signal mainly due to multiple paths during the downlink transmission; (3) correlation induced by filters and emission and fading models.
An equalization method is proposed for this model within a Bayesian estimation framework [27]. It consists in considering the posterior distribution of the sampled symbols conditionally to the sequence of the received samples: p e jT s 1≤ j≤Nr r jT s 1≤ j≤Nr .
An estimation procedure is then implemented by computing the MAP estimate of distribution (26). Monte Carlo estimation methods developed in the previous paragraphs can be slightly modified in order to take into account the parameters of the complete transmission chain (cf. points (1) and (2) above). The correlation of the samples at the receiver stage mainly comes from the linear filters in the channel. In fact, this problem yields to the estimation of parameter p: the number of received samples per symbol rate p = T/T s , as parameters of Chebyshev filters at the emission and multiplexing stages depend on its value.
Computing the correlation of the received samples makes it possible to give an estimate of p [31] in the case where this quantity is an integer, and thus to estimate the parameters of the filters. In the sequel, we consider that this is the case, assuming that a proper synchronization processing has been performed at the receiver stage. This task can also be achieved via Monte Carlo simulation methods [24]. This parameter p will be used in an explicit manner in the recursive equalization algorithm introduced in Section 4.3. An MCMC simulation scheme [29,30], for the batch processing of received data, was studied in [31]. A sequential simulation method for sampling the distribution (26) is now introduced, as many applications in the field of telecommunication require in-line processing methods when data is available sequentially.

Sequential simulation method
A sequential method for sampling distribution (26) can be implemented via particle filtering techniques [15,16]. The wide scope of this approach, originally developed for the recursive estimation of nonlinear state space models [17,18,20,21], is well suited for the sampling task of this equalization problem. The basic idea of particle filtering is to generate iteratively sequences of the variables of interest, each of them denoted as a "particle," here written as x 0 (i), x 1 (i), . . . , x t (i) 1≤i≤M (27) such that particles (x t (1), . . . , x t (M)) at time t are distributed from the desired distribution, denoted as p t (x). This goal can be reached with the use of two "tuning factors" in the algorithm: (i) the way the particles are propagated or diffused, x t (i) → x t+1 (i), in the sampling space, namely, the choice of a proposal or candidate distribution; (ii) the way the distribution of particles (x t (1), . . . , x t (M)) approximates the target distribution p t (x): by affecting a weight w t (i) to each particle depending on the proposal distribution, and updating these weights with an appropriate recursive scheme.
Algorithm 1: Equalization algorithm. explained in the next paragraph. Such Monte Carlo simulation scheme is now proposed to tackle the sequential sampling of distribution (26).

Equalization algorithm
In the present case, phase samples φ j = φ( jT s ) of the emitted signal are directly sampled. The simulation scheme which is considered is the bootstrap filter [15,16,17,18] and is given in Algorithm 1.
The important sampling and computation steps (28) and (29) are detailed hereinafter.
The information brought by parameter p, number of received samples per symbol duration, is taken into account via the proposal distribution (28). Indeed, candidates for particles φ j (i) can be naturally sampled from the following scheme: This sampling scheme is very simple and can easily be improved by considering φ kp (i) ∼ U {π/4,3π/4,5π/4,7π/4} and φ kp+s (i) = φ kp (i) for 1 ≤ k ≤ p − 1, for instance. However, the scheme above gives sufficiently accurate results as a first approximation, due to its flexibility (if a false symbol is chosen, there is probability to switch to other symbols again) and its ability to deal with possible uncertainty on the value of parameter p. This scheme is also efficient to limit the negative effect of sample impoverishment due to the resampling step, as detailed hereinafter. The proposed scheme, however, does not take into account completely the information coming from emission and received signals and if some coding techniques are used to generate the symbols, this knowledge should be introduced in the sampling scheme (28) if possible.
The computation of weights (29) is realized by using similar Monte Carlo approaches to the ones introduced previously, including filters and their parameters in expressions (14) and (18). In this respect, Algorithm 1 can be seen as a Rao-Blackwellized particle filter [20,21] where the parameters of the channel, considered here as nuisance parameters, are integrated out. This generally helps to lead to more robust estimates in practice [16,33,34].
A crucial point in the implementation of particle filtering techniques lies in the resampling stage, step (4) in Algorithm 1. As the computations for sampling the candidates (28) and updating the weights (29) can be performed in parallel, the resampling step gives the main contribution in the computing time of the algorithm as its achievement needs the interaction of all the particles. This stage is compulsory in practice if one wants the sampler to work efficiently. This is mainly due to the fact that the sequential importance sampling algorithm without resampling increases naturally the variance of the weights (w j (i)) 1≤i≤M with time [20,22,35]. In such case, only a few particles are affected nonnegligible weights after several iterations, implying a poor approximation of the target distribution and a waste of computation.
To limit this effect, several approaches can be considered [15]. One consists in using very large numbers of particles M and/or in performing the resampling step for each iteration [17,18]. However, resampling too many times often leads to severe sample impoverishment [16,20,21]. Other methods, also aiming at minimizing computational and memory costs, consist in using efficient sampling schemes for diffusing the particles [20] and performing occasionally the resampling stage when it seems to be needed [15,21]. When to perform resampling can be decided by measuring the variance of weights via the computation of the effective sample size M/(1 + var( w j (i))), whose one estimate is given by [15,16,21,35]. In this case, the resampling stage can be performed each time the estimated effective sample size is small, measuring how the propagation of the particles in the sampling space is efficient. This quantity equals M for uniformly weighted particles and equals 1 for degenerated cases where all the particles have zero weights except one.
It is also possible to compute the entropy of the weights, describing "how far" the distribution of the weights is from the uniform distribution. Indeed, the entropy is maximized for uniform weights and minimized for the degenerated configurations as mentioned above. In this sense, the entropy of weights quantifies the information of the samples and measures the efficiency of representation for a given population of particles. This approach is adopted in [24,36], for instance, and also in our algorithm as follows.
Step (4) of Algorithm 1 is therefore replaced by the computing of entropy of the weights and a resampling/selection step is processed only if the con- holds, assuming that λ is a threshold value set by the user.
To show that the estimated effective sample size and entropy lead to similar results for the resampling task, their values for one realization of the algorithm are depicted in Figure 11. Also, the resampling step can be performed via different techniques [18,21]. In the sequel, we use the general multinomial sampling procedure [16,18].
As the variable of interest φ is distributed from a set of discrete values, using large numbers of particles M happened to be unnecessary. Simulations have shown that several dozens are sufficient for the problem considered here. This is also the case for other Monte Carlo simulation methods used for estimating telecommunication models [37]. Moreover, the degeneracy phenomenon is not really a nuisance in the present case as the simulation algorithm Algorithm 1 is only used in a MAP estimation purpose and not for computing mean integral estimates as (18) and (A.4), for instance.
Some simulations concerning performance of the equalization algorithm with this sequential sampling scheme are now presented.

NUMERICAL EXPERIMENTS
The equalization method was run for 100 realizations of sequences of 1 000 samples for each value of SNR at the receiver stage and for SNR equal to 12 dB at the emission stage. The number of received samples per symbol rate is p = 8. The channel model is depicted in Figure 1.
(a) The amplifier model is described in Section 2.2 where parameters are set as the first row of Table 1. For each realization, the emitted symbol sequence is estimated by considering the MAP trajectory computed from a Monte Carlo approximation of the distribution (26) with M = 50 particles (27). Weights (29) are computed with the help of Monte Carlo estimation techniques introduced in Section 3.2, considering sequences of 100 samples. In our simulations, the value for threshold (31) λ = 0.1 gave acceptable estimation results. The mean values (straight lines) and their associated variances (dotted lines) of the BER of MAP estimates are depicted in Figure 12.
The equalization algorithm was also run for different values of uplink SNR: 10 dB and 15 dB. The mean values of the BER computed from the estimated phases are depicted in Figure 13. Curves from the bottom to the top are associated to a decreasing SNR e .
One of the advantages of the proposed equalization method is its robustness with respect to nonstationarities of the transmission channel. This property comes from the MAP estimation procedure considering the distribution (18) marginalized with respect to channel parameters. Simulations including perturbations of the parameters of the transmission chain lead to similar results to those presented in Figures 12 and 13. On the other hand, in case of dysfunc- tion in devices of the amplifier and/or of the filters or if sudden change of noise intensities happens during the transmission, the estimation method remains almost insensitive to these changes. The approach currently developed cannot thus be used for diagnostic purposes, as it is the case for certain methods based on neural networks [4]. An interesting hybrid approach is proposed in the conclusion to cope with this task. In order to compare the performance of the equalization method, BERs computed from the MAP estimates of the symbols are compared with ones obtained with signals estimated by an equalizer built from a 2-10-4 multilayer neural network, using hyperbolic tangents as activation functions [3]. The mean values (straight lines) and standard deviations (dotted lines) of BER computed from Monte Carlo MAP estimates are represented in Figures 14 and 15. The mean values of BER computed from signals estimated with the neural network method are depicted in dashed lines.
The two methods give similar results for this configuration. Nevertheless, an important and interesting characteristic of the sequential Monte Carlo estimation method is that it does not require any learning sequence for equalizing the transmission chain, contrary to approaches based on neural networks. The proposed method is thus efficient in the context of blind communication. A calibration step is at least necessary in order to estimate the number of received samples per symbol rate. This tuning can be realized in a simple manner by computing the correlation of the received samples.

CONCLUSION
The particle filtering equalization method proposed in this paper makes it possible to estimate sequentially digital signals transmitted through a satellite communication channel. The approach takes explicitly into account the nonlinearities induced by the amplification stage in the satellite thanks to a Bayesian estimation framework implemented via Monte Carlo simulation methods. This approach enables to estimate recursively the distribution of interest, namely, the posterior distribution of the variables, marginalized with respect to the parameters of the transmission chain. An advantage of this approach is that it is robust to nonstationarities of the channel. On the contrary, the method cannot detect these nonstationarities and thus be employed to predict some dysfunction of transmission devices, as it is the case for certain neural networks approach. Therefore, it seems interesting to use Markov chain and/or sequential Monte Carlo methods to train appropriate neural networks models. This approach was successfully applied in the field of statistical learning, for instance [16,Section 17].
As for many particle-filtering-based methods, the implementation of the sampling algorithm is generally demanding in terms of computing time, if compared to classical nonsimulation approaches, especially for tackling problems arising in the field of digital communication [19].
However, an advantage of the proposed Monte Carlo equalization method is that it does not require the knowledge of any learning input sequence in order to update the equalizer parameters. To the best of our knowledge, it seems at the moment the only method for achieving directly the blind equalization task of such transmission channel. This blind property can be of premium importance in the framework of communication in noncooperative context. This is the case in passive listening, for instance, or when transmission of learning sequences cannot be completed correctly because of intense noises during the transmission. As stated in Section 3.1, the problem is to compute likelihood (13). A solution consists in integrating this expression with respect to y, the amplified signal (cf. Figure 7): The right probability density function in integral (A.2) can be computed by marginalizing it with respect to x, the signal to amplify (cf. Figure 7):

ACKNOWLEDGMENTS
This work was partly supported by the research program BLISS (IST-1999-14190) from the European Commission. The first author is grateful to the Japan Society for the Promotion of Science and the Grant-in-Aid for Scientific Research in Japan for their funding. The authors thank P. Comon and C. Jutten for helpful comments and are grateful to the anonymous reviewers for their helpful suggestions which have greatly improved the presentation of this paper.