On equivalence of major relaxation methods for minimum ellipsoid covering intersection of ellipsoids

Automatica 103 (2019) 337–345

Contents lists available at ScienceDirect

Automatica journal homepage: www.elsevier.com/locate/automatica

Brief paper

On equivalence of major relaxation methods for minimum ellipsoid covering intersection of ellipsoids✩ ∗

Zhiguo Wang, Xiaojing Shen , Yunmin Zhu Department of Mathematics, Sichuan University, Chengdu, Sichuan 610064, China

article

info

Article history: Received 27 May 2018 Received in revised form 7 November 2018 Accepted 15 January 2019 Available online xxxx Keywords: Set membership filter The intersection of ellipsoids S-procedure Bounding ellipsoid Multisensor fusion

a b s t r a c t This paper investigates the problem on the minimum trace or volume ellipsoid covering intersection of ellipsoids. This problem can be solved by three major relaxation methods involving semi-definite programming relaxation, S-procedure relaxation and parameterized bounding ellipsoid relaxation. They are proposed from the different ideas or viewpoints, so that it is difficult to judge which relaxation method is the most suitable for the different practical implementations. To the best of our knowledge, the problem of proving the convexity of the optimization problem derived by the third relaxation method remains open. By rigorous analysis, we reveal the equivalence among the three relaxation methods, and address the open problem, i.e., the optimization problem derived by the third relaxation method can be equivalently transformed to a convex optimization problem. The computation complexity of these relaxation methods are also discussed, and the equivalent result is helpful to set membership filter and multisensor fusion obtaining a tighter ellipsoid to contain the state with lower computation complexity. Finally, the performance of each method is evaluated by several typical numerical examples in information fusion and filtering. © 2019 Elsevier Ltd. All rights reserved.

1. Introduction Problems involving the intersection of ellipsoids occur in system control, information fusion, and target tracking, etc. Boyd, El Ghaoui, Feron, and Balakrishnan (1994), Schweppe (1968), Shen, Zhu, Song, and Luo (2011), Shi, Chen, and Lin (2015) and Wang and Rong Li (2012). In system control, when the uncertain system noise is unknown but bounded, Schweppe proposed the set membership filter (Schweppe, 1968), which is to obtain a parameterized bounding ellipsoid containing the true state. The key of set membership filter is determining a minimum ellipsoid containing the intersection of the predicted ellipsoid and the measurement ellipsoid (Becis-Aubry, Boutayeb, & Darouach, 2008; Joachim & Deller, 2006; Polyak, Nazin, Durieu, & Walter, 2004; Ros, Sabater, & Thomas, 2002; Wang, Shen, Zhu, & Pan, 2018). Recently, Sinyakov (2015) considers the problem of constructing exterior and interior approximations to the reachability sets for applied control. ✩ This work was supported in part by the open research funds of BACC-STAFDL of China under Grant No. 2013afdl011, the special funds of NEDD of China under Grant No. 201314 the NSFC No. 61673282 and the PCSIRT16R53. The material in this paper was not presented at any conference. This paper was recommended for publication in revised form by Associate Editor Andrea Garulli under the direction of Editor Torsten Söderström. ∗ Corresponding author. E-mail addresses: [email protected] (Z. Wang), [email protected] (X. Shen), [email protected] (Y. Zhu). https://doi.org/10.1016/j.automatica.2019.02.001 0005-1098/© 2019 Elsevier Ltd. All rights reserved.

In information fusion, when the cross-correlation of local sensor estimation errors is unknown or impractical (Chen, Ho, Zhang, & Yu, 2017; Gu, He, & Han, 2015; Noack, Sijs, Reinhardt, & Hanebeck, 2017; Wang, Shen, & Zhu, 2017), the covariance intersection (CI) algorithm is derived to fuse the local estimates (He, Xue, & Fang, 2018; Julier & Uhlmann, 2009; Uhlmann, 1996; Zhang & Shi, 2018). The geometric interpretation of CI approach is that finding a minimum ellipsoid to enclose the intersection region of local estimated ellipsoids (Uhlmann, 1996). In fact, it is a nonconvex optimization problem for determining the minimum trace or volume ellipsoid containing the intersection of several ellipsoids, since we need to maximize a convex function over a convex set. To the best of our knowledge, this problem has been considered in many papers and books (Boyd et al., 1994; Calafiore & El Ghaoui, 2004; Chernousko, 1993; Durieu, Walter, & Polyak, 2001; Kurzhanski & Mesyats, 2013; Kurzhanski & Vâlyi, 1997; Luo, Ma, So, Ye, & Zhang, 2010), and references therein. They provide different approximate solutions by different ideas or viewpoints. In Boyd et al. (1994), the authors apply the S-procedure to approximate the original problem to a semi-definite programming (SDP) which can be solved by the interior point methods (Nesterov & Nemirovski, 1994). Since the ellipsoidal approximation of the intersection of the ellipsoids involves solving nonconvex quadratically constrained quadratic programs, the semi-definite relaxation (SDR) technique (Luo et al., 2010) can be used to deal with this

338

Z. Wang, X. Shen and Y. Zhu / Automatica 103 (2019) 337–345

problem. In Henrion, Tarbouriech, and Arzelier (2001), the authors survey various linear matrix inequality relaxation techniques for evaluating the maximum norm vector within the intersection of several ellipsoids. However, it is different from the problem of minimum trace or volume ellipsoid containing the intersection of multiple ellipsoids. In fact, the described problem in Henrion et al. (2001) is similar to the Chebyshev center of a convex set (Eldar, Beck, & Teboulle, 2008; Wu, Zhou, & Hu, 2013). In Durieu et al. (2001), Polyak et al. (2004) and Ros et al. (2002), the authors employ the parameterized bounding ellipsoid relaxation method (PBER) to contain the intersection of ellipsoids, but the problem of proving the convexity of the optimization problem by the PBER remains open (Polyak et al., 2004). Although these relaxation methods can obtain the outer ellipsoids to contain the intersection, respectively, it is unclear for the interrelationships among these methods on the tightness and computational complexity. The goal of this paper is to clarify the relationship among the three typical relaxation methods involving SDP relaxation, S-procedure relaxation and PBER, and insight to the pros and cons of these methods. Our contribution mainly lies in the following three aspects.

