Non-coherent Sensor Fusion via entropy regularized optimal mass transport

Non-coherent Sensor Fusion via
entropy regularized optimal mass transport


This work presents a method for information fusion in source localization applications. The method utilizes the concept of optimal mass transport in order to construct estimates of the spatial spectrum using a convex barycenter formulation. We introduce an entropy regularization term to the convex objective, which allows for low-complexity iterations of the solution algorithm and thus makes the proposed method applicable also to higher-dimensional problems. We illustrate the proposed method’s inherent robustness to misalignment and miscalibration of the sensor arrays using numerical examples of localization in two dimensions.

Non-coherent Sensor Fusion via

entropy regularized optimal mass transport

Filip Elvander, Isabel Haasler, Andreas Jakobsson, and Johan Karlssonthanks: This work was supported in part by the Swedish Research Council, Carl Trygger’s foundation, and the Royal Physiographic Society in Lund.
Div. of Mathematical Statistics, Lund University, Sweden
Dept. of Mathematics, KTH Royal Institute of Technology, Sweden
emails: {filipelv, aj}, {haasler, johan.karlsson}

Index Terms—  Optimal mass transport, Entropy regularization, Target localization, Sensor fusion, Non-coherent processing.

1 Introduction

Source localization using an array of sensors is a well established and still highly active field of research in signal processing, with applications in areas such as radar, sonar, and audio processing [1, 2, 3, 4]. The localization problem is often formulated as a spatial spectral estimation problem, where signal sources are identified with peaks of power in the spectral representation. In contrast to the closely related problem of direction of arrival (DoA) estimation, where one aims to find the directions of incoming planar waves from far field sources, the problem of localization also includes the estimation of the range (i.e., distance to the sources) and thus needs to utilize that the incoming wavefronts are non-planar. Classical methods, originally used mainly for DoA estimation, such as beamforming techniques [5] and subspace methods like MUSIC [6] and ESPRIT [7], may be readily extended to the localization problem, with more recent efforts including works exploiting spectral sparsity [8, 9, 10]. Many of these methods rely on the assumption of well-calibrated arrays, i.e., known sensor positions, in order to produce reliable estimates. However, this assumption may in some applications be violated due to, e.g., calibration errors, perturbations of the sensor localizations, or synchronization errors. Therefore, considerable effort has been directed towards formulating estimators robust to such calibration errors (see, e.g., [11, 12, 13, 14, 15]). In some of the approaches to address this problems, techniques have been developed that divide the available sensors into subarrays, each assumed to be perfectly calibrated, although relative errors between the different subarrays are allowed for [16]. The DoA estimation and/or localization is then performed by minimizing measures of discrepancy between the measurements, i.e., the estimated joint covariance matrix of the complete set of sensors, and the assumed array geometry [15].

In this work, we instead consider localization by fusing information from separate sensor arrays, where we only assume to have access to measurements in the form of covariance matrices for the individual arrays, so-called non-coherent processing (see, e.g., [17, 18, 19, 20]), i.e., we assume no information about cross-correlation between different arrays. Such a scenario may arise in distributed sensor networks, wherein sensor nodes transmit information to a central processing unit that performs the estimation [21, 22]. In such settings, it is common to have restrictions on the communication bandwidth, which requires transmitting only aggregated data [21], e.g., in the form of covariance matrices. Also, as the separate sensor arrays may not be jointly synchronized, cross-correlations cannot be estimated from the raw data.

Extending on the work presented in [23], which proposes a robust method for DoA sensor fusion using optimal mass transport (OMT), we here generalize this approach to the source localization problem. One of the main ideas is to use OMT as a notion of distance between spectra (see, e.g., [24, 25, 26]), and in particular, the proposed method performs sensor fusion by computing the corresponding barycenter solutions based on incomplete measurements. In addition, to make the scheme computationally efficient, we introduce an entropy regularization term in the problem formulation and show how to solve the resulting optimization problem iteratively, similar to the popular Sinkhorn iteration scheme [27] (see also [28]). The robustness of the proposed method is demonstrated in numerical examples.

2 Background

