On SDP and ESDP Relaxation of Sensor Network Localization

Localization

Paul Tseng

Mathematics, University of Washington Seattle

7th Int. Conf. Num. Optim. Num. Lin. Algeb., Lijiang August 9, 2009

(joint work with Ting Kei Pong)

Talk Outline

• Sensor network localization

• SDP, ESDP relaxations: properties and soln accuracy certificate

Talk Outline

• Sensor network localization

• SDP, ESDP relaxations: properties and soln accuracy certificate

• A robust version of ESDP to handle noises

• Log-barrier penalty CGD method

Talk Outline

• Sensor network localization

• SDP, ESDP relaxations: properties and soln accuracy certificate

• A robust version of ESDP to handle noises

• Log-barrier penalty CGD method

• Numerical simulations

• Conclusion & Ongoing work

Sensor Network Localization

Basic Problem

• n pts in <².

• Know last n − m pts (‘anchors’) xm+1, ..., x_n and Eucl. dist. estimate for pairs of ‘neighboring’ pts

d_ij ≥ 0 ∀(i, j) ∈ A with A ⊆ {(i, j) : 1 ≤ i, j ≤ n}.

• Estimate first m pts (‘sensors’).

Sensor Network Localization

Basic Problem

• n pts in <².

• Know last n − m pts (‘anchors’) xm+1, ..., x_n and Eucl. dist. estimate for pairs of ‘neighboring’ pts

d_ij ≥ 0 ∀(i, j) ∈ A with A ⊆ {(i, j) : 1 ≤ i, j ≤ n}.

• Estimate first m pts (‘sensors’).

History? Graph realization/rigidty, Euclidean matrix completion, position estimation in wireless sensor network, ...

Optimization Problem Formulation

υ_opt := min

x₁,...,x_m

X

(i,j)∈A

kx_i − x_jk² − d²_ij 

Optimization Problem Formulation

υ_opt := min

x₁,...,x_m

X

(i,j)∈A

kx_i − x_jk² − d²_ij 

• Objective function is nonconvex. m can be large (m ≥ 1000).

6. .

_

• Problem is NP-hard (reduction from PARTITION).

6. .

_

• Local improvement heuristics can fail badly.

6. .

_

Optimization Problem Formulation

υ_opt := min

x₁,...,x_m

X

(i,j)∈A

kx_i − x_jk² − d²_ij 

• Objective function is nonconvex. m can be large (m ≥ 1000).

6. .

_

• Problem is NP-hard (reduction from PARTITION).

6. .

_

• Local improvement heuristics can fail badly.

6. .

_

• Use a convex (SDP, SOCP) relaxation (& local improvement).

Low soln accuracy OK. Distributed computation preferred.

SDP Relaxation

Let X := [x₁ · · · x_m]. Y = X^TX ⇐⇒ Z =  Y X^T

X I



 0, rankZ = 2

SDP Relaxation

Let X := [x₁ · · · x_m]. Y = X^TX ⇐⇒ Z =  Y X^T

X I



 0, rankZ = 2 SDP relaxation (Biswas,Ye ’03):

υ_sdp := min

Z

X

(i,j)∈A,i≤m<j

y_ii − 2x^T_j x_i + kx_jk² − d²_ij 

+ X

(i,j)∈A,i<j≤m

y_ii − 2y_ij + y_jj − d²_ij  s.t. Z =  Y X^T

X I



 0

Adding the nonconvex constraint rankZ = 2 yields original problem.

But SDP relaxation is still expensive to solve for m large..

ESDP Relaxation

ESDP relaxation (Wang, Zheng, Boyd, Ye ’06):

υ_esdp := min

Z

X

(i,j)∈A,i≤m<j

y_ii − 2x^T_j x_i + kx_jk² − d²_ij 

+ X

(i,j)∈A,i<j≤m

y_ii − 2y_ij + y_jj − d²_ij  s.t. Z =  Y X^T

X I







y_ii y_ij x^T_i y_ij y_jj x^T_j x_i x_j I



  0 ∀(i, j) ∈ A, i < j ≤ m

0 ≤ υ_esdp ≤ υ_sdp ≤ υ_opt. In simulation, ESDP is nearly as strong as SDP, and solvable much faster by IP method.

