Prony's method under an almost sharp multivariate Ingham inequality

TL;DR

This paper extends Prony's method to multivariate signals by establishing a sharper Ingham inequality, improving parameter identifiability conditions and refining the theoretical understanding of multivariate moment-based reconstruction.

Contribution

It introduces a new multivariate Ingham inequality with a logarithmic constant and refines identifiability conditions for multivariate Prony methods.

Findings

Improved dimension-dependent constant from square-root to logarithmic in Ingham inequality.

Enhanced understanding of parameter identifiability in multivariate Prony's method.

Connections to flat extension principle and Hilbert function stagnation.

Abstract

The parameter reconstruction problem in a sum of Dirac measures from its low frequency trigonometric moments is well understood in the univariate case and has a sharp transition of identifiability with respect to the ratio of the separation distance of the parameters and the order of moments. Towards a similar statement in the multivariate case, we present an Ingham inequality which improves the previously best known dimension-dependent constant from square-root growth to a logarithmic one. Secondly, we refine an argument that an Ingham inequality implies identifiability in multivariate Prony methods to the case of commonly used max-degree by a short linear algebra argument, closely related to a flat extension principle and the stagnation of a generalized Hilbert function.

Tables

Table 1. Table 2.1: Explicit constants c d = 2 ​ min p ∈ 2 ​ ℕ ⁡ C p ​ d p subscript 𝑐 𝑑 2 subscript 𝑝 2 ℕ subscript 𝐶 𝑝 𝑝 𝑑 c_{d}=2\min_{p\in 2\mathbb{N}}C_{p}\sqrt[p]{d} , see Equation ( 2.1 ), or its minor improvements via Remark 2.2 and 2.3 .

$d$	$1$	$2$	$3$	$4$	$10$	$16$	$20$	$64$	$100$	$256$
$c_d$	$1.0$	$1.4$	$1.7$	$2.0$	$2.7$	$3.0$	$3.2$	$4.0$	$4.3$	$5.0$

Equations

$$q:=\underset{r\in {\mathbb{Z}}^d,j\ne \ell }{\min }\Vert t_j-t_{\ell }+r{\Vert }_{\infty }.$$

$$f:{\mathbb{Z}}^d\to \mathbb{C},k\mapsto {\int }_{[0,1{)}^d}{\text{e}}^{-2\pi \text{i}kt}\mathrm{d}\tau (t)=\sum _{j=1}^M{\hat{f}}_j{z}_j^k,$$

$$V(K):=\{z\in {\mathbb{C}}^d:p(z)=0\text{ for all }p\in K\}.$$

$$T:=T_n:={\left(f(k-\ell )\right)}_{k,\ell \in {\mathbb{N}}_0^d,|k|,|\ell |\le n}\in {\mathbb{C}}^{(n+1{)}^d\times (n+1{)}^d}.$$

$$A:=A_n:={\left(z_j^k\right)}_{\begin{array}{c}j=1,\dots ,M\\ k\in {\mathbb{N}}_0^d,\Vert k{\Vert }_{\infty }\le n\end{array}}\in {\mathbb{C}}^{M\times (n+1{)}^d},z_j:={\text{e}}^{-2\pi \text{i}{t}_j}\in {\mathbb{C}}^d,$$