2.1 Near field localization

Consider an array consisting of sensors located at , for , and assume that point-like sources, located in a space , emit spherical waves impinging on the sensors. Then, the spatial distribution of the energy of the source signals may be representated by a non-negative function, or measure, on , referred to as the power spectrum. We will denote the space of such functions by . Letting denote the wavelength of the impinging waves, the covariance matrix of the sensor measurements is , where the linear operator is given by (see, e.g., [29])


where is the array manifold vector [30], i.e.,


with denoting the Lebesgue measure. In other words, a positive semi-definite matrix is a valid covariance for a sensor array if and only if , for some . From this, it may be noted that performing source localization based on sensor measurements corresponds to inferring a power spectrum from observations of the array covariance matrix .

We will here consider a scenario wherein several sensor arrays measure the superposition of the impinging waves, i.e., when covariance matrices , for , are available, corresponding to separate sensor arrays. For the case when all arrays are calibrated, i.e., when the position of the sensors are known perfectly, source localization thus corresponds to finding a single power spectrum that is consistent with the observed matrices, i.e., , for , where is the linear operator corresponding to array . However, for the case of miscalibration, the (assumed) operators contain errors, which may result in there being no single consistent with all covariance matrices. We will here adress this problem by utilizing OMT in order to induce robustness to the sensor fusion problem.

2.2 Optimal mass transport

Let and be two non-negative mass distributions on the space . A transport plan for such a pair is a non-negative distribution on , with , where describes the amount of mass transported from location to location . By defining a cost function describing the cost of moving mass, we may define the total cost associated with a transport plan as


The OMT problem is then to find the feasible transport plan from to incurring minimal transportation cost (see, e.g., [31]), i.e., to solve

subject to

where the constraints ensure that the transport plan is valid, i.e., matches the marginals and . The suitability of utilizing OMT costs as measures of distance has previously been considered for imposing metric structure on the space of power spectra [24], as well as to induce distances and smooth interpolants for covariance matrices [25]. In this work, the optimal transport cost will be used for interpolating between mass distributions and . It is worth noting that the OMT cost quantifies location errors based on the distances on the underlying space, in this case , thereby allowing for offsets due to, e.g., calibration errors. This is in contrast to standard metrics (e.g., , ) which only compare functions point by point without taking the geometry of the underlying space, i.e., relative location of the points, into account.

3 Proposed method

Consider the localization scenario with the source signals impinging on a set of arrays, with array manifold vectors , generating the covariance matrices , . Given the array geometries and covariance matrix estimates of each individual array, we aim to estimate the spatial spectrum . That is, we assume that each sensor array is calibrated within itself, but there might be global calibration errors resulting in inaccurate knowledge of each sensor array’s position and orientation. For the far field case, i.e., DoA estimation, an OMT formulation for fusing information from different sensor arrays has been proposed in [23]. There, the estimate was formed as the spectrum closest, in the OMT sense, to a set of spectra, each consistent with a corresponding covariance matrix. Correspondingly, in the localization setting the spatial spectrum may be found as the solution to the convex barycenter problem

subject to

where denotes the linear operator (1) of array . In order to find practical solutions to the barycenter problem in (5), the space may be discretized as to be represented by points . Accordingly, the marginal spatial spectra are represented by the vectors , and the cost function and transport plans are represented by matrices and , respectively, where , and denotes the mass transported from to . Although such a discretization yields finite-dimensional linear programs, finding an optimal transport plan can be computationally demanding due to the large number of variables.

For discretizations of the standard OMT problem (4), this was in [27] addressed by regularizing the problem by introducing an entropy term to the objective function, yielding the strictly convex program

subject to

where is an entropy term, a small constant, and denotes a vector of ones. In [27], it was shown that this problem may be solved in a computationally efficient way, using so-called Sinkhorn iterations.

3.1 Efficient sensor fusion

In order to arrive at a computationally efficient scheme for computing spectral barycenters, we propose to regularize the problem in (5) according to

subject to