Example 1

n = 3, m = 1, d₁₂ = d₁₃ = 2

Problem

0 = min

x₁∈<²

|kx₁ − ₁

0 k² − 4| + |kx₁ − ₋₁

0  k² − 4|

SDP/ESDP Relaxation

0 = min

x1=[^αβ]^∈<2

y11∈<

|y₁₁ − 2α − 3| + |y₁₁ + 2α − 3|

s.t.





y₁₁ α β

α 1 0

β 0 1



  0

If solve SDP/ESDP by IP method, then likely get analy. center y₁₁ = 3, x₁ = h

0 0

i

Example 2

n = 4, m = 1, d₁₂ = d₁₃ = 2, d₁₄ = 1

Problem

0 = min

x₁∈<²

|kx₁ − ₁

0 k² − 4| + |kx₁ − ₋₁

0  k² − 4| + |kx₁ − h

√1 3

ik² − 1|

SDP/ESDP Relaxation

0 = min

x1=[^αβ]^∈<2

y11∈<

|y₁₁ − 2α − 3| + |y₁₁ + 2α − 3| + |y₁₁ − 2α − 2√

3β + 3|

s.t.





y₁₁ α β

α 1 0

β 0 1



  0

SDP/ESDP has unique soln y11 = 3, x₁ = h

√0 3

i

Properties of SDP & ESDP Relaxations

Assume each i ≤ m is conn. to some j > m in the graph ({1, ..., n}, A).

Fact 0

• Sol(SDP) and Sol(ESDP) are nonempty, closed, convex.

• If

d_ij = kx^true_i − x^true_j k ∀ (i, j) ∈ A “noiseless case”

(x^true_i = x_i ∀ i > m), then

υ_opt = υ_sdp = υ_esdp = 0 and

Z^true :=  X^true I ^T

 X^true I 

is a soln of SDP and ESDP (i.e., Z^true ∈ Sol(SDP) ⊆ Sol(ESDP)).

Let tr_i[Z] := y_ii − kx_ik², i = 1, ..., m. “ith trace”

Fact 1

tr_i[Z] = 0 ∃Z ∈ ri(Sol(ESDP)) =⇒ x_i is invariant over Sol(ESDP) (so x_i = x^true_i in noiseless case) Still true with “ESDP” changed to “SDP”.

Let tr_i[Z] := y_ii − kx_ik², i = 1, ..., m. “ith trace”

Fact 1

tr_i[Z] = 0 ∃Z ∈ ri(Sol(ESDP)) =⇒ x_i is invariant over Sol(ESDP) (so x_i = x^true_i in noiseless case) Still true with “ESDP” changed to “SDP”.

Fact 2

tr_i[Z] = 0 ∀Z ∈ Sol(ESDP) ⇐= x_i is invariant over Sol(ESDP).

Proof is by induction, starting from sensors that neighbor anchors.

(Q: True for SDP?)

Proof idea

• If (i, j) ∈ A and x_i, x_j are invar. over Sol(ESDP), then tr_i[Z] = tr_j[Z]

∀Z ∈ Sol(ESDP).

• Suppose ∃i ≤ m such that x_i is invar. over Sol(ESDP) but tr_i[ ¯Z] > 0 for some ¯Z ∈ Sol(ESDP). Consider maximal ¯I ⊂ {1, . . . , m} such that x_i is invar. over Sol(ESDP) and tri[ ¯Z] > 0 ∀i ∈ ¯I.

• Then x_i is not invar. over Sol(ESDP) ∀i ∈ N (¯I).

So ∃Z ∈ ri(Sol(ESDP)) with x_i 6= ¯x_i ∀i ∈ N (¯I).

• Let Z^α = α ¯Z + (1 − α)Z with α > 0 suff. small.

Can rotate x^α_i ∀i ∈ ¯I and Z^α still remains in Sol(ESDP). ⇒⇐

In practice, there are measurement noises:

d²_ij = kx^true_i − x^true_j k² + δ_ij ∀(i, j) ∈ A.

When δ := (δij)_(i,j)∈A ≈ 0, does tr_i[Z] = 0 (with Z ∈ ri(Sol(ESDP))) imply x_i ≈ x^true_i ?