• By rigorous analysis, we reveal the equivalence among these relaxation methods through the proposed decoupled SDP relaxation method, in the sense that the three methods provide the same outer-bounding ellipsoid containing the intersection of ellipsoids. • The open problem on the convexity of the optimization problem of the PBER is addressed, i.e., it can be equivalently transformed to a convex optimization problem. • The equivalence result is helpful to set-membership filter and sensor fusion for obtaining a tighter bounding ellipsoid with lower computation complexity. The structure of the paper is as follows. In Section 2, we describe the problem on the minimum ‘‘size’’ ellipsoid containing the intersection of multiple ellipsoids. Section 3 presents the SDP relaxation method, the S-procedure relaxation method and the parameterized bounding ellipsoid relaxation method. The equivalence among the three relaxation methods is established in Section 4. Section 5 shows the significance of the equivalence. Simulations are given in Section 6. Conclusions and potential research are drawn in Section 7. Notations: For a matrix X , X ⪰ 0 means X is positive semidefinite matrix. The superscript ‘‘T ’’ means the transpose. X + denotes the (Moore–Penrose) pseudo-inverse of X . ∥X ∥ denotes the standard Euclidean norm of X . diag(·) stands for a block diagonal matrix. Ellipsoid is described as E = {x ∈ Rn : (x − xˆ )T P −1 (x − xˆ ) ≤ 1}, where xˆ and P are the center and shape matrix of the ellipsoid E . The ‘‘size’’ of the ellipsoid is the function of shape matrix P, and it is denoted by f (P). Throughout this paper, the ‘‘size’’ of the ellipsoid is either trace(P), which means the sum of squares of semi-axes lengths of the ellipsoid E , or logdet(P), which corresponds to the volume of the ellipsoid E . 2. Problem statement

P0 ,x0

s.t. (x − x0 )T P0−1 (x − x0 ) ≤ 1, ∀x ∈ F ,

F =

(1) (2)

where the vector x0 ∈ Rn and the symmetric positive-definite matrix P0 ∈ Rn×n are the decision variables. The set F denotes

m ⋂

Ei ,

(3)

i=1

−1

Ei = {x ∈ Rn : (x − xi )T Pi

(x − xi ) ≤ 1},

(4)

i = 1, . . . , m. Actually, the constraint (2) has infinite inequalities, which enforces that the intersection of these ellipsoids is contained in the optimized ellipsoid. Thus, the optimization problem (1)–(2) is also called semi-infinite optimization problem (Boyd & Vandenberghe, 2004), and the computational complexity of the optimization problem ⋂m (1)–(2) is NP-hard. In addition, simply verifying that E0 ⊃ i=1 Ei holds, given E1 , . . . , Em , is NP-complete (Boyd et al., 1994). In next sections, we do not expect to get the optimal solution of the problem (1)–(2) in polynomial time, but provide efficient relaxation methods for approximation solutions. We concentrate on clarifying the relationship among the following relaxation methods, and analyze the pros and cons of these methods. 3. Three major relaxation methods For solving the problem (1)–(2), there are three major relaxation methods involving SDP relaxation, S-procedure relaxation and parameterized bounding ellipsoid relaxation, which are derived from different ideas or viewpoints. 3.1. SDP relaxation Let X = xxT , then (4) can be written equivalently as max f0 (X , x)

(5)

s.t. (X , x) ∈ F ,

(6)

where F = {(X , x) : fi (X , x) ≤ 0, 1 ≤ i ≤ m, X = xxT },

and, for 0 ≤ i ≤ m, fi (X , x) = trace(Pi−1 X ) − 2xTi Pi−1 x + xTi Pi−1 xi − 1. The objective function (5) is linear in decision variable (X , x), but the set F is not convex due to X = xxT . In order to obtain a relaxation of the optimization problem (5)–(6), usually F is replaced by the following convex set T = {(X , x) : fi (X , x) ≤ 0, 1 ≤ i ≤ m, xxT ⪯ X }.

Based on the relaxed convex set T , a suboptimal ellipsoid can be obtained to contain the intersection F by the SDP relaxation method as follows. Lemma 1 (SDP Relaxation Boyd & Vandenberghe, 2004). The optimization problem (1)–(2) can be relaxed to a convex optimization problem as follows min f (P0 )

(7)

s.t. λi ≥ 0, i = 1, . . . , m,

(8)

⎡

Consider the optimization problem of computing the minimum ellipsoid containing the intersection of multiple ellipsoids min f (P0 )

the intersection of m ellipsoids, which is defined as follows

P0−1 −

m ⎢ ∑ ⎢ λi Pi−1 ⎢ ⎢ i=1 ⎢ ⎢ ⎢ −˜xT + ⎢ 0 ⎢ ⎢ m ⎢∑ ⎢ λi xTi Pi−1 ⎢ ⎣ i=1

0

−˜x0 + m ∑ λi Pi−1 xi i=1

−1 − −

m ∑

m ∑

λi

0

⎤

x˜ T0

⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⪯0 ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

i=1

λi xTi Pi−1 xi

i=1

x˜ 0

−P0−1

(9)

Z. Wang, X. Shen and Y. Zhu / Automatica 103 (2019) 337–345

where x˜ 0 = P0−1 x0 . Moreover, the relaxed optimization problem (7)–(9) is convex in the variables P0−1 , x˜ 0 , λ1 , . . . , λm . Algorithm 1 SDP relaxation (1) Solve the optimization problem (7)–(9) to obtain the optimal solutions P0−1 and x˜ 0 . (2) Compute the shape matrix P0 and the center x0 = P0 x˜ 0 .

339

⎤

⎡

−1 −P0−1 x0 ⎦ ⎣ P0 −xT0 P0−1 xT0 P0−1 x0 − 1 ⎡ ⎤ m −1 −1 ∑ P −Pi xi ⎦ ⪯ τi ⎣ i T −1 T −1 −xi Pi xi Pi xi − 1 i=1