Here, , for , are the vectorizations of the covariances matrices, and are matrix representations of the discretized operators . Also, the vectors allow for noisy covariance matrix estimates, penalized by the standard -norm, as controlled by the positive penalty parameter . The addition of the entropy term, , renders the problem strictly convex and allows for efficiently computing the barycenter , using the observation in [28], that Sinkhorn iterations are equivalent to a block coordinate ascent in a Lagrange dual of the entropy regularized transport problem. This then allows for increasing the dimensionality of the sensor fusion problem, i.e., making it feasible to include also range estimation or considering higher spatial dimensions. The extended Sinkhorn iterations, constituting a globally convergent solution algorithm for (7), are detailed in the following proposition.

  Require: Initial guess and let for
  while Sinkhorn not converged do
     for j=1:J do
        while Newton not converged do
           , with determined by a linesearch
        end while
     end for
     for j=1:J do
     end for
  end while
Algorithm 1 Sinkhorn-Newton method for the optimal mass transport barycenter problem with partial information.
Proposition 1

A block coordinate ascent scheme for the dual of the barycenter problem (7) is to iteratively perform

  1. For , let be the solution to



  2. Update


Here, , and , , and denotes the elementwise exponential function, elementwise division, and elementwise product, respectively. The barycenter solution is then given by , where is the limit point of the iterations. Note that we here convert all complex quantities to an equivalent real-valued representation such that the vectorized covariances, , are mapped according to

where and denote the real and imaginary parts, respectively, and where and the discretized operators are converted as to be consistent with this.

A proof of the proposition may be found in Appendix A. It may be noted that for the noise-free case with complete observations of the underlying spectra, i.e., and , as well as , one obtains the standard Sinkhorn iterations. The roots of (8) are determined by the means of a Newton method, using the Jacobian

The complete estimation method for computing the spectral barycenter is summarized in Algorithm 1.

Fig. 1: Spectral estimate as given by the proposed estimator in (7). The alignment error is 6.7 degrees.
Fig. 2: Pseudo-spectrum as given by non-coherent MUSIC. The alignment error is 6.7 degrees.
Fig. 3: Average distance from estimates to signal sources, as function of the misalignment angle.

4 Numerical results

In this section, we illustrate the proposed method’s behavior in sensor fusion scenarios for localization in 2D, demonstrating the robustness to array alignment errors that the OMT criterium induces. Specifically, we consider two uncorrelated sources impinging on two asynchronous sensor arrays; one ellipsoidal shaped array consisting of 8 sensors, and one linear array consisting of 7 sensors. To simulate misalignment between the two arrays, we introduce an unknown rotation to the ellipsoidal array. The scenario, as shown in Figure 1, is then used in a simulation study. Varying the misalignment, i.e., the rotation angle, between 0 and 10 degrees, we for each considered angle perform 100 Monte Carlo simulations, were the locations of the signal sources are randomized uniformly on the square . We let the source signals be uncorrelated circularly symmetric Gaussian white noises and add an uncorrelated Gaussian noise to the sensors. The source signals have variance 100 and the variance of the sensor noise is 1. The wavelength of the impinging waves is twice that of the smallest sensor spacing in the linear array. Estimating the covariance matrices of the separate sensor arrays using the sample covariance matrix and 500 signal snapshots, we apply the proposed estimator in (7) to the covariance matrix estimates and record the distance from the location estimates to the actual target positions. Throughout, we use the parameters and , and the cost function , for grid points .

As comparison, we also obtain estimates using the non-coherent MUSIC and MVDR estimators, as described in [18], as well as the least-squares (LS) estimator from [32], and the non-coherent SPICE estimator from [19]. Figures 1 and 2 show the estimates obtained in one simulation using the proposed and MUSIC estimators, respectively. The misalignment angle is 6.7 degrees. As can be seen, for the proposed estimator, the misalignment only causes a small perturbation on the estimated source locations, demonstrating the ability of the OMT formulation to perform smooth interpolation on the underlying domain. In contrast, the MUSIC (pseudo-) spectral estimate contains spurious peaks, as well as larger deviations from the true source locations. The robustness of the proposed estimator is further demonstrated in Figure 3, displaying the results from the simulations. From this, it may be noted that the proposed estimator deviates slightly more from the ground truth for the non-perturbed case as compared to the other methods. However, the robustness to misalignment is significantly higher; the rate of error increase is considerably lower than for the other methods.