$$\phi (x):=\{\begin{array}{ll}{\left(1-{\left(\frac{2x}{q}\right)}^2\right)}^r,& |x|<\frac{q}{2},\\ 0,& \text{otherwise},\end{array}$$

$$\psi =\left({\left(2\pi n\right)}^p-(-1{)}^r\sum _{s=1}^d\frac{{\partial }^p}{\partial x_s^p}\right)\underset{\ell =1}{\overset{d}{⨂}}\phi \ast \phi ,$$

$$\hat{\psi }(v)\{\begin{array}{ll}\ge 0,& \Vert v{\Vert }_p\le n\\ \le 0,& \Vert v{\Vert }_p\ge n,\end{array}$$

$${\phi }^{(r)}(x)=\{\begin{array}{ll}(-1{)}^r{q}^{-r}{4}^{r}r!P_r\left(\frac{2x}{q}\right),& |x|<\frac{q}{2},\\ 0,& \text{otherwise},\end{array}$$

$$\hat{\psi }(v)=\left({\left(2\pi n\right)}^p-\sum _{s=1}^d{\left(2\pi v_s\right)}^p\right)\prod _{\ell =1}^d{\left(\hat{\phi }(v_{\ell })\right)}^2.$$

$$\phi \ast \phi (0)=\frac{q}{2}{\int }_{-1}^1{\left(1-x^2\right)}^p\mathrm{d}x=\frac{q\sqrt{\pi }p!}{2\Gamma (p+\frac{3}{2})}$$

$${\phi }^{(r)}\ast {\phi }^{(r)}(0)={\left(\frac{4^{r}r!}{q^r}\right)}^2\frac{(-1{)}^{r}q}{2}{\int }_{-1}^1(P_r(x){)}^2\mathrm{d}x=\frac{4^p(r!{)}^2(-1{)}^r}{(p+1)q^{p-1}}.$$

$$\psi (0)$$

$$={\left(\phi \ast \phi (0)\right)}^{d-1}\left({\left(2\pi n\right)}^p\frac{q\sqrt{\pi }p!}{2\Gamma (p+\frac{3}{2})}-\frac{d{4}^p(r!{)}^2}{(p+1)q^{p-1}}\right)>0,$$

$$C_p:={\left(\frac{\Gamma \left(\frac{p}{2}+1\right)\Gamma \left(p+\frac{3}{2}\right)}{{\pi }^p\Gamma \left(\frac{p+3}{2}\right)}\right)}^{\frac{1}{p}}\le \frac{{\left(\Gamma \left(p+\frac{3}{2}\right)\right)}^{\frac{1}{p}}}{\pi }\le \frac{2p+3}{\text{e}\pi },$$

$$\phi (x):=\{\begin{array}{ll}1-{\left(\frac{2x}{q}\right)}^2,& |x|<\frac{q}{2},\\ 0,& \text{otherwise},\end{array}$$

$$\phi \ast \phi (0)=\frac{q}{2}{\int }_{-1}^1(1-x^2{)}^2\mathrm{d}x=\frac{8q}{15},{\phi }^{\prime }\ast {\phi }^{\prime }(0)=-\frac{4q}{3}\cdot \frac{4}{q^2},$$

$$\phi (x):=\{\begin{array}{ll}\cos \frac{\pi x}{q},& |x|<\frac{q}{2},\\ 0,& \text{otherwise},\end{array}$$

$${\phi }^{\prime }\ast {\phi }^{\prime }(0)=\phi \ast {\phi }^{\prime \prime }(0)=-\frac{{\pi }^2}{q^2}\phi \ast \phi (0).$$

$$\phi (x):=\{\begin{array}{ll}{\left(1-{\left(\frac{2x}{q}\right)}^2\right)}^2,& |x|<\frac{q}{2},\\ 0,& \text{otherwise},\end{array}$$

$$\phi \ast \phi (0)=\frac{q}{2}{\int }_{-1}^1(1-x^2{)}^4\mathrm{d}x=\frac{128q}{315},{\phi }^{\prime \prime }\ast {\phi }^{\prime \prime }(0)=\frac{64q}{5}\cdot \frac{16}{q^4},$$

$$\phi (x):=\{\begin{array}{ll}1+\cos \frac{2\pi x}{q},& |x|<\frac{q}{2},\\ 0,& \text{otherwise},\end{array}$$

$${\phi }^{\prime \prime }(x)=\frac{4{\pi }^2}{q^2}\left({\chi }_{\left(-\frac{q}{2},\frac{q}{2}\right)}(x)-\phi (x)\right)$$

$${\phi }^{\prime \prime }\ast {\phi }^{\prime \prime }=\frac{16{\pi }^4}{q^4}\left(\phi \ast \phi -2\phi \ast {\chi }_{\left(-\frac{q}{2},\frac{q}{2}\right)}+{\chi }_{\left(-\frac{q}{2},\frac{q}{2}\right)}\ast {\chi }_{\left(-\frac{q}{2},\frac{q}{2}\right)}\right),$$

$$\psi (0)$$

$$=3\cdot 8{\pi }^4{\left(\phi \ast \phi (0)\right)}^{d-1}q\left(n^4-\frac{d}{3}{q}^{-4}\right)>0$$

$$\phi (x):=\{\begin{array}{ll}\frac{\cosh \sigma }{\cosh \sigma -\cos \sigma }\cos (2\sigma x/q)-\frac{\cos \sigma }{\cosh \sigma -\cos \sigma }\cosh (2\sigma x/q),& |x|<\frac{q}{2},\\ 0,& \text{otherwise},\end{array}$$

$$\frac{{\mathrm{d}}^4}{\mathrm{d}{x}^4}\phi (x)=\frac{16{\sigma }^4}{q^4}\phi (x),|x|<\frac{q}{2},\phi (\pm \frac{q}{2})={\phi }^{\prime }(\pm \frac{q}{2})=0.$$

$${\phi }^{\prime \prime }\ast {\phi }^{\prime \prime }(0)=\phi \ast {\phi }^{(4)}(0)=\frac{16{\sigma }^4}{q^4}\phi \ast \phi (0).$$

$$\underset{r\in {\mathbb{Z}}^d,j\ne \ell }{\min }\Vert t_j-t_{\ell }+r{\Vert }_{\infty }>q$$

$$\sum _{\begin{array}{c}k\in {\mathbb{Z}}^d\\ \Vert k{\Vert }_p\le n\end{array}}{\left|\sum _{j=1}^M{\hat{f}}_j{\text{e}}^{2\pi \text{i}k{t}_j}\right|}^2\ge c\sum _{j=1}^M{\left|{\hat{f}}_j\right|}^2,$$

Full text

Key words and phrases : frequency analysis, exponential sum, moment problem, super-resolution, Ingham inequality.

2010 AMS Mathematics Subject Classification : 65T40, 42C15, 30E05, 65F30

1 Introduction

The reconstruction of a Dirac ensemble from its low frequency trigonometric moments or equivalently the parameter reconstruction problem in exponential sums has been studied since [13] as Prony’s method, see e.g. [37, 35, 16, 41, 32, 15] and the survey [30]. Popular different approaches to this reconstruction problem include, without being exhaustive, convex optimization [8, 7, 4, 3], where the problem has been termed ‘super-resolution’ and can be solved via a lifting technique by semidefinite optimization, and maximum likelihood techniques [39, 5, 6, 10, 34, 27, 28, 24]. Multivariate variants of Prony-like methods have been considered e.g. in [18, 38, 1, 33, 31, 14, 21, 36, 20] using generic arguments, projections to univariate problems, or Gröbner basis techniques.

The previously best known condition for an Ingham inequality to hold is $nq>\sqrt{d}/2$, where $n$ denotes the order of the used moments, $q$ the separation distance of the parameters, and $d$ the dimension, respectively. Together with a Gröbner basis construction which asks for a total-degree setting and thus introduces a factor $d$, this results in a deterministic a-priori bound $(n-d-1)q>d^{3/2}$ for the success of Prony’s method, see also [20].

In this paper we refine an appropriate Ingham inequality as well as the algebra techniques to improve this a-priori condition to $(n-1)q>3+2\log d$. After fixing notation, Section 2.2 improves the Ingham inequality [17, 19, 33] by constructing a function with compact support in space domain and zero crossing at an $\ell^{p}$-ball in frequency domain. In Section 2.3, we refine our previous analysis [20] of the vanishing ideal of the parameter set and prove that an interpolation condition, equivalent to some full rank condition, implies that the parameters can be identified by polynomials of certain max-degree. In total, the multivariate Prony method correctly identifies the parameters under a deterministic condition which is sharp up to a logarithmic factor in the space dimension.

2 Main results

In what follows, we establish a deterministic condition for the success of Prony’s method relying on two ingredients - a new Ingham inequality and an analysis of the parameters’ vanishing ideal. Ingham’s inequality generalizes Parseval’s identity in the sense that a norm equality between a space and a frequency representation of a vector is replaced by a certain norm equivalence. With respect to the identifiability in Prony’s method, we are interested in establishing conditions on the parameters such that a non trivial lower bound in this inequality exists. Theorem 2.4 and Corollary 2.5 give such a result under a separation condition on the parameters and improve upon [17, 19, 33] by a smaller space dimension dependent constant. Subsequently, Theorem 2.8 and Corollary 2.10 apply the established Ingham inequality to prove success of Prony’s method by adapting standard arguments from algebraic geometry [40, App. B] to our so-called max-degree setting.

2.1 Preliminaries

Throughout the paper, $d\in\mathbb{N}$ always denotes the dimension and $\mathbb{T}^{d}:=\{z\in\mathbb{C}:|z|=1\}^{d}$ the torus with parameterization $\mathbb{T}^{d}\ni z=\textnormal{e}^{-2\pi\textnormal{i}{t}}$ for a unique $t\in[0,1)^{d}$. Now let $M\in\mathbb{N}$, pairwise distinct parameters $t_{j}\in[0,1)^{d}$, $j=1,\ldots,M$, be given, and define their (wrap-around) separation distance by

[TABLE]

For given coefficients $\hat{f}_{j}\in\mathbb{C}\setminus\{0\}$, the trigonometric moment sequence of the complex Dirac ensemble $\tau:\mathcal{P}([0,1)^{d})\rightarrow\mathbb{C}$, $\tau=\sum_{j=1}^{M}\hat{f}_{j}\delta_{t_{j}}$, is the exponential sum

[TABLE]

with parameters $z_{j}:=\textnormal{e}^{-2\pi\textnormal{i}{t_{j}}}=(\textnormal{e}^{-2\pi\textnormal{i}{t_{j,1}}},\ldots,\textnormal{e}^{-2\pi\textnormal{i}{t_{j,d}}})\in\mathbb{T}^{d}$. A convenient choice for the truncation of this sequence is $\|k\|_{\infty}:=\max\{|k_{1}|,\ldots,|k_{d}|\}\leq n$ and our aim is to reconstruct the parameters $t_{j}\in[0,1)^{d}$ and coefficients $\hat{f}_{j}\in\mathbb{C}\setminus\{0\}$ from the trigonometric moments $f(k)$, $k\in\mathbb{Z}^{d}$, $\|k\|_{\infty}\leq n$.

Prony’s method tries to realize the parameters $z_{j}$ as common zeros of certain polynomials as follows. We identify $\mathbb{C}^{(n+1)^{d}}$ with the space of polynomials of max-degree $n$, i.e., $\Pi_{n}:=\{p:p(z)=\sum_{k\in\mathbb{N}_{0}^{d},\|k\|_{\infty}\leq n}p_{k}z^{k}\}$, and for $K\subset\mathbb{C}^{(n+1)^{d}}$ we let

[TABLE]

Given the above moments, it turns out that an appropriate set of polynomials is the kernel of the multilevel Toeplitz matrix

[TABLE]

Direct computation easily shows the factorization $T=A^{*}DA$ with diagonal matrix $D:=\operatorname*{diag}(\hat{f}_{j})\in\mathbb{C}^{M\times M}$ and Vandermonde matrix

[TABLE]

which is our entry point for the subsequent analysis of the problem.

2.2 Ingham’s inequality

The univariate case is well understood and the popular paper [26] used a so-called Beurling-Selberg majorant and minorant to prove a discrete Ingham inequality (which has been generalized to the unit disk recently [2]). Tensorizing the majorant easily gives an upper bound in a multivariate discrete Ingham inequality but simple attempts to provide a minorant and thus a lower bound failed, see e.g. [25]. A construction of a valid minorant by linear combinations of univariate majorants and minorants can be found e.g. in [9], these seem to become effective only if $nq>C\cdot d$, i.e., show a linear dependence in the space dimension. A similar bound follows from [22] using localized trigonometric polynomials and a packing argument. Moreover note that the Fourier analytic approach to sphere packing problems [11, 42] asks for a similar but not quite the same construction of a function. The classical and currently best known construction [17, 19, 33] uses an eigenfunction of the Laplace operator on the cube and becomes effective if $nq>\frac{1}{2}\sqrt{d}$. Subsequently, we replace the Laplace operator by the sum of higher order derivatives, improve the previous bound for $d>5$, and show a logarithmic growth of the space dependent constant.

Lemma 2.1**.**

Let $d,r\in\mathbb{N}$, $p:=2r$, $n,q>0$, and the functions $\varphi:\mathbb{R}\rightarrow\mathbb{R}$,

[TABLE]

and $\psi:\mathbb{R}^{d}\rightarrow\mathbb{R}$,

[TABLE]

be given, see Figures 2.1 and 2.2 for examples, then

1. i)

the Fourier transform $\hat{\psi}(v):=\int_{\mathbb{R}^{d}}\psi(x)\textnormal{e}^{-2\pi\textnormal{i}{vx}}\mathrm{d}x$ is bounded and obeys

[TABLE] 2. ii)

$\textnormal{supp}\psi=[-q,q]^{d}$, 3. iii)

and $\psi(0)>0$ if $nq>C_{p}\sqrt[p]{d}$ with $C_{p}\leq(2p+3)/(\textnormal{e}\pi)$.

Proof.

We have the weak $r$-th derivative

[TABLE]

where $P_{r}$ denotes the $r$-th Legendre polynomial with normalization $P_{r}(1)=1$. This derivative is of bounded variation and thus, the Fourier transform of the function $\varphi$ obeys $|\hat{\varphi}(v)|\leq C(1+|v|)^{-r-1}$. Now the first assertion follows from

[TABLE]

The second claim easily follows from $\textnormal{supp}\varphi=[-\frac{q}{2},\frac{q}{2}]$. Moreover, we have

[TABLE]

and noting $\varphi^{(r)}$ being odd for $r$ odd, we get

[TABLE]

Finally, this implies

[TABLE]

provided the term in brackets is positive. Using the Legendre duplication formula for $\Gamma(p+2)$, this is equivalent to $nq>C_{p}\sqrt[p]{d}$ with a constant

[TABLE]

where the first inequality follows from monotonicity and the second by Stirling’s approximation. ∎

Remark 2.2**.**

The construction of Lemma 2.1 for $p=2$ reads as

[TABLE]

leading to

[TABLE]

and we get $\psi(0)>0$ for $nq>\frac{\sqrt{5/2}}{\pi}\cdot\sqrt{d}\approx 0.5033\sqrt{d}$.

A minor improvement is possible by choosing the function $\varphi:\mathbb{R}\rightarrow\mathbb{R}$,

[TABLE]

which has a distributional $2$-nd derivative $\varphi^{\prime\prime}=-\pi^{2}q^{-2}\varphi+C\delta_{\pm\frac{q}{2}}$ and thus

[TABLE]

Proceeding as in Lemma 2.1, we get $\psi(0)>0$ for $nq>\frac{1}{2}\cdot\sqrt{d}$, see also [19, 33].

Remark 2.3**.**

The construction of Lemma 2.1 for $p=4$ reads as

[TABLE]

leading to

[TABLE]

and we get $\psi(0)>0$ for $nq>\frac{\sqrt[4]{63/2}}{\pi}\cdot\sqrt[4]{d}\approx 0.7541\cdot\sqrt[4]{d}$.

An alternative choice is the function $\varphi:\mathbb{R}\rightarrow\mathbb{R}$,

[TABLE]

which has a weak $2$-nd derivative of bounded variation

[TABLE]

and the Fourier transforms obey the bounds $|\hat{\varphi}(v)|\leq C(1+|v|)^{-3}$ and $|\hat{\psi}(v)|\leq C^{\prime}$. Moreover, we have

[TABLE]

$\varphi*\varphi(0)=\frac{3}{2}q$, and $\varphi*\chi_{\left(-\frac{q}{2},\frac{q}{2}\right)}(0)=q$, which yields

[TABLE]

if $nq>\sqrt[4]{d/3}\approx 0.7598\cdot\sqrt[4]{d}$.

Finally, a minor improvement is possible as follows. Let $\sigma\approx 2.365$ be the first positive root of $\cos t\sinh t+\cosh t\sin t$ and the function $\varphi:\mathbb{R}\rightarrow\mathbb{R}$,

[TABLE]

be given, then $\varphi$ solves the biharmonic eigenvalue problem

[TABLE]

Globally, we have a $1$-st derivative, a weak $2$-nd derivative, and distributional $3$-rd and $4$-th derivatives. In particular, we have $\varphi^{(4)}=16\sigma^{4}q^{-4}\varphi^{(4)}+C\delta_{\pm\frac{q}{2}}^{\prime}$ and thus

[TABLE]

Proceeding as in Lemma 2.1, we get $\psi(0)>0$ for $nq>\frac{\sigma}{\pi}\cdot\sqrt[4]{d}\approx 0.7528\cdot\sqrt[4]{d}$.

Theorem 2.4**.**

Let $n>0$, $d,p,M\in\mathbb{N}$, $p$ even, $t_{j}\in[0,1)^{d}$, $j=1,\ldots,M$, with

[TABLE]

and $nq>\frac{2p+3}{\textnormal{e}\pi}\sqrt[p]{d}$, then there exists a constant $c>0$, depending only on $d,p,n,q$ such that

[TABLE]

for all choices $\hat{f}_{j}\in\mathbb{C}$, $j=1,\ldots,M$.

Proof.

Let the function $\psi$ be as in Lemma 2.1, then

[TABLE]

where the inequality is implied by Lemma 2.1 i), Poisson’s summation formula yields the first equality, and Lemma 2.1 ii) and the separation of the nodes $t_{j}$ the last equality. Finally, the assertion follows since Lemma 2.1 iii) assures $\psi(0)>0$. ∎

Corollary 2.5**.**

With the notation of Section 2.1, the condition $nq>c_{d}$, see Table 2.1, or $nq>3+2\log d$ implies $\operatorname{rank}A_{n}=M$.

Proof.

First note that

[TABLE]

and thus $A_{n}$ has rank $M$ if and only if the second matrix on the right hand side has this rank. Dropping the last column for odd $n$, this is assured by Theorem 2.4 if $nq>2\frac{2p+3}{\textnormal{e}\pi}\sqrt[p]{d}$, since

[TABLE]

for $\hat{f}\in\mathbb{C}^{M}\setminus\{0\}$. Choosing $p:=2\lceil\log d\rceil$ yields the assertion since

[TABLE]

∎

Remark 2.6**.**

Corollary 2.5 can also be found in [22, Cor. 4.7] and [33, Lem. 3.1] under the stronger conditions $(n+1)q>2d$ and $nq>\sqrt{d}$, respectively. Regarding a trivial lower bound of the constant $c_{d}$, a $d$-fold Cartesian product of equispaced parameters in each dimension consists of $M=q^{-d}$ parameters in total and thus $\operatorname{rank}A\leq(n+1)^{d}<M$ if $(n+1)q<1$. Moreover, a non-trivial lower bound in [29, Section 5.7] indicates that functions as in Lemma 2.1 necessarily imply that the constant $c_{d}$ tends to infinity for $d\rightarrow\infty$. Going beyond the rank of the Vandermonde matrix $A_{n}$ and considering its condition number or equivalently upper and lower bounds in the Ingham inequality, [22, Cor. 4.7] implies a uniform bound

[TABLE]

for each $q$-separated parameter set with $(n+1)q\geq 4d$ and a certain diagonal preconditioner $W$. On the other hand, every weaker condition $(n+1)q\leq c\cdot d^{1-\varepsilon}$, $c>0$, $\varepsilon\in(0,1)$, implies at least a subexponential growing condition number, since [22, Cor. 4.11] leads to

[TABLE]

for equispaced nodes with $(n+1)q\not\in\mathbb{N}$. For a numerical example, let $d=256$ and consider equispaced nodes with $\mathbb{N}\not\ni(n+1)q\approx 4.98\approx c_{256}$ - while Corollary 2.5 assures full rank of $A_{n}$, we have $\kappa(A_{n}A_{n}^{*})=(5/4)^{256}\geq 10^{24}$.

2.3 Prony’s method

As outlined in Section 2.1, Prony’s method tries to realize the unknown parameters $z_{j}$ as common roots of $d$-variate polynomials belonging to the kernel of the multilevel Toeplitz matrix $T_{n}$. In [21], we provided the well-known a-priori condition $n>M$ which however implies the need of $\mathcal{O}(M^{d})$ trigonometric moments for the reconstruction of $M$ parameters. We improved this in [20] to the a-priori condition $(n-d-1)q>d^{3/2}$ using Gröbner basis arguments and a variant of the flat extension principle [12, 23]. This allows for the reconstruction of $M$ well distributed parameters from $\mathcal{O}(M)$ trigonometric moments but the constant deteriorates with larger space dimension. Subsequently, we refine a variant of the flat extension principle to our max-degree setting and give a simple linear algebra proof. Together with the improved Ingham inequality, Prony’s method succeeds under the weakened condition $(n-1)q>3+2\log d$.

Lemma 2.7**.**

With the notation of Section 2.1, if $\operatorname{rank}A_{n}=\operatorname{rank}A_{n+1}$ for some $n\in\mathbb{N}_{0}$, then $\operatorname{rank}A_{\ell}=M$ for all $\ell\geq n$.

Proof.

Let $r:=\operatorname{rank}A_{n}$ and pick $k_{1},\dots,k_{r}\in\mathbb{N}_{0}^{d}$ with $\left\|k_{i}\right\|_{\infty}\leq n$ such that the matrix

[TABLE]

has rank $r$. Subsequently, we let $m\in\mathbb{N}_{0}^{d}$, $m\notin\left\{k_{1},\dots,k_{r}\right\}$, $\|m\|_{\infty}\leq n+2$ and

[TABLE]

If $\left\|m\right\|_{\infty}\leq n+1$, then $B_{m}$ is a submatrix of $A_{n+1}$ and $r=\operatorname{rank}B\leq\operatorname{rank}B_{m}\leq\operatorname{rank}A_{n+1}=r$ by assumption. If $\left\|m\right\|_{\infty}=n+2$, then $m=s+k$, $\|s\|_{\infty}=1$, $\|k\|_{\infty}=n+1$, implies

[TABLE]

and thus $(z_{j}^{m})_{j=1,\ldots,M}$ is linear dependent on the columns of $A_{n+1}$ and thus on the columns of $B$, i.e., $\operatorname{rank}B_{m}=r$. This shows $\operatorname{rank}A_{n+2}=r$ and inductively $\operatorname{rank}A_{\ell}=r$ for all $\ell\geq n$ Finally, the space $\Pi_{M}$ of polynomials of max-degree at most $M$ interpolates any set of $M$ points and thus $\operatorname{rank}A_{M}=M$ and $r=M$. ∎

Theorem 2.8**.**

With the notation of Section 2.1, if $\operatorname{rank}A_{n}=M$, then $V(\ker A_{n+1})=\{z_{1},\ldots,z_{M}\}$.

Proof.

Clearly, “$\supset$” holds. To prove the reverse inclusion, assume that there is a

[TABLE]

For $\ell\in\mathbb{N}_{0}$ consider the matrix

[TABLE]

i.e. $A_{\ell}$ extended by the row $(z_{0}^{k})_{\left\|k\right\|_{\infty}\leq\ell}$. We show that $\ker C_{\ell}=\ker A_{\ell}$ for all $\ell\leq n+1$. The inclusion “$\subset$” is clear. Let $p\in\ker A_{\ell}$. Since $\ell\leq n+1$, we can consider $p$ as an element of $\ker A_{n+1}$, which, by assumption, yields $p(z_{0})=0$. By the definition of $C_{\ell}$ we have $C_{\ell}p=0$ (considering $p$ as an element of $\ker A_{\ell}$ again). Thus $\ker C_{\ell}=\ker A_{\ell}$ for all $\ell\leq n+1$, and therefore $\operatorname{rank}C_{\ell}=\operatorname{rank}A_{\ell}$. In particular, $\operatorname{rank}C_{n}=\operatorname{rank}A_{n}=M$ and $\operatorname{rank}C_{n+1}=\operatorname{rank}A_{n+1}=M$, thus $\operatorname{rank}C_{n}=\operatorname{rank}C_{n+1}$. By Lemma 2.7 we have $\operatorname{rank}C_{n}=\lvert\left\{z_{0},\dots,z_{M}\right\}\rvert=M+1$, a contradiction. ∎

Remark 2.9**.**

Inspection of the proofs readily shows that the results of Lemma 2.7 and Theorem 2.8 hold for arbitrary pairwise different $z_{j}\in\mathbb{C}^{d}$. The Vandermonde matrix $A_{n}$ is the representation matrix of the evaluation homomorphism which maps a polynomial to its values at the prescribed points $z_{j}$. In algebraic geometry and with the max-degree replaced by the total degree, the rank of such maps is called Hilbert function. It is well known that such a function is strictly increasing with $n$ up to stagnation at the value $M$, see e.g. [40, App. B]. The previous Lemma 2.7 shows by elementary arguments that this is true also for the max-degree.

Corollary 2.10**.**

With the notation of Section 2.1, the condition $(n-1)q>c_{d}$, see Table 2.1, or $(n-1)q>3+2\log d$ implies $V(\ker T_{n})=\{z_{1},\ldots,z_{M}\}$.

Proof.

Corollary 2.5 implies $\operatorname{rank}A_{n-1}=M$ and thus $V(\ker A_{n})=\{z_{1},\ldots,z_{M}\}$ follows from Theorem 2.8. Finally note that always $\ker A_{n}\subset\ker T_{n}$ and equality follows from $\operatorname{rank}A_{n}=M$ since $p\in\ker T_{n}$ implies $A_{n}p\in\ker(A_{n}^{*}D)=\operatorname{im}(D^{*}A_{n})^{\bot}=\{0\}$ and thus $p\in\ker A_{n}$. ∎

3 Summary

We have shown that Prony’s method identifies a Dirac ensemble from its first moments provided an associated Vandermonde matrix has full rank. This is indeed fulfilled if the number of given coefficients exceeds a constant divided by the separation distance of the unknown parameters. Explicit bounds for the involved constant grow logarithmically in the space dimension and improve over previous constructions also for small space dimensions.

Acknowledgment. The authors thank the referees for their valuable suggestions and additional pointers to related literature. Moreover, we gratefully acknowledge support by the DFG within the research training group 1916: Combinatorial structures in geometry and by the DAAD-PRIME-program.