Sensor Array Design Through Submodular
Optimization
Gal Shulkind, Student Member, IEEE, Stefanie Jegelka, Gregory W. Wornell, Fellow, IEEE

arXiv:1705.06616v2 [] 28 Dec 2017

Abstract
We consider the problem of far-field sensing by means of a sensor array. Traditional array geometry design techniques are
agnostic to prior information about the far-field scene. However, in many applications such priors are available and may be utilized
to design more efficient array topologies.
We formulate the problem of array geometry design with scene prior as one of finding a sampling configuration that enables
efficient inference, which turns out to be a combinatorial optimization problem. While generic combinatorial optimization problems
are NP-hard and resist efficient solvers, we show how for array design problems the theory of submodular optimization may be
utilized to obtain efficient algorithms that are guaranteed to achieve solutions within a constant approximation factor from the
optimum. We leverage the connection between array design problems and submodular optimization and port several results of
interest. We demonstrate efficient methods for designing arrays with constraints on the sensing aperture, as well as arrays respecting
combinatorial placement constraints.
This novel connection between array design and submodularity suggests the possibility for utilizing other insights and
techniques from the growing body of literature on submodular optimization in the field of array design.
Index Terms
Array design, Submodular optimization, Far field sensing

I. I NTRODUCTION
Sensor arrays for spatial sensing are widely deployed in a wide range of applications including radar, sonar, medical imaging
and radio astronomy and there is a vast literature on the topics of array design and array processing from the last century [1],
[2].
A major goal in designing arrays is efficiently meeting sensing specifications with a limited budget of sensing elements,
which is often a main determinant of system cost in terms of dimensions, weight, complexity and manufacturing cost. However,
traditionally most studies assume a uniform and linear design, often a truncated half-wavelength array, or restrictions thereof.
Indeed, the problem of designing non-uniform arrays hints at combinatorial optimization and is computationally intractable as
we discuss later.
In beamforming arrays attaining a desired resolution level is often a primary concern, achieved by means of manipulating
beam pattern parameters such as lobe widths and positions. The problem of designing efficient array geometries that facilitate
desired beam patterns has been widely studied in the past. Techniques involving array thinning start with a dense uniform
geometry, removing elements while maintaining performance within specified bounds [3]. Other approaches consider various
methods such as swarm optimization [4], dynamic programming [5], genetic algorithms [6], inversion [7] and Bayesian
compressive sampling [8].
In direction-of-arrival estimation applications other specialized techniques have emerged for finding efficient array designs,
such as optimizing the corresponding Cramer-Rao error bound [9], or designing according to the nested array methodology
for increasing the available number of degrees of freedom [10], and theoretical results for estimation performance in these
settings have been developed [24], [25], [26]. Alternatively, other studies have addressed array design from the point of view
of performing efficient compressive sampling of a Wide Sense Stationary (WSS) process, e.g. for a setting with a parametrized
covariance matrix [11], for wide-band power spectrum estimation [12] and evaluating performance in the presence of modeling
errors [13]. However, these samplers are optimized for reconstructing the second order statistics of the scene and not the actual
scene realization.
In virtually all design methodologies some assumptions are made on the scene of interest. These represent beliefs, constraints
or knowledge that hold over the unobserved scene. A favorable design is one meeting requirements while taking into account
these assumptions. For example, in direction of arrival applications we may assume some limit on the number of point targets
present in the scene [9]. In beamforming we assume some separation level between objects of interest that implies a certain
resolution level is required [1], or some scene sparsity structure [14], [15].
In this paper we study the problem of inference on a scene of interest through measurements collected by an array of sensors.
The approach we take for designing array geometries is novel in that we consider settings where some Bayesian prior on the
This work was supported, in part, by NSF CAREER Award 1553284.
This work was supported, in part, by DARPA under Contract No. HR0011-16-C-0030, and by NSF under Grant No. CCF-1319828.
The authors are with the Department of Electrical Engineering and Computer Science, Massachusetts Institute of Technology, Cambridge, MA 02139 USA
(e-mail: {shulkind,stefje,gww}@mit.edu).

1

scene is available at the time of design and propose exploiting this knowledge by adapting the array geometry accordingly to
achieve efficient inference with a limited budget for sensors.
Prior knowledge is in fact not a rarity at all. Frequently, the same device is used to sense multiple similar scenes, where past
examples are indicative. A medical imaging device, for instance, is typically used to image the same organ across different
patients. In other cases, we may have prior knowledge in the form of scene properties such as smoothness or adherence to
spatial constraints. We incorporate such knowledge as a prior in a Bayesian model, such that sensing the scene is equivalent
to performing inference in this model, and the problem of array geometry design asks to select a geometry that optimizes the
quality of inference. This naturally turns out to be a combinatorial optimization problem, which can be extremely hard to solve
even approximately without exploiting additional structure.
We show how to choose a cost function for the array geometry design problem that holds the property of submodularity
[16]. A submodular set function is one that exhibits diminishing marginal gains, i.e. adding additional elements results in
diminishing benefit. Recently, there has been significant progress on the theory of optimizing submodular functions [17], [18],
[19], [20]. In particular, these results state that, while NP-hard, submodular maximization admits variants of greedy algorithms
that are guaranteed to achieve near-optimal solutions, i.e. within a constant factor. Submodularity has been used in connection
with sensor placement problems, for example in [21], [22], [23], however these works are not tailored to the far field scenes
and models that we focus on here.
Our connections and formulations open avenues for leveraging those results for efficient array design with strong guarantees.
We demonstrate this by showcasing the design of array geometries in settings with arbitrary apertures and settings where sensor
placement respects matroid combinatorial constraints by porting results from the growing body of literature on submodular
optimization. We demonstrate that together with our new adaptive formulations, the exploitation of prior knowledge leads to
higher quality inference at lower cost in terms of the number of sensing elements.
The paper is structured as follows: In Section II, we review far-field sensing and introduce a statistical framework for
modeling a-priori scene distributions. In Section III, we define a cost-function to quantify the inference loss associated with
array configurations and derive a corresponding array design optimization problem. We show that this optimization problem
belongs to the class of submodular optimization problems and adapt established algorithms and guarantees to our setting. In
Section IV, we review the theory of matroids and showcase its use in designing array geometries with combinatorial placement
constraints. Finally, in Section V we summarize numerical experiments exemplifying and validating our theory. The appendix
provides proofs for some of our claims.
II. P ROBLEM FORMULATION
In this section we formulate the array design problem. We focus on far-field sensing applications, although, as will become
apparent, the same techniques could also be generalized to other settings where the measurement process is linear. For simplicity
we consider a one dimensional setting as the extension to multiple dimensions is straightforward. We begin with a review of
the far-field sensing model to establish notation, and then pose the sensing problem as a Bayesian inference problem.
A. Far-Field Sensing
The far-field sensing setup is depicted in Fig. 1. At a distance from an observation axis x we have a scene of interest,
characterized by an illumination function β̃(θ) that depends on the angle θ∈[− π2 , + π2 ] between the direction of observation
and the axis of measurement broadside. We are to place N sensors along the observation axis at positions S ≡ {x1 , . . . , xN }.

Fig. 1. Far-field sensing. Antenna elements are depicted as blue dots.