Appendix A Proof of Proposition

Proof 1

In the following, we denote the set by , and similarly for , and . The Lagrangian of (7) with the dual variables and is


Minimizing the Lagrangian with respect to gives the expression for the mass transport matrix solutions,


with , and , for . Further, minimizing (9) with respect to the variables and , one gets a dual problem to (7), formulated as

subject to


A block coordinate ascent method is to maximize (11) iteratively with respect to and . Maximization with respect to requires

which gives the first step (8). To maximize (11) with respect to , consider the corresponding Lagrangian with dual variable , given by

Maximization with respect to requires

Further, due to the constraint in (11), it follows

This completes the second step of the method. Finally, with (10), we get the expression for the barycenter marginal


  • [1] H. Krim and M. Viberg, “Two Decades of Array Signal Processing Research,” IEEE Signal Process. Mag., pp. 67–94, July 1996.
  • [2] J. C. Chen, Kung Yao, and R. E. Hudson, “Source localization and beamforming,” IEEE Signal Process. Mag., vol. 19, no. 2, pp. 30–39, Aug. 2002.
  • [3] P. Stoica and R. Moses, Spectral Analysis of Signals, Prentice Hall, Upper Saddle River, N.J., 2005.
  • [4] N. Dey and A. S. Ashour, Direction of Arrival Estimation and Localization of Multi-Speech Sources, Springer, 2018.
  • [5] J. Capon, “High Resolution Frequency Wave Number Spectrum Analysis,” Proc. IEEE, vol. 57, pp. 1408–1418, 1969.
  • [6] R. Schmidt, “Multiple emitter location and signal parameter estimation,” in Proceedings of RADC Spectrum Estimation Workshop, 1979, pp. 243–258.
  • [7] R. Roy and T. Kailath, “ESPRIT – Estimation of Signal Parameters via Rotational Invariance Techniques,” IEEE Trans. Acoust., Speech, Signal Process., vol. 37, no. 7, pp. 984–995, July 1989.
  • [8] P. Stoica, P. Babu, and J. Li, “SPICE : a novel covariance-based sparse estimation method for array processing,” IEEE Trans. Signal Process., vol. 59, no. 2, pp. 629 –638, Feb. 2011.
  • [9] S. I. Adalbjörnsson, T. Kronvall, S. Burgess, K. Åström, and A. Jakobsson, “Sparse Localization of Harmonic Audio Sources,” IEEE Transactions on Audio, Speech, and Language Processing, vol. 24, no. 1, pp. 117–129, January 2016.
  • [10] Z. Yang and L. Xie, “Enhancing Sparsity and Resolution via Reweighted Atomic Norm Minimization,” IEEE Trans. Signal Process., vol. 64, no. 4, pp. 995–1006, Feb 2016.
  • [11] J. Yang and A. Swindlehurst, “The Effects of Array Calibration Errors on DF-Based Signal Copy Performance,” IEEE Trans. Signal Process., vol. 43, no. 11, pp. 2724–2732, November 1995.
  • [12] S. D. Somasundaram and A. Jakobsson, “Degradation of Covariance Reconstruction-Based Robust Adaptive Beamformers,” in Sensor Signal Processing for Defense, Edinburgh, Sept. 8-9 2014.
  • [13] M. Rübsamen and M. Pesavento, “Maximally Robust Capon Beamformer,” IEEE Trans. Signal Process., vol. 61, no. 8, pp. 2030–2041, 2013.
  • [14] R. G. Lorenz and S. Boyd, “Robust Minimum Variance Beamforming,” IEEE Trans. Signal Process., vol. 53, no. 5, pp. 1684–1696, May 2005.
  • [15] M. Pesavento, A. B. Gershman, and K. M. Wong, “Direction finding in partly calibrated sensor arrays composed of multiple subarrays,” IEEE Trans. Signal Process., vol. 50, no. 9, pp. 2103–2115, Sep. 2002.
  • [16] C. M. See and A. B. Gershman, “Direction-of-arrival estimation in partly calibrated subarray-based sensor arrays,” IEEE Trans. Signal Process., 2004.
  • [17] M. Wax and T. Kailath, “Detection of Signals by Information Theoretic Criteria,” IEEE Trans. Acoust., Speech, Signal Process., vol. 33, no. 2, pp. 387–392, April 1985.
  • [18] D. W. Rieken and D. R. Fuhrmann, “Generalizing MUSIC and MVDR for Multiple Noncoherent Arrays,” IEEE Trans. Signal Process., vol. 52, no. 9, pp. 2396–2406, Sep. 2004.
  • [19] W. Suleiman, P. Parvazi, M. Pesavento, and A. M. Zoubir, “Non-Coherent Direction-of-Arrival Estimation Using Partly Calibrated Arrays,” IEEE Trans. Signal Process., vol. 66, no. 21, pp. 5776–5788, Nov. 2018.
  • [20] H. Kim, A. M. Haimovich, and Y. C. Eldar, “Non-Coherent Direction of Arrival Estimation from Magnitude-Only Measurements,” IEEE Signal Process. L, vol. 22, no. 7, pp. 925–929, Jul. 2015.
  • [21] L. M. Kaplan, “Global node selection for localization in a distributed sensor network,” IEEE Trans. Aerosp. Electron. Syst., vol. 42, no. 1, pp. 113–135, Jan. 2006.
  • [22] D. Li, K. D. Wong, Y. H. Hu, and A. M. Sayeed, “Detection, classification, and tracking of targets,” IEEE Signal Process. Mag., vol. 19, no. 2, pp. 17–29, Mar. 2002.
  • [23] F. Elvander, I. Haasler, A. Jakobsson, and J. Karlsson, “Tracking and Sensor Fusion in Direction of Arrival Estimation Using Optimal Mass Transport,” in 26th European Signal Processing Conference, Rome, Italy, Sep. 3-7 2018.
  • [24] T. T. Georgiou, J. Karlsson, and M. S. Takyar, “Metrics for power spectra: an axiomatic approach,” IEEE Trans. Signal, vol. 57, no. 3, pp. 859–867, Mar. 2009.
  • [25] F. Elvander, A. Jakobsson, and J. Karlsson, “Interpolation and Extrapolation of Toeplitz Matrices via Optimal Mass Transport,” IEEE Trans. Signal. Process, vol. 66, no. 20, pp. 5285 – 5298, Oct. 2018.
  • [26] F. Elvander, A. Jakobsson, and J. Karlsson, “Using Optimal Mass Transport for Tracking and Interpolation of Toeplitz Covariance Matrices,” in 43rd IEEE Int. Conf. on Acoustics, Speech, and Signal Processing, Calgary, Canada, April 15-20 2018, pp. 4469–4473.
  • [27] M. Cuturi, “Sinkhorn distances: Lightspeed computation of optimal transport,” Proc. Adv. Neural Inf. Process. Syst., pp. 2292–2300, 2013.
  • [28] J. Karlsson and A. Ringh, “Generalized Sinkhorn iterations for regularizing inverse problems using optimal mass transport,” SIAM Journal on Imaging Sciences, vol. 10, no. 4, pp. 1935–1962, 2017.
  • [29] T. T. Georgiou, “Solution of the general moment problem via a one-parameter imbedding,” IEEE Trans. Autom. Control, vol. 50, no. 6, pp. 811–826, 2005.
  • [30] D. H. Johnson and D. E. Dudgeon, Array Signal Processing: Concepts and Techniques, Prentice Hall, Englewood Cliffs, N.J., 1993.
  • [31] C. Villani, Optimal transport: old and new, Springer Science & Business Media, 2008.
  • [32] J-A. Luo, K. Yu, Z. Wang, and Y-H. Hu, “Passive source localization from array covariance matrices via joint sparse representations,” Neurocomputing, vol. 270, pp. 82–90, 2017.
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description