(18)

(19)

with variables P0 , x0 , τi , i = 1, . . . , m. 3.3. Parameterized bounding ellipsoid relaxation

As we know, this optimization problem (7)–(9) can be efficiently solved by interior point methods (Vandenberghe & Boyd, 1996). In addition, we can use the convex optimization toolbox CVX (Grant, Boyd, & Ye, 2009) to solve (7)–(9) in Matlab, and Algorithm 1 gives us a suboptimal ellipsoid that contains the intersection of ellipsoids E1 , . . . , Em . 3.2. S-procedure relaxation In system and control theory, one often encounters the constraint that a quadratic form is negative when other quadratic forms are all negative. In some cases, this constraint can be relaxed as a linear matrix inequality (LMI) based on the S-procedure method. It has long been known that the ellipsoids can be written as quadratic forms, therefore, Boyd et al. (1994) use the S-procedure relaxation technique to deal with the constraint (2). Now, a suboptimal ellipsoid for the intersection F can be derived by the S-procedure relaxation technique. Let ξ = [xT 1]T , then (4) can be equivalent to fi (X , x) = ξ T Ai ξ ≤ 0, 0 ≤ i ≤ m,

(10)

⎡

−1

P Ai = ⎣ i −xTi Pi−1

−1

(11)

xTi Pi−1 xi

f0 (X , x) = ξ A0 ξ ≤ 0 T

i=1

Let D be the convex set ∑ of all vectors t = (t1 , . . . , tm ) ∈ Rm , with m ti ≥ 0, i = 1, . . . , m, and i=1 ti = 1. If x ∈ F , then m ∑

ti (x − xi )T Pi−1 (x − xi ) ≤ 1, ∀ t ∈ D + .

(20)

i=1

After simple calculations and transformations, the following equation can be obtained ti (x − xi )T Pi−1 (x − xi )

i=1

= (x − xt )T Pt−1 (x − xt ) + δt ,

(21)

where (12)

whenever

Pt−1 =

m ∑

ti Pi−1

(22)

i=1

fi (X , x) = ξ T Ai ξ ≤ 0, 1 ≤ i ≤ m.

(13)

By S-procedure (Yakubovich, 1971), we can obtain E0 containing ⋂ m i=1 Ei if there exist τ1 , . . . , τm ≥ 0 such that

τi Ai .

(14)

xt = Pt

m ∑

ti Pi−1 xi

(23)

i=1

δt =

m ∑

ti xTi Pi−1 xi − xTt Pt−1 xt .

(24)

i=1

i=1

Using the definition of Ai in (11), we can write Eq. (14) as an LMI

⎡

m ⋂ {x : (x − xi )T Pi−1 (x − xi ) ≤ 1, i = 1, . . . , m}.

+

⎤

−Pi xi ⎦ , 0 ≤ i ≤ m. −1 ⋂m Thus, the problem E0 ⊃ i=1 Ei is equivalent to

A0 ⪯

F =

m ∑

where

m ∑

The authors (Durieu et al., 2001) divide the problem of computing the minimum trace or volume ellipsoid of the intersection of ellipsoids into two steps. First, one needs to find a parameterized ellipsoid to contain the intersection. Second, the optimal ellipsoid can be derived by optimizing the parameter. Assume that the intersection F of m ellipsoids is a nonempty bounded region, which is defined by (3), i.e.,

⎤

−1 −P0−1 x0 ⎦ ⎣ P0 −xT0 P0−1 xT0 P0−1 x0 − 1 ⎡ ⎤ m −1 −1 ∑ P − P x i i ⎦. ⪯ τi ⎣ i −xTi Pi−1 xTi Pi−1 xi − 1 i=1

(15)

Substituting (21) into (20) obtains F ⊆ Et , where the ellipsoid Et is denoted by Et = {x : (x − xt )T ((1 − δt )Pt )−1 (x − xt ) ≤ 1}. Therefore, we have found an ellipsoid to contain the intersection of multiple ellipsoids, but it depends on the parameter ti , i = 1, . . . , m. In order to derive an optimal ellipsoid Et , we need to minimize the function of the shape matrix (1 − δt )Pt . The optimization problem can be written as follows min f ((1 − δt )Pt )

Then, based on the S-procedure relaxation technique, the best such outer-bounding ellipsoid of the intersection F can be derived by solving the following optimization problem min f (P0 )

(16)

s.t. τi ≥ 0, i = 1, . . . , m,

(17)

(25)

s.t. (t1 , . . . , tm ) ∈ D , +

(26)

where the vector t = (t1 , . . . , tm ) is the optimization variable of this problem. Note that, in Polyak et al. (2004), Polyak et al. proposed that proving the convexity of the scalar functions trace((1 − δt )Pt ) and logdet((1 − δt )Pt ) (25) over the constraint (26) is an open problem.

340

Z. Wang, X. Shen and Y. Zhu / Automatica 103 (2019) 337–345

Based on the parameterized bounding ellipsoid relaxation method, Algorithm 2 derives an outer suboptimal ellipsoid E0 to contain the intersection F . Its shape matrix and the center are P0 and x0 , respectively. Algorithm 2 Parameterized bounding ellipsoid relaxation (1) Solve the optimization problem (25)–(26) to obtain the optimal solutions t ∗ . ∗

(2) Compute the center x0 = Pt

m ∑

− 1− ⎢ ⎢ ⎢ m ⎢∑ ⎢ λi (xTi Pi−1 xi − 1) ⎢ ⎢ i=1 ⎢ m ⎢ ∑ ⎢ λi xTi Pi−1 ⎢ ⎢ i=1 ⎢ ⎢ ⎣ −xT P −1 0 0

P0

4.1. The relationship between parameterized bounding ellipsoid relaxation and SDP relaxation In this subsection, to show the relationship between parameterized bounding ellipsoid relaxation and SDP relaxation, we first use a decoupled technique (Calafiore & El Ghaoui, 2004) to obtain an analytic formula of P0 and x0 in the convex optimization problem (7)–(9), respectively. Proposition 1. When the objective function is trace or logdet function, the optimization problem (7)–(9) based on SDP relaxation technique can be equivalently decoupled to