The passive sensors pick up the radiation emitted from the scene, which is assumed to be narrow-band around a fixed
wavelength λ. Recall that in an ideal noiseless setting, r(x), the reading collected at position x, is given according to [1]:
1

π

Z+ 2
r(x)=

β̃(θ)ej

2π
λ x sin(θ)

Z+ 2
cos(θ)dθ=

−π
2

β(ψ)ej

4π
λ xψ

dψ,

(1)

− 21

where we have introduced ψ ≡ sin2 θ such that ψ ∈ [− 12 , + 21 ] and β(ψ) ≡ β̃(sin−1 (2ψ)).
Eq. (1) bears close resemblance to the definition of the Fourier transform. Namely, defining f (t) to be the Fourier transform
pair of β(ψ)
Z
Z
F
j2πtψ
f (t) = β(ψ)e
dψ  β(ψ)= f (t)e−j2πtψ dt
(2)
F −1

2

2
we immediately identify, comparing (1)
 2 and (2) that r(x) = f ( λ x). That is, collecting measurements at x ∈ 1{xn }1 is equivalent
to sampling the function f (t) at t ∈ λ xn . Notice that since the support of β(ψ) is restricted to ψ ∈ [− 2 , + 2 ] we observe
that f (t) is band-limited with bandwidth 1. Thus, as a consequence of the Whittaker-Kotelnikov-Shannon sampling theorem
we have that perfect reconstruction of f (t), and subsequently the scene β(ψ), is possible from an infinite set of samples taken
at t ∈ {. . . , −1, 0, +1, . . .} which are exactly the samples collected by an infinite λ2 -spaced array.
Finally, we consider the effect of noise by introducing the vector of N noisy measurements f˜ = [f˜1 , . . . , f˜N ]> modelled
according to:


2
˜
xn + wn
n = 1, . . . , N
(3)
fn ≡ r(xn ) + wn = f
λ

where wn are additive noise components. Stacked into a vector w = [w1 , . . . , wN ]> , we assume throughout that the noise is
a complex, circular, Gaussian random vector w ∼ CN (0, Σww ) [27]. We assume the noise to be independent across different
2
sensors, i.e., Σww = σw
IN , where IN is an N ×N identity matrix such that the noise is i.i.d. across different sensors.
B. Bayesian Formulation
The sensing problem we consider in this paper entails the estimation of the illumination function β(ψ) or equivalently its
inverse Fourier transform f (t) from the set of noisy measurements f˜. Even in the noiseless setting this problem is gravely
ill-posed as infinitely many wildly varying scenes map to any given finite set of observations1 .
To cope with this ill-posedness, some prior belief or knowledge pertaining to β(ψ) (or equivalently f (t)) must hence be
incorporated into the model, and this could be achieved in several ways. Wingham [29] proposed selecting one specific f (t) of
the multiple such functions consistent with the samples, namely the minimum norm solution. Alternatively, some constraints
or other preferences may be imposed on the solution by penalizing the inversion. For example one may require the solution
to satisfy some constraints, e.g. lie in some pre-specified sub-space of the function space, or regularize the inversion, e.g. to
promote smoothness or low total variation [30].
In this study we take a Bayesian approach and impose a prior distribution on the scene β(ψ). Subsequently, sensing is
equivalent to performing inference in this model. The prior may be assigned based on past observations over the distribution
of scenes or based on a-priori knowledge of scene properties.
We conceptually think of our scene as a natural image, and as such, our Bayesian formulation entails posing a prior that is
consistent with its expected characteristics, tailored to the specific application that is of interest to the user. Utilizing Gaussian
distributions as priors for natural scenes and images has been successfully applied in the past [31], [32], [33] and our work
uses similar notions, with the Gaussian prior specifically expressed in the frequency domain (e.g., as in [9]). We note that this
Gaussian description naturally captures the statistics of spread-out scenes and is not immediately suitable for serving as a prior
for point-source target models. The formulations and methods that we develop in subsequent sections can be adapted for use
with other scene priors in lieu of the Gaussian one we focus on for this study. In particular, our formulations may account for
settings with point-source scenes by posing a point-process prior (e.g. a determinantal point process, DPP [34]) that is better
suited for capturing priors of point-targets in specific environments. With this alternative prior we could follow similar steps to
develop efficient array designs for those applications. Studying such alternative priors and resulting formulations is a possible
interesting extension to our work.
In what follows, to be concrete we detail one specific choice for the form of the scene prior, expressed as a Gaussian
distribution in the frequency domain:
1) Discrete Representation: Assigning a prior on a continuous function β(ψ) may seem like a daunting task. However, in
many scenarios this can be considerably simplified if the scene may be efficiently expanded in a countable basis of functions
such that the prior may be imposed in thediscrete domain of expansion coefficients. In what follows, we consider expanding
β(ψ) by means of Fourier basis functions ej2πmψ |m = . . . , −1, 0, +1, . . . in ψ ∈ [− 21 , + 12 ]. Alternatively, other expansions
may be used with similar results. We write:
X
1 1
β(ψ) =
βm ej2πmψ ,
ψ ∈ [− , + ]
(4)
2 2
m
1

Z+ 2
βm ≡

β(ψ)e−j2πmψ dψ

(5)

− 12

where {βm } are the Fourier expansion coefficients. The usual Parseval relation holds:
Z
X
|β(ψ)|2 dψ =
|βm |2
m
1 The

mapping between a band-limited function f (t) and a finite set of its (generally non-uniform) samples {f (tn )} is not bijective [28].

(6)

3

In lieu of the prior on β(ψ) we will impose a prior distribution over {βm }. As will become evident such a description is
especially suited for applications involving real world smooth scenes β(ψ) as suggested by the following properties of the
Fourier series expansion [35]:
|m|→∞

(a) (Riemann-Lebesgue lemma) Let β(ψ) be any integrable function. Then |βm | −→ 0.
α
(b) Let β(ψ) ∈ C r where C r is the space of r-times continuously differentiable functions over some domain. Then |βm | ≤ |m|
r
∂r
with α = supψ | ∂ψr β(ψ)| [35].
Thus, for any nicely behaved scene β(ψ) the high frequency Fourier expansion coefficients diminish to zero with an
asymptotically polynomial rate determined by the level of smoothness.
2) Prior: In the sequel we use Gaussian priors on the coefficients {βm }:
2
βm ∼ CN (0, σm
)

(7)

2
where βm are independent, complex, circular and Gaussian, and σm
are the corresponding variances, which are assumed known
in advance. Using (6) we define the expected scene power:
Z
X
2
P ≡ E |β(ψ)|2 dψ =
σm
(8)
m

 2
The σm
can be set following some initial measurements of sample functions β(ψ) or taking into account prior knowledge.
For example if we have a-priori knowledge that β(ψ) ∈ C r we may use

1
m=0
2
(9)
σm
=
1
m
6= 0
2r
m
which is a very simple distribution respecting the polynomial variance decay. For the rest of the paper we adopt the prior in
(9) and take r = 1 to promote continuously differentiable functions. Also notice that the resulting prior distribution on β(ψ) is
2
Wide Sense Stationary (WSS), with the coefficients σm
determining the shape of the power spectrum, based on the available
prior information we have on the scene.
C. Observation Model
With the prior stated in the discrete {βm } domain as described above, our next goal is to circumvent β(ψ) by directly
stating the problem in terms of the measurements f˜n and the coefficients βm , replacing the continuous representation with
manageable discrete counterparts. Substituting the sum (4) into Equation (1) we have:
1

