Prony Method

Prony method extracts damped complex exponential functions (or sinusoids) from a given signal by solving a set of linear equations (Lobos et al., 2003; Mitrofanov and Priimenko, 2015; Peter and Plonka, 2013; Prony, 1795). The Prony method allows for estimation of frequencies, amplitudes and phases of a signal (For details see Appendix). Assume we want to solve the problem:

$$x[n] = \sum_{k=1}^{M}A_k e^{(\alpha_k + j\omega_k)(n-1)\Delta t + j\phi_k},$$ (2)

if let $h_k = A_k e^{j\phi_k}, z_k = e^{(\alpha_k + j\omega_k)\Delta t}$, we derive the concise form

$$x[n] = \sum_{k=1}^{M}h_k z_k^{n-1}.$$ (3)

The above M equations can be written in a matrix form:

$$\left[ \begin{array}{cccc} z_1^0 & z_2^0 & \cdots & z_M^0\\ z_1^... ...ht] = \left[ \begin{array}{c} x[1] x[2] \vdots x[M] \end{array} \right].$$ (4)

The above $z_k, k=1,2, \cdots , M$ of equation 4 can be computed by solving a polynomial of the form:

$$\mathbf{P}(z) = \prod_{k=1}^{M}(z-z_k),$$ (5)

equation 5 can also be written in a form:

$$\mathbf{P}(z) = a_0 z^M + a_1 z^{M-1} + \cdots + a_{M-1}z + a_M.$$ (6)

The coefficients $a_k$ of the polynomial can be computed by solving the following equation:

$$\sum_{m=0}^{M}a_mx[n-m] = 0.$$ (7)

We use the method proposed by Toh and Trefethen (1994) to compute the roots $z_k$ of equation 6. If the roots are solved, the $h_k$ can be computed using equation 3. Finally, the components are computed based on the equation below(For details see Appendix):

$$c_k[n] = h_kz_k^{n-1} = h_kz_k^{n-1}, k=1,2,\cdots,M.$$ (8)