State dimension reduction and analysis of quantized estimation systems

Signal Processing 105 (2014) 363–375

Contents lists available at ScienceDirect

Signal Processing journal homepage: www.elsevier.com/locate/sigpro

State dimension reduction and analysis of quantized estimation systems$ Hui Zhang a,n, Ying Shen b, Li Chai c, Li-Xin Gao d a State Key Laboratory of Industrial Control Technology, Institute of Industrial Process Control, Department of Control Science and Engineering, Zhejiang University, Hangzhou 310027, China b School of Engineering, Huzhou Teachers College, Huzhou 31300, China c Engineering Research Center of Metallurgical Automation and Measurement Technology, Wuhan University of Science and Technology, Wuhan 430081, China d Institute of Intelligent Systems and Decision, Wenzhou University, Zhejiang 325027, China

a r t i c l e i n f o

abstract

Article history: Received 28 August 2013 Received in revised form 14 April 2014 Accepted 10 June 2014 Available online 23 June 2014

The problem of state dimension reduction and quantizer design under communication constraints is discussed for state estimation in quantized linear systems. Subject to the limited signal power, number and bandwidth of the parallel channels, a differential pulse code modulation (DPCM)-like structure is adopted to generate the quantized innovations as the transmitted signals, and the multi-level quantized Kalman filter (MLQ-KF) is used to serve as the pre- and post-filters. The dimension reduction matrix and quantizer are designed jointly under the MMSE criterion of estimation at the channel receiver. To demonstrate the validity of state estimation under the adopted framework, the state estimability based on quantized innovations is analyzed by using information theoretic method. This leads to a sufficient and necessary condition of a certain estimability Gramian matrix having full rank. The quantized Gramian is proved to converge to that of the original unquantized system when the quantization intervals turn to zero. Our work also provides an auxiliary analytic support for the estimation under 1-bit quantization. Simulations show that under communication constraints, the estimation performance is satisfactory when the designed dimension reduction method and quantizer are applied. The analytic conclusion of estimability is also verified by illustrative simulations. & 2014 Elsevier B.V. All rights reserved.

Keywords: Communication constraint Quantized linear system Quantizer design State dimension reduction Estimability Mutual information

1. Introduction Recently there have been plenty of works concerning the quantized estimation problem, e.g. [1–7]. A modified Kalman filter called SOI-KF was proposed by using the sign of

☆ This work is supported by the Natural Science Foundation of China under Grants 61174063 & 61171160, and Zhejiang Provincial Natural Science Function of China under Grant LY13F030048. n Corresponding author. Tel./fax: þ86 571 87952520. E-mail addresses: [email protected] (H. Zhang), [email protected] (Y. Shen), [email protected] (L. Chai), [email protected] (L.-X. Gao).

http://dx.doi.org/10.1016/j.sigpro.2014.06.007 0165-1684/& 2014 Elsevier B.V. All rights reserved.

innovations or the 1-bit quantization of innovations [1]. The filtering problem based on multi-level quantizers was studied in [2–4], where the optimal filters and the multilevel quantization schemes are designed simultaneously. Lu et al. [5] and Liu et al. [6] discussed the problem of state estimation under quantization, while Farias and Brossier [7] considered the quantizer design for estimation. Another important methodology [40–42] dealing with quantized estimation adopts the so-called expectation maximization (EM) algorithm involving the Monte Carlo (MC) techniques. In this paper, we address the problem of joint state dimension reduction and quantizer design for state estimation in linear dynamical systems with a differential pulse code modulation

364

H. Zhang et al. / Signal Processing 105 (2014) 363–375

(DPCM, [8])-like structure, where the multi-level quantized Kalman filter (MLQ-KF) proposed in [3] is adopted as the estimator. To deal with the case of limited number of simultaneous connections, this paper adopts the method of reducing the dimension of transmitted data of system states. Dimension reduction has been extensively used in various applications in fields of communication and control, such as the system model reduction [9], reduced-dimension coding [10,11], and speech data compression [12]. The main difference between model reduction [9] and the state dimension reduction here is, the former is looking for a projection matrix which aggregates the states while making the reduced-order model approximate the input–output property of the original system as precise as possible, whereas the later is looking for a compression matrix which reduces the quantity of transmission data while making the state estimation at the receiver as precise as possible. Contributions of Zhu et al. [10] and Goela and Gastpar [11] consider the communication in network setting by using reduced-dimension transformation to deal with the bandwidth limitation, while our work focuses on the state estimation under reduced-dimension transmission where the bandwidth limitation is further dealt with by quantization. DPCM is an effective way of reducing the redundancy of transmission data [12–15]. Our work bases on a DPCM-like scheme where the quantized innovations are transmitted. There are a pre-filter which generates innovations to be transmitted and a post-filter which estimates the states at the receiver of noiseless channel. The present structure is different from the traditional DPCM in that, in our system the estimation performance is pursued partly by designing the dimension reduction matrix (which can also be considered as a designable system output matrix). Our design work is also enlightened by but different from Msechu et al. [2], You et al. [3], Fu and de Souza [4], Curry [13], Ramabadran and Sinha [12] and Arildsen et al. [15]. Msechu et al. [2], You et al. [3], Fu and de Souza [4] and Curry [13] design the quantized estimator and quantizer without considering state dimension reduction, while Ramabadran and Sinha [12] and Arildsen et al. [15] consider the coder-reduction (measurements selection) and decoder (filter) design with a representation of linear additive quantization noise model (LAQNM). Specially, in the classic book [13] Curry provided the methods of designing the predictive quantizer under the Gaussian assumption (GA) [13, Section 3.3], and the joint optimization of quantizer and filter under LAQNM [13, Section 4.5]. In fact, the dimension reduction, pre-filter and quantizer play jointly the role of preprocessing or coding, and the post-filter is the decoder. In this sense, the design problem in the present paper is similar to the precoding design for fixed receiver [16], joint source-channel coding [17], and joint design of precoder and decoder [18,19]. Different from Wiesel et al. [16], Persson et al. [17], Sampath et al. [18] and Khaled et al. [19] the dimension reduction matrix and quantizer are involved here for design, and the MLQ-KF [3] under GA is adopted in our system as both of the pre- and post-filters. Our scheme is also similar to the recently proposed approaches involving predictor in the pre-processing part, the generalized quantization

scheme (GQS) [42] and its special form of the noise shaping quantizer (NSQ) [43], but different from them in two points. One point is that in our scheme the adopted MLQ-KF [3] (see Section 2.3) provides analytic formula of the estimator, while in [42] the EM algorithm is adopted based on the MC technique. The other point is in our scheme the dimension reduction matrix is designed jointly with the quantizer levels. The reduction matrix can also be considered as a static filter under the framework of the GQS [42]. Although many works concerning the quantized estimation algorithm and performance have been published, a problem is slightly ignored: the validity or effectiveness of the estimation may be destroyed by improper quantization scheme. Few contributions concerning the property of state reconstruction under quantization were found. In [45] the conditions of input and quantization design were proposed for guaranteeing the active observability under binary valued observations. Here, we focus on the estimability property of quantized systems which is important from both the system theory point of view and the engineering point of view. This is another topic of this paper. The concept of estimability was proposed by Baram and Kailath to evaluate the validity of state estimation under the minimum mean square error (MMSE) criterion in [20], where it is stated that the system is estimable if and only if the posterior estimation error is strictly less than the prior estimation error. The estimability was also discussed from information theoretic viewpoint [21–23] based on the minimum maximum error entropy (MMEE) estimation [24–26]. Referring to this methodology, the present paper discusses the estimability given quantized innovations, which is different from that in [23] concerning quantized outputs. Our discussion involves the ‘true’ quantizer. This is different from that based on the assumption of additive uniform noise for approximating quantization [43,46]. In short, this paper consists of two main works: the optimal design of state dimension reduction and quantizer for state estimation in linear systems under limited channel (or subchannel) number, bandwidth and power of transmitted signal; the estimability analysis given quantized innovations. Under the MMSE criterion, the reduction matrix and quantizer are designed jointly based on a DPCM-like structure and the estimator of MLQ-KF. Simulation results show that the method presented in this paper has satisfactory estimation performance. The analysis work leads to an algebraic condition of estimability for quantized linear Gaussian systems, i.e. certain Gramian matrix having full rank, and shows the quantized Gramian converges to that of the original unquantized system when the quantization intervals turn to zero. The relation between the present result and that in [23] is discussed. Analytical results are also illustrated by simulations. The design work was presented primarily in the conference paper [27], and the estimability analysis is an extension of [28,23]. The rest of the paper is organized as follows: Section 2 introduces the system and the problem of state preprocessing, reviews the quantizer and the MLQ-KF; The state dimension reduction and quantizer are designed jointly under MMSE criterion in Section 3; Section 4 introduces