r(xn ) =

Z+ 2 X
− 12

m

βm ej2πmψ ej

4π
λ xn ψ

dψ =

X

Knm βm

(10)

m

where
1

Z+ 2
Knm ≡

2
2
ej2π( λ xn +m)ψ dψ = sinc(m+ xn )
λ

(11)

− 21

with sinc(x) ≡

sin(πx)
πx .

Plugging this into (3) we finally have the full observation model
X
n = 0, . . . , N − 1
f˜n =
Knm βm + wn

(12)

m

and the sensing problem amounts to estimating the coefficients {βm } given the observation vector f˜. As we have assumed a
Gaussian distribution for the noise vector w as well as for the prior over {βm } the posterior distribution p({βm } |f˜) turns out
to be Gaussian with a convenient analytic expression, as we detail in Section III-D.
III. A RRAY D ESIGN W ITH AN A PERTURE C ONSTRAINT
In the previous section we formulated the problem of far-field sensing in a Bayesian setting with a prior on the distribution
of the underlying scene. In this section we design corresponding array geometries to facilitate efficient sensing that exploit the
prior knowledge and model. Initially we consider simple constraints on the total number of sensors we can place and on the
aperture where these can be situated. Specifically, we assume that up to N sensors are free to be placed over some aperture
A on the real line, e.g. for simplicity we can consider A ≡ {x| − a ≤ x ≤ a} , a ∈ R+ . In Section IV we consider more
sophisticated combinatorial placement constraints and show that the same formulations we develop here may be adapted to
those more challenging use cases.

4

A. Formulating a Cost Function
In order to make the problem of optimal array design well-posed in our Bayesian setting we need to specify a performance
measure that will be used as a cost function to compare different designs and choose the best one. Revisiting our observation
model (12) we notice that the array design determines the coefficients Kmn through the set of sensor positions S, as defined in
Section II-A. With any given design the sensing problem amounts to inferring the posterior p({βm } |f˜). A natural performance
measure in this setting quantifies the quality of inference, i.e., the information gained by performing the sensing experiments
which results in updating our beliefs about the coefficients {βm } from the prior p({βm }) to the posterior p({βm } |f˜). This
problem has been extensively studied in the context of statistical inference and experimental design [36], [37]. If we are
interested in making general inference to learn the state of the world (represented by the posterior distribution of the Fourier
coefficients) then the natural performance measure, which is referred to as Bayes D-optimality is the mutual information
between the expansion coefficients and the collected measurements [36], [37]:
I(f˜S ; {βm }) = H({βm }) − H({βm } |f˜S )

(13)

where I(·; ·) is the mutual information and H(·) the Shannon entropy. The subscript S explicitly emphasizes the dependence
of the measurements on the set of sensor positions S. Notice that this cost function measures the reduction in uncertainty of
the scene expansion coefficients before and after measurements are made. The larger the mutual information the more we trust
the values of the coefficients {βm } |f˜S .
With the cost function in place the array design problem becomes
S ? = argmax I(f˜S ; {βm }),

(14)

S⊆A,|S|≤N

which is an NP-hard combinatorial optimization problem. However, we will show later that we can obtain a constant-factor
approximation for (14) using efficient computational techniques.

B. Finite Dimensional Approximation
Solving (14) under our model (12) involves manipulations of the infinite sequence {βm } which may not be amenable to
computer representation. To make our formulation tractable we approximate the infinite sequence {βm } with a finite truncated
set of coefficients {βm |m ∈ M}, where M is some finite set. Stacked in vector form β, we consider the simplified finite
dimensional approximation of (12):
fˆS = KS β + w

(15)

where KS is an N ×|M| matrix formed by restricting Knm on m ∈ M and the dependency on the sampling set S again
made explicit. Notice that hat notation is replacing the previous tilde.
We show that the approximate finite-dimensional model (15) is a good proxy for the original infinite-dimensional model
(12) for a suitably selected M. Indeed, if M is chosen to only exclude those βm coefficients that in expectation contribute a
marginally small part of the energy of the full infinite sequence {βm } then the mutual information derived from the approximate
model (15) will closely track that of the infinite dimensional model in (12). More precisely we have the following result:
P 2
3
2
2
Lemma III.1. Let the prior on {βm } be i.i.d. according to βm ∼ CN (0, σm
), M fixed, and  ≡
σm satisfies  < σw
N−2 .
m∈M
/

We have:
3

−N log(1+

3

N 2
N 2
)≤I(f˜S ; {βm })−I(fˆS ; β)≤−N log(1− 2 )
2
σw
σw

Proof. See Appendix A.
3

→0
2
ˆ
˜
By virtue of the last lemma and N log(1± N
2 ) −→ 0 we have that I(fS ; β) is an arbitrarily accurate proxy for I(fS ; {βm })
σw
for  small enough, such that in lieu of problem (14) we now continue with the simplified finite dimensional approximation:

S ? = argmax I(fˆS ; β)
S⊆A,|S|≤N

and the results will be accurate to within the approximation bounds from Lemma III.1.

(16)

5

C. Grid Discretization
Next we turn to discretizing the aperture A to cast the array design problem in the form of a generic finite selection problem.
We thus restrict the choice of sampling positions to the finite set
o
n
(17)
V ≡ x1 , . . . , x|V| ⊂ A
For the sequel we take V to be a uniform δ-spaced grid of positions in A. We adapt the array design problem (16) accordingly
as:
Sd? = argmax I(fˆS ; β)

(18)

S⊆V,|S|≤N

with the subscript d implying discretization. The next result can be used to quantify the level of discretization δ necessary to
guarantee performance close to optimum within some specified error bound:
Lemma III.2. With V a uniform grid of sampling positions with distance δ between adjacent positions we have:
3

I(fˆSd? ; β) ≤ I(fˆS ? ; β) ≤ I(fˆSd? ; β)+N log(1+

4δP (1+δ)N 2
)
2
λσw

Proof. See Appendix B.
With this last lemma in place the array design problem (16) may be further approximated in the more convenient finite
combinatorial problem form of (18) with guarantees on the accuracy of the resulting designs. In the sequel we assume that
δ is chosen such as to meet desired accuracy levels, as prescribed in Lemma III.2, and work with the simplified formulation
(18).
D. Sensing
With the observation model of (15), coupled with a Gaussian distribution for the noise vector w and a Gaussian prior for the
coefficients vector β, calculation of the posterior distribution β|fˆS is particularly simple and can be performed analytically.
Concretely, we have as a result of all random variables being Gaussian β|fˆS ∼ CN (µ̂, Σ̂) and the parameters are given
according to the conventional Gaussian conditional parameters:
fˆ
µ̂ = Σβ fˆΣ−1
fˆfˆ S
Σ†
Σ̂ = Σββ − Σβ fˆΣ−1
fˆfˆ β fˆ
where Σβ fˆ = Σββ KS† ,

(19)