λi Pi−1

)−1

) (27)

i=1

s.t. λi ≥ 0, i = 1, . . . , m,

⎡ ⎢1−

m ∑

λi +

m ∑

⎢ i=1 ⎢ ⎢∑ m ⎢ ⎢ λi xTi Pi−1 xi ⎢ ⎢ i=1 ⎢ ⎢ ∑ m ⎣ λi Pi−1 xi i=1

⎤ λi xTi Pi−1 ⎥ ⎥ i=1 ⎥ ⎥ ⎥ ⎥ ⪰ 0, ⎥ ⎥ ⎥ ⎥ m ∑ −1 ⎦ λ i Pi

(28)

i=1

where λi , i = 1, . . . , m, are the optimization variables of this problem. If the above problem is feasible, then there exists an optimal solution for optimization problem (7)–(9). In this case, calling λ∗i optimal values of the variable λi , i = 1, . . . , m, then the optimal shape matrix and the center of problem (7)–(9) satisfy P0−1 =

m ∑

λ∗i Pi−1 ,

⎢ ⎢ xT ⎢ 0 ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣I

m ∑

I

λi xTi Pi−1 xi

i=1

1+ m ∑

m ∑

λi

i=1

λi xTi Pi−1

x0 = P 0

λ∗i Pi−1 xi .

⎤

m ∑

⎥ λi Pi−1 xi ⎥ ⎥ ⎥ i=1 ⎥ ⎥ ⎥ ⪰ 0, ⎥ ⎥ ⎥ m ⎥ ∑ −1 ⎦ λ i Pi

(31)

i=1

where I is the identity matrix with appropriate dimension. Based on Corollary A.1 (Decoupling) in Calafiore and El Ghaoui (2004) to (31), the optimization problem (7)–(9) can be equivalent to (27)–(30). Algorithm 3 Decoupled SDP relaxation (1) Solve the optimization problem (27)–(28) to obtain the optimal solutions λi , i = 1, . . . , m. (2) Compute the shape matrix P0 and the center x0 by (29)–(30).

It is interesting to see that the analytic expressions (29)–(30) by the decoupled SDP relaxation method are similar to Eqs. (22)–(23), which encourages us to prove the equivalence between the optimization problems (25)–(26) and (27)–(28). The optimization problem (27)–(30) is also a significant bridge for deriving the equivalence among three typical relaxation methods. Proposition 2. When the objective function is logdet or trace function, the optimization problem (25)–(26) based on the parameterized bounding ellipsoid relaxation technique is equivalent to the convex optimization problem (27)–(28) based on the decoupled SDP relaxation method, i.e., they have the same optimal outer-bounding ellipsoid to contain the intersection F . Proof. The proof has three steps. Firstly, when the objective function is logdet function, we ∑ continue to simplify the optimization m λi problem (27)–(28). Let η = i=1 λi and ti = η , i = 1, . . . , m, then ∑m i=1 ti = 1. Since η > 0 (Boyd et al., 1994), if the constraint (28) is scaled by a factor of η1 , then the optimization problem (27)–(28) can be rewritten as min − log det

m (∑

)

ti Pi−1 − n · log(η)

(32)

i=1

(29)

i=1 m ∑

i i

x0

i=1