H. Zhang et al. / Signal Processing 105 (2014) 363–375

the concept of estimability, gives analytic results based on quantized innovations, and proves the convergence of quantized estimability Gramian; Sections 5 and 6 are illustrative simulations and conclusion, respectively. Throughout this paper, we use Y k to denote a vector ½yT ðk0 Þ yT ðk0 þ 1Þ ⋯ yT ðkÞT which consists of a sequence fyðkÞ; k Zk0 g of vectors (or variables), where k denotes the time instants. The initial instant can be set as k0 ¼ 0 for time-invariant systems without losing generality. We use x^ ðjjkÞ to denote the MMSE or MMEE estimation of system state xðjÞ at instant j given data Y k . Other notations are as follows. ℝ, ℝn denote the set of all real numbers and the n-dimensional real space, respectively. AT , tr(A) and rankfAg denote respectively the transposition, trace and rank of a matrix (or vector) A. For random variables or vectors x and y, pðxÞ and pðxjyÞ denote the probability density of x and that of x given y, respectively; Efxg and Efxjyg denote respectively the mathematical expectation of x and the conditional expectation of x given y; HðxÞ, HðxjyÞ and Iðx; yÞ denote the entropy of x, the conditional entropy of x given y and the mutual information between x, y, respectively. 2. System design framework 2.1. System state The system state is modeled as xðk þ 1Þ ¼ AxðkÞ þBwðkÞ

ð1Þ

where xðkÞ A ℝn , wðkÞ A ℝm ; A, B are constant matrices with appropriate dimensions; wðkÞ is a zero-mean white Gaussian random vector with constant covariance matrix Q, and is uncorrelated with the initial state x(0) which is Gaussian with zero mean and bounded covariance. The state xðkÞ is measured and preprocessed before transmission, and reconstructed at the receiver. 2.2. State preprocessing The channel under discussion is assumed to be the parallel channel [18,19,29] which consists of p independent noiseless subchannels with limited bandwidth and a power constraint on the transmitted signal, where p is generally assumed to be smaller than n (the dimension of system state). The number of parallel data streams is then required to be not greater than p. This will be met by state dimension reduction. We reduce the state dimension linearly, namely generate a linear combination of state measurements by aggregation. In general, the sensing is affected by thermal and other unobservable tiny noises, and the reduction error is inevitable. In areas of communication and data compression, independent white Gaussian noise is commonly used to approximate kinds of noise, uncertainty and disturbance [30,31]. Then in our problem, various uncertainties caused by sensing and dimension reduction are also approximated by an additive white Gaussian process, according to the central limit theorem [32]. In this way, the state dimension reduction can be described as yðkÞ ¼ CxðkÞ þvðkÞ

ð2Þ

365

where yðkÞ A ℝp is the output of state dimension reduction; C is the reduction matrix to be designed; vðkÞ A ℝp is a Gaussian white noise with zero-mean and constant covariance R, and is uncorrelated with w(k) and x(0). In the design problem, we consider the case of p ¼ 1, i.e. yðkÞ; vðkÞ A ℝ. Due to the bandwidth constraint, quantization is involved. We adopt a DPCM-like framework, where the innovations are generated and quantized for transmission: ~ yðkÞ ¼ yðkÞ  y^ ðkjk  1Þ

ð3Þ

k1 where y^ ðkjk  1Þ ¼ EfyðkÞjY~ q g is the prediction of yðkÞ k1 given Y~ q ¼ ½y~ q ð0Þ: y~ q ð1Þ⋯y~ q ðk  1ÞT which contains the

quantized innovations from instant 0 to k  1. The designable quantizer is a symmetrical one, by which the innovation is partitioned into 2L þ1 contiguous and nonoverlapping intervals: 8 ~ aL ; aL o yðkÞ > > > > > ⋮ > <⋮ ~ r a2 a1 o yðkÞ ~ ð4Þ y~ q ðkÞ ¼ Q s ðyðkÞÞ ¼ a1 ; > > > ~ a0 ; 0 r yðkÞ ra1 > > > : Q ð  yðkÞÞ; ~ ~ o0 yðkÞ s ~ where y~ q ðkÞ denotes the quantized yðkÞ, Q s ð U Þ is the symmetrical quantizer, fa  ðL þ 1Þ ; a  L ; …; a  1 ; a0 ; a1 ; …; aL ; aL þ 1 g are the thresholds, 0 ¼ a0 o a1 o ⋯ o aL o þ1 ¼ aL þ 1 , a  i ¼ ai . We use aL to denote the set of thresholds fa1 ; a2 ; …; aL g. In contrast to the case of transmitting states directly, data redundancy is reduced when innovations are transmitted, and the signal power is cut down. The power of innovation is

σ 2y~ ðkÞ ¼ Efy~ 2 ðkÞg:

ð5Þ

The channel is assumed to be lossless. The receiver k estimates xðkÞ by using quantized innovations Y~ q . The predictor and the filter will be introduced in the following context. The system with a DPCM-like structure is shown in Fig. 1. 2.3. MLQ-KF For system (1)–(4), the MLQ-KF [3] is adopted as both the pre- and post-filters (i.e., the predictor and the filter in Fig. 1). The MLQ-KF is an MMSE estimator which seeks k x^ ðkjkÞ ¼ EfxðkÞjY~ g under the GA [1,3,13] that the state q

distribution conditioned on the quantized innovations is Gaussian, and possesses the property of convergence to standard Kalman filter when the number of quantization

Fig. 1. System structure.

366

H. Zhang et al. / Signal Processing 105 (2014) 363–375

levels turns to infinite. Equations of the MLQ-KF algorithm used here are as follows: One-step-ahead predictions of the state and output at time instant k: ^ xðkjk  1Þ ¼ Ax^ ðk  1jk  1Þ;

ð6Þ

^ y^ ðkjk 1Þ ¼ C xðkjk  1Þ:

ð7Þ

Estimation of the state at time instant k:

σ 2y~ ðkÞ rW:

f ðaL ; y~ q ðkÞÞPðkjk 1ÞC T ^ ^ xðkjkÞ ¼ xðkjk  1Þ þ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi : CPðkjk  1ÞC T þR

ð8Þ

Covariance of the state prediction error at time instant k: Pðkjk  1Þ ¼ APðk  1jk  1ÞAT þQ :

ð9Þ

Covariance of the state estimation error at time instant k: T

Pðkjk  1ÞC CPðkjk  1Þ : PðkjkÞ ¼ Pðkjk  1Þ  FðaL Þ CPðkjk  1ÞCT þ R

ð10Þ

In these equations: y~ q ðkÞ is produced by (3) and (4):

ϕðal =σ y~ ðkÞÞ  ϕðal þ 1 =σ y~ ðkÞÞ ; f ða ; y~ q ðkÞÞ ¼ ∑ I fal g ðy~ q ðkÞÞ Tðal þ 1 =σ y~ ðkÞÞ Tðal =σ y~ ðkÞÞ l ¼ L L

L

ð11Þ ðϕðal =σ y~ ðkÞÞ  ϕðal þ 1 =σ y~ ðkÞÞÞ2 ; l ¼ 0 Tðal þ 1 =σ y~ ðkÞÞ Tðal =σ y~ ðkÞÞ L

FðaL Þ ¼ 2 ∑