2
Σfˆfˆ = KS Σββ KS† + Σww and Σββ = diag[σ12 , . . . , σM
].

E. Submodular Optimization
The next step is to prescribe an efficient algorithm for the solution of (18). As we will show shortly (18) is an instance of
a known NP-hard problem such that it is widely believed that no efficient algorithm for its solution exists. However, due to
the structure of the cost function an efficient approximation algorithm is known to exist with strong theoretical guarantees. In
this subsection we survey the relevant results and adapt them to our needs.
1) Submodularity: We begin by invoking the submodularity property of set functions [16]:
Definition III.1. Let G : 2V → R be a set function.
(a) G is submodular if it satisfies the property of decreasing marginals: ∀S, T ⊆ V such that S ⊆ T and x ∈ V\T it holds
that G(S ∪ {x}) − G(S) ≥ G(T ∪ {x}) − G(T ).
(b) G is monotonic (increasing) if ∀S, T ⊆ V s.t. S ⊆ T we have G(S) ≤ G(T ).
As it turns out, our cost function is monotonic and submodular as the next result shows (this is almost identical to corollary
4 in [38]. We reproduce it here with adaptations to our setting for completeness):
Lemma III.3. Let V be defined as before, and define the set function G : 2V → R according to G(S) = I(fˆS ; β). Then G is
submodular and monotonic (increasing).
Proof. Expanding the mutual information according to I(x; y) = H(x)−H(x|y) we have:
G(S∪ {x})−G(S) = H(fˆS∪{x} )−H(fˆS )
−[H(fˆS∪{x} |β)−H(fˆS |β)] = H(fˆ{x} |fˆS )−H(fˆ{x} |β)

(20)

6

where in the last equality we used the conditional independence of the components of fˆS∪{x} given β. Substituting T for S
we immediately get:
[G(S∪ {x})−G(S)]−[G(T ∪ {x})−G(T )]
= H(fˆ{x} |fˆS )−H(fˆ{x} |fˆT )

(21)

Using S ⊆ T we have H(fˆ{x} |fˆS ) ≥ H(fˆ{x} |fˆT ) such that G(S∪ {x})−G(S) ≥ G(T ∪ {x})−G(T ) and G is submodular.
To prove monotonicity it is enough to show G(S∪ {x})−G(S) ≥ 0. This time expand the mutual information according to
I(x; y) = H(y)−H(y|x):
G(S∪ {x})−G(S) = H(β|fˆS )−H(β|fˆS∪{x} )

(22)

Conditioning can never increase entropy so H(β|fˆS )≥H(β|fˆS∪{x} ) and the result follows.
2) Efficient Solvers: Our optimization problem (18) is the maximization of a monotonic submodular function. A greedy
algorithm, shown as Algorithm 1, solves this problem to within the best possible approximation factor, as the following
fundamental theorem states.
Algorithm 1 Greedy Submodular Maximization
1: function G REEDY M AX (G(·), V, N )
2:
S←∅
3:
for i = 1 to N do
4:
x? = argmaxx∈V\S G(S ∪ {x})
5:
S ← S ∪ {x? }
6:
Return S
The following theorem guarantees that the greedy algorithm achieves approximately optimal performance:
Theorem III.1. (Nemhauser [16]) Let G be a monotonic, submodular set function and S ? = argmax G(S).
S⊆V,|S|≤N

Let S gr be the set retrieved by the greedy maximization Algorithm 1. We have the following guarantee for the performance of
the greedy algorithm:
1
1
G(S gr ) ≥ (1 − (1 − )N )G(S ? ) ≥ (1 − )G(S ? )
N
e
Moreover, no polynomial time algorithm can provide a better approximation guarantee unless P=NP [39].
Combining the guarantees of Theorem III.1, Lemma III.1 and III.2 we derive an approximation bound on the original
problem (14):
Corollary III.1.
"
#
3
2
2
1
λσ
+4δP
(1+δ)N
w
(1− ) I(f˜S ? ; {βm })−N log
≤I(fˆS gr ; β)
2 −λN 23
e
λσw

(23)

We hence apply Algorithm 1 to solve our optimization problem (18). The algorithm runs in time O(|V|N ), linear in the
size of the set V and the number of selected elements N [40] such that it is easily implementable for problems of large size.
However, as it turns out even more efficient variants of the algorithm have been introduced and studied. One such variant
commonly referred to as the ’lazy greedy’ algorithm was studied in [40] and it was shown to offer substantial running-time
improvements in practice (with an unlikely worst-case theoretical performance upper bounded by that of the conventional
greedy algorithm). Our numerical experiments described in Section V implement this more efficient variant to reduce running
time.
3) Improved Approximation Bounds: While Theorem III.1 guarantees an approximation bound of (1 − 1e ) ≈ 63% for the
efficient greedy algorithm this guarantee is not tight. It is possible to derive a tighter data-dependent online bound on the gap
between the cost of the greedy solution G(S gr ) and that of the optimal solution G(S ? ) [41]. Specifically, we have:
X
G(S gr )≥G(S ? )− max
(G(S gr ∪ {e})−G(S gr ))
(24)
B:|B|≤N

e∈B

which takes O(|V| log |V|) evaluations of G(S) to compute and sort. We use (24) to improve the distance from optimality
bound in some of our numerical solutions in Section V.

7

F. Design Example: A Simple Ideal Setting
In the previous subsections we formulated the array design problem in a setting with constraints on the aperture A and
the number of sensors N and showed how a greedy algorithm (Algorithm 1) is guaranteed to efficiently find an approximate
solution.
Here, we study a particular instance of that problem, where the Signal to Noise Ratio (SNR) is high, and the aperture is
effectively unconstrained (the N sensors may be placed anywhere on the real line). Under these conditions a truncated λ2 -spaced
array is traditionally considered the design of choice in the conventional non-Bayesian setting (this is a truncated version of
the infinite λ2 -spaced design mentioned in Section II-A). We show next that the truncated λ2 -spaced design naturally emerges as
the approximately optimal solution as retrieved by our schemes in Bayesian settings where the a-priori β distribution satisfies
some conditions. Specifically, we have the following result:
Theorem III.2. Consider the high SNR regime σP2 → ∞ and assume the prior from (7) takes a symmetric, monotonically
w
2
2
2
2
decreasing form, i.e. σm
= σ−m
and σm
≥ σm
whenever 0 ≤ m1 < m2 . In addition, take V as an arbitrarily dense set of
1
2
sampling points on R, and M = −M, . . . , M with M → ∞.
We then have that a greedy solver on (18) will return a length N , λ2 -spaced truncated uniform array centered around x = 0.
Proof. See Appendix C.
The last theorem studies one class of simple idealized problems where the greedy solution is reminiscent of generic nonBayesian array designs. However, notice that our formalism is also useful in more challenging design problems such as when
the aperture A takes on arbitrary forms, and the effects of noise and application-tailored priors are considered. Furthermore, in
Section IV we demonstrate that additional combinatorial constraints may be naturally incorporated into our formulation such
that even more challenging problems become manageable.
IV. A RRAY D ESIGN WITH M ATROID C ONSTRAINTS
In Section III we formalized the array design problem in a setting where we imposed constraints on the aperture and the
number of sensors. Specifically, the constraints were S ⊆ V, |S| ≤ N . In many practical scenarios these may be too simplistic
to accurately represent real world design constraints. For example in applications where sensors are heavy and mounted on
support beams we may want to restrict the number of sensors in specific sections of the aperture.
In this section we briefly review key elements from matroid theory which is a branch in combinatorics [42] and survey
results from submodular optimization with matroid constraints guaranteeing the existence of efficient approximate solvers for
this class of problems. We continue to show that these mathematical structures may be utilized to impose constraints of interest
in array design enriching the set of problems our Bayesian formulation can describe and solve.
A. Submodular Optimization with Matroid Constraints
We begin by defining matroids and their corresponding independent sets [16]:
Definition IV.1. A finite matroid M is a pair (V, I) where V is a ground set and I is a collection of subsets of V (the
independent sets) that satisfies the following properties:
1) The empty set is independent: ∅ ∈ I
2) A subset of an independent set is independent:
X ⊂ Y, Y ∈ I ⇒ X ∈ I
3) If X is an independent set and Y is a larger indepedent set, X can be augmented to a larger independent set by adding
an element from Y \X:
X, Y ∈I, |X|<|Y | ⇒ ∃e ∈ Y \X s.t. X∪ {e} ∈ I
A matroid structure may be used to classify subsets of a ground set V into permissible subsets which belong to I and
non permissible subsets which do not belong to I. In the next subsection we show that using this formalism we can express
interesting array design constraints.
From the theory of submodular optimization we have the following results for submodular optimization with matroid
constraints [43], [44]. Let M = (V, I) be a matroid and G(S) a monotonic, submodular set function. There exists an efficient
approximate solver for the problem argmaxS∈I G(S). Specifically, a greedy solver (maximizing the immediate marginal benefit
at each step) taking the form
(
)
S gr ←S gr ∪