In practice, there are measurement noises:

d²_ij = kx^true_i − x^true_j k² + δ_ij ∀(i, j) ∈ A.

When δ := (δij)_(i,j)∈A ≈ 0, does tr_i[Z] = 0 (with Z ∈ ri(Sol(ESDP))) imply x_i ≈ x^true_i ? No!

6. .

_

Fact 3

tr_i[Z] = 0 ∃Z ∈ ri(Sol(ESDP)) 6=⇒ x_i ≈ x^true_i . Still true with “ESDP” changed to “SDP”.

Proof is by counter-example.

An example of sensitivity of ESDP solns to measurement noise:

Problem data: m = 2, n = 6;

d₁₂ = p4 + (1 − )², d₁₃ = 1 + , d₁₄ = 1 − , d₂₅ = d₂₆ = √

2 ( > 0)

Thus, even when Z ∈ Sol(ESDP) is unique, tri[Z] = 0 fails to certify accuracy of xi in the noisy case!

Robust ESDP

Fix any ρij > |δ_ij| ∀(i, j) ∈ A (ρ > |δ|).

Let Sol(ρESDP) denote the set of Z =  Y X^T

X I



satisfying

|y_ii − 2x^T_j x_i + kx_jk² − d²_ij| ≤ ρ_ij ∀(i, j) ∈ A, i ≤ m < j

|y_ii − 2y_ij + y_jj − d²_ij| ≤ ρ_ij ∀(i, j) ∈ A, i < j ≤ m





y_ii y_ij x^T_i y_ij y_jj x^T_j x_i x_j I



  0 ∀(i, j) ∈ A, i < j ≤ m

Note: Z^true =  X^true I ^T

 X^true I  ∈ Sol(ρESDP).

Let

Z^ρ,δ := arg min

Z∈Sol(ρESDP)

X

(i,j)∈A,i<j≤m

− ln det





y_ii y_ij x^T_i y_ij y_jj x^T_j x_i x_j I





Let

Z^ρ,δ := arg min

Z∈Sol(ρESDP)

X

(i,j)∈A,i<j≤m

− ln det





y_ii y_ij x^T_i y_ij y_jj x^T_j x_i x_j I





Fact 4

|δ|<ρ→0x^ρ,δ_i = x^true_i

tr_i[Z^ρ,δ] > ₁₀^η ∃|δ| < ρ ≤ ¯ρe =⇒ x_i not invar. over Sol(ESDP) when δ = 0 Moreover,

kx^ρ,δ_i − x^true_i k ≤ p

2|A| + m q

tr_i[Z^ρ,δ] ∀ |δ| < ρ.

Log-barrier Penalty CGD Method

Efficiently compute Z^ρ,δ? Let h_a(t) := 1

2(t − a)²₊ + 1

2(−t − a)²₊ (|t| ≤ a ⇐⇒ h_a(t) = 0) and

f_µ(Z) := X

(i,j)∈A,i≤m<j

h_ρij(y_ii − 2x^T_j x_i + kx_jk² − d²_ij)

+ X

(i,j)∈A,i<j≤m

h_ρij(y_ii − 2y_ij + y_jj − d²_ij)

+µ X

(i,j)∈A,i<j≤m

− ln det





y_ii y_ij x^T_i y_ij y_jj x^T_j x_i x_j I





• f_µ is partially separable, strictly convex & diff. on its domain.

• For each fixed ρ > |δ|, argminfµ → Z^ρ,δ as µ → 0.

• f_µ is partially separable, strictly convex & diff. on its domain.

• For each fixed ρ > |δ|, argminfµ → Z^ρ,δ as µ → 0.

Idea

Yun ’06)

Log-barrier Penalty CGD Method

Given Z in domfµ, compute gradient ∇Z_if_µ of fµ w.r.t.

Z_i := {x_i, y_ii, y_ij : (i, j) ∈ A} for each i.

• If k∇Z_if_µk ≥ max{µ, 10⁻⁷} for some i, update Zi by moving along the Newton direction −

∂_Z²

iZ_if_µ⁻¹

∇_Zif_µ with Armijo stepsize rule.