where I fal g ð:Þ is a standard sign function: ( 1; x ¼ al ; I fal g ðxÞ ¼ 0; otherwise:

be minimal. The effect of noise can be attenuated by increasing the power of the useful signal. However the signal power is always limited by measurement equipment and other practical factors [8]. Many works characterized this limitation in terms of signal-to-noise ratio (SNR) constraints [16,33,34]. Instead of this approach, here we assume a power constraint W on the transmitted innovations as in [18]:

ð12Þ

ð13Þ

ϕð U Þ denotes the probability density function of a standard

Thus our problem boils down to the problem as follows: mintrðPðkjkÞÞ C;aL

s:t: σ 2y~ ðkÞ r W:

CðkÞPðkjk  1ÞCðkÞT þ R rW

ð17Þ

T

where CðkÞPðkjk 1ÞCðkÞ is the power of the useful signal, R is the power of the noise. Hence, the SNR is maximal when CðkÞPðkjk  1ÞCðkÞT þ R ¼ W:

ð18Þ

The optimal estimation performance is expected when (18) is satisfied. This will be implied in Section 3.2. We employ the iterative method to solve this problem: compute the dimension reduction matrix C(k) and quantizer aL ðkÞ ¼ fa1 ðkÞ; a2 ðkÞ; …; aL ðkÞg at every time instant k, which are expected to approach to a constant matrix C and a threshold set aL of a time-invariant symmetrical quantizer, respectively. From (9), trðPðkjkÞÞ ¼ trðPðkjk  1ÞÞ L

 Fða ðkÞÞtr

Pðkjk  1ÞCðkÞT CðkÞPðkjk 1Þ

!

CðkÞPðkjk  1ÞCðkÞT þ R

CðkÞP 2 ðkjk  1ÞCðkÞT ¼ trðPðkjk  1ÞÞ  FðaL ðkÞÞ CðkÞPðkjk  1ÞCðkÞT þR ð19Þ At time instant k, Pðkjk  1Þ is known, then

There are three factors affecting the state reconstruction performance: the dimension reduction, the quantizer and the quantized filter. When the filter is defined, the design problem is then to look for the reduction matrix C and the quantization algorithm under preferred optimal principle of state estimation at the channel receiver. There is no reason to think uncritically that independently designed reduction and quantizer will be optimal for state estimation. This leads to a joint optimization problem stated as follows. 3.1. Design method

CðkÞP 2 ðkjk 1ÞCðkÞT min trðPðkjkÞÞ 3 max FðaL ðkÞÞ CðkÞPðkjk 1ÞCðkÞT þ R

ð14Þ

ð20Þ

where “ 3 ” denotes “be equivalent to”. Let y~ 0 ðkÞ ¼ ~ yðkÞ= σ y~ ðkÞ. For y~ 0 ðkÞ is a standard Gaussian variable with constant covariance 1, the quantizer for y~ 0 ðkÞ can be designed as aL ¼ fa1 ; a2 ; …; aL g, then the thresholds aL ðkÞ at time instant k can be obtained by al ðkÞ ¼ al σ y~ ðkÞ; l ¼ 1; 2; …; L. Thus the value of FðaL ðkÞÞ formulated as (12) ~ is not affected by the covariance of innovation yðkÞ. As a result, min trðPðkjkÞÞ 3 max FðaL Þ

It is expected to recover the states from the quantized innovations. Under MMSE principle, we look for a linear dimension reduction matrix and a quantizer which make trðPðkjkÞÞ

ð16Þ

Constraint (15) is equivalent to

Gaussian variable with zero mean and unit covariance [3], and T(.) denotes its probability distribution function; σ 2y~ ðkÞ ¼ CPðkjk  1ÞC T þ R is the innovation covariance, the ~ innovation yðkÞ is defined as (3). 3. Joint design of reduction matrix and quantizer

ð15Þ

CðkÞP 2 ðkjk 1ÞCðkÞT CðkÞPðkjk  1ÞCðkÞT þ R

where ðϕðal Þ  ϕðal þ 1 ÞÞ2 l ¼ 0 Tðal þ 1 Þ  Tðal Þ L

FðaL Þ ¼ 2 ∑

ð21Þ

H. Zhang et al. / Signal Processing 105 (2014) 363–375

is related with the quantizer only and irrelevant to C(k). While, CðkÞP 2 ðkjk 1ÞCðkÞT CðkÞPðkjk 1ÞCðkÞT þ R is related with C(k) only and irrelevant to the quantizer aL ðkÞ. Thus the quantizer and the dimension reduction matrix can be designed separately at each time instant. Consequently, problem (16) is decomposed into the following two sub-problems: Sub-problem I: Design of the optimal quantizer, i.e. 2 L max FðaL Þ ¼ 2 ∑ ðϕðal Þ  ϕðal þ 1 ÞÞ L a Tða Þ Tða Þ lþ1 l l¼0

s:t: 0 o a1 o ⋯ o aL o þ 1:

ð22Þ

Sub-problem II: Design of the optimal dimension reduction matrix, i.e. ! CðkÞP 2 ðkjk 1ÞCðkÞT max CðkÞ CðkÞPðkjk 1ÞCðkÞT þ R s:t: CðkÞPðkjk  1ÞCðkÞT þ R rW:

ð23Þ

3.2. Solutions In the inequality constrained optimization problem of Sub-problem I, the objective function FðaL Þ is bounded and infinitely differentiable, while the constraints define a convex set of the standard Gaussian variables al ; l ¼ 1; 2; …; L. Due to these properties and the special form of the constraints matrices, the gradient and Hessian conditions of the optimality [47] can be checked analytically for simple cases, e.g. L ¼ 1. Numerical methods are practical consideration for seeking solutions to the cases of high levels. From which we prefer the sequential quadratic programming (SQP) method [48] based on expanding quadratically the objective function and checking the Karush–Kuhn–Tucker conditions. By checking the feasible and bounded conditions [48], it is seen that the SQP is applicable to solve the Sub-problem I, and the local optimal solution is expected to be obtained by evoking the following Matlab command: ½x; y ¼ fminconðf un; x0 ; A; b; Aeq; beq; lb; ubÞ: The “fmincon” attempts to find a constrained minimum of a scalar function of several variables starting at an initial estimate, and it is a gradient-based approach to problems where the objective and constraint functions are both continuous and have continuous first derivatives, so its use in this sub-problem is reasonable. This command starts at an initial value x0 and attempts to find a minimizer x of the function described by fun subject to the linear inequalities A Ux rb and linear equalities Aeq Ux ¼ beq. The solution x is set in lb rx rub, and the minimal value of fun is returned by y. In our problem, only the inequalities A U x rb are required and formulated by using constraints on aL in (22); f un ¼ FðaL Þ. The “fmincon” uses one of three algorithms: active-set, interiorpoint, or trust-region-reflective, and the optimization result is local minimum possible. In our procedure, the

367

active-set algorithm which is implemented by SQP method will be selected and illustrated by Monte Carlo (MC) simulations. The solution aL is then used to get aL ðkÞ ¼ pffiffiffiffiffi ffi aL U σ y~ ðkÞ ¼ aL U W . This means aL ðkÞ ¼ aL is constant and can be computed independently of the solution to Subproblem II. The Lagrange multiplier method is employed to solve Sub-problem II. The object function is then written as JðCðkÞÞ ¼

CðkÞP 2 ðkjk 1ÞCðkÞT CðkÞPðkjk  1ÞCðkÞT þ R

þ ρðCðkÞPðkjk 1ÞCðkÞT þR  WÞ

ð24Þ where ρ is the Lagrange multiplier. The optimal C(k) can be obtained by solving dJðCðkÞÞ=dCðkÞ ¼ 0:

ð25Þ

Let λ1 ðkÞ Z λ2 ðkÞ Z ⋯ Z λn ðkÞ be the eigenvalues of Pðkjk  1Þ and ε1 ðkÞ; ε2 ðkÞ; …; εn ðkÞ be the corresponding unit orthogonal eigenvectors, then the optimal reduction matrix at time k is qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi CðkÞ ¼ ðW  RÞ=λ1 ðkÞ:εT1 ðkÞ: ð26Þ Condition (18) is satisfied by using (26). Let the difference value between the traces of estimation error covariance matrices of the adjacent iteration times be dðiÞ ¼ trðPði 1ji 1ÞÞ  trðPðijiÞÞ, the expected precision of d(i) be delta. The sub-problem II is then solved by the following algorithm: Step 0: Set quantization level L, power constraint W, and compute the optimal quantizer (the solution to sub-problem I); Step 1 (initiation): Let iteration times i¼0, set prediction error covariance Pð0j  1Þ, delta, and dð0Þ ¼ 1, compute C(0); Step 2: Let i ¼ i þ1, compute Pðiji  1Þ, then C(i) by (26) and d(i) by using (19); Step 3: If dðiÞ r delta, then stop; otherwise, go to step 2. When the algorithm stops, the obtained matrix CðiÞ ¼ C is the constant dimension reduction matrix we are looking for. The parameter delta represents the convergence degree of the algorithm. The smaller the delta is, the closer to the optimal point the solution is. The precision is at the expense of computation time. For time invariant systems, we can choose to get the dimension reduction matrix offline. In this case, we pay more attention to precision than computation, and hence tend to set delta as a very small value. 4. Estimability analysis In order to evaluate the validity of state estimation in systems involving the DPCM-like structure as Fig. 1, here we discuss the concept of estimability proposed for the time-varying linear system (27) which takes system (1), (2) as a special case: ( xðk þ 1Þ ¼ AðkÞxðkÞ þ BðkÞwðkÞ ð27Þ yðkÞ ¼ CðkÞxðkÞ þ vðkÞ

368

H. Zhang et al. / Signal Processing 105 (2014) 363–375

where xðkÞ A ℝn , wðkÞ A ℝm , yðkÞ; vðkÞ A ℝp ; A(k), B(k), C(k) are time-varying matrices with appropriate dimensions; the initial state possesses zero mean and bounded covariance; Efxðk0 ÞwT ðkÞg ¼ 0, Efxðk0 ÞvT ðkÞg ¼ 0, EfwðkÞvT ðjÞg ¼ 0, EfwðkÞwT ðjÞg ¼ Q ðkÞδk;j , EfvðkÞvT ðjÞg ¼ RðkÞδk;j , where δk;j ¼ 1 when k ¼ j, otherwise δk;j ¼ 0. It was proposed in [20] that the system (27) is estimable if the optimal posterior mean-square error of state estimation is strictly smaller than the prior estimation error. This leads to a necessary and sufficient condition that the following estimability Gramian has full-rank: k

W e ½k0 ; k ¼ ∑ Φðk; jÞΠ ðjÞC T ðjÞCðjÞΠ ðjÞΦ ðk; jÞ; T

ð28Þ

j ¼ k0

Assumption 1. In system (27) the initial state xðk0 Þ, noises wðkÞ and vðkÞ are mutually independent Gaussian variables; The outputs are scalar, i.e. yðkÞ; vðkÞ A ℝ; The outputs are preprocessed to generate quantized innovations by (29) before transmission, y~ q ðkÞ are decoded by (30) at the lossless channel receiver and then used to estimate the states. It was shown in [23] that the quantized estimability means the mutual information between state and quantized outputs is larger than zero. While, it can be seen [24–26] that when n the optimal posterior errors are achieved, Hðg T x~ ðkÞÞ ¼ k k 8 g A ℝn ; g a0, Hðg T xðkÞjY~ Þ ¼ Hðg T xðkÞÞ Iðg T xðkÞ; Y~ Þ, d

d

x^ ðkÞ ¼ EfxðkÞjY~ d g, Y~ d ¼ where x~ ðkÞ ¼ xðkÞ  x^ ðkÞ, ½y~ d ðk0 Þy~ d ðk0 þ1Þ⋯y~ d ðkÞT . Then referring to [22–26], we prefer the following condition for the analysis of estimation given quantized innovations: The system (27), (29) and (30) is estimable if 8 k Zk0 þn 1, 8 g A ℝn ; g a 0, n

n

where Φðk; jÞ ¼ Aðk 1ÞAðk  2Þ⋯AðjÞ; jo k, Φðk; kÞ ¼ I, and Π ðkÞ ¼ EfxðkÞxT ðkÞg. Eq. (28) can be interpreted by using the mutual information between state and the data received by the estimator [21–23]. In (27), estimation given outputs is equivalent to that given innovations [35]. However, the relation between estimability based on quantized innovations and that based on quantized outputs [23] needs to be discussed.

k

n

Iðg T xðkÞ; Y~ d Þ 4 0; k

ð31Þ

or equivalently, Hðg T xðkÞÞ 4 Hðg T xðkÞjY~ d Þ: k

4.1. Estimability based on quantized innovations We discuss the estimability of system (27) (with p ¼ 1) under the framework of Fig. 1 with a slight difference that the measurement output matrix is not designable. In order to get general result, here we do not focus on a certain quantizer, but on the general time-varying N-level quantization [13,14] described by (29) which takes (4) as a special case. ~ ~ A κ l ðkÞ; y~ q ðkÞ ¼ Q k ðyðkÞÞ ¼ zl ðkÞ for yðkÞ

l ¼ 1; 2; …; N ð29Þ

~ where yðkÞ ¼ yðkÞ  y^ ðkjk 1Þ are innovations, y^ ðkjk  1Þ ¼ k1 EfyðkÞjY~ q g is the prediction of yðkÞ in (27) given quank1 Y~ q

tized innovations ¼ ½y~ q ðk0 Þy~ q ðk0 þ1Þ⋯y~ q ðk 1Þ ; y~ q ðkÞ A fz1 ðkÞ; z2 ðkÞ; …; zN ðkÞg are quantizer outputs; κ l ðkÞ ¼  N þ 1 ai ðkÞ i ¼ 1 , ðal ðkÞ; al þ 1 ðkÞ are quantization intervals,  1 ¼ a1 ðkÞ oa2 ðkÞ o⋯ oaN þ 1 ðkÞ ¼ þ1g are thresholds. Also for the consideration of generality, we assume that the quantized innovations y~ q ðkÞ are decoded by y~ d ðkÞ ¼ Dk ðy~ q ðkÞÞ

T

ð30Þ

ð32Þ

As information theoretic quantities [29], entropy describes the information amount of a variable, while mutual information measures the information amount commonly contained in and the statistic dependence between two variables. So, (31) and (32) imply that the quantized system is estimable if the received data contain the information of any element of the state, which can be extracted by proper estimator. It is known that Ið U ; U Þ Z 0 with the equality holds if and only if these two variables are mutual independent. So (31) also means that any direction of state space is not orthogonal to all the past decoded data, i.e. 8 g A ℝn ; g a0; g T EfxðkÞy~ d ðjÞg a0; (j r k; 8 k Zk0 þn 1:

The one-to-one property of decoding implies that state estimation given quantized innovations is equivalent to that given the decoded signals. The quantized estimability is then discussed under the following assumption.

ð33Þ

From (33) the following result is gotten. Theorem 1. Under Assumption 1, the quantized Gaussian system (27), (29), (30) is estimable if and only if q

rank fW e ½k0 ; kg ¼ n; 8 k Z k0 þn  1;

ð34Þ

where the quantized estimability Gramian is k

at the channel receiver before used to estimate, where Dk ð UÞ is a one-to-one decoding mapping. A common decoding model in accord with (30) is [13] Dk ðy~ q ðkÞÞ ¼ ~ EfyðkÞj y~ q ðkÞ ¼ zl ðkÞg; l ¼ 1; 2; …; N. In our system stated in Section 2, the one-to-one mapping is implied in the MLQKF [3]: qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Dk ðy~ q ðkÞÞ ¼ f ðaL ; y~ q ðkÞÞ CPðkjk  1ÞC T þ R:

k

W qe ½k0 ; k ¼ ∑ ψ 2 ðjÞΦðk; jÞΠ ðjÞC T ðjÞCðjÞΠ ðjÞΦ ðk; jÞ; T

ð35Þ

j ¼ k0

with 0 0 1 Dj ðzl ðjÞÞ B B al ðjÞ þ CðjÞx^ ðjjj 1Þ C ψ ðjÞ ¼ ∑ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi@ϕ@qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiA l¼1 CðjÞΠ ðjÞC T ðjÞ þ RðjÞ CðjÞΠ ðjÞC T ðjÞ þ RðjÞ N

0

11

^ 1ÞCC Ba 1 ðjÞ þ CðjÞxðjjj  ϕ@ql þffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi AA; T CðjÞΠ ðjÞC ðjÞ þ RðjÞ

ð36Þ

ϕð U Þ is the probability density function of a standard Gaussian variable as in Section 2.3.

H. Zhang et al. / Signal Processing 105 (2014) 363–375

Proof. From Eqs. (27), (29) and (30) we have EfxðkÞy~ d ðjÞg ¼ EfðΦðk; jÞxðjÞ kj1

þ ∑ Φðk; jþ i þ 1ÞBðj þ iÞwðj þ iÞÞy~ d ðjÞg i¼0

¼ Φðk; jÞEfxðjÞy~ d ðjÞg:

ð37Þ

Let x ¼ xðjÞ, x^ ¼ x^ ðjjj 1Þ, y~ d ¼ y~ d ðjÞ, y^ ¼ y^ ðjjj  1Þ, R ¼ RðjÞ for notation simplicity, then Z N EfxðjÞy~ d ðjÞg ¼ ∑ x UDj ðzl Þ U pðx; y~ d ¼ Dj ðzl ÞÞdx l¼1

N

¼ ∑ Dj ðzl Þ l¼1

Z

x A ℝn

x A ℝn

x Upðx; y~ d ¼ Dj ðzl ÞÞdx:

ð38Þ

Based on the Bayesian law pðx; y~ d ¼ Dj ðzl ÞÞ ¼ pðy~ d ¼ Dj ðzl ÞjxÞpðxÞ ¼ pðal o y~ r al þ 1 jxÞpðxÞ ¼ pðal o Cx þv C x^ ral þ 1 jxÞpðxÞ ¼ pððal þ C x^  CxÞ o v r ðal þ 1 þ C x^ CxÞjxÞpðxÞ ¼ pððal þ C x^  CxÞ o v r ðal þ 1 þ C x^ CxÞÞpðxÞ

x A ℝn

o v rðal þ 1 þ C x^  CxÞÞpðxÞdx      Z ^ Cx a þC x a þ C x^ Cx T l pffiffiffi Gð0; Π Þdx ¼ x T l þ 1 pffiffiffi R R x A ℝn ð39Þ

where Tð UÞ is the probability distribution function value of a standard Gaussian variable with zero mean and covariance 1 as in Section 2.3, Gð0; Π Þ denotes the probability density function which means that the stochastic vector x is Gaussian with zero-mean and covariance matrix Π (Π here denotes Π ðjÞ for short). Compute the integral in (39) by variable replacements, from (38) we have EfxðjÞy~ d ðjÞg

0 0

1

N Dj ðzl Þ B B al þ CðjÞx^ ðjjj  1Þ C ¼ ∑ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi@ϕ@qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi A T l¼1 CðjÞΠ ðjÞC T ðjÞ þRðjÞ CðjÞΠ ðjÞC ðjÞ þ RðjÞ 0 11 ^ a þ CðjÞ x ðjjj  1Þ B lþ1 CC T  ϕ@qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi AAΠ ðjÞC ðjÞ: CðjÞΠ ðjÞC T ðjÞ þRðjÞ

Then (37) shows that (33) is equivalent to 0 0 1 N ^ D ðz ðkÞÞ a ðkÞ þ CðjÞ x ðjjj  1Þ j B B C l l g T ∑ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi@ϕ@qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiA T T l¼1 CðjÞΠ ðjÞC ðjÞ þRðjÞ CðjÞΠ ðjÞC ðjÞ þRðjÞ 0 11 ^ 1ÞCC Bal þ 1 ðkÞ þCðjÞxðjjj ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi AAΦðk; jÞΠ ðjÞC T ðjÞ a0;  ϕ@ q T CðjÞΠ ðjÞC ðjÞ þ RðjÞ ( jr k; 8 k Zk0 þ n 1; 8 g A ℝn ; g a0: Thus (34) is obtained.

can be preserved under quantization, even if the quantizer is 1-bit [1,44,45]. To see this, suppose a time-varying quantizer is designed as qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Dj ðzl ðjÞÞ ¼ cl CðjÞΠ ðjÞC T ðjÞ þ RðjÞ; l ¼ 1; 2; …; N qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ^ 1Þ; l ¼ 1; 2; …; N þ 1 al ðjÞ ¼ dl CðjÞΠ ðjÞC T ðjÞ þRðjÞ CðjÞxðjjj where cl and dl are chosen to make ψ ðjÞ nonzero constant, N i.e. ψ ðjÞ ¼ ψ ¼ Σ l ¼ 1 cl ðϕðdl Þ  ϕðdl þ 1 ÞÞ a 0; j ¼ k0 ; k0 þ1; …. Then comparing (35) with (28) we know such a quantizer does not change the estimability of the original system. Especially, when N ¼ 2 in the above equations, let c1 p ¼ffiffiffiffiffiffiffiffi  1,ffi c2 ¼ 1, d1 ¼  1, d2 ¼ 0, d3 ¼ þ 1, then ψ ðjÞ  2=π . Namely a 1-bit coarse quantizer can preserve the estimability. k ^ Remark 2. If the estimation xðkÞ given Y~ d is optimal, k T T~ ~ ~ ^ then the Hðg xðkÞjY d Þ ¼ Hðg xðkÞÞ, where xðkÞ ¼ xðkÞ  xðkÞ, quantized system is estimable if the posterior error entropy is strictly less than the prior entropy:

where the last equality holds because xðjÞ and vðjÞ are mutually independent. Then Z Z x U pðx; y~ d ¼ Dj ðzl ÞÞdx ¼ x U pððal þ C x^  CxÞ x A ℝn

369

□

Remark 1. Comparing (35) with (28) shows that the quantized estimability is defined by the quantizer and the intrinsic property of the original system jointly. When ψ ðjÞ a0; 8 j Zk0 , the estimability of the original system

~ o Hðg T xðkÞÞ; 8 k Z k0 þ n  1; g A ℝn ; g a0: Hðg T xðkÞÞ

ð40Þ

This condition is an information theoretic counterpart of that in [20], i.e. the optimal posterior mean-square estimation error is strictly smaller than the prior error. If k ~ x^ ðkÞ is not optimal, Hðg T xðkÞÞ 4Hðg T xðkÞjY~ d Þ, then (40) is a sufficient condition for practical consideration. Remark 3. There are three supplementary explanations to the applicability of Theorem 1 to the system in Sections 2 and 3 where the design method is proposed: (1) Because of the one-to-one property of the decoder (30), k k Iðg T xðkÞ; Y~ Þ ¼ Iðg T xðkÞ; Y~ Þ, so (31) is equivalent to d k

q

Iðg T xðkÞ; Y~ q Þ 40 which can be applied to the analysis for k the case that the state is estimated by using Y~ q (as in Sections 2 and 3) without decoding (30); (2) The designed procedure in Section 3 provides innovations under a power level (15), while Theorem 1 is derived for any innovation power; (3) The reduction matrix in (2) of the system in Sections 2 and 3 is designable, while the output matrix in (27) of the framework in this section is not, this does not affect the estimability analysis. In fact the Gramian (35) involves two key items playing roles in quantized estimability, the C concerning the dimension reduction, and the ψ concerning the quantization. The estimability yields constraints on the design of C and the quantizer. If no quantization is involved, the constraints will be reduced to constraint on the dimension of C, because when C is designed under the index of (16) or (23), the information loss is minimal in the procedure of dimension reduction.

Remark 4. Gramian (35) based on quantized innovations is slightly different from that based on quantized outputs [23]: k

W qe n ½k0 ; k ¼ ∑ ψ n 2 ðjÞΦðk; jÞΠ ðjÞC T ðjÞCðjÞΠ ðjÞΦ ðk; jÞ j ¼ k0

T

ð41Þ

370

H. Zhang et al. / Signal Processing 105 (2014) 363–375

with D0 ðz0 ðjÞÞ

N

0 0

1 a0 ðjÞ

j l B B C l ψ n ðjÞ ¼ ∑ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi @ϕ@qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi A T T l¼1 CðjÞΠ ðjÞC ðjÞ þ RðjÞ CðjÞΠ ðjÞC ðjÞ þ RðjÞ

0

11 a0l þ 1 ðjÞ B CC  ϕ@qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiAA T CðjÞΠ ðjÞC ðjÞ þRðjÞ where D0j , z0l ðjÞ, a0l ðjÞ are respectively the decoding, quantizer outputs and levels corresponding to system outputs. However, if the quantizer (29) and decoder (30) are applied for the both cases, a0l ðjÞ can be replaced by ^ al ðjÞ þ CðjÞxðjjj 1Þ, so can D0j , z0l ðjÞ. It can be seen from the following equivalent form of ψ ðjÞ (42) that the Gramian (35) is equivalent to (41). So under the same quantization and decoding, the estimability based on quantized innovations is equivalent to that based on quantized outputs. 4.2. Convergence Let Δ ¼

sup 1 r l r N;k Z k0

Δl ðkÞ, where Δl ðkÞ ¼ al þ 1 ðkÞ  al ðkÞ.

Consider the Gramian in Theorem 1 when Δ-0. Theorem 2. Under Assumption 1, the estimability Gramian (35) of quantized system (27), (29) and (30) converges to the estimability Gramian (28) of unquantized system (27) when Δ-0. Proof. We have 0 0 1 ^ ^ CðjÞxðjjj  1Þ B B al ðjÞ þ CðjÞxðjjj  1Þ C ∑ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi@ϕ@qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi A l¼1 CðjÞΠ ðjÞC T ðjÞ þ RðjÞ CðjÞΠ ðjÞC T ðjÞ þ RðjÞ 0 11 ^ a ðjÞ þ CðjÞ x ðjjj 1Þ B CC 1  ϕ@ql þffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi AA CðjÞΠ ðjÞC T ðjÞ þ RðjÞ 0 0 1 ^ ^ CðjÞxðjjj 1Þ a ðjÞ þ CðjÞ x ðjjj 1Þ B B 1 C ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi@ϕ@qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiA CðjÞΠ ðjÞC T ðjÞ þ RðjÞ CðjÞΠ ðjÞC T ðjÞ þRðjÞ 0 11 ^ a ðjÞ þCðjÞ xðjjj 1Þ B Nþ1 CC  ϕ@ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi AA ¼ 0 CðjÞΠ ðjÞC T ðjÞ þRðjÞ N

11

^ 1ÞCC Ba 1 ðjÞ þ CðjÞxðjjj  ϕ@ql þffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi AA T CðjÞΠ ðjÞC ðjÞ þ RðjÞ

Δl -0

     znl ðjÞ a~ ðjÞ a~ ðjÞ þ Δ ffiffiffiffi ϕ pl ffiffiffiffi  ϕ l pffiffiffiffi l ¼ ∑ lim p U U U l ¼ 1Δl -0    pffiffiffiffi  pffiffiffiffi  n 1 zl ðjÞ ϕ a~ l ðjÞ= U  ϕ a~ l ðjÞ þ Δl = U Δ  ffiffiffiffi :  pffiffiffil ffi ¼ ∑ lim p pffiffiffiffi U U  Δl = U l ¼ 1Δl -0   n 1 zl ðjÞ d Δ ffiffiffiffi ϕðsÞjs ¼ a~ l ðjÞ=pffiffiffiU  pffiffiffil ffi ¼ ∑ lim p ds U U Δ -0 l l¼1 n   n pffiffiffiffi2 o   1 znl ðjÞ 2 exp  1=2 zl ðjÞ= U Δ pffiffiffiffiffiffi pffiffiffil ffi ¼ ∑ lim pffiffiffiffi U U 2π l ¼ 1Δl -0 n  o    p ffiffiffi ffi 2 Z  n 2 exp  1=2 zn ðjÞ= U zl ðjÞ dznl ðjÞ l pffiffiffiffiffiffi pffiffiffiffi pffiffiffiffi : ¼ U U 2π R pffiffiffiffi pffiffiffiffi n n Denote r ¼ zl ðjÞ= U , dr ¼ dzl ðjÞ= U . Then lim ψ ðjÞ ¼ Δ-0 R 2  r2 =2 pffiffiffiffiffiffi = 2π Þdr ¼ 1, which is in fact an Euler integral R r ðe of Gamma function [36]. Thus complete the proof. □ 1

Theorem 2 coincides in some sense with the conclusion in [3] that the MLQ-KF convergences to the standard Kalman filter when the quantization intervals go to zero. 5. Simulation Simulations of different models show that the estimation performance of the proposed method is satisfactory. Analytic results on estimability given quantized innovations are also verified by simulations. In this section, only the simulation results of 3 models are stated. The design method of dimension reduction and quantizer is illustrated mainly by Model 1 and Model 2, and the analysis of estimability is illustrated further by Model 3. These models are simulated with different power constraints (W) and different quantization levels (L for the symmetrical quantizer (4) or N for the general quantizer (29)). 5.1. Model 1 and Model 2

because ϕð 1Þ-0, ϕð þ 1Þ-0 for a1 ¼  1, aN þ 1 ¼ þ 1. Then ψ ðjÞ (36) can be re-expressed as 0 0 1 N D ðz ðjÞÞ þ CðjÞx ^ ^ ðjjj  1Þ B Bqaffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi l l ðjÞ þCðjÞxðjjj 1Þ C ψ ðjÞ ¼ ∑ qj ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi @ϕ@ A l¼1 CðjÞΠ ðjÞC T ðjÞ þ RðjÞ CðjÞΠ ðjÞC T ðjÞ þ RðjÞ 0

Noting N-1, Δl -0 when Δ-0, from (42) we have      N zn ðjÞ a~ ðjÞ a~ þ 1 ðjÞ l ffiffiffiffi ffiffiffiffi ϕ pl ffiffiffiffi  ϕ lp lim ψ ðjÞ ¼ lim ∑ p N-1 l ¼ 1 U U U Δ-0

ð42Þ

Denote a~ l ðjÞ ¼ al ðjÞ þ CðjÞx^ ðjjj  1Þ, then Δl ðjÞ ¼ a~ l þ 1 ðjÞ  a~ l ðjÞ. qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffi Let U ¼ CðjÞΠ ðjÞC T ðjÞ þ RðjÞ, znl ðjÞ ¼ Dj ðzl ðjÞÞ þ CðjÞx^ ðjjj  1Þ.

Model 1 are: 2  1:1 6 0:01 6 A¼6 4 0:275

described by (1) is from [37]. The parameters

0:06

0

3

1

0

0 0

1 0

07 7 7; 15

0

0

0

2 3 1 607 6 7 B ¼ 6 7; 405

Q ¼ 1;

R ¼ 1:

0

Model 2 from [38] is as xðk þ1Þ ¼ AxðkÞ þ BuðkÞ þ Γ wðkÞ; where u(k) is the control input which does not affect the state estimation. We then let uðkÞ ¼ 0 and Γ ¼ B in accord with (1). The parameters of Model 2 are: 2 3 0:3718  0:3111 0:1357  0:0557 6 0:9575 0:1524 0:0685 0:1043 7 6 7 A¼6 7; 4 0:3408  0:3287 0:3923  0:2060 5 0:0529

 0:2922

0:9320

0:4835

H. Zhang et al. / Signal Processing 105 (2014) 363–375

2

1:2800 6 1:6400 6 B¼6 4 1:0400

1:0400

0:8400

2:6000

3

1:4000 7 7 7; 2:4000 5

Q ¼ I;

R ¼ 1:

These two models are mainly used to illustrate the algorithms in Section 3.2 which takes the symmetrical model (4) as the quantizer and the MLQ-KF as the predictor and filter. MC simulation method is employed to illustrate the design method which involves numerical technique, the SQP. The active-set algorithm is selected here when using the Matlab function “fmincon”. Set delta¼0.00001 in the algorithm. By using “fmincon”, the optimization depends on the initial setting of the normalized quantizer thresholds al ; l ¼ 1; 2; …; L, in (22). We set the initial values of these standard Gaussian variables randomly under the constraints in (22). Simulation results illustrate our design method. For the case of L ¼ 3, 200 optimization results of quantizer thresholds are shown in Fig. 2(a), where ani denote the optimization results of ai , the star points are defined by values of triple {a1 ; a2 ; a3 }. We evaluate the corresponding estimation performance by using the trace of the steady estimation error covariance P s . For Model 1, 200 simulation results of trðP s Þ in the case of W ¼ 3 and

371

L ¼ 3, corresponding to the quantizer thresholds in Fig. 2(a), are shown in Fig. 2(b). Results of one simulation realization for different power constraints and quantization levels are also shown respectively in Tables 1 and 2, and Figs. 3 and 4. Along with the trðP s Þ, time average square errors defined as follows are also given: ðΔxi Þ2 ¼

1 T ∑ ðxi ðkÞ  x^ i ðkÞÞ2 Tk¼1

where xi ðkÞ and x^ i ðkÞ are respectively the ith elements of the actual states and estimates, T denotes simulation time. Data in Tables 1 and 2 are Δx2i and trðP s Þ in the cases of different W and L (T¼100). It is seen that the square error decreases when signal power or quantization level increases. This is consistent with the intuition. The data show that the performance of L ¼ 3 is close to that of L ¼ 1 (i.e. the signals are not quantized) when the transmitted power remains invariant. This illustrates the advantage of the method. Performances of cases of W ¼ 3, L ¼ 1 and W ¼ 6, L ¼ 1 are very close. This suggests that the performance decrease caused by lower quantization level can be compensated by higher power. Only the 1st state elements of these models in the case of W¼ 3 and L¼1 are shown by figures for saving text space. Figs. 3(a) and 4(a) show respectively the actual states by solid lines and the state estimates by dotted lines, while Figs. 3(b) and 4(b) are estimation errors, where the “estimation error of x1 ” denotes x1 ðkÞ  x^ 1 ðkÞ. Simulations of Models 1 and 2 also illustrate our conclusion of estimability. In Model 1, for example, the obtained reduction matrix is C ¼ 0:9556 0:0486 0:1174 0:0282 when W¼ 3 and L ¼1. All the eigenvalues of A are inside the unit circle. For the stationary system {A, B, C, Q, R}, computation shows that rankð½Π C T AΠ C T A2 Π C T A3 Π C T Þ ¼ 4, where Π is the static state covariance, which means that this unquantized system is estimable [20]. We get ψ ðjÞ a0; j ¼ 0; 1; … by computing (36), so theoretically this quantized system is Table 1 Simulation results of model 1.

W ¼ 3, W ¼ 3, W ¼ 3, W ¼ 6, W ¼ 6, W ¼ 6,

L¼1 L¼3 L¼1 L¼1 L¼3 L¼1

Δx21

Δx22

Δx23

Δx24

trðP s Þ

0.8487 0.5876 0.4147 0.4631 0.2212 0.1734

0.0255 0.0235 0.0242 0.0260 0.0139 0.0112

0.0324 0.0263 0.0255 0.0257 0.0142 0.0113

0.0017 0.0015 0.0013 0.0014 0.0007 0.0005

1.0566 0.6976 0.6099 0.5869 0.3040 0.2368

Table 2 Simulation results of Model 2.

Fig. 2. (a) Thresholds, L ¼ 3 and (b) trðP s Þ, Model 1, W ¼3, L ¼3.

W ¼ 3, W ¼ 3, W ¼ 3, W ¼ 6, W ¼ 6, W ¼ 6,

L¼1 L¼3 L¼1 L¼1 L¼3 L¼1

Δx21

Δx22

Δx23

Δx24

trðP s Þ

1.5055 1.4906 1.0252 1.1475 0.9749 0.9009

2.9696 2.5241 2.3152 2.4792 2.3826 1.8176

3.7679 3.1872 2.1876 2.5080 2.4776 1.8104

6.4590 5.2531 3.0971 4.1069 3.3500 2.4057

17.8801 13.9071 12.4909 12.4848 8.2446 7.0180

372

H. Zhang et al. / Signal Processing 105 (2014) 363–375

actual value and estimate of x1

6 actual value estimate

4 2 0 −2 −4 −6 −8

0

20

40

60

80

100

[29] Hðχ Þ ¼ ð1=2Þln 2π e þð1=2Þln Χ . In quantized systems, ~ the probability distribution of g T xðkÞ is undefined but can be supposed to make its entropy maximal according to the well known Jaynes’ “maximum entropy principle” [39,26]. The maximum entropy distribution is Gaussian when the ~ covariance is given, the posterior error entropy Hðg T xðkÞÞ in simulation can also be calculated by using that as the

above formula. When g T ¼ 1 0 0 0 , condition (40) means Hðx1 Þ 4 Hðx~ 1 Þ which is shown in Fig. 3(c), where the “entropy difference of x1 ” is Hðx1 Þ  Hðx~ 1 Þ. This implies that the estimability of the original system (Model 1) is preserved under the designed framework proposed in Sections 2 and 3, even if a coarse quantizer (L¼1) is used. For Model 2, we have also verified estimability results for different cases. For example, Fig. 4(c) shows that Hðx1 Þ 4 Hðx~ 1 Þ in the case of W¼3 and L ¼1.

k

5.2. Model 3

estimation error of x1

3

Model 3 from [4] is simulated to verify the analysis of quantized estimability. Its parameters (as in (1) and (2)) are: 2 3 2:4744  2:811 1:7038  0:5444 0:0723 6 1 0 0 0 0 7 6 7 6 7 6 1 0 0 0 7 A¼6 0 7; 6 7 0 0 0 1 5 4 0

2

1

0

−1

2

−2

−3

0

20

40

60

80

100

0 3

0

1 607 6 7 6 7 07 B¼6 6 7; 607 4 5 0

C¼

Q ¼ 1;

1 : 16

h

0:245

1

0:236

0

0:384

0

0:146

i 0:035 ;

k 1.4

entropy difference of x1

1.2 1 0.8 0.6 0.4 0.2 0

0

20

40

60

80

100

k Fig. 3. (a) Model 1, W ¼ 3, L ¼ 1, (b) Model 1, W ¼3, L ¼1, and (c) Model 1, W ¼ 3, L ¼1.

also estimable according to Remark 1. To illustrate this computing conclusion, we use the condition (40). The entropy of a Gaussian variable χ with covariance Χ is

R¼

Different from that in the designed algorithm, here we choose the filter of Gaussian fit algorithm [13] and the Max–Lloyd quantizer [8] to illustrate the general conclusion. Cases of quantization level N ¼3 and 2 in (28) are shown. The quantizer outputs are respectively σ y~ ðkÞ {  1.2240, 0, 1.2240} and σ y~ ðkÞ{  0.7979, 0.7979} for ~ N ¼3 and 2, where σ yðkÞ is the standard deviation of yðkÞ. ~ By [20] we know Model 3 is estimable without quantization. Computation shows ψ ðjÞ a0; j ¼ 0; 1; 2; …, for both cases of N¼ 3 and 2 (Figs. 5(a) and 6(a)). So, the estimability is not changed by quantization according to Remark 1. Figs. 5(b) and 6(b) show respectively the actual values of x1 ðkÞ by solid lines and its estimates x^ 1 ðkÞ by dotted lines. Figs. 5(c) and 6(c) are differences between the prior and posterior entropies. The “entropy difference of x1 ”¼ Hðx1 ðkÞÞ Hðx~ 1 ðkÞÞ, where Hðx~ 1 ðkÞÞ of estimation errors x~ 1 ðkÞ are also computed under the “maximum entropy principle” as above. We observe that all of the entropy differences are greater than zero, i.e. the posterior error entropies are strictly smaller than prior entropies. This indicates the computation conclusion that the quantized systems with the 3-level quantizer and the 1-bit quantizer are both

H. Zhang et al. / Signal Processing 105 (2014) 363–375

0.7

5 actual value estimate

0.65

3

0.6

2

0.55 0.5

1

psi

actual value and estimate of x1

4

0

0.45

−1

0.4

−2

0.35

−3

0.3

−4

0.25

−5

373

0.2 0

20

40

60

80

100

0

20

40

80

100

15

3

actual value estimate

actual value and estimate of x1

2

estimation error of x1

60

k

k

1 0 −1 −2

10

5

0

−5

−3 −4

−10 0

20

40

60

80

100

0

20

40

k

80

100

60

80

100

1.4

0.4

1.2

entropy difference of x1

0.35

entropy difference of x1

60

k

0.3 0.25 0.2 0.15

1 0.8 0.6 0.4

0.1 0.2 0.05 0

0 0

20

40

60

80

100

0

40

k

k Fig. 4. (a) Model 2, W ¼ 3, L ¼ 1, (b) Model 2, W¼ 3, L ¼ 1, and (c) Model 2, W ¼3, L ¼ 1.

20

Fig. 5. (a) Model 3, N¼ 3, ψ, (b) Model 3, N ¼3, and (c) Model 3, N ¼3.

6. Conclusion estimable, according to (40), and also illustrated the Remark 4 concerning the estimability equivalence between quantizing innovation and quantizing outputs.

The present paper discussed the problem of joint design of state dimension reduction and quantizer for

374

H. Zhang et al. / Signal Processing 105 (2014) 363–375

0.6 0.55 0.5 0.45

psi

0.4 0.35 0.3 0.25 0.2 0.15 0.1

0

20

40

60

80

100

k

15

actual value and estimate of x1

actual value estimate

10

5

0

−5

−10

−15

0

20

40

60

80

100

k

1.4

Acknowledgment

1.2

entropy difference of x1

MMSE criterion and the MLQ-KF. To illustrate the validity of state estimation in the designed system, the state estimability given quantized innovations was analyzed based on the quantity of mutual information, along with the design. The quantized estimability Gramian was obtained, and its convergence property was also analyzed. It was shown that the estimability of the original system can be preserved even if a 1-bit quantizatier is involved. This provides an auxiliary support for the results on 1-bit quantized estimation [1,44,45]. The estimability equivalence between cases of quantized innovations and quantized outputs [23] was also discussed. Monte Carlo simulation results show that the performance of state estimation under communication constraints is satisfactory when the derived dimension reduction and quantizers are applied. The analytic conclusion on estimability was also verified by illustrative simulations. In the present paper, the design method was given and verified for the case of one subchannel. For the case of more than one parallel subchannels, it is suspected by the authors that the performance of estimation can be improved partly by considering the signal power allocation. On the other hand, as stated in Section 1, the scheme of the pre-processing part is similar to the GQS [42] and its special form of the NSQ [43]. A very interesting result is that the NSQ improves the accuracy of the estimates. Also, while the estimability equivalence between quantizing innovations and quantizing output observations was gotten in this paper, an important fact was exploited in [4,40] that quantizing innovations requires fewer bits than quantizing observations directly. These contributions provide valuable hints for further consideration of estimation performance of predictive quantizer or the DPCM-like scheme.

The authors thank the anonymous reviewers whose comments and suggestions are very valuable and helpful for improving the quality of our work.

1 0.8

References 0.6 0.4 0.2 0

0

20

40

60

80

100

k Fig. 6. (a) Model 3, N¼ 2, ψ,. (b) Model 3, N¼ 2, and (c) Model 3, N ¼ 2.

state estimation in quantized linear dynamic systems with a DPCM-like structure. For systems under constraints of limited signal power, number and bandwidth of the communication (sub) channels, the dimension reduction matrix and quantizer thresholds were obtained based on

[1] A. Ribeiro, G.B. Giannakis, S.I. Roumeliotis, SOI-KF: distributed Kalman filtering with low-cost communications using the sign of innovations, IEEE Trans. Signal Process. 54 (12) (2006) 4782–4795. [2] E.J. Msechu, S.I. Roumeliotis, A. Ribeiro, G.B. Giannakis, Decentralized quantized Kalman filtering with scable communication cost, IEEE Trans. Signal Process. 56 (8) (2008) 3727–3741. [3] K.-Y. You, L.-H. Xie, S.-L. Sun, W.-D. Xiao, Quantized filtering of linear stochastic systems, Trans. Inst. Meas. Control 33 (6) (2011) 683–698. [4] M.-Y. Fu, C.E. de Souza, State estimation for linear discrete-time systems using quantized measurements, Automatica 45 (12) (2009) 2937–2945. [5] R.-Q. Lu, Y. Xu, A.-K. Xue, H1 filtering for singular systems with communication delays, Signal Process. 90 (2010) 1240–1248. [6] M. Liu, J. You, X.-C. Ma, H1 filtering for sampled-data stochastic systems with limited capacity channel, Signal Process. 91 (2011) 1826–1837. [7] R.C. Farias, J.-M. Brossier, Adaptive quantizers for estimation, Signal Process. 93 (2013) 3076–3087. [8] J.G. Proakis, Digital Communications, 4th edition, The McGraw-Hill Companies, Singapore, 2001, 4 (115, 127).

H. Zhang et al. / Signal Processing 105 (2014) 363–375

[9] S. Gugercin, A.C. Antoulas, A survey of model reduction by balanced truncation and some new results, Int. J. Control 77 (8) (2004) 748–766. [10] Y. Zhu, E. Song, J. Zhou, Z. You, Optimal dimensionality reduction of sensor data in multisensory estimation fusion, IEEE Trans. Signal Process. 53 (5) (2005) 1631–1639. [11] N. Goela, M. Gastpar, Reduced-dimension linear transform coding of correlated signals in networks, IEEE Trans. Signal Process. 60 (6) (2012) 3174–3187. [12] T.V. Ramabadran, D. Sinha, Speech data compression through sparse coding of innovations, IEEE Trans. Speech Audio Process. 2 (2) (1994) 274–284. [13] R.E. Curry, Estimation and Control with Quantized Measurements, The M.I.T Press, London, 1970, 23–25. [14] R.M. Gray, D.L. Neuhoff, Quantization, IEEE Trans. Inf. Theory 44 (6) (1998) 2325–2383. [15] T. Arildsen, M.N. Murthi, S.V. Andersen, S.H. Jensen, On predictive coding for erasure channels using a Kalman framework, IEEE Trans. Signal Process. 57 (11) (2009) 4456–4466. [16] A. Wiesel, Y.C. Eldar, S. Shamai, Linear precoding via conic optimization for fixed MIMO receivers, IEEE Trans. Signal Process. 54 (1) (2006) 161–176. [17] D. Persson, J. Kron, M. Skoglund, E.G. Larsson, Joint source-channel coding for the MIMO broadcast channel, IEEE Trans. Signal Process. 60 (4) (2012) 2085–2090. [18] H. Sampath, P. Stoica, A. Paulraj, Generalized linear precoder and decoder design for MIMO channels using the weighted MMSE criterion, IEEE Trans. Commun. 49 (12) (2001) 2198–2206. [19] N. Khaled, S. Thoem, L. Deneire, Optimizing the joint transmit and receive MMSE design using mode selection, IEEE Trans. Commun. 53 (4) (2005) 730–737. [20] Y. Baram, T. Kailath, Estimability and regulability of linear systems, IEEE Trans. Autom. Control 33 (12) (1988) 1116–1121. [21] H. Zhang, Y.-X. Sun, Estimability of stochastic systems – an information theory approach, Control Theory Appl. 13 (5) (1996) 567–572. [22] A.R. Liu, R.R. Bitmead, Stochastic observability in network state estimation and control, Automatica 47 (2011) 65–78. [23] H. Zhang, Y. Shen, Estimability and stochastic observability of quantized linear systems, Int. J. Syst. Sci. (2014), http://dx.doi.org/ 10.1080/00207721.2014.919713. in press. [24] H.L. Weidemann, E.B. Stear, Entropy analysis of estimating systems, IEEE Trans. Inf. Theory IT-16 (3) (1970) 264–270. [25] Y. Tomita, S. Ohmatsu, T. Soeda, An application of information theory to estimation problems, Inf. Control 32 (2) (1976) 101–111. [26] P. Kalata, R. Priemer, Linear prediction, filtering, and smoothing: an information-theoretic approach, Inf. Sci. 17 (1) (1979) 1–14. [27] Y. Shen, H. Zhang, Joint states dimension reduction and quantizer design for estimation under communication constraints, in: Proceedings of IEEE 2nd International Conference on Computational Problem - Soving, Chengdu, China, October 21–23, 2011. [28] Y. Shen, H. Zhang, Estimability of quantized linear systems, in: Proceedings of the 2011 International Conference on Wireless

[29] [30] [31] [32] [33] [34]

[35] [36]

[37]

[38] [39]

[40]

[41] [42] [43]

[44] [45]

[46]

[47] [48]

375

Communications and Signal Processing, Nanjing, China, November 9–11, 2011. T.M. Cover, J.A. Thomas, Elements of Information Theory, John Wiley & Sons, Inc., New Jersey, 2003. J.A. Betts, Signal Processing, Modulation and Noise, American Elsevier Pub. Co, New York, 1971. M. Akçakaya, V. Tarokh, Shannon-theoretic limits on noisy compressive sampling, IEEE Trans. Inf. Theory 56 (1) (2010) 492–504. S.M. Ross, A First Course in Probability, 5th edition, Prentice-Hall, New Jersey, US, 1998. E.I. Silva, G.C. Goodwin, D.E. Quevedo, Control system design subject to SNR constraints, Automatica 46 (2010) 428–436. J.S. Freudenberg, R.H. Middleton, V. Solo, Stabilization and disturbance attenuation over a Gaussian communication channel, IEEE Trans. Autom. Control 55 (3) (2010) 795–799. B.D.O. Anderson, J.B. Moore, Optimal Filtering, Prentice-Hall, Inc., New Jersey, 1979. F.W.J. Olver, D.W. Lozier, R.F. Boisvert, C.W. Clark, NIST Handbook of Mathematical Functions, Cambridge University Press, New York, 2010. D.A. Wagie, R.S. Skelton, A projection approach to covariance equivalent realizations of discrete systems, IEEE Trans. Autom. Control 31 (12) (1986) 1114–1120. M.S. Ahmed, A New Algorithm for State Estimation of stochastic discrete systems, IEEE Trans. Autom. Control 39 (8) (1994) 1652–1656. (a) E.T. Jaynes, Information theory and statistical mechanics I, Phys. Rev. 106 (4) (1957) 620–630; (b) E.T. Jaynes, Information theory and statistical mechanics II, Phys. Rev. 108 (2) (1957) 171–190. D. Marelli, K.-Y. You, M.-Y. Fu, Identification of ARMA models using intermittent and quantised output observations, Automatica 49 (2) (2013) 360–369. B.I. Godoy, et al., On the identification of FIR systems having quantised output data, Automatica 47 (9) (2011) 1905–1915. B.I. Godoy, et al., Identification of sparse FIR systems using a general quantisation scheme, Int. J. Control 87 (4) (2014) 874–888. J.C. Aguero, G.C. Goodwin, J.I. Yuz, System identification using quantised data, in: Proceedings of 46th Conference on Decision & Control, New Orleans, LA, December 12–14, 2007, pp. 4263–4268. L.Y. Wang, System Identification with Quantised Observations, Birkhauser, New York, US, 2010. L.Y. Wang, G.H. Xu, G. Yin, State reconstruction for linear timeinvariant systems with binary-valued output observations, Syst. Control Lett. 57 (2008) 958–963. F. Gustafsson, R. Karlsson, Statistical results for system identification based on quantized observations, Automatica 45 (12) (2009) 2794–2801. I. Griva, S.G. Nash, A. Sofer, Linear and Nonlinear Optimization, 2nd edition, SIAM, Philadelphia, US, 2009. J.F. Bonnans, et al., Numerical Optimization: Theoretical and Practical Aspects, 2nd edition, Springer, Berlin, 2008.