argmax

[G(S gr ∪ {e})−G(S gr )]

(25)

e:e∈S
/ gr ,S gr ∪{e}∈I

and stopping when no more elements e can be added is guaranteed to achieve a half-approximation bound:
1
G(S gr ) ≥ max G(S).
2 S∈I

1
The constant factor may be tightened to 1 − e by utilizing specialized randomized algorithms [44].

(26)

8

B. Combinatorial Constraints in Array Design
Here we invoke one specific well known matroid structure and show its application in expressing useful array design
constraints. Let V be a ground
S set of gridTpoints where sensors are allowed to be placed as before. Let V1 , . . . , VK be a
partition of the set V, i.e.
Vk = V, Vi Vj = ∅, ∀i 6= j, and let N, n1 , . . . , nK be a set of integers.
k

We define the (cardinality constrained) partition matroid [44] M = (V, I) with the following definition for the collection
of independent sets: a subset S⊆V is an independent subset S∈I if it holds |S∩V| ≤ N, |S∩Vj | ≤ nj , ∀j.
In the context of array design the partition matroid may be useful in expressing practical constraints over sensor placement
configurations. For example if the subsets Vi represent closed line sections, e.g. a physical partitioning of the aperture into
zones, and I represents the collection of all permissible designs then the structure of the matroid limits the number of sensors
that may be placed in the ith zone to ni which may be an important engineering constraint coupled with some specific
application. We solve:
S ? = argmaxS∈I I(fˆS ; β)

(27)

Applying the results from the previous subsection we immediately have an efficient approximate solver for the array design
problem coupled with a partition matroid constraint. In Section V we detail such a design for one numerical example.
V. N UMERICAL E XPERIMENTS
In this section we perform numerical experiments validating our theoretical results and exemplifying them. We showcase
an array design with cardinality and aperture constraints as prescribed in Section III and design arrays adhering to matroid
constraints as prescribed in Section IV.
Our initial setting is as follow. We fix λ = 1 throughout as the wavelength only serves to scale the x axis. The aperture is set
as A = {x| − 3.5 ≤ x ≤ 3.5} and the selection set V is chosen as a uniform grid of 113 positions from A spaced δ = 0.0625
apart. We set out to design an array consisting of N = 11 sensor locations. The prior for {βm } is set as per (9) with r = 1
and normalized to sum to P = 1. For the simulations we consider the truncated vector β formed when restricting the set of m
coefficients to 901 consecutive elements centered around the origin, i.e., we set M = {−450, . . . , +450}. For the preliminary
design we implement the lazy greedy algorithm of Section III and plot the results in the left column of Fig. 2 as a function
of the SNR which we define here as NPσ2 . Blue markers denote the full selection set V and red markers delineate the active
w
S selected by the algorithm.
In the high SNR regime (SNR = 10dB or higher values) the resulting design is a truncated λ2 uniform array as predicted
according to Theorem III.2. As the SNR decreases the reliability of the measurements deteriorates and the algorithm prefers
locating samplers right next to each other on expense of widening the array as this serves to average out the noise.
The performance in terms of mutual information I(fˆS ; β) for the selected locations S appears in the title of the plots (in
natural units). Notice for example that for the 5dB SNR design the achieved mutual information is 12.54. Using Theorem
III.1 we have that the optimal design cannot achieve mutual information better than 1−1 1 12.54 = 19.83. This bound can be
e
improved using the improved bounding method briefly described following Theorem III.1 to show that the optimal performance
is not greater than 17.45.
The truncation level dictated by our choice of M translates to  ∼
= 1e−4 and the truncation bounds from Lemma III.1 read
(for the lower, extreme SNR case): −0.45≤I(f˜S ; {βm })−I(fˆS ; β)≤0.47. We find empirically that these bounds are extremely
loose and M can be shrunk considerably without substantially compromising accuracy. To achieve a similar upper bound on
I(fˆS ? ; β) − I(fˆSd? ; β) as per Lemma III.2 a discretization level of δ ∼
= 2e−4 is needed. However, we empirically find that
our choice of δ = 0.0625 is accurate enough as further refining the grid does not significantly change the design. Our lemmas
prove to be extremely pessimistic as is expected given that the proofs take into account worst-case scenarios.
The array geometries above, derived according to the formulations of Section III, are designed to optimize the quality
of inference between the measurements and the scene expansion coefficients β. Many sensing applications of interest
specifically involve imaging the scene, that is reconstructing β(ψ) from the measurements. Our next experiment was designed
to empirically evaluate the Mean Square Error (MSE) performance in scene reconstruction from measurements collected using
the prescribed designs. First, we designed five array geometries as described above, optimized for several target SNR levels
{30dB, 12dB, 10dB, 5dB, 0dB}. We set up a Monte-Carlo experiment where 1000 scenes were randomly drawn from the
distribution of Section II-B2. For each scene, noisy measurements were collected by each of the five optimized arrays. The
measurements were repeated with five different synthetic noise levels corresponding to the five target SNR levels.
We repeatedly performed maximum likelihood estimation [45] of the expansion coefficients β, and synthesized an estimated
scene β̂(ψ) according to (4). The MSE discrepancy between β̂(ψ) and the true scene is depicted in Fig. 3. It is evident that
the quality of inference criterion is indicative of MSE performance, as each of the five geometries yielded the best MSE
performance at its specified target SNR level.
Next we solve a corresponding set of design problems with matroid constraints installed to limit the number of sensors in given
aperture segments. Specifically, we use the (cardinality constrained) partition matroid from Section IV-B with N = 11, ni =
1, ∀i and Vi spanning consecutive line segments of length 0.5: Vi =[−0.25+0.5·i, +0.25+0.5·i). The matroid constraints limits

