- Email: [email protected]
Quadratic estimation from uncertain observations with white plus coloured noises using covariance information
Applied Mathematics and Computation 155 (2004) 65–79 www.elsevier.com/locate/amc
Quadratic estimation from uncertain observations with white plus coloured noises using covariance information b S. Nakamori a,*, R. Caballero-Aguila , c A. Hermoso-Carazo , J. Linares-Perez c a
Department of Technology, Faculty of Education, Kagoshima University, 1-20-6 Kohrimoto, Kagoshima 890-0065, Japan b Departamento de Estadıstica e Investigaci on Operativa, Universidad de Ja en, Paraje Las Lagunillas, s/n, 23071 Ja en, Spain c Departamento de Estadıstica e Investigaci on Operativa, Universidad de Granada, Campus Fuentenueva, s/n, 18071 Granada, Spain
Abstract In this paper recursive least mean-squared error quadratic filtering and fixed-point smoothing algorithms to estimate signals from uncertain observations are obtained for the case of white plus coloured observation noises. It is assumed that the state-space model of the signal is not known and only the information on the moments, up to the fourth one, of the involved processes and the probability that the signal exists in the observations are available. The estimators require the covariance functions of the signal and coloured noise, as well as the covariance functions of its second-order powers in a semi-degenerate kernel form. Ó 2003 Elsevier Inc. All rights reserved. Keywords: Uncertain observation; Quadratic estimation; Covariance information; False alarm probability; Coloured noise
*
Corresponding author. E-mail addresses: [email protected] (S. Nakamori), [email protected] (R. Caballero-Aguila), [email protected] (A. Hermoso-Carazo), [email protected] (J. Linares-Perez). 0096-3003/$ - see front matter Ó 2003 Elsevier Inc. All rights reserved. doi:10.1016/S0096-3003(03)00760-4
66
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
1. Introduction This paper studies the least mean-squared error (LMSE) quadratic (or second-order polynomial) estimation problem from uncertain observations when the observation equation, besides the multiplicative noise component, is also affected by white plus coloured noises. Systems with uncertain observations include a multiplicative noise component in the observation equation, described by a sequence of Bernoulli random variables, whose values––one or zero––indicate the presence or absence of the signal in the observation, respectively. Systems with uncertain observations have been widely studied since they constitute an appropriate model for the analysis of situations in which the observation may not contain the signal to be estimated and, hence, it consists only of noise; for example, situations of fading or reflection of the transmitted signals from the ionosphere and, also, certain cases involving sampling, gating or amplitude modulation. Because of the multiplicative noise component in the observation equation, these systems are non-Gaussian and hence, the LMSE estimator is not a linear function of the observations and, generally, it is not easy to be obtained. The difficulty in obtaining the optimal estimator has motivated the necessity of looking for suboptimal estimators which are easier to obtain, such as the linear estimators or, even, polynomial estimators which improve the widely used linear ones. Besides these cases, there exists a considerable number of situations in which the widely used assumption of Gaussian noises must be removed in order to obtain a more realistic statistical description of the random processes involved. De Santis et al. [1] obtained the optimal quadratic filter to estimate signals in non-Gaussian systems using a state-space approach; this filter provides the optimal quadratic solution for the non-Gaussian finite-horizon regulator problem [2]. Carravetta et al. [3,4] generalized the study of De Santis et al. [1], considering arbitrary-order polynomial estimators; this theory was also used by Dalla Mora et al. [5] for the restoration of images corrupted by additive non-Gaussian noise. Also, in systems with uncertain observations, considering different hypotheses on the noises, lineal and polynomial filtering algorithms have been derived by some authors [6–8]. All the above mentioned papers on polynomial estimation in non-Gaussian systems use the full knowledge of the state-space model. However, as it is well known, this model can be unavailable in some practical situations. The linear estimation problem in systems with uncertain observations using as information the covariance functions of the processes involved has been considered for the case of white observation noise [9] and white plus coloured observation noises [10]. Recently, in [11], the results of [9] are generalized by obtaining algorithms for the quadratic estimation problem.
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
67
In this paper, on the one hand, we generalize the results in [10] by obtaining algorithms for the LMSE quadratic estimators and, on the other, the study of [11] is extended for the case of white plus coloured noises. The problem is approached by using the same technique as that in [11], which consists of augmenting the signal and observation vectors by aggregating them with their second-order powers. Thus, the LMSE linear estimator of the augmented signal based on the augmented observations provides the LMSE quadratic estimator for the original signal. The problem is then reduced to find the linear estimator for the augmented signal. It is proven that the augmented observation satisfies the necessary conditions to apply the linear algorithms proposed in [10]; hence, these algorithms allow us to obtain the linear estimators of the augmented signal and, from them, the quadratic estimators for the original signal are deduced, as we have already indicated. The LMSE second-order polynomial estimation problem and the hypotheses on the model are formulated in Section 2. The augmented observation equation and the statistical properties of the augmented processes are analyzed in Section 3. Also, in this section, the filtering and fixed-point smoothing algorithms are established. Finally, in Section 4, these algorithms are applied to the estimation of a signal generated by a first-order autoregressive model, showing the effectiveness of the quadratic estimators in comparison with the linear ones.
2. System model and problem formulation Let zðkÞ and yðkÞ be n 1 vectors which describe the signal that we wish to estimate and the observation of this signal at time k, respectively. Let us suppose that the observations of the signal are affected by a multiplicative noise, fuðkÞ; k P 0g, and by additive white and coloured noises, fvðkÞ; k P 0g and fv0 ðkÞ; k P 0g, respectively; so, the observation equation is given by yðkÞ ¼ uðkÞzðkÞ þ vðkÞ þ v0 ðkÞ:
ð1Þ
For each instant of time k, uðkÞ is a random variable which takes on the values 0 or 1 with P ½uðkÞ ¼ 1 ¼ pðkÞ; hence the observation yðkÞ may not contain the signal ðuðkÞ ¼ 0Þ, in which case it will only consist of noise. The probability 1 pðkÞ that the observation at time k is only noise is named false alarm probability. Our aim is to obtain the least mean-squared error quadratic (or second-order polynomial) estimator, ^zðk; LÞ, of the signal zðkÞ based on the observations yð1Þ; . . . ; yðLÞ (L P k). By defining the random vectors y ½2 ðiÞ ¼ yðiÞ yðiÞ, where denotes the Kronecker product [12], and by assuming that E½y ½2 T ðiÞy ½2 ðiÞ < 1, the required estimator is the orthogonal projection of zðkÞ
68
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
on the space of n-dimensional linear transformations of yð1Þ; . . . ; yðLÞ and their second-order powers y ½2 ð1Þ; . . . ; y ½2 ðLÞ. To treat the LMSE quadratic estimation problem, we assume the following hypotheses on the signal and the noise processes involved in Eq. (1): (H.1) The signal process fzðkÞ; k P 0g has zero mean and its autocovariance function, Kz ðk; sÞ ¼ E½zðkÞzT ðsÞ , as well as the autocovariance function of their second-order powers, Kz½2 ðk; sÞ ¼ E½ðz½2 ðkÞ E½z½2 ðkÞ Þðz½2 ðsÞ
T E½z½2 ðsÞ Þ , are expressed in a semi-degenerate kernel form, namely, AðkÞBT ðsÞ; 0 6 s 6 k Kz ðk; sÞ ¼ BðkÞAT ðsÞ; 0 6 k 6 s aðkÞbT ðsÞ; 0 6 s 6 k Kz½2 ðk; sÞ ¼ bðkÞaT ðsÞ; 0 6 k 6 s where the n M 0 matrix functions A, B and the n2 L0 matrix functions a, b are known. Moreover, let us suppose that the covariance function of the signal and their second-order powers, Kzz½2 ðk; sÞ ¼ E½zðkÞz½2 T ðsÞ , can be also expressed as c1 ðkÞcT2 ðsÞ; 0 6 s 6 k Kzz½2 ðk; sÞ ¼ d1 ðkÞd2T ðsÞ; 0 6 k 6 s where c1 , c2 , d1 and d2 are n N 0 , n2 N 0 , n P 0 and n2 P 0 known matrix functions, respectively. (H.2) The noise process fvðkÞ; k P 0g is a zero-mean white sequence and its moments, up to the fourth one, are also known, and denoted as follows Rv ðkÞ ¼ E½vðkÞvT ðkÞ ;
Rvv½2 ðkÞ ¼ E½vðkÞv½2 T ðkÞ ;
Rv½2 ðkÞ ¼ E½ðv½2 ðkÞ E½v½2 ðkÞ Þðv½2 ðkÞ E½v½2 ðkÞ ÞT : (H.3) The coloured noise fv0 ðkÞ; k P 0g has zero mean and its autocovariance function, Kv0 ðk; sÞ ¼ E½v0 ðkÞvT0 ðsÞ , as well as the autocovariance function ½2 ½2 ½2 of their second-order powers, Kv½2 ðk; sÞ ¼ E½ðv0 ðkÞ E½v0 ðkÞ Þðv0 ðsÞ
0 ½2 E½v0 ðsÞ ÞT , are expressed in a semi-degenerate kernel form, aðkÞbT ðsÞ; 0 6 s 6 k Kv0 ðk; sÞ ¼ bðkÞaT ðsÞ; 0 6 k 6 s cðkÞdT ðsÞ; 0 6 s 6 k Kv½2 ðk; sÞ ¼ 0 dðkÞcT ðsÞ; 0 6 k 6 s where the n M 00 matrix functions a, b and the n2 L00 matrix functions c, d are known.
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
69
Also, let us suppose that the covariance function of the coloured noise ½2 T and their second-order powers, Kv v½2 ðk; sÞ ¼ E½v0 ðkÞv0 ðsÞ , can be ex0 0 pressed in a similar way, namely, 1 ðkÞT2 ðsÞ; 0 6 s 6 k Kv v½2 ðk; sÞ ¼ 0 0 q1 ðkÞqT2 ðsÞ; 0 6 k 6 s where 1 , 2 , q1 and q2 are n N 00 , n2 N 00 , n P 00 and n2 P 00 known matrix functions, respectively. (H.4) The multiplicative noise fuðkÞ; k P 0g is a sequence of independent Bernoulli random variables and let us suppose that the probability of existence of the signal into each observation, P ½uðkÞ ¼ 1 ¼ pðkÞ, is available. (H.5) The signal process, fzðkÞ; k P 0g, and the noise processes, fuðkÞ; k P 0g, fvðkÞ; k P 0g and fv0 ðkÞ; k P 0g, are mutually independent.
3. LMSE quadratic estimation problem In order to treat the LMSE quadratic estimation problem, let us define the augmented signal and observation vectors by aggregating to the original vectors zðkÞ and yðkÞ their second-order powers z½2 ðkÞ and y ½2 ðkÞ, that is, zðkÞ yðkÞ ZðkÞ ¼ ; YðkÞ ¼ : z½2 ðkÞ y ½2 ðkÞ Then, the vector constituted by the first n entries of the LMSE linear estimator of the augmented signal, ZðkÞ, based on the augmented observations, Yð1Þ; . . . ; YðLÞ, provides the LMSE quadratic estimator of the original signal zðkÞ. To obtain that estimator, we first analyze the properties of the random vectors ZðkÞ and YðkÞ. Previously, taking into account the Kronecker product properties and since uðkÞ ¼ u2 ðkÞ, we obtain that the vector YðkÞ can be expressed as a function of ZðkÞ similar to that Eq. (1), namely, YðkÞ ¼ uðkÞZðkÞ þ VðkÞ þ V0 ðkÞ where VðkÞ ¼
vðkÞ ; f ðkÞ
V0 ðkÞ ¼
v0 ðkÞ f0 ðkÞ
with f ðkÞ ¼ uðkÞðIn2 þ Kn2 ÞðzðkÞ vðkÞÞ þ ðIn2 þ Kn2 Þðv0 ðkÞ vðkÞÞ þ ðuðkÞ pðkÞÞðIn2 þ Kn2 ÞðzðkÞ v0 ðkÞÞ þ v½2 ðkÞ
ð2Þ
70
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
and ½2
f0 ðkÞ ¼ pðkÞðIn2 þ Kn2 ÞðzðkÞ v0 ðkÞÞ þ v0 ðkÞ
ð3Þ
being In2 the n2 n2 identity matrix and Kn2 the n2 n2 commutation matrix satisfying Kn2 ðzðkÞ vðkÞÞ ¼ vðkÞ zðkÞ. Since the augmented signal, ZðkÞ, and the noises, VðkÞ and V0 ðkÞ, have no zero mean, taking into account that E½uðkÞZðkÞ ¼ pðkÞE½ZðkÞ and by denoting ZðkÞ ¼ ZðkÞ E½ZðkÞ , V ðkÞ ¼ VðkÞ E½VðkÞ þ ½uðkÞ pðkÞ E½ZðkÞ and V0 ðkÞ ¼ V0 ðkÞ E½V0 ðkÞ , it follows that Y ðkÞ ¼ YðkÞ E½YðkÞ satisfy the following equation Y ðkÞ ¼ uðkÞZðkÞ þ V ðkÞ þ V0 ðkÞ
ð4Þ
in which the signal, ZðkÞ, and the noises, V ðkÞ and V0 ðkÞ, have zero mean. The remaining statistical properties of these processes are established in the following propositions. Proposition 1. Let us suppose that the hypotheses (H.1)–(H.5) are satisfied. Then the process fZðkÞ; k P 0g has zero mean and its covariance function is expressed in a semi-degenerate kernel form, specifically AðkÞBT ðsÞ; 0 6 s 6 k T KZ ðk; sÞ ¼ E½ZðkÞZ ðsÞ ¼ BðkÞAT ðsÞ; 0 6 k 6 s being AðkÞ ¼ BðkÞ ¼
AðkÞ 0n2 M 0
c1 ðkÞ 0n2 N 0
BðkÞ 0n2 M 0
0nN 0 c2 ðkÞ
0nP 0 d2 ðkÞ d1 ðkÞ 0n2 P 0
0nL0 : aðkÞ 0nL0 : bðkÞ
Moreover, the process fZðkÞ; k P 0g is independent of the noise fuðkÞ; k P 0g. Proof. Nakamori et al. [11].
h
Proposition 2. Under the hypotheses (H.1)–(H.5), the noise process fV ðkÞ; k P 0g of Eq. (4) is a sequence of zero-mean, mutually uncorrelated random vectors with covariance matrices given by RV ðkÞ ¼ E½V ðkÞV T ðkÞ ¼ pðkÞð1 pðkÞÞE½ZðkÞ E½ZT ðkÞ þ
Rv ðkÞ RTvv½2 ðkÞ
Rvv½2 ðkÞ R22 ðkÞ
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
71
being E½ZðkÞ ¼
0 vecð AðkÞBT ðkÞÞ
ð5Þ
and R22 ðkÞ ¼ pðkÞðIn2 þ Kn2 ÞðAðkÞBT ðkÞ Rv ðkÞÞðIn2 þ Kn2 Þ þ ðIn2 þ Kn2 ÞðaðkÞbT ðkÞ Rv ðkÞÞðIn2 þ Kn2 Þ þ pðkÞð1 pðkÞÞ ðIn2 þ Kn2 ÞðAðkÞBT ðkÞ aðkÞbT ðkÞÞðIn2 þ Kn2 Þ þ Rv½2 ðkÞ:
ð6Þ
Moreover, fV ðkÞ; k P 0g is uncorrelated with the processes fZðkÞ; k P 0g and fuðkÞZðkÞ; k P 0g. Remark. vecðÞ denote the vec or stack operator [12]. Proof. Using the independence hypotheses on the model it is easy to see that, for k 6¼ s, E½V ðkÞV T ðsÞ ¼ 0, and, for k ¼ s, Rvv½2 ðkÞ Rv ðkÞ RV ðkÞ ¼ pðkÞð1 pðkÞÞE½ZðkÞ E½ZT ðkÞ þ RTvv½2 ðkÞ R22 ðkÞ T
with R22 ðkÞ ¼ E½ðf ðkÞ E½f ðkÞ Þðf ðkÞ E½f ðkÞ Þ . Since u v ¼ vecðvuT Þ, using hypothesis (H.1), it is immediate that E½ZðkÞ verifies (5). On the other hand, using the Kronecker product properties, the hypotheses on the model and Eq. (2), we obtain that R22 ðkÞ verifies (6). Finally, the uncorrelation between fV ðkÞ; k P 0g and the processes fZðkÞ; k P 0g, fuðkÞZðkÞ; k P 0g is derived in a similar way, by using hypotheses (H.1)–(H.5) and employing the Kronecker product properties. h Proposition 3. Let us suppose that the hypotheses (H.1)–(H.5) are satisfied. Then the process fV0 ðkÞ; k P 0g of Eq. (4) is a zero-mean coloured noise and its covariance function is expressed in a semi-degenerate kernel form, specifically CðkÞDT ðsÞ; 0 6 s 6 k KV0 ðk; sÞ ¼ E½V0 ðkÞV0T ðsÞ ¼ ð7Þ DðkÞCT ðsÞ; 0 6 k 6 s being CðkÞ ¼ DðkÞ ¼
aðkÞ 0n2 M 00
1 ðkÞ 0n2 N 00
0nP 00 0nL00 q2 ðkÞ cðkÞ
bðkÞ 0n2 M 00
0nN 00 2 ðkÞ
q1 ðkÞ 0n2 P 00
0nL00 dðkÞ
0nM 0 M 00 pðkÞðIn2 þ Kn2 ÞðAðkÞ aðkÞÞ
0nM 0 M 00 : pðkÞðIn2 þ Kn2 ÞðBðkÞ bðkÞÞ
Moreover, fV0 ðkÞ; k P 0g is uncorrelated with the processes fV ðkÞ; k P 0g, fZðkÞ; k P 0g and fuðkÞZðkÞ; k P 0g.
72
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
Proof. Expression (7) for KV0 ðk; sÞ is obtained taking into account that, from hypothesis (H.3), for 0 6 s 6 k, KV0 ðk; sÞ aðkÞbT ðsÞ ¼ q1 ðkÞqT2 ðsÞ
1 ðkÞT2 ðsÞ 2 p ðkÞðIn2 þ Kn2 ÞðAðkÞBT ðsÞ aðkÞbT ðsÞÞðIn2 þ Kn2 Þ
and, for 0 6 k 6 s, KV0 ðk; sÞ bðkÞaT ðsÞ ¼ q2 ðkÞqT1 ðsÞ
2 ðkÞT1 ðsÞ : p2 ðkÞðIn2 þ Kn2 ÞðBðkÞAT ðsÞ bðkÞaT ðsÞÞðIn2 þ Kn2 Þ
Finally, as in Proposition 2, the uncorrelation between fV0 ðkÞ; k P 0g and the processes fV ðkÞ; k P 0g, fZðkÞ; k P 0g and fuðkÞZðkÞ; k P 0g is derived by using hypotheses (H.1)–(H.5) and employing the Kronecker product properties. h In view of the properties of the processes involved in Eq. (4), which have been established in Propositions 1–3, the recursive algorithms given in [10] can be applied to obtain the linear filtering and fixed-point smoothing estimators, b ðk; LÞ, L P k, of the signal ZðkÞ based on the observations Y ð1Þ; . . . ; Y ðLÞ. Z These algorithms, which are presented in the following theorem, allow us to obtain the required quadratic filtering and fixed-point smoothing estimators of b ðk; LÞ. the original signal zðkÞ, just by extracting the first n entries of Z Theorem 1. The filtering and fixed-point smoothing algorithm of the augmented signal ZðkÞ based on the augmented observations Y ð1Þ; . . . ; Y ðLÞ, L P k, is given by b ðk; LÞ ¼ Z b ðk; L 1Þ þ hðk; L; LÞmðLÞ Z where mðLÞ, the innovation, is mðLÞ ¼ Y ðLÞ pðLÞAðLÞOðL 1Þ CðLÞQðL 1Þ: The M 0 1 and N 0 1 vectors OðLÞ and QðLÞ, respectively, are recursively calculated by OðLÞ ¼ OðL 1Þ þ J ðL; LÞmðLÞ;
Oð0Þ ¼ 0
QðLÞ ¼ QðL 1Þ þ IðL; LÞmðLÞ;
Qð0Þ ¼ 0
being J ðL; LÞ ¼ ½pðLÞðBT ðLÞ rðL 1ÞAT ðLÞÞ cðL 1ÞCT ðLÞ P 1 ðLÞ
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
73
IðL; LÞ ¼ ½DT ðLÞ pðLÞcT ðL 1ÞAT ðLÞ dðL 1ÞCT ðLÞ P 1 ðLÞ where PðLÞ, the covariance matrix of the innovation, is given by PðLÞ ¼ RV ðLÞ þ pðLÞ½BðLÞ pðLÞAðLÞrðL 1Þ CðLÞcT ðL 1Þ AT ðLÞ þ ½DðLÞ pðLÞAðLÞcðL 1Þ CðLÞdðL 1Þ CT ðLÞ: The functions r, c and d are M 0 M 0 , M 0 N 0 and N 0 N 0 matrices, respectively, verifying rðLÞ ¼ rðL 1Þ þ J ðL; LÞ½pðLÞðBðLÞ AðLÞrðL 1ÞÞ
CðLÞcT ðL 1Þ ;
rð0Þ ¼ 0
cðLÞ ¼ cðL 1Þ þ J ðL; LÞ½DðLÞ pðLÞAðLÞcðL 1Þ
CðLÞdðL 1Þ ;
cð0Þ ¼ 0
dðLÞ ¼ dðL 1Þ þ IðL; LÞ½DðLÞ pðLÞAðLÞcðL 1Þ
CðLÞdðL 1Þ ;
dð0Þ ¼ 0:
The smoothing gain, hðk; L; LÞ, is given by hðk; L; LÞ ¼ ½pðLÞðBðkÞAT ðLÞ Eðk; L 1ÞAT ðLÞÞ
F ðk; L 1ÞCT ðLÞ P 1 ðLÞ where Eðk; LÞ and F ðk; LÞ are n M 0 and n N 0 matrices satisfying Eðk; LÞ ¼ Eðk; L 1Þ þ hðk; L; LÞ½pðLÞðBðLÞ AðLÞrðL 1ÞÞ
CðLÞcT ðL 1Þ Eðk; kÞ ¼ AðkÞrðkÞ F ðk; LÞ ¼ F ðk; L 1Þ þ hðk; L; LÞ½DðLÞ pðLÞAðLÞcðL 1Þ
CðLÞdðL 1Þ F ðk; kÞ ¼ AðkÞcðkÞ: b ðk; kÞ, which provides the initial condition for the fixedThe filtering estimate, Z b ðk; kÞ ¼ AðkÞOðkÞ. point smoothing algorithm, is given by Z The fixed-point smoothing and filtering error covariance matrices, P ðk; LÞ, L P k, satisfy P ðk; LÞ ¼ P ðk; L 1Þ hðk; L; LÞPðLÞhT ðk; L; LÞ; P ðk; kÞ ¼ AðkÞ½BT ðkÞ rðkÞAT ðkÞ :
L>k
74
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
4. A numerical simulation example In order to show the effectiveness of the proposed quadratic estimators, we have performed a program in MATLAB, which simulates the signal value at each iteration, and provides the linear and quadratic estimates, as well as the corresponding error covariance matrices. This program has been applied to a scalar signal fzðkÞ; k P 0g which, as in [11], is generated by the following first-order autoregressive model zðk þ 1Þ ¼ 0:95zðkÞ þ wðkÞ where fwðkÞ; k P 0g is a zero-mean white Gaussian noise with Var½wðkÞ ¼ 0:1; 8k P 0. The autocovariance functions of the signal and their second-order powers are given in a semi-degenerate kernel form, specifically, Kz ðk; sÞ ¼ 1:025641 0:95k s ; Kz2 ðk; sÞ ¼ 2:1038795 0:952ðk sÞ ;
06s6k
and the covariance function of the signal and their second-order powers is given by Kzz2 ðk; sÞ ¼ 0;
8s; k
According to hypothesis (H.1), the functions which constitute these covariance functions are as follows AðkÞ ¼ 1:025641 0:95k ; bðkÞ ¼ 0:95 2k ;
BðkÞ ¼ 0:95 k ;
aðkÞ ¼ 2:1038795 0:952k ;
c1 ðkÞ ¼ c2 ðkÞ ¼ d1 ðkÞ ¼ d2 ðkÞ ¼ 0:
Let us consider that the observations of the signal are given by yðkÞ ¼ uðkÞzðkÞ þ vðkÞ þ v0 ðkÞ where fuðkÞ; k P 0g is a sequence of independent Bernoulli random variables with P ½uðkÞ ¼ 1 ¼ p;
8k P 0
The noise fvðkÞ; k P 0g is a sequence of independent random variables with P ½vðkÞ ¼ 8 ¼ 18; P vðkÞ ¼ 87 ¼ 78; 8k P 0 and so, E½vðkÞ ¼ 0;
Rv ðkÞ ¼ 9:142857;
Rvv2 ðkÞ ¼ 62:693878;
Rv2 ðkÞ ¼ 429:900875: The coloured noise fv0 ðkÞ; k P 0g is generated by the following first-order autoregressive model
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
75
v0 ðk þ 1Þ ¼ 0:5v0 ðkÞ þ v1 ðkÞ where fv1 ðkÞ; k P 0g is a zero-mean white Gaussian noise with Var½v1 ðkÞ ¼ 0:075; 8k P 0. The autocovariance functions of the coloured noise and their second-order powers are given by in a semi-degenerate kernel form, namely, Kv0 ðk; sÞ ¼ 0:1 0:5k s ;
Kv2 ðk; sÞ ¼ 0:02 0:52ðk sÞ ; 0
06s6k
and the covariance function of the signal and their second-order powers is given by Kv0 v2 ðk; sÞ ¼ 0; 0
8s; k:
Hence, according to hypothesis (H.3), the functions which constitute these covariance functions are as follows aðkÞ ¼ 0:1 0:5k ; dðkÞ ¼ 0:5 2k ;
bðkÞ ¼ 0:5 k ;
cðkÞ ¼ 0:02 0:52k ;
1 ðkÞ ¼ 2 ðkÞ ¼ q1 ðkÞ ¼ q2 ðkÞ ¼ 0:
Finally, in accordance with the hypothesis (H.5) imposed in the theoretic study that we have developed, it is also assumed that the processes fzðkÞ; k P 0g, fuðkÞ; k P 0g, fvðkÞ; k P 0g and fv0 ðkÞ; k P 0g are mutually independent. In order to compare the linear (see [10]) and quadratic estimates, we have performed 200 iterations of the respective algorithms, considering different values of the false alarm probability 1 p, specifically, 1 p ¼ 0:5, 1 p ¼ 0:25 and 1 p ¼ 0, case in which the signal is always present in the observations. Fig. 1 displays the linear and quadratic filtering error variances, which measure the performance of the filters; this figure shows that, for all the values of the false alarm probability, the error variances corresponding to the quadratic filter are less than the linear filter ones; also we observe that as the false alarm probability is smaller, the error variances are also smaller and, consequently, the performance of the filter is better, as one could expect intuitively. The simulated signal, the linear and the quadratic filtering estimates for the values 1 p ¼ 0:5, 1 p ¼ 0:25 and 1 p ¼ 0, are displayed in Figs. 2–4, respectively. These figures show that the quadratic filtering estimate follows the signal evolution better and, also, that as the false alarm probability, 1 p, is smaller both the linear and the quadratic filtering estimates follow the signal evolution better, according to the results displayed in Fig. 1. Finally, Tables 1 and 2 summarize the mean-square values (MSVs) of the filtering and fixed-point smoothing errors, respectively, for the different values P of the false alarm probability. The MSVs are calculated by 200 zði; iÞÞ2 = i¼1 ðzðiÞ ^ P200 P10 2 200, for the filtering errors, and by i¼1 j¼1 ðzðiÞ ^zði; i þ jÞÞ =2000, for the fixed-point smoothing errors. Tables 1 and 2 show, on the one hand, that the
76
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79 1
0.9 (a)
Variances of Filtering Error
0.8 (b)
0.7 (c)
0.6
0.5
(a’)
0.4 (b’)
0.3
(c’)
0.2 0
5
10
15
20
25
30
Time k
Fig. 1. Linear filtering error variances (solid line) and quadratic filtering error variances (dashed line) when (a), (a0 ): 1 p ¼ 0:5, (b), (b0 ): 1 p ¼ 0:25, (c), (c0 ): 1 p ¼ 0.
4 Signal Linear Filtering Estimate Quadratic Filtering Estimate
Signal and Filtering Estimates
3
2
1
0
–1
–2
0
20
40
60
80
100
120
140
160
180
200
Time k
Fig. 2. Signal (solid line), linear filtering estimate (dashed line) and quadratic filtering estimate (dotted line) for 1 p ¼ 0:5.
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
77
4 Signal Linear Filtering Estimate Quadratic Filtering Estimate
Signal and Filtering Estimates
3
2
1
0
–1
–2 0
20
40
60
80
100
120
140
160
180
200
Time k
Fig. 3. Signal (solid line), linear filtering estimate (dashed line) and quadratic filtering estimate (dotted line) for 1 p ¼ 0:25.
5 Signal Linear Filtering Estimate Quadratic Filtering Estimate
Signal and Filtering Estimates
4
3
2
1
0
–1
–2
0
20
40
60
80
100
120
140
160
180
200
Time k
Fig. 4. Signal (solid line), linear filtering estimate (dashed line) and quadratic filtering estimate (dotted line) for 1 p ¼ 0.
78
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
Table 1 Mean-square values of the linear filtering errors and the quadratic filtering errors for different values of the false alarm probability False alarm probability
MSVs of linear filtering errors
MSVs of quadratic filtering errors
1 p ¼ 0:50 1 p ¼ 0:25 1 p ¼0
0.6967 0.5046 0.3939
0.3069 0.2315 0.1527
Table 2 Mean-square values of the linear fixed-point smoothing errors and the quadratic fixed-point smoothing errors for different values of the false alarm probability False alarm probability
MSVs of linear smoothing errors
MSVs of quadratic smoothing errors
1 p ¼ 0:50 1 p ¼ 0:25 1 p ¼0
0.5832 0.4421 0.3503
0.2543 0.1697 0.1180
estimation accuracy of the smoother is superior to that of the filter and, on the other, not only that the MSVs of the quadratic filtering and smoothing errors are less than the MSVs of the linear ones, but also that as the false alarm probability is smaller, the MSVs are smaller too, agreeing again with the theoretical results showed in Fig. 1.
5. Conclusions In this paper, using covariance information, recursive LMSE algorithms for the quadratic filter and fixed-point smoother are derived from uncertain observations, for the case of white plus coloured observation noises. The proposed estimators do not require the knowledge of state-space model of the signal; the available information is only the autocovariance and crosscovariance functions of the signal and their second-order powers, defined by the Kronecker product, as well as the corresponding functions of the additive noises. It is also assumed that the probability that the signal exists in the observed values is known. By defining a suitably augmented observation equation, the LMSE linear estimator of the augmented signal based on the augmented observations provides the LMSE quadratic estimator for the original signal. The effectiveness of the quadratic estimators in contrast to the linear ones is shown by applying the filtering and fixed-point algorithms to the estimation of
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
79
a signal generated by a first-order autoregressive model, from uncertain observations with different values of the false alarm probability.
Acknowledgements This work has been partially supported by the ‘‘Ministerio de Ciencia y Tecnologıa’’ under contract BFM2002-00932.
References [1] A. De Santis, A. Germani, M. Raimondi, Optimal quadratic filtering of linear discrete-time non-Gaussian systems, IEEE Transactions on Automatic Control AC-40 (1995) 1274–1278. [2] A. Germani, G. Mavelli, Optimal quadratic solution for the non-Gaussian finite-horizon regulator problem, Systems and Control Letters 38 (1999) 321–331. [3] F. Carravetta, A. Germani, M. Raimondi, Polynomial filtering for linear discrete time nonGaussian systems, SIAM Journal of Control and Optimization 34 (1996) 1666–1690. [4] F. Carravetta, A. Germani, M. Raimondi, Polynomial filtering of discrete-time stochastic linear systems with multiplicative state noise, IEEE Transactions on Automatic Control AC-42 (1997) 1106–1126. [5] M. Dalla Mora, A. Germani, A. Nardecchia, Restoration of images corrupted by additive nonGaussian noise, IEEE Transactions on Circuits and Systems I: Fundamental Theory and Applications 48 (2001) 859–875. [6] A. Hermoso, J. Linares, Linear smoothing for discrete-time systems in the presence of correlated disturbances and uncertain observations, IEEE Transactions on Automatic Control 40 (1995) 1486–1488. [7] W. NaNakara, E.E. Yaz, Recursive estimator for linear and nonlinear systems with uncertain observations, Signal Processing 62 (1997) 215–228. [8] R. Caballero, A. Hermoso, J. Linares, Polynomial filtering with uncertain observations in stochastic linear systems, International Journal of Modelling and Simulation 23 (2003) 22–28. [9] S. Nakamori, R. Caballero, A. Hermoso, J. Linares, New design of estimators using covariance information with uncertain observations in linear discrete-time systems, Applied Mathematics and Computation 135 (2003) 429–441. [10] S. Nakamori, R. Caballero, A. Hermoso, J. Linares, Linear estimation from uncertain observations with white plus coloured noises using covariance information, Digital Signal Processing 13 (3) (2003) 552–568. [11] S. Nakamori, R. Caballero, A. Hermoso, J. Linares, Second-order polynomial estimators from uncertain observations using covariance information, Applied Mathematics and Computation 143 (2003) 319–338. [12] J.R. Magnus, H. Neudecker, Matrix Differential Calculus with Applications in Statistics and Econometrics, Wiley & Sons, New York, 1999.
Quadratic estimation from uncertain observations with white plus coloured noises using covariance information b S. Nakamori a,*, R. Caballero-Aguila , c A. Hermoso-Carazo , J. Linares-Perez c a
Department of Technology, Faculty of Education, Kagoshima University, 1-20-6 Kohrimoto, Kagoshima 890-0065, Japan b Departamento de Estadıstica e Investigaci on Operativa, Universidad de Ja en, Paraje Las Lagunillas, s/n, 23071 Ja en, Spain c Departamento de Estadıstica e Investigaci on Operativa, Universidad de Granada, Campus Fuentenueva, s/n, 18071 Granada, Spain
Abstract In this paper recursive least mean-squared error quadratic filtering and fixed-point smoothing algorithms to estimate signals from uncertain observations are obtained for the case of white plus coloured observation noises. It is assumed that the state-space model of the signal is not known and only the information on the moments, up to the fourth one, of the involved processes and the probability that the signal exists in the observations are available. The estimators require the covariance functions of the signal and coloured noise, as well as the covariance functions of its second-order powers in a semi-degenerate kernel form. Ó 2003 Elsevier Inc. All rights reserved. Keywords: Uncertain observation; Quadratic estimation; Covariance information; False alarm probability; Coloured noise
*
Corresponding author. E-mail addresses: [email protected] (S. Nakamori), [email protected] (R. Caballero-Aguila), [email protected] (A. Hermoso-Carazo), [email protected] (J. Linares-Perez). 0096-3003/$ - see front matter Ó 2003 Elsevier Inc. All rights reserved. doi:10.1016/S0096-3003(03)00760-4
66
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
1. Introduction This paper studies the least mean-squared error (LMSE) quadratic (or second-order polynomial) estimation problem from uncertain observations when the observation equation, besides the multiplicative noise component, is also affected by white plus coloured noises. Systems with uncertain observations include a multiplicative noise component in the observation equation, described by a sequence of Bernoulli random variables, whose values––one or zero––indicate the presence or absence of the signal in the observation, respectively. Systems with uncertain observations have been widely studied since they constitute an appropriate model for the analysis of situations in which the observation may not contain the signal to be estimated and, hence, it consists only of noise; for example, situations of fading or reflection of the transmitted signals from the ionosphere and, also, certain cases involving sampling, gating or amplitude modulation. Because of the multiplicative noise component in the observation equation, these systems are non-Gaussian and hence, the LMSE estimator is not a linear function of the observations and, generally, it is not easy to be obtained. The difficulty in obtaining the optimal estimator has motivated the necessity of looking for suboptimal estimators which are easier to obtain, such as the linear estimators or, even, polynomial estimators which improve the widely used linear ones. Besides these cases, there exists a considerable number of situations in which the widely used assumption of Gaussian noises must be removed in order to obtain a more realistic statistical description of the random processes involved. De Santis et al. [1] obtained the optimal quadratic filter to estimate signals in non-Gaussian systems using a state-space approach; this filter provides the optimal quadratic solution for the non-Gaussian finite-horizon regulator problem [2]. Carravetta et al. [3,4] generalized the study of De Santis et al. [1], considering arbitrary-order polynomial estimators; this theory was also used by Dalla Mora et al. [5] for the restoration of images corrupted by additive non-Gaussian noise. Also, in systems with uncertain observations, considering different hypotheses on the noises, lineal and polynomial filtering algorithms have been derived by some authors [6–8]. All the above mentioned papers on polynomial estimation in non-Gaussian systems use the full knowledge of the state-space model. However, as it is well known, this model can be unavailable in some practical situations. The linear estimation problem in systems with uncertain observations using as information the covariance functions of the processes involved has been considered for the case of white observation noise [9] and white plus coloured observation noises [10]. Recently, in [11], the results of [9] are generalized by obtaining algorithms for the quadratic estimation problem.
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
67
In this paper, on the one hand, we generalize the results in [10] by obtaining algorithms for the LMSE quadratic estimators and, on the other, the study of [11] is extended for the case of white plus coloured noises. The problem is approached by using the same technique as that in [11], which consists of augmenting the signal and observation vectors by aggregating them with their second-order powers. Thus, the LMSE linear estimator of the augmented signal based on the augmented observations provides the LMSE quadratic estimator for the original signal. The problem is then reduced to find the linear estimator for the augmented signal. It is proven that the augmented observation satisfies the necessary conditions to apply the linear algorithms proposed in [10]; hence, these algorithms allow us to obtain the linear estimators of the augmented signal and, from them, the quadratic estimators for the original signal are deduced, as we have already indicated. The LMSE second-order polynomial estimation problem and the hypotheses on the model are formulated in Section 2. The augmented observation equation and the statistical properties of the augmented processes are analyzed in Section 3. Also, in this section, the filtering and fixed-point smoothing algorithms are established. Finally, in Section 4, these algorithms are applied to the estimation of a signal generated by a first-order autoregressive model, showing the effectiveness of the quadratic estimators in comparison with the linear ones.
2. System model and problem formulation Let zðkÞ and yðkÞ be n 1 vectors which describe the signal that we wish to estimate and the observation of this signal at time k, respectively. Let us suppose that the observations of the signal are affected by a multiplicative noise, fuðkÞ; k P 0g, and by additive white and coloured noises, fvðkÞ; k P 0g and fv0 ðkÞ; k P 0g, respectively; so, the observation equation is given by yðkÞ ¼ uðkÞzðkÞ þ vðkÞ þ v0 ðkÞ:
ð1Þ
For each instant of time k, uðkÞ is a random variable which takes on the values 0 or 1 with P ½uðkÞ ¼ 1 ¼ pðkÞ; hence the observation yðkÞ may not contain the signal ðuðkÞ ¼ 0Þ, in which case it will only consist of noise. The probability 1 pðkÞ that the observation at time k is only noise is named false alarm probability. Our aim is to obtain the least mean-squared error quadratic (or second-order polynomial) estimator, ^zðk; LÞ, of the signal zðkÞ based on the observations yð1Þ; . . . ; yðLÞ (L P k). By defining the random vectors y ½2 ðiÞ ¼ yðiÞ yðiÞ, where denotes the Kronecker product [12], and by assuming that E½y ½2 T ðiÞy ½2 ðiÞ < 1, the required estimator is the orthogonal projection of zðkÞ
68
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
on the space of n-dimensional linear transformations of yð1Þ; . . . ; yðLÞ and their second-order powers y ½2 ð1Þ; . . . ; y ½2 ðLÞ. To treat the LMSE quadratic estimation problem, we assume the following hypotheses on the signal and the noise processes involved in Eq. (1): (H.1) The signal process fzðkÞ; k P 0g has zero mean and its autocovariance function, Kz ðk; sÞ ¼ E½zðkÞzT ðsÞ , as well as the autocovariance function of their second-order powers, Kz½2 ðk; sÞ ¼ E½ðz½2 ðkÞ E½z½2 ðkÞ Þðz½2 ðsÞ
T E½z½2 ðsÞ Þ , are expressed in a semi-degenerate kernel form, namely, AðkÞBT ðsÞ; 0 6 s 6 k Kz ðk; sÞ ¼ BðkÞAT ðsÞ; 0 6 k 6 s aðkÞbT ðsÞ; 0 6 s 6 k Kz½2 ðk; sÞ ¼ bðkÞaT ðsÞ; 0 6 k 6 s where the n M 0 matrix functions A, B and the n2 L0 matrix functions a, b are known. Moreover, let us suppose that the covariance function of the signal and their second-order powers, Kzz½2 ðk; sÞ ¼ E½zðkÞz½2 T ðsÞ , can be also expressed as c1 ðkÞcT2 ðsÞ; 0 6 s 6 k Kzz½2 ðk; sÞ ¼ d1 ðkÞd2T ðsÞ; 0 6 k 6 s where c1 , c2 , d1 and d2 are n N 0 , n2 N 0 , n P 0 and n2 P 0 known matrix functions, respectively. (H.2) The noise process fvðkÞ; k P 0g is a zero-mean white sequence and its moments, up to the fourth one, are also known, and denoted as follows Rv ðkÞ ¼ E½vðkÞvT ðkÞ ;
Rvv½2 ðkÞ ¼ E½vðkÞv½2 T ðkÞ ;
Rv½2 ðkÞ ¼ E½ðv½2 ðkÞ E½v½2 ðkÞ Þðv½2 ðkÞ E½v½2 ðkÞ ÞT : (H.3) The coloured noise fv0 ðkÞ; k P 0g has zero mean and its autocovariance function, Kv0 ðk; sÞ ¼ E½v0 ðkÞvT0 ðsÞ , as well as the autocovariance function ½2 ½2 ½2 of their second-order powers, Kv½2 ðk; sÞ ¼ E½ðv0 ðkÞ E½v0 ðkÞ Þðv0 ðsÞ
0 ½2 E½v0 ðsÞ ÞT , are expressed in a semi-degenerate kernel form, aðkÞbT ðsÞ; 0 6 s 6 k Kv0 ðk; sÞ ¼ bðkÞaT ðsÞ; 0 6 k 6 s cðkÞdT ðsÞ; 0 6 s 6 k Kv½2 ðk; sÞ ¼ 0 dðkÞcT ðsÞ; 0 6 k 6 s where the n M 00 matrix functions a, b and the n2 L00 matrix functions c, d are known.
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
69
Also, let us suppose that the covariance function of the coloured noise ½2 T and their second-order powers, Kv v½2 ðk; sÞ ¼ E½v0 ðkÞv0 ðsÞ , can be ex0 0 pressed in a similar way, namely, 1 ðkÞT2 ðsÞ; 0 6 s 6 k Kv v½2 ðk; sÞ ¼ 0 0 q1 ðkÞqT2 ðsÞ; 0 6 k 6 s where 1 , 2 , q1 and q2 are n N 00 , n2 N 00 , n P 00 and n2 P 00 known matrix functions, respectively. (H.4) The multiplicative noise fuðkÞ; k P 0g is a sequence of independent Bernoulli random variables and let us suppose that the probability of existence of the signal into each observation, P ½uðkÞ ¼ 1 ¼ pðkÞ, is available. (H.5) The signal process, fzðkÞ; k P 0g, and the noise processes, fuðkÞ; k P 0g, fvðkÞ; k P 0g and fv0 ðkÞ; k P 0g, are mutually independent.
3. LMSE quadratic estimation problem In order to treat the LMSE quadratic estimation problem, let us define the augmented signal and observation vectors by aggregating to the original vectors zðkÞ and yðkÞ their second-order powers z½2 ðkÞ and y ½2 ðkÞ, that is, zðkÞ yðkÞ ZðkÞ ¼ ; YðkÞ ¼ : z½2 ðkÞ y ½2 ðkÞ Then, the vector constituted by the first n entries of the LMSE linear estimator of the augmented signal, ZðkÞ, based on the augmented observations, Yð1Þ; . . . ; YðLÞ, provides the LMSE quadratic estimator of the original signal zðkÞ. To obtain that estimator, we first analyze the properties of the random vectors ZðkÞ and YðkÞ. Previously, taking into account the Kronecker product properties and since uðkÞ ¼ u2 ðkÞ, we obtain that the vector YðkÞ can be expressed as a function of ZðkÞ similar to that Eq. (1), namely, YðkÞ ¼ uðkÞZðkÞ þ VðkÞ þ V0 ðkÞ where VðkÞ ¼
vðkÞ ; f ðkÞ
V0 ðkÞ ¼
v0 ðkÞ f0 ðkÞ
with f ðkÞ ¼ uðkÞðIn2 þ Kn2 ÞðzðkÞ vðkÞÞ þ ðIn2 þ Kn2 Þðv0 ðkÞ vðkÞÞ þ ðuðkÞ pðkÞÞðIn2 þ Kn2 ÞðzðkÞ v0 ðkÞÞ þ v½2 ðkÞ
ð2Þ
70
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
and ½2
f0 ðkÞ ¼ pðkÞðIn2 þ Kn2 ÞðzðkÞ v0 ðkÞÞ þ v0 ðkÞ
ð3Þ
being In2 the n2 n2 identity matrix and Kn2 the n2 n2 commutation matrix satisfying Kn2 ðzðkÞ vðkÞÞ ¼ vðkÞ zðkÞ. Since the augmented signal, ZðkÞ, and the noises, VðkÞ and V0 ðkÞ, have no zero mean, taking into account that E½uðkÞZðkÞ ¼ pðkÞE½ZðkÞ and by denoting ZðkÞ ¼ ZðkÞ E½ZðkÞ , V ðkÞ ¼ VðkÞ E½VðkÞ þ ½uðkÞ pðkÞ E½ZðkÞ and V0 ðkÞ ¼ V0 ðkÞ E½V0 ðkÞ , it follows that Y ðkÞ ¼ YðkÞ E½YðkÞ satisfy the following equation Y ðkÞ ¼ uðkÞZðkÞ þ V ðkÞ þ V0 ðkÞ
ð4Þ
in which the signal, ZðkÞ, and the noises, V ðkÞ and V0 ðkÞ, have zero mean. The remaining statistical properties of these processes are established in the following propositions. Proposition 1. Let us suppose that the hypotheses (H.1)–(H.5) are satisfied. Then the process fZðkÞ; k P 0g has zero mean and its covariance function is expressed in a semi-degenerate kernel form, specifically AðkÞBT ðsÞ; 0 6 s 6 k T KZ ðk; sÞ ¼ E½ZðkÞZ ðsÞ ¼ BðkÞAT ðsÞ; 0 6 k 6 s being AðkÞ ¼ BðkÞ ¼
AðkÞ 0n2 M 0
c1 ðkÞ 0n2 N 0
BðkÞ 0n2 M 0
0nN 0 c2 ðkÞ
0nP 0 d2 ðkÞ d1 ðkÞ 0n2 P 0
0nL0 : aðkÞ 0nL0 : bðkÞ
Moreover, the process fZðkÞ; k P 0g is independent of the noise fuðkÞ; k P 0g. Proof. Nakamori et al. [11].
h
Proposition 2. Under the hypotheses (H.1)–(H.5), the noise process fV ðkÞ; k P 0g of Eq. (4) is a sequence of zero-mean, mutually uncorrelated random vectors with covariance matrices given by RV ðkÞ ¼ E½V ðkÞV T ðkÞ ¼ pðkÞð1 pðkÞÞE½ZðkÞ E½ZT ðkÞ þ
Rv ðkÞ RTvv½2 ðkÞ
Rvv½2 ðkÞ R22 ðkÞ
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
71
being E½ZðkÞ ¼
0 vecð AðkÞBT ðkÞÞ
ð5Þ
and R22 ðkÞ ¼ pðkÞðIn2 þ Kn2 ÞðAðkÞBT ðkÞ Rv ðkÞÞðIn2 þ Kn2 Þ þ ðIn2 þ Kn2 ÞðaðkÞbT ðkÞ Rv ðkÞÞðIn2 þ Kn2 Þ þ pðkÞð1 pðkÞÞ ðIn2 þ Kn2 ÞðAðkÞBT ðkÞ aðkÞbT ðkÞÞðIn2 þ Kn2 Þ þ Rv½2 ðkÞ:
ð6Þ
Moreover, fV ðkÞ; k P 0g is uncorrelated with the processes fZðkÞ; k P 0g and fuðkÞZðkÞ; k P 0g. Remark. vecðÞ denote the vec or stack operator [12]. Proof. Using the independence hypotheses on the model it is easy to see that, for k 6¼ s, E½V ðkÞV T ðsÞ ¼ 0, and, for k ¼ s, Rvv½2 ðkÞ Rv ðkÞ RV ðkÞ ¼ pðkÞð1 pðkÞÞE½ZðkÞ E½ZT ðkÞ þ RTvv½2 ðkÞ R22 ðkÞ T
with R22 ðkÞ ¼ E½ðf ðkÞ E½f ðkÞ Þðf ðkÞ E½f ðkÞ Þ . Since u v ¼ vecðvuT Þ, using hypothesis (H.1), it is immediate that E½ZðkÞ verifies (5). On the other hand, using the Kronecker product properties, the hypotheses on the model and Eq. (2), we obtain that R22 ðkÞ verifies (6). Finally, the uncorrelation between fV ðkÞ; k P 0g and the processes fZðkÞ; k P 0g, fuðkÞZðkÞ; k P 0g is derived in a similar way, by using hypotheses (H.1)–(H.5) and employing the Kronecker product properties. h Proposition 3. Let us suppose that the hypotheses (H.1)–(H.5) are satisfied. Then the process fV0 ðkÞ; k P 0g of Eq. (4) is a zero-mean coloured noise and its covariance function is expressed in a semi-degenerate kernel form, specifically CðkÞDT ðsÞ; 0 6 s 6 k KV0 ðk; sÞ ¼ E½V0 ðkÞV0T ðsÞ ¼ ð7Þ DðkÞCT ðsÞ; 0 6 k 6 s being CðkÞ ¼ DðkÞ ¼
aðkÞ 0n2 M 00
1 ðkÞ 0n2 N 00
0nP 00 0nL00 q2 ðkÞ cðkÞ
bðkÞ 0n2 M 00
0nN 00 2 ðkÞ
q1 ðkÞ 0n2 P 00
0nL00 dðkÞ
0nM 0 M 00 pðkÞðIn2 þ Kn2 ÞðAðkÞ aðkÞÞ
0nM 0 M 00 : pðkÞðIn2 þ Kn2 ÞðBðkÞ bðkÞÞ
Moreover, fV0 ðkÞ; k P 0g is uncorrelated with the processes fV ðkÞ; k P 0g, fZðkÞ; k P 0g and fuðkÞZðkÞ; k P 0g.
72
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
Proof. Expression (7) for KV0 ðk; sÞ is obtained taking into account that, from hypothesis (H.3), for 0 6 s 6 k, KV0 ðk; sÞ aðkÞbT ðsÞ ¼ q1 ðkÞqT2 ðsÞ
1 ðkÞT2 ðsÞ 2 p ðkÞðIn2 þ Kn2 ÞðAðkÞBT ðsÞ aðkÞbT ðsÞÞðIn2 þ Kn2 Þ
and, for 0 6 k 6 s, KV0 ðk; sÞ bðkÞaT ðsÞ ¼ q2 ðkÞqT1 ðsÞ
2 ðkÞT1 ðsÞ : p2 ðkÞðIn2 þ Kn2 ÞðBðkÞAT ðsÞ bðkÞaT ðsÞÞðIn2 þ Kn2 Þ
Finally, as in Proposition 2, the uncorrelation between fV0 ðkÞ; k P 0g and the processes fV ðkÞ; k P 0g, fZðkÞ; k P 0g and fuðkÞZðkÞ; k P 0g is derived by using hypotheses (H.1)–(H.5) and employing the Kronecker product properties. h In view of the properties of the processes involved in Eq. (4), which have been established in Propositions 1–3, the recursive algorithms given in [10] can be applied to obtain the linear filtering and fixed-point smoothing estimators, b ðk; LÞ, L P k, of the signal ZðkÞ based on the observations Y ð1Þ; . . . ; Y ðLÞ. Z These algorithms, which are presented in the following theorem, allow us to obtain the required quadratic filtering and fixed-point smoothing estimators of b ðk; LÞ. the original signal zðkÞ, just by extracting the first n entries of Z Theorem 1. The filtering and fixed-point smoothing algorithm of the augmented signal ZðkÞ based on the augmented observations Y ð1Þ; . . . ; Y ðLÞ, L P k, is given by b ðk; LÞ ¼ Z b ðk; L 1Þ þ hðk; L; LÞmðLÞ Z where mðLÞ, the innovation, is mðLÞ ¼ Y ðLÞ pðLÞAðLÞOðL 1Þ CðLÞQðL 1Þ: The M 0 1 and N 0 1 vectors OðLÞ and QðLÞ, respectively, are recursively calculated by OðLÞ ¼ OðL 1Þ þ J ðL; LÞmðLÞ;
Oð0Þ ¼ 0
QðLÞ ¼ QðL 1Þ þ IðL; LÞmðLÞ;
Qð0Þ ¼ 0
being J ðL; LÞ ¼ ½pðLÞðBT ðLÞ rðL 1ÞAT ðLÞÞ cðL 1ÞCT ðLÞ P 1 ðLÞ
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
73
IðL; LÞ ¼ ½DT ðLÞ pðLÞcT ðL 1ÞAT ðLÞ dðL 1ÞCT ðLÞ P 1 ðLÞ where PðLÞ, the covariance matrix of the innovation, is given by PðLÞ ¼ RV ðLÞ þ pðLÞ½BðLÞ pðLÞAðLÞrðL 1Þ CðLÞcT ðL 1Þ AT ðLÞ þ ½DðLÞ pðLÞAðLÞcðL 1Þ CðLÞdðL 1Þ CT ðLÞ: The functions r, c and d are M 0 M 0 , M 0 N 0 and N 0 N 0 matrices, respectively, verifying rðLÞ ¼ rðL 1Þ þ J ðL; LÞ½pðLÞðBðLÞ AðLÞrðL 1ÞÞ
CðLÞcT ðL 1Þ ;
rð0Þ ¼ 0
cðLÞ ¼ cðL 1Þ þ J ðL; LÞ½DðLÞ pðLÞAðLÞcðL 1Þ
CðLÞdðL 1Þ ;
cð0Þ ¼ 0
dðLÞ ¼ dðL 1Þ þ IðL; LÞ½DðLÞ pðLÞAðLÞcðL 1Þ
CðLÞdðL 1Þ ;
dð0Þ ¼ 0:
The smoothing gain, hðk; L; LÞ, is given by hðk; L; LÞ ¼ ½pðLÞðBðkÞAT ðLÞ Eðk; L 1ÞAT ðLÞÞ
F ðk; L 1ÞCT ðLÞ P 1 ðLÞ where Eðk; LÞ and F ðk; LÞ are n M 0 and n N 0 matrices satisfying Eðk; LÞ ¼ Eðk; L 1Þ þ hðk; L; LÞ½pðLÞðBðLÞ AðLÞrðL 1ÞÞ
CðLÞcT ðL 1Þ Eðk; kÞ ¼ AðkÞrðkÞ F ðk; LÞ ¼ F ðk; L 1Þ þ hðk; L; LÞ½DðLÞ pðLÞAðLÞcðL 1Þ
CðLÞdðL 1Þ F ðk; kÞ ¼ AðkÞcðkÞ: b ðk; kÞ, which provides the initial condition for the fixedThe filtering estimate, Z b ðk; kÞ ¼ AðkÞOðkÞ. point smoothing algorithm, is given by Z The fixed-point smoothing and filtering error covariance matrices, P ðk; LÞ, L P k, satisfy P ðk; LÞ ¼ P ðk; L 1Þ hðk; L; LÞPðLÞhT ðk; L; LÞ; P ðk; kÞ ¼ AðkÞ½BT ðkÞ rðkÞAT ðkÞ :
L>k
74
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
4. A numerical simulation example In order to show the effectiveness of the proposed quadratic estimators, we have performed a program in MATLAB, which simulates the signal value at each iteration, and provides the linear and quadratic estimates, as well as the corresponding error covariance matrices. This program has been applied to a scalar signal fzðkÞ; k P 0g which, as in [11], is generated by the following first-order autoregressive model zðk þ 1Þ ¼ 0:95zðkÞ þ wðkÞ where fwðkÞ; k P 0g is a zero-mean white Gaussian noise with Var½wðkÞ ¼ 0:1; 8k P 0. The autocovariance functions of the signal and their second-order powers are given in a semi-degenerate kernel form, specifically, Kz ðk; sÞ ¼ 1:025641 0:95k s ; Kz2 ðk; sÞ ¼ 2:1038795 0:952ðk sÞ ;
06s6k
and the covariance function of the signal and their second-order powers is given by Kzz2 ðk; sÞ ¼ 0;
8s; k
According to hypothesis (H.1), the functions which constitute these covariance functions are as follows AðkÞ ¼ 1:025641 0:95k ; bðkÞ ¼ 0:95 2k ;
BðkÞ ¼ 0:95 k ;
aðkÞ ¼ 2:1038795 0:952k ;
c1 ðkÞ ¼ c2 ðkÞ ¼ d1 ðkÞ ¼ d2 ðkÞ ¼ 0:
Let us consider that the observations of the signal are given by yðkÞ ¼ uðkÞzðkÞ þ vðkÞ þ v0 ðkÞ where fuðkÞ; k P 0g is a sequence of independent Bernoulli random variables with P ½uðkÞ ¼ 1 ¼ p;
8k P 0
The noise fvðkÞ; k P 0g is a sequence of independent random variables with P ½vðkÞ ¼ 8 ¼ 18; P vðkÞ ¼ 87 ¼ 78; 8k P 0 and so, E½vðkÞ ¼ 0;
Rv ðkÞ ¼ 9:142857;
Rvv2 ðkÞ ¼ 62:693878;
Rv2 ðkÞ ¼ 429:900875: The coloured noise fv0 ðkÞ; k P 0g is generated by the following first-order autoregressive model
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
75
v0 ðk þ 1Þ ¼ 0:5v0 ðkÞ þ v1 ðkÞ where fv1 ðkÞ; k P 0g is a zero-mean white Gaussian noise with Var½v1 ðkÞ ¼ 0:075; 8k P 0. The autocovariance functions of the coloured noise and their second-order powers are given by in a semi-degenerate kernel form, namely, Kv0 ðk; sÞ ¼ 0:1 0:5k s ;
Kv2 ðk; sÞ ¼ 0:02 0:52ðk sÞ ; 0
06s6k
and the covariance function of the signal and their second-order powers is given by Kv0 v2 ðk; sÞ ¼ 0; 0
8s; k:
Hence, according to hypothesis (H.3), the functions which constitute these covariance functions are as follows aðkÞ ¼ 0:1 0:5k ; dðkÞ ¼ 0:5 2k ;
bðkÞ ¼ 0:5 k ;
cðkÞ ¼ 0:02 0:52k ;
1 ðkÞ ¼ 2 ðkÞ ¼ q1 ðkÞ ¼ q2 ðkÞ ¼ 0:
Finally, in accordance with the hypothesis (H.5) imposed in the theoretic study that we have developed, it is also assumed that the processes fzðkÞ; k P 0g, fuðkÞ; k P 0g, fvðkÞ; k P 0g and fv0 ðkÞ; k P 0g are mutually independent. In order to compare the linear (see [10]) and quadratic estimates, we have performed 200 iterations of the respective algorithms, considering different values of the false alarm probability 1 p, specifically, 1 p ¼ 0:5, 1 p ¼ 0:25 and 1 p ¼ 0, case in which the signal is always present in the observations. Fig. 1 displays the linear and quadratic filtering error variances, which measure the performance of the filters; this figure shows that, for all the values of the false alarm probability, the error variances corresponding to the quadratic filter are less than the linear filter ones; also we observe that as the false alarm probability is smaller, the error variances are also smaller and, consequently, the performance of the filter is better, as one could expect intuitively. The simulated signal, the linear and the quadratic filtering estimates for the values 1 p ¼ 0:5, 1 p ¼ 0:25 and 1 p ¼ 0, are displayed in Figs. 2–4, respectively. These figures show that the quadratic filtering estimate follows the signal evolution better and, also, that as the false alarm probability, 1 p, is smaller both the linear and the quadratic filtering estimates follow the signal evolution better, according to the results displayed in Fig. 1. Finally, Tables 1 and 2 summarize the mean-square values (MSVs) of the filtering and fixed-point smoothing errors, respectively, for the different values P of the false alarm probability. The MSVs are calculated by 200 zði; iÞÞ2 = i¼1 ðzðiÞ ^ P200 P10 2 200, for the filtering errors, and by i¼1 j¼1 ðzðiÞ ^zði; i þ jÞÞ =2000, for the fixed-point smoothing errors. Tables 1 and 2 show, on the one hand, that the
76
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79 1
0.9 (a)
Variances of Filtering Error
0.8 (b)
0.7 (c)
0.6
0.5
(a’)
0.4 (b’)
0.3
(c’)
0.2 0
5
10
15
20
25
30
Time k
Fig. 1. Linear filtering error variances (solid line) and quadratic filtering error variances (dashed line) when (a), (a0 ): 1 p ¼ 0:5, (b), (b0 ): 1 p ¼ 0:25, (c), (c0 ): 1 p ¼ 0.
4 Signal Linear Filtering Estimate Quadratic Filtering Estimate
Signal and Filtering Estimates
3
2
1
0
–1
–2
0
20
40
60
80
100
120
140
160
180
200
Time k
Fig. 2. Signal (solid line), linear filtering estimate (dashed line) and quadratic filtering estimate (dotted line) for 1 p ¼ 0:5.
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
77
4 Signal Linear Filtering Estimate Quadratic Filtering Estimate
Signal and Filtering Estimates
3
2
1
0
–1
–2 0
20
40
60
80
100
120
140
160
180
200
Time k
Fig. 3. Signal (solid line), linear filtering estimate (dashed line) and quadratic filtering estimate (dotted line) for 1 p ¼ 0:25.
5 Signal Linear Filtering Estimate Quadratic Filtering Estimate
Signal and Filtering Estimates
4
3
2
1
0
–1
–2
0
20
40
60
80
100
120
140
160
180
200
Time k
Fig. 4. Signal (solid line), linear filtering estimate (dashed line) and quadratic filtering estimate (dotted line) for 1 p ¼ 0.
78
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
Table 1 Mean-square values of the linear filtering errors and the quadratic filtering errors for different values of the false alarm probability False alarm probability
MSVs of linear filtering errors
MSVs of quadratic filtering errors
1 p ¼ 0:50 1 p ¼ 0:25 1 p ¼0
0.6967 0.5046 0.3939
0.3069 0.2315 0.1527
Table 2 Mean-square values of the linear fixed-point smoothing errors and the quadratic fixed-point smoothing errors for different values of the false alarm probability False alarm probability
MSVs of linear smoothing errors
MSVs of quadratic smoothing errors
1 p ¼ 0:50 1 p ¼ 0:25 1 p ¼0
0.5832 0.4421 0.3503
0.2543 0.1697 0.1180
estimation accuracy of the smoother is superior to that of the filter and, on the other, not only that the MSVs of the quadratic filtering and smoothing errors are less than the MSVs of the linear ones, but also that as the false alarm probability is smaller, the MSVs are smaller too, agreeing again with the theoretical results showed in Fig. 1.
5. Conclusions In this paper, using covariance information, recursive LMSE algorithms for the quadratic filter and fixed-point smoother are derived from uncertain observations, for the case of white plus coloured observation noises. The proposed estimators do not require the knowledge of state-space model of the signal; the available information is only the autocovariance and crosscovariance functions of the signal and their second-order powers, defined by the Kronecker product, as well as the corresponding functions of the additive noises. It is also assumed that the probability that the signal exists in the observed values is known. By defining a suitably augmented observation equation, the LMSE linear estimator of the augmented signal based on the augmented observations provides the LMSE quadratic estimator for the original signal. The effectiveness of the quadratic estimators in contrast to the linear ones is shown by applying the filtering and fixed-point algorithms to the estimation of
S. Nakamori et al. / Appl. Math. Comput. 155 (2004) 65–79
79
a signal generated by a first-order autoregressive model, from uncertain observations with different values of the false alarm probability.
Acknowledgements This work has been partially supported by the ‘‘Ministerio de Ciencia y Tecnologıa’’ under contract BFM2002-00932.
References [1] A. De Santis, A. Germani, M. Raimondi, Optimal quadratic filtering of linear discrete-time non-Gaussian systems, IEEE Transactions on Automatic Control AC-40 (1995) 1274–1278. [2] A. Germani, G. Mavelli, Optimal quadratic solution for the non-Gaussian finite-horizon regulator problem, Systems and Control Letters 38 (1999) 321–331. [3] F. Carravetta, A. Germani, M. Raimondi, Polynomial filtering for linear discrete time nonGaussian systems, SIAM Journal of Control and Optimization 34 (1996) 1666–1690. [4] F. Carravetta, A. Germani, M. Raimondi, Polynomial filtering of discrete-time stochastic linear systems with multiplicative state noise, IEEE Transactions on Automatic Control AC-42 (1997) 1106–1126. [5] M. Dalla Mora, A. Germani, A. Nardecchia, Restoration of images corrupted by additive nonGaussian noise, IEEE Transactions on Circuits and Systems I: Fundamental Theory and Applications 48 (2001) 859–875. [6] A. Hermoso, J. Linares, Linear smoothing for discrete-time systems in the presence of correlated disturbances and uncertain observations, IEEE Transactions on Automatic Control 40 (1995) 1486–1488. [7] W. NaNakara, E.E. Yaz, Recursive estimator for linear and nonlinear systems with uncertain observations, Signal Processing 62 (1997) 215–228. [8] R. Caballero, A. Hermoso, J. Linares, Polynomial filtering with uncertain observations in stochastic linear systems, International Journal of Modelling and Simulation 23 (2003) 22–28. [9] S. Nakamori, R. Caballero, A. Hermoso, J. Linares, New design of estimators using covariance information with uncertain observations in linear discrete-time systems, Applied Mathematics and Computation 135 (2003) 429–441. [10] S. Nakamori, R. Caballero, A. Hermoso, J. Linares, Linear estimation from uncertain observations with white plus coloured noises using covariance information, Digital Signal Processing 13 (3) (2003) 552–568. [11] S. Nakamori, R. Caballero, A. Hermoso, J. Linares, Second-order polynomial estimators from uncertain observations using covariance information, Applied Mathematics and Computation 143 (2003) 319–338. [12] J.R. Magnus, H. Neudecker, Matrix Differential Calculus with Applications in Statistics and Econometrics, Wiley & Sons, New York, 1999.