m (∑

λ i P i xi ⎥ ⎥ ⎥ ⎥ −1 −P0 x0 ⎥ ⎥ ⎥ ⎥ ⪯ 0. ⎥ ⎥ −1 P0 ⎥ ⎥ ⎥ m ⎥ ∑ −1 ⎦ − λP i=1

i=1

⎡

Although the optimization problems (7)–(9), (16)–(19), and (25)–(26) have different forms, in this section, we clarify the relationship among the typical three relaxation methods.

min f

⎤ −1

By Schur complement again, it can be equivalent to

4. The equivalence of three relaxation methods

(

m ∑

ti Pi xi and the shape matrix

i=1

∗

xT0 P0−1 x0

∗ −1

P0 = (1 − δt )Pt by (22)–(24). ∗

⎡

(30)

i=1

Proof. According to Schur complement (Boyd et al., 1994) and the definition x˜ 0 = P0−1 x0 , the inequality (9) is equivalent to

⎡

1

m ∑

− ti ⎢ ⎢ η ⎢ m i=1 ⎢ ∑ ⎢ ti xT P −1 xi s.t. ⎢+ ⎢ i=1 i i ⎢ m ⎢ ∑ ⎣ ti Pi−1 xi

i=1

m ∑

⎤ ti xTi Pi−1

i=1

m ∑ i=1

ti Pi−1

⎥ ⎥ ⎥ ⎥ ⎥ ⎥⪰0 ⎥ ⎥ ⎥ ⎦

(33)

Z. Wang, X. Shen and Y. Zhu / Automatica 103 (2019) 337–345

341

Proposition 2, the problem (25)–(26) is equivalent to a convex optimization problem (27)–(28), which can be solved efficiently by the modern convex optimization approach. Thus this open problem can be addressed by Proposition 2. 4.2. The relationship between S-procedure relaxation and SDP relaxation Replacing the variable x˜ 0 by x˜ 0 = P0−1 x0 , (19) is equivalent to (9) by Schur complement. Therefore, we can obtain the following proposition. Fig. 1. The equivalence among these relaxation techniques. Table 1 Comparisons of the different relaxation algorithms. Algorithm

logdet(P)

Computing time

Algorithm 1 Algorithm 2 Algorithm 3 Algorithm 4

0.5394 0.5394 0.5394 0.5409

1.7156 1.9265 1.6119 1.3943

m ∑

ti = 1, ti ≥ 0, i = 1, . . . , m.

Proposition 3. The optimization problem (7)–(9) based on the SDP relaxation technique is equivalent to the optimization problem (16)–(19) by the S-procedure relaxation method, i.e., they can derive the same optimal outer-bounding ellipsoid to contain the intersection of these ellipsoids. Remark 2. Based on Propositions 1–3, we can derive that the SDP relaxation method, S-procedure method and parameterized bounding ellipsoid relaxation method obtain the same solution for the optimization problem (1)–(2). The equivalence among these important relaxation techniques is shown in Fig. 1.

(34)

5. The significance of the equivalence

i=1

Let ξ = η1 , and using Schur complement, the above optimization problem (32)–(34) is equivalent to the following problem min − log det

s.t. ξ −

m ∑

m (∑

∑

ti +

∑

5.1. The computational complexity ti xTi Pi−1 xi

(36)

i=1

ti xTi Pi−1

i=1 m ∑

(35)

i=1 m

i=1 m

−

)

ti Pi−1 + n · log(ξ )

In this section, we present the significance of the equivalence from three aspects involving the computational complexity, the application for set membership filter and the application for multisensor fusion.

m (∑

ti Pi−1

m )−1 ∑

ti Pi−1 xi ≥ 0

i=1

i=1

ti = 1, ti ≥ 0, i = 1, . . . , m.

(37)

i=1

Secondly, because of the monotonicity of the log function, we can introduce an optimization variable ξ , then the optimization problem (25)–(26) min n · log(1 − δt ) − log det

m (∑

ti Pi−1

)

(38)

i=1

s.t. (t1 , . . . , tm ) ∈ D +

(39)

is equivalent to min n · log(ξ ) − log det

m (∑

ti Pi−1

)

(40)

i=1

s.t. ξ − 1 + δt ≥ 0, (t1 , . . . , tm ) ∈ D + .

(41) (42)

Thirdly, according to the definition of δt in (24), it is shown that the optimization problems (40)–(42) and (35)–(37) have the same optimal solution. When the objective function is trace function, the proof is similar. Thus the proof of Proposition 2 is finished. Remark 1. It is an open problem to prove the convexity of the scalar functions trace((1 − δt )Pt ) and logdet((1 − δt )Pt ) (25) over the constraint (26) (Polyak et al., 2004). However, according to

From Proposition 1, we can use the decoupled SDP relaxation method to calculate an ellipsoid to contain the intersection F by (27)–(30). Comparing with the optimization problem (7)–(9), we can see that the dimension of the constraint matrix (28) is n + 1, which is smaller than that of matrix (9) based on the SDP relaxation method. In addition, the number of the decision variables of the optimization problem (27)–(28) is m, but that of the optimization n(n+1) problem (7)–(9) is m + 2 + n. Therefore, if we use the decoupled SDP relaxation method to calculate the ellipsoid containing the intersection of the ellipsoids, compared with the SDP relaxation method, the computing time may be reduced much more, which can be clearly seen in Table 1 in the numerical section. The optimization problem (16)–(19) by S-procedure relaxation cannot be solved directly, since the constraint matrix involves the inverse of matrix P0 . In practice, it is accustomed to translate the problem (16)–(19) into the convex optimization problem (7)–(9). Thus, the computation complexity of S-procedure relaxation and SDP relaxation are assumed to be the same. If the intersection F is produced by two ellipsoids, i.e., m = 2, according to the result of equivalence, we can reduce the search of the optimal outer-bounding ellipsoid to a one-dimensional optimization (25)–(26). It does not need to solve the SDP problem (7)–(9). However, with the increase of the number of the intersected ellipsoid, Ros et al. (2002) show that the effective computation of the optimum parameterized bounding ellipsoid involves the conceptually simple, but algebraically tedious. Durieu et al. (2001) also recognized that the direct acquisition of the optimal value of the problem (25)–(26) becomes marginal in most examples treated so far, since it is at cost of more computation. All in all, if we hope to find a tighter ellipsoid to contain the intersection of ellipsoids with lower computational complexity, the decoupled SDP relaxation method is the best choice so far.

342

Z. Wang, X. Shen and Y. Zhu / Automatica 103 (2019) 337–345

∑m

5.2. The application for set membership filter Since the measurement updated step of set-membership filter involves to solve the problem, which is determining a minimum updated ellipsoid to cover the intersection of predicted ellipsoid and measurement ellipsoid. The equivalence among three relaxation methods is helpful for finding a tighter updated ellipsoid to enclose the state. Actually, in Durieu et al. (2001), the authors have proved δt ≥ 0, then the ellipsoid Et′ = {x : (x − xt )T Pt−1 (x − xt ) ≤ 1} ⊇ Et . Since the optimization problem (25)–(26) may lead to high computational complexity, they advise to find the optimal value of vector t by deleting (1 − δt ) in (25), which yields the following optimization problem (Durieu et al., 2001), min f (Pt )

(43)

s.t. (t1 , . . . , tm ) ∈ D . +

(44)

Here the vector t = (t1 , . . . , tm ) is the optimization variable. If the optimal value t ∗ of t takes place within the ellipsoid Et′ , the final ellipsoidal approximation is improved by taking Et ∗ = {x : (x − xt ∗ )T ((1 − δt ∗ )Pt ∗ )−1 (x − xt ∗ ) ≤ 1}.

(45)

Finally, we can use Algorithm 4 to obtain a suboptimal ellipsoid containing the intersection of these ellipsoids E1 , . . . , Em . Algorithm 4 Parameterized bounding ellipsoid relaxation without optimizing 1 − δt (1) Solve the optimization problem (43)–(44) to obtain the optimal solutions ti , i = 1, . . . , m. (2) Compute the shape matrix Pt and the center xt of the outer ellipsoid by (22)–(23); and calculate δt in (24). (3) Obtain the suboptimal ellipsoid with the shape matrix P0 = (1 − δt )Pt and the center x0 = xt .

Corollary 1. When the objective function is logdet or trace function, Algorithm 3 based on the decoupled SDP relaxation method can derive a tighter ellipsoid to contain the intersection than Algorithm 4. Proof. Since the optimal solution of the optimization problem (43)–(45) is feasible for the problem (25)–(26), according to Proposition 2, we can obtain the result of Corollary 1. Remark 3. Algorithm 4 is widely applied to the measurement updated step of set-membership filter for target tracking (Schweppe, 1968), since it has less computing time, which can be seen in Table 1. However, based on Corollary 1, the minimum trace or volume ellipsoid calculated by Algorithm 3 can be tighter than that by Algorithm 4. Therefore, if Algorithm 3 is used in set-membership filter, it can obtain a better estimation performance (see Figs. 3–4). 5.3. The application for multisensor fusion In the distributed estimation fusion setting, when the crosscorrelation of local sensor estimation errors is unknown or impractical, the covariance intersection (CI) algorithm (He et al., 2018; Julier & Uhlmann, 2009; Uhlmann, 1996) is widely used to deal with this problem. CI algorithm can be expressed as Pc−1 =

m ∑

wi Pi−1 ,

(46)

i=1

xˆ c = Pc

m (∑ i=1

wi Pi

−1

) xˆ i ,

where i=1 wi = 1 and wi ≥ 0, i = 1, . . . , m. Let Ec denote the fused ellipsoid in the fusion center, xˆ c and Pc are the center and the shape matrix of Ec , respectively. The local ellipsoid is denoted by Ei , the center and the shape matrix of Ei are xˆ i and Pi , i = 1, . . . , m, respectively. Interestingly, the definitions of Pt and xt in (22)–(23) are same as the expressions of Pc and xˆ c in (46)–(47), then the steps 1–2 of Algorithm 4 are same as the CI algorithm. Corollary 2. Algorithm 3 based on the decoupled SDP relaxation method can derive a tighter ellipsoid to cover the intersection than CI algorithm. Proof. Due to 0 ≤ δt ≤ 1, the shape matrix P0 = (1 − δt )Pt should be smaller than the fused covariance Pc , thus the ellipsoid derived by Algorithm 4 is tighter than that by CI algorithm. According to Corollary 1, we can obtain the result. Therefore, if we use the decoupled SDP relaxation method to solve the distributed estimation fusion with unavailable crosscorrelation among multiple sensors, the estimation performance may be further improved (see Fig. 5–6). 6. Simulation results In this section, we consider two cases: static case and dynamic case. The following simulation results are under Matlab R2015a with CVX (Grant et al., 2009). (1) Static Case: Suppose that there exist two local estimated ellipsoids E1 , E2 , the shape matrix and the center are denoted as follows 2.5 −1

[ P1 =

] [ −1 0.8 , P2 = 1.2 −0.5

] −0.5 4

,

x1 = [0.5 1]T , x2 = [2 1]T . The goal is determining an ellipsoid E containing the intersection of the two ellipsoids. The matrix P is denoted by the shape matrix ellipsoid E , and the objective function f (P) is logdet(P). We use the Algorithms 1–4 to compute the ellipsoid E containing the intersection of the two ellipsoids. Table 1 shows the performance of the Algorithms 1–4. The values of logdet(P) and the computing time are calculated by the average of 100 Monte Carlo runs. It is easy to see that Algorithms 1, 2 and 3 do derive the same ellipsoid to contain the intersection, which is consistent with the result of Propositions 1 and 2. The objective function in (25) is very complex, so that the CVX solver cannot be used, then we use the brute-force search method. From Table 1, the computing time of Algorithm 3 is less than that of Algorithms 1 and 2, the reason is that Algorithm 3 has less optimization variables. The minimum volume ellipsoids are derived by Algorithm 3 and the method in Kurzhanski and Vâlyi (1997) in Fig. 2, respectively. It presents that the method in Kurzhanski and Vâlyi (1997) is not equivalent to the three major relaxation methods, and the size of the ellipsoid derived by Algorithm 3 is smaller than that of the method in Kurzhanski and Vâlyi (1997). The reason may be that the method in Kurzhanski and Vâlyi (1997) is based on two-step relaxation, one is that the Minkowski sum of the two ellipsoids contains the intersection, the other is that find an external ellipsoid to contain the Minkowski sum of two ellipsoids, which leads to lose more information. According to Fig. 1, Table 1, Remark 3, Algorithms 1–3 may provide the smallest volume ellipsoid compared to the other two algorithms. If we consider the computing time, the Algorithm 4 might be the appropriate choice. (2) Dynamic Case: Consider a constant-velocity moving target with three local sensors, and the dynamic model is as follows:

(47)

xk+1 = Fk xk + wk

Z. Wang, X. Shen and Y. Zhu / Automatica 103 (2019) 337–345

Fig. 2. Comparing Algorithm 3 with the method in Kurzhanski and Vâlyi (1997).

343

Fig. 4. Comparison of the volume of updated ellipsoid for sensor 1.

Fig. 5. Comparison of the RMSE of state estimation in the fusion center. Fig. 3. Comparison of the RMSE of state estimation for sensor 1.

yik+1 = xk+1 + vki , i = 1, 2, 3,

[

]

1 T and sampling time T = 1. The process noise 0 1 and measurement noise are assumed to be confined to specified ellipsoidal sets

where Fk =

Wk = {wk : wkT Qk−1 wk ≤ 1} Vk = {vk : vkT Rik

−1

vk ≤ 1}, i = 1, 2, 3,

where

⎡

T3

⎢ Qk = ⎣ 32 T

R2k

[ 2 =

18 0

T2

⎤ [

20 1 2 ⎥ ⎦ , Rk = 0 T

]

[

0 22 , R3k = 22 0

]

0 , 20

]

0 . 18

The goal of this example is to find an optimal outer-bounding ellipsoid to contain the true state by the set-membership filter and the distributed estimation fusion, respectively. In the predicted step of the set-membership filter, the shape matrix Pki +1|k and the center xik+1|k of the predicted ellipsoid Eki +1|k are calculated by Durieu et al. (2001). The updated step is to determine an updated ellipsoid Eki +1|k+1 to contain the intersection

Fig. 6. Comparison of the volume of fused ellipsoid in the fusion center.

of the predicted ellipsoid Eki +1|k and the measurement ellipsoid by each sensor, where the measurement ellipsoid is defined by −1 i Ekv+1 = {xk+1 : (yik+1 − xk+1 )T Rik+1 (yik+1 − xk+1 ) ≤ 1}. In this example, the target starts with x0 = [1 1]. Assume that the center and the shape matrix of the initial bounding ellipsoid are xˆ 0 = [2 2]T and P0 = diag(5, 5), respectively.

344

Z. Wang, X. Shen and Y. Zhu / Automatica 103 (2019) 337–345

Firstly, we compare the estimation precision of Algorithm 3 with that of Algorithm 4 by a single sensor. Since Algorithm 4 is widely applied to set-membership filter in updated step, we use the measurement of sensor 1 to calculate the updated ellipsoid Ek1+1|k+1 by Algorithms 3 and 4, respectively. Figs. 3–4 show the root mean square error (RMSE) of the state estimation and the volume of the updated ellipsoid Ek1+1|k+1 . It is clear to see that Algorithm 3 based on the decoupled SDP method performs better than Algorithm 4, which is consistent with the result of Corollary 1. Therefore, if we want to obtain the higher estimation precision, Algorithm 3 is a better choice. Secondly, in the fusion center, we compare the estimation precision of Algorithm 3 with that of CI method. Assume that each sensor can obtain an updated ellipsoid Eki +1|k+1 by the set-membership filter with Algorithm 3, i = 1, 2, 3, and these local updated ellipsoids are sent to the fusion center. The goal of the fusion center is to calculate an optimal ellipsoid containing the intersection of Eki +1|k+1 , i = 1, 2, 3. Based on Algorithm 3 and CI method, we can derive the fused ellipsoid, respectively. The results are plotted in Figs. 5–6. From Figs. 5–6, we can see that the estimation precisions of Algorithm 3 and CI method are higher than that of each sensor. The Algorithm 3 performs the best, which is also consistent with the result of Corollary 2. 7. Conclusions and future work This paper has investigated the problem on the minimum trace or volume ellipsoid covering the intersection of ellipsoids. In fact, it is a non-convex optimization problem, which needs to compute the maximum of the quadratic function. There are three major relaxation methods involving SDP relaxation, S-procedure relaxation and PBER, but it is unclear for the interrelationships among these methods. By rigorous analysis, we have revealed the equivalence among the three relaxation methods, and solved the open problem, i.e., the optimization problem derived by PBER can be equivalently transformed to a convex optimization problem. In addition, the equivalent result is helpful to the research of the set membership filter and multisensor fusion. Since S-procedure relaxation in complex linear space is lossless for the case of two constraints (Yakubovich, 1971), the future research may extend the equivalent result into complex space. In addition, we hope to apply the equivalent result for maneuvering target tracking based on multisensor networks. References Becis-Aubry, Y., Boutayeb, M., & Darouach, M. (2008). State estimation in the presence of bounded disturbances. Automatica, 44, 1867–1873. Boyd, S., El Ghaoui, L., Feron, E., & Balakrishnan, V. (1994). Linear matrix inequalities in system and control theory. Philadelphia, PA: SIAM (Studies in Applied Mathematics). Boyd, S., & Vandenberghe, L. (2004). Convex optimization. Cambridge University Press. Calafiore, G., & El Ghaoui, L. (2004). Ellipsoidal bounds for uncertain equations and dynamical systems. Automatica, 40, 773–787. Chen, Bo, Ho, Daniel W. C., Zhang, Wen-An, & Yu, Li (2017). Networked fusion estimation with bounded noises. IEEE Transactions on Automatic Control, 62(10), 5415–5421. Chernousko, F. L. (1993). State estimation of dynamic systems. Boca Raton, FL: CRC Press. Durieu, C., Walter, E., & Polyak, B. T. (2001). Multi-input multi-output ellipsoidal state bounding. Journal of Optimization Theory and Applications, 111(2), 273–303. Eldar, Yonina C., Beck, Amir, & Teboulle, Marc (2008). A minimax chebyshev estimator for bounded error estimation. IEEE Transactions on Signal Processing, 56(4), 1388–1397. Grant, Michael, Boyd, Stephen, & Ye, Yinyu (2009). CVX: Matlab software for disciplined convex programming [online]. Available: http://www.stanford.edu/ boyd/cvx. Gu, Feng, He, Yuqing, & Han, Jianda (2015). Active persistent localization of a threedimensional moving target under set-membership uncertainty description through cooperation of multiple mobile robots. IEEE Transactions on Industrial Electronics, 62, 4958–4971.

He, Xingkang, Xue, Wenchao, & Fang, Haitao (2018). Consistent distributed state estimation with global observability over sensor network. Automatica, 92, 162–172. Henrion, Didie, Tarbouriech, Sophie, & Arzelier, Denis (2001). LMI Approximations for the radius of the intersection of ellipsoids: Survey. Journal of Optimization Theory and Applications, 108(1), 1–28. Joachim, Dale, & Deller, Jr. John R. (2006). Multiweight optimization in optimal bounding ellipsoid algorithms. IEEE Transactions on Signal Processing, 54(2), 679–690. Julier, Simon, & Uhlmann, Jeffrey K. (2009). Handbook of multisensor data fusion: theory and practice, General decentralized data fusion with covariance intersection (CI). New York: CRC Press. Kurzhanski, Alexander B., & Mesyats, Aleksei I. (2013). Ellipsoidal motions for applied control: From theory to computation. In IEEE conference on decision and control. Kurzhanski, Alexander, & Vâlyi, Istvân (1997). Ellipsoidal calculus for estimation and control. systems and control: Foundations and applications. Boston, MA: Birkhauser. Luo, Zhiquan, Ma, Wingkin, So, Anthony Man-Cho, Ye, Yinyu, & Zhang, Shuzhong (2010). Semidefinite relaxation of quadratic optimization problems. IEEE Signal Processing Magazine, 27(3), 20–34. Nesterov, Y., & Nemirovski, A. (1994). Interior point polynomial methods in convex programming: Theroy and applications. Philadelphia, PA: SIAM. Noack, Benjamin, Sijs, Joris, Reinhardt, Marc, & Hanebeck, Uwe D. (2017). Decentralized data fusion with inverse covariance intersection. Automatica, 79, 35–41. Polyak, Boris T., Nazin, Sergey A., Durieu, Cêcile, & Walter, Eric (2004). Ellipsoidal parameter or state estimation under model uncertainty. Automatica, 40, 1171–1179. Ros, Llus, Sabater, Assumpta, & Thomas, Federico (2002). An ellipsoidal calculus based on propagation and fusion. IEEE Transactions on Systems, Man and Cybernetics, Part B (Cybernetics), 32(4), 430–442. Schweppe, Fred C. (1968). Recursive state estimation: Unknown but bounded errors and system inputs. IEEE Transactions on Automatic Control, AC-13(1), 22–28. Shen, Xiaojing, Zhu, Yunmin, Song, Enbin, & Luo, Yingting (2011). Minimizing euclidian state estimation error for linear uncertain dynamic systems based on multisensor and multi-algorithm fusion. IEEE Transactions on Information Theroy, 57(10), 7131–7146. Shi, Ke, Chen, Hongsheng, & Lin, Yao (2015). Probabilistic coverage based sensor scheduling for target tracking sensor networks. Information Sciences, 292, 95–110. Sinyakov, V. V. (2015). Method for computing exterior and interior approximations to the reachability sets of bilinear differential systems. Differential Equations, 51, 1097–1111. Uhlmann, Jeffrey K. (1996). General data fusion for estimates with unknown cross covariances. Proceedings of the SPIE Aerosense Conference, 2755(14), 536–547. Vandenberghe, L., & Boyd, S. (1996). Semidefinite programming. SIAM Review, 38(1), 49–95. Wang, Yimin, & Rong Li, X. (2012). Distributed estimation fusion with unavailable cross-correlation. IEEE Transactions on Aerospace and Electronic Systems, 48(1), 259–278. Wang, Zhiguo, Shen, Xiaojing, & Zhu, Yunmin (2017). Set-membership information fusion for multisensor nonlinear dynamic systems. In Proceedings of the 20th international conference of information fusion, July Xi’an. Wang, Zhiguo, Shen, Xiaojing, Zhu, Yunmin, & Pan, Jianxin (2018). A tighter set-membership filter for some nonlinear dynamic systems. IEEE Access, 6, 25351–25362. Wu, Duzhi, Zhou, Jie, & Hu, Aiping (2013). A new approximate algorithm for the chebyshev center. Automatica, 49(8), 2483–2488. Yakubovich, V. A. (1971). S-procedure in nonlinear control theory (pp. 62–77). Vestnik Leningrad University, English translation in Vestnik Leningrad University, (1977) 73–93. Zhang, Wen-An, & Shi, Ling (2018). Sequential fusion estimation for clustered sensor networks. Automatica, 89, 358–363.

Zhiguo Wang received the Ph.D. degree from Sichuan University, Chengdu, China, in 2018. He is currently working as a post-Doctor in school of science and engineering, The Chinese University of Hong Kong, Shenzhen. His research interests are in the fields of information fusion, target tracking, signal and data processing, machine learning, nonconvex optimization, and analysis and control of uncertain systems.

Z. Wang, X. Shen and Y. Zhu / Automatica 103 (2019) 337–345 Xiaojing Shen graduated from Sichuan University in 2003, obtained his M.S. at Jilin University in 2006 and his Ph.D. at Sichuan University in 2009. He is the recipient of the China National Excellent Doctoral Dissertation Award (the highest honor of Ph.D. degree). During 2009–2011, he worked as a post-Doctor and assistant professor in school of computer science, Sichuan University. Since 2013, he has been working as a professor in school of Mathematics, Sichuan University. He was with Syracuse University as a postdoctoral research associate during 2012–2013. His current research interests are in the fields of multisensor decision and estimation fusion, multitarget data association and tracking, optimization theory with applications to information fusion, statistical theory with applications to information fusion. He is one of the authors of two books—Networked Multisensor Decision and Estimation Fusion: Based on Advanced Mathematical Methods (CRC Press, 2012) and Nonlinear Estimation and Applications to Industrial Systems Control (G. Rigatos, ed., ch. 3, pp. 61–88, Nova Science, 2012). He has also served as a cooperator to several major research institutes.

345

Yunmin Zhu received the B.S. degree from the Department of Mathematics and Mechanics, Beijing University, Beijing, China, in 1968. From 1968 to 1978, he was with Luoyang Tractor Factory, Luoyang, Henan, China, as a steel worker and a machine engineer. From 1981 to 1994, he was with Institute of Mathematical Sciences, Chengdu Institute of Computer Applications, Chengdu Branch, Academia Sinica. Since 1995 he has been with Department of Mathematics, Sichuan University as a Professor. From 1986 to 1987, 1989–1990, 1993–1996, 1998–1999, 2001–2002, he was a Visiting Associate or Visiting Professor, at the Lefschetz Centre for Dynamical Systems and Division of Applied Mathematics, Brown University; Department of Electrical Engineering, McGill University; Communications Research Laboratory, McMaster University; and the Department of Electrical Engineering, University of New Orleans. He is the author or coauthor of over 80 papers in international and Chinese journals. He is the author of Multisensor Decision and Estimation Fusion (Kluwer Academic Publishers, 2002), Networked Multisensor Decision and Estimation Fusion: Based on Advanced Mathematical Methods (CRC Press, 2012), and Multisensor Distributed Statistical Decision (Science Press, Chinese Academy of Science, Beijing, 2000), and coauthor of Stochastic Approximations (Shanghai Scientific & Technical Publishers, 1996). He is on the editorial boards of the Journal of Systems Science and Complexity and Systems Science and Mathematics. His research interests include stochastic approximations, adaptive filtering, other stochastic recursive algorithms and their applications in estimations, optimizations, and decisions for dynamic system as well as for signal processing, information compression. In particular, his present major interest is multisensor distributed estimation and decision fusion.