9

Fig. 2. (left) Array designs with an aperture constraint and various SNR levels. (right) Array designs with combinatorial placement constraints and various
SNR levels.

Fig. 3. Reconstruction MSE for various arrays designed for fixed SNR levels, with actual SNR levels swept.

the proximity between sensor elements, which may be useful in real world use cases. We plot the results for the matroid
constrained designs in the right column of Fig. 2. Notice that while the theoretical guarantees pertaining to the greedy matroid
optimization scheme of Section IV is 0.5 compared to 1− 1e with the cardinality constraints of Section III, the actual performance
achieved in the constrained design instances is not far from those achieved with the simple cardinality constraints.
VI. C ONCLUDING REMARKS
We introduced a novel framework for designing sensor arrays. Our setting is Bayesian in the sense that our model incorporates
the notion of prior beliefs about the scene of interest. Our goal is designing arrays that are adapted to the scene and perform
efficient inference. We showed that this NP-hard combinatorial selection problem may be efficiently approximated by porting
results from the theory of submodular set function optimization.

10

Initially, we showed how to apply our formalism in design problems with straightforward cardinality and aperture constraints.
Later we showed that more challenging combinatorial constraints may be enforced by utilizing results and concepts from the
field of matroids and submodular optimization with matroid constraints.
Our formulation connecting the array design problem and submodularity may further be extended to reap additional benefits,
namely tackling problems such as robust array design or adaptive design that evolve as the scene is learned. We leave the
treatment of some of these and other problems to separate study [46].
A PPENDIX A
P ROOF OF L EMMA III.1
For simplicity of notation, we suppress the subscript S throughout. Begin by expanding the mutual information expressions:
I(f˜; {βm }) = H(f˜) − H(f˜| {βm })
I(fˆ; β) = H(fˆ) − H(fˆ|β)

(28)

Examining (12) and (15) we have H(f˜| {βm }) = H(fˆ|β) = H(w) and so:
I(f˜; {βm }) − I(fˆ; β) = H(f˜) − H(fˆ)

(29)

Next, notice that f˜ and fˆ are both circular, complex, Gaussian random N -length vectors, such that their entropies are given
according to [27]:
H(f˜) = log((πe)N det(Σ̃))
H(fˆ) = log((πe)N det(Σ̂))

Σ̃ij = E[f˜i f˜j∗ ]
Σ̂ij = E[fˆi fˆ∗ ]
j

(30)

∗
2
Further expand using the independence between w and βm and E[βm βm
0 ] = δmm0 σm :
X
X
∗
∗
∗
Σ̃ij = E[(
Kim βm + wi )(
Kjm
0 βm0 + wj )]
m

= [Σww ]ij +

X

m0
∗
2
Kim Kjm σm

(31)

m

Σ̂ij = E[(

X

Kim βm + wi )(

X

∗
∗
∗
Kjm
0 βm0 + wj )]

m0 ∈M

m∈M

X

= [Σww ]ij +

∗
2
Kim Kjm
σm

(32)

m∈M

Comparing (31) and (32) and using:
X

∗
2
Kim Kjm
σm
≤

m∈M
/

X

∗
2
Kim Kjm
σm
≤

X

2
σm
=

(33)

m∈M
/

m∈M
/

we have |Σ̃ij − Σ̂ij | ≤  such that we may write:
Σ̂ = Σ̃ + X

(34)

for some N × N matrix X satisfying |Xij | ≤ 1. We use (34) to bound the determinants. Σ̃ is positive-definite and invertible
such that we can write:
det(Σ̂) = det(Σ̃ + X) = det(Σ̃)det(IN + Σ̃−1 X)

(35)

H(fˆ)−H(f˜) = log(det(IN +Σ̃−1 X))= log(det(X̃))

(36)

Substituting (35) in (30) we have:

with X̃ ≡ IN + Σ̃−1 X. We turn next to bounding the term log(det(X̃)). First notice:
[Σ̃−1 X]ij = 

X

Σ̃−1
im Xmj ≤ 

m

≤

X

Σ̃−1
im Xmj

m

√
−1
Σ̃−1
k∞ ≤  N kΣ̃−1 k2
im ≤ kΣ̃

m

√

X

√
=  N σmax (Σ̃−1 ) =  N

√
1
 N
≤ 2
σw
σmin (Σ̃)

(37)

11

√
Where we have used the matrix norm equivalence kAk∞ ≤ N kAk2 (for N ×N matrices) and σmax (·) (σmin (·)) is the
2
maximal (minimal)
singular value such that σmin (Σ̃) ≥ σw
. Thus we have that X̃ has
√ diagonal elements centered around 1:
√
P
−1)
 N