• Decrease µ when k∇_Zif_µk < max{µ, 10⁻⁶} ∀ i.

µ_initial = 10, µ_final = 10⁻¹⁴. Decrease µ by a factor of 10 each time.

Coded in Fortran. Compute Newton direc. by sparse Cholesky.

Computation easily distributes.

Simulation Results

• Compare ρESDP as solved by LPCGD method with ESDP as solved by Sedumi 1.05 Sturm (with the interface to Sedumi coded by Wang et al).

Simulation Results

• Compare ρESDP as solved by LPCGD method with ESDP as solved by Sedumi 1.05 Sturm (with the interface to Sedumi coded by Wang et al).

• Anchors and sensors x^true₁ , ..., x^true_n uniformly distributed in [−.5, .5]², m = .9n. (i, j) ∈ A whenever kx^true_i − x^true_j k ≤ rr. Set

d_ij = kx^true_i − x^true_j k · |1 + σ · _ij|, where _ij ∼ N (0, 1).

Simulation Results

• Compare ρESDP as solved by LPCGD method with ESDP as solved by Sedumi 1.05 Sturm (with the interface to Sedumi coded by Wang et al).

• Anchors and sensors x^true₁ , ..., x^true_n uniformly distributed in [−.5, .5]², m = .9n. (i, j) ∈ A whenever kx^true_i − x^true_j k ≤ rr. Set

d_ij = kx^true_i − x^true_j k · |1 + σ · _ij|, where _ij ∼ N (0, 1).

• Sensor i is judged as “accurately positioned” if

tr_i[Z^found] < (.01 + 30σ)d^avg_ij .

ρESDP_LPCGD ESDP_Sedumi

n m σ rr cpu/m_ap/err_ap cpu(cpus)/m_ap/err_ap 1000 900 0 .06 7/662/1.7e-3 182(104)/669/2.1e-3 1000 900 .01 .06 5/660/2.2e-2 119(42)/720/3.1e-2 2000 1800 0 .06 26/1762/3.1e-4 1157(397)/1742/3.9e-4 2000 1800 .01 .06 20/1699/1.4e-2 966(233)/1746/2.4e-2 10000 9000 0 .02 77/7844/2.3e-3 16411(1297)/6481/2.5e-3 10000 9000 .01 .02 63/8336/1.0e-2 16368(1264)/8593/8.7e-3

• cpu(sec) times are on a HP DL360 workstation, running Linux 3.5. ESDP is solved by Sedumi; cpus:= run time for Sedumi.

• Set ρij = d²_ij · ((1 − 2σ)⁻² − 1).

• m_ap := # accurately positioned sensors.

err_ap := max_iaccurate. pos. kx_i − x^true_i k.

900 sensors, 100 anchors, rr = 0.06, σ = 0.01, solve ρESDP by LPCGD method. x^true_i (shown as ∗) and x^ρ,δ_i (shown as •) are joined by blue line segment; anchors are shown as ◦.

60 sensors, 4 anchors at corners, rr = 0.3, σ = 0.1. x_i (shown as ∗) and x^ρ,δ_i (shown as

•) are joined by blue line segment; anchors are shown as ◦. Left: Soln of ρESDP found by LPCGD method. Right: After local gradient improvement.

−0.5 0 0.5

−0.5

−0.4

−0.3

−0.2

−0.1 0 0.1 0.2 0.3 0.4 0.5

−0.5 0 0.5

−0.5

−0.4

−0.3

−0.2

−0.1 0 0.1 0.2 0.3 0.4 0.5

Conclusion & Ongoing work

• SDP and ESDP solns are sensitive to measurement noise. Has soln accuracy certificate under no noise only (though it works well enough in simulation).

• ρESDP solns are more stable. Has soln accuracy certificate under low noise (which works well enough in simulation). Needs to estimate the noise level δ to set ρ. Can ρ > |δ| be

relaxed?

• SDP, ESDP, ρESDP solns can be further refined by local improvement. This improves the rmsd when noise level is high (e.g., σ = 0.1).

• Approximation bounds? Extensions to handle lower bounds on distances (e.g., (i, j) 6∈ A imply kx^true_i − x^true_j k > rr)?

Thanks for coming!

6. .

^