X̃im ≤  Nσ(N
.
X̃ii − 1 ≤ σ2 and the row-sums over non-diagonal entries satisfy
2
w

w

m6=i

Applying the Gershgorin circle theorem we have for the eigenvalues of X̃:
√
√
 NN
 NN
1−
≤
|λ
|
≤
1
+
i
2
2
σw
σw
det(X̃) is a positive real number as the quotient of the determinants of two positive definite matrices det(X̃) =
Q
Q
that we can write det(X̃) = i λi = i |λi | and consequently:
3

N log(1 −

det(Σ̂)
det(Σ̃)

such

3

N 2
N 2
) ≤ log(det(X̃)) ≤ N log(1 + 2 )
2
σw
σw

(38)

This finally leads to:
3

−N log(1+

3

N 2
N 2
) ≤ I(f˜; {βm })−I(fˆ; β) ≤ −N log(1− 2 )
2
σw
σw

(39)

A PPENDIX B
P ROOF OF L EMMA III.2
The left inequality is trivial. We have I(fˆSd? ; β) ≤ I(fˆS ? ; β) as the second optimization is over a larger set.
To prove the right inequality we show that for every S ⊆ A there is Sd ⊆ V such that
3

I(fˆS ; β) ≤ I(fˆSd ; β)+N log(1+

4δP (1+δ)N 2
)
2
λσw

we will show this for |S|=N but the proof for other cardinalities is identical.

With distance δ between adjacent elements of V, for every S= {x1 , . . . , xN } ⊂ A there is a set Sd = xd1 , . . . , xdN
that Sd ⊆ V and |xi − xdi | ≤ 2δ for all i. We have, similarly to Appendix A:

(40)

such

I(fˆS ; β) − I(fˆSd ; β) = H(fˆS ) − H(fˆSd )

(41)

H(fˆS ) = log((πe)N det(Σ̂S ))
H(fˆS ) = log((πe)N det(Σ̂Sd ))

(42)

Using the model (15):

d

where:
Σ̂S ≡ E[fˆS fˆS† ] = KS Σββ KS† + Σww
Σ̂Sd ≡ E[fˆS fˆ† ] = KS Σββ K † + Σww
d

Sd

d

Sd

†

Σββ ≡ E[ββ ]

(43)

and H(fˆS ) − H(fˆSd ) = log(det(Σ̂S )) − log(det(Σ̂Sd )).
Both KS and KSd are size N × |M| matrices defined as per the definition in (11), such that:
2
2
|[KS ]nm −[KSd ]nm | = |sinc(m+ xn )−sinc(m+ xdn )|
λ
λ
2
2
≤ 2 |xq −xdq | ≤ δ
λ
λ

(44)

where for the first inequality we have used the fact that sinc(·) is Lipschitz with constant smaller than 2. We can thus define
∆ ≡ KS − KSd and we have |∆nm | ≤ λ2 δ. Substitution in (43) yields:
Σ̂S = Σ̂Sd + ∆Σββ ∆† + ∆Σββ KS† d + KSd Σββ ∆†

(45)

12

We bound the perturbation terms by noticing:
X

[∆Σββ ∆† ]ij =

∆im [Σββ ]mm ∆jm

m

X
2
4
≤ ( δ)2
[Σββ ]mm ≤ 2 δ 2 P
λ
λ
m
[∆Σββ KS† d ]ij =

X

∆im [Σββ ]mm [KSd ]jm

m

≤

X
2
2
δ·1·
[Σββ ]mm ≤ δP
λ
λ
m

(46)

and overall we have:
|Σ̂Sij − Σ̂Sijd | ≤

4
4
4
δP + 2 δ 2 P = δP (1 + δ)
λ
λ
λ

(47)

define: 0 ≡ λ4 δP (1 + δ) and we have
Σ̂S = Σ̂Sd + 0 X

(48)

with N × N matrix X satisfying |Xij | ≤ 1 which is akin to (34). We thus port the results from Appendix A here (we only
need the lower bound):
3

0 N 2
−N log(1 + 2 ) ≤ I(fˆSd ; β) − I(fˆS ; β)
σw

(49)

which, upon substitution of 0 is equivalent to (40).
A PPENDIX C
P ROOF OF T HEOREM III.2
The greedy algorithm sequentially selects elements according to the rule x? = argmaxx∈V\S I(fˆS∪{x} ; β) where S is the
set of elements selected so far. We recursively show that the added elements can be selected on a λ2 -spaced grid centered
around x = 0. Expanding the mutual information as in Appendix A we have:
argmaxx∈V\S I(fˆS∪{x} ; β)=argmaxx∈V\S H(fˆS∪{x} )

(50)

We begin by showing that the first selected element is x1 = 0. Indeed, using the results from Appendix A:
argmax H(fˆ{x} ) = argmax det(Σ̂11 ) = argmax Σ̂11
x∈V

x∈V

where again using Appendix A and under the assumptions of the theorem (high SNR):
X
X
2
∗
2
2
Σ̂11 =
K1m K1m
σm
=
sinc2 (m+ x)σm
λ
m
m
X
2
2
2
≤ σ0
sinc (m+ x) = σ02
λ
m
2
where we used σm
≤ σ02 ,

(51)

x∈V

∀m 6= 0 and the identity:
X
sinc(m+a)sinc(m+b) = sinc(b−a)

(52)

(53)

m

It easy to see that in (52) equality is achieved for the choice x1 = 0 which
 is the claim.
Next, assume that the greedy algorithm has already picked a set S = x1 , . . . , x|S| of adjacent elements on a λ2 -spaced
grid centered around x = 0 and show that the next element to be added is an adjacent location on the same λ2 -spaced grid.
We have using H(x, y)=H(x)+H(y|x):
2
argmax H(fˆS∪{x} )=argmax H(fˆ{x} |fˆS )=argmax σx|S
x∈V\S

x∈V\S

(54)

x∈V\S

2
where σx|S
is the conditional variance of the Gaussian observation collected at x given the Gaussian observations made at the
set S:
†
2
σx|S
= σx2 − ΣxS Σ−1
SS ΣxS

(55)

13

with the usual definitions for the covariance matrices (as in Appendix A):
X
2
2
2
2
2
[ΣxS ]1i =
sinc(m+ x)sinc(m+ xi )σm
=sinc(m(i)+ x)σm(i)
λ
λ
λ
m
X
2
2
2
2
[ΣSS ]ij =
sinc(m+ xi )sinc(m+ xj )σm
=δii σm(i)
λ
λ
m
where in the last equations almost all sinc(·) functions nulled out as the xi ’s are situated on a
m(i)≡− λ2 xi , such that I ≡ {m(i)} is a set of consecutive integers. Additionally, we have:
X
2
2
σx2 =
sinc2 (m+ x)σm
λ
m

(56)
λ
2

grid, and we have defined

(57)

Substituting back into (55) we have:
2
σx|S
=

X

=

X

m

m∈I
/

X
2
2
2
2
sinc2 (m+ x)σm
−
sinc2 (m(i)+ x)σm(i)
λ
λ
i∈I

2
2
sinc (m+ x)σm
λ
2

(58)

and this is similar to the optimization over the selection of the first location in (52) with the optimum achieved by selecting
/ I which is the next adjacent location on the λ2 grid which completes the proof.
x such that λ2 x = m for the first m ∈
R EFERENCES
[1] H. L. Van Trees, Detection, estimation, and modulation theory, optimum array processing. John Wiley & Sons, 2004.
[2] H. Krim and M. Viberg, “Two decades of array signal processing research: the parametric approach,” Signal Processing Magazine, IEEE, vol. 13, no. 4,
pp. 67–94, 1996.
[3] R. L. Haupt, “Thinned arrays using genetic algorithms,” Antennas and Propagation, IEEE Transactions on, vol. 42, no. 7, pp. 993–999, 1994.
[4] M. M. Khodier and C. G. Christodoulou, “Linear array geometry synthesis with minimum sidelobe level and null control using particle swarm
optimization,” Antennas and Propagation, IEEE Transactions on, vol. 53, no. 8, pp. 2674–2679, 2005.
[5] M. Skolnik, G. Nemhauser, and J. Sherman, “Dynamic programming applied to unequally spaced arrays,” IEEE Transactions on Antennas and
Propagation, vol. 12, no. 1, pp. 35–43, 1964.
[6] R. L. Haupt and D. H. Werner, Genetic algorithms in electromagnetics. John Wiley & Sons, 2007.
[7] B. P. Kumar and G. Branner, “Design of unequally spaced arrays for performance improvement,” IEEE Transactions on Antennas and Propagation,
vol. 47, no. 3, pp. 511–523, 1999.
[8] G. Oliveri, M. Carlin, and A. Massa, “Complex-weight sparse linear array synthesis by bayesian compressive sampling,” Antennas and Propagation,
IEEE Transactions on, vol. 60, no. 5, pp. 2309–2326, 2012.
[9] H. Gazzah and S. Marcos, “Cramer-rao bounds for antenna array design,” IEEE Transactions on Signal Processing, vol. 54, no. 1, pp. 336–345, 2006.
[10] P. P. Vaidyanathan and P. Pal, “Sparse sensing with co-prime samplers and arrays,” IEEE Transactions on Signal Processing, vol. 59, no. 2, pp. 573–586,
2011.
[11] D. Romero, R. Lopez-Valcarce, and G. Leus, “Compression limits for random vectors with linearly parameterized second-order statistics,” IEEE
Transactions on Information Theory, vol. 61, no. 3, pp. 1410–1425, 2015.
[12] D. D. Ariananda and G. Leus, “Compressive wideband power spectrum estimation,” IEEE Transactions on signal processing, vol. 60, no. 9, pp.
4775–4789, 2012.
[13] A. Koochakzadeh and P. Pal, “Performance of uniform and sparse non-uniform samplers in presence of modeling errors: A cramér-rao bound based
study,” IEEE Transactions on Signal Processing, vol. 65, no. 6, pp. 1607–1621, 2017.
[14] J. D. Krieger, Y. Kochman, and G. W. Wornell, “Design and analysis of multi-coset arrays,” in 2013 IEEE International Conference on Acoustics, Speech
and Signal Processing. IEEE, 2013, pp. 3781–3785.
[15] ——, “Multi-coset sparse imaging arrays,” IEEE Transactions on Antennas and Propagation, vol. 62, no. 4, pp. 1701–1715, 2014.
[16] G. L. Nemhauser, L. A. Wolsey, and M. L. Fisher, “An analysis of approximations for maximizing submodular set functions,” Mathematical Programming,
vol. 14, no. 1, pp. 265–294, 1978.
[17] S. Fujishige, Submodular functions and optimization. Elsevier, 2005, vol. 58.
[18] N. Buchbinder, M. Feldman, J. Seffi, and R. Schwartz, “A tight linear time (1/2)-approximation for unconstrained submodular maximization,” SIAM
Journal on Computing, vol. 44, no. 5, pp. 1384–1402, 2015.
[19] J. Vondrák, “Symmetry and approximability of submodular maximization problems,” SIAM Journal on Computing, vol. 42, no. 1, pp. 265–304, 2013.
[20] B. Mirzasoleiman, A. Karbasi, R. Sarkar, and A. Krause, “Distributed submodular maximization: Identifying representative elements in massive data,”
in Advances in Neural Information Processing Systems, 2013, pp. 2049–2057.
[21] A. Krause, A. Singh, and C. Guestrin, “Near-optimal sensor placements in gaussian processes: Theory, efficient algorithms and empirical studies,” Journal
of Machine Learning Research, vol. 9, no. Feb, pp. 235–284, 2008.
[22] M. Shamaiah, S. Banerjee, and H. Vikalo, “Greedy sensor selection: Leveraging submodularity,” in 49th IEEE conference on decision and control (CDC).
IEEE, 2010, pp. 2572–2577.
[23] A. Krause and C. Guestrin, “Submodularity and its applications in optimized information gathering,” ACM Transactions on Intelligent Systems and
Technology (TIST), vol. 2, no. 4, p. 32, 2011.
[24] P. Stoica and A. Nehorai, “Performance study of conditional and unconditional direction-of-arrival estimation,” IEEE Transactions on Acoustics, Speech,
and Signal Processing, vol. 38, no. 10, pp. 1783–1795, 1990.
[25] M. Wang and A. Nehorai, “Coarrays, music, and the cramér–rao bound,” IEEE Transactions on Signal Processing, vol. 65, no. 4, pp. 933–946, 2017.
[26] A. Koochakzadeh and P. Pal, “Cramér–rao bounds for underdetermined source localization,” IEEE Signal Processing Letters, vol. 23, no. 7, pp. 919–923,
2016.
[27] F. D. Neeser and J. L. Massey, “Proper complex random processes with applications to information theory,” Information Theory, IEEE Transactions on,
vol. 39, no. 4, pp. 1293–1302, 1993.
[28] J. L. Yen, “On nonuniform sampling of bandwidth-limited signals,” Circuit Theory, IRE Transactions on, vol. 3, no. 4, pp. 251–257, 1956.

14

[29] D. J. Wingham, “The reconstruction of a band-limited function and its fourier transform from a finite number of samples at arbitrary locations by singular
value decomposition,” Signal Processing, IEEE Transactions on, vol. 40, no. 3, pp. 559–570, 1992.
[30] Y. C. Eldar, Sampling Theory: Beyond Bandlimited Systems. Cambridge University Press, 2015.
[31] M. K. Mihcak, I. Kozintsev, K. Ramchandran, and P. Moulin, “Low-complexity image denoising based on statistical modeling of wavelet coefficients,”
IEEE Signal Processing Letters, vol. 6, no. 12, pp. 300–303, 1999.
[32] A. Levin, R. Fergus, F. Durand, and W. T. Freeman, “Image and depth from a conventional camera with a coded aperture,” ACM transactions on graphics
(TOG), vol. 26, no. 3, p. 70, 2007.
[33] A. Levin, Y. Weiss, F. Durand, and W. T. Freeman, “Understanding and evaluating blind deconvolution algorithms,” in Computer Vision and Pattern
Recognition, 2009. CVPR 2009. IEEE Conference on. IEEE, 2009, pp. 1964–1971.
[34] A. Kulesza, B. Taskar et al., “Determinantal point processes for machine learning,” Foundations and Trends R in Machine Learning, vol. 5, no. 2–3,
pp. 123–286, 2012.
[35] G. Friesecke, “Lecture notes in fourier analysis [MA4064],” 2013.
[36] K. Chaloner and I. Verdinelli, “Bayesian experimental design: A review,” Statistical Science, pp. 273–304, 1995.
[37] J. M. Bernardo, “Expected information as expected utility,” The Annals of Statistics, pp. 686–690, 1979.
[38] A. Krause and C. E. Guestrin, “Near-optimal nonmyopic value of information in graphical models,” Proceedings of the Twenty-First Conference on
Uncertainty in Artificial Intelligence, 2012.
[39] U. Feige, “A threshold of ln n for approximating set cover,” Journal of the ACM (JACM), vol. 45, no. 4, pp. 634–652, 1998.
[40] M. Minoux, “Accelerated greedy algorithms for maximizing submodular set functions,” in Optimization Techniques. Springer, 1978, pp. 234–243.
[41] J. Leskovec, A. Krause, C. Guestrin, C. Faloutsos, J. VanBriesen, and N. Glance, “Cost-effective outbreak detection in networks,” in Proceedings of the
13th ACM SIGKDD international conference on Knowledge discovery and data mining. ACM, 2007, pp. 420–429.
[42] J. G. Oxley, Matroid theory. Oxford University Press, USA, 2006, vol. 3.
[43] A. Krause and D. Golovin, “Submodular function maximization,” Tractability: Practical Approaches to Hard Problems, vol. 3, p. 19, 2012.
[44] G. Calinescu, C. Chekuri, M. Pál, and J. Vondrák, “Maximizing a monotone submodular function subject to a matroid constraint,” SIAM Journal on
Computing, vol. 40, no. 6, pp. 1740–1766, 2011.
[45] T. Adali, P. J. Schreier, and L. L. Scharf, “Complex-valued signal processing: The proper way to deal with impropriety,” IEEE Transactions on Signal
Processing, vol. 59, no. 11, pp. 5101–5125, 2011.
[46] G. Shulkind, S. Jegelka, and G. W. Wornell, “Multiple wavelength sensing array design,” in Acoustics, Speech and Signal Processing (ICASSP), 2017
IEEE International Conference on. IEEE, 2017, pp. 3424–3428.

