Prony analysis (Prony's method) was developed by Gaspard Riche de Prony in 1795. However, practical use of the method awaited the digital computer. Similar to the Fourier transform, Prony's method extracts valuable information from a uniformly sampled signal and builds a series of damped complex exponentials or sinusoids. This allows for the estimation of frequency, amplitude, phase and damping components of a signal.
Let                     f        (        t        )                 be a signal consisting of                     N                 evenly spaced samples. Prony's method fits a function
                                                        f              ^                                      (        t        )        =                  ∑                      i            =            1                                M                                    A                      i                                    e                                    σ                              i                                      t                          cos                (        2        π                  f                      i                          t        +                  ϕ                      i                          )                to the observed                     f        (        t        )                . After some manipulation utilizing Euler's formula, the following result is obtained. This allows more direct computation of terms.
                                                                                                                                    f                      ^                                                                      (                t                )                                                            =                                  ∑                                      i                    =                    1                                                        M                                                                    A                                      i                                                                    e                                                            σ                                              i                                                              t                                                  cos                                (                2                π                                  f                                      i                                                  t                +                                  ϕ                                      i                                                  )                                                                                                  =                                  ∑                                      i                    =                    1                                                        M                                                                                        1                    2                                                                    A                                      i                                                                    e                                      ±                    j                                          ϕ                                              i                                                                                                              e                                                            λ                                              i                                                              t                                                                                              where:
                              λ                      i                          =                  σ                      i                          ±        j                  ω                      i                                   are the eigenvalues of the system,                              σ                      i                                   are the damping components,                              ϕ                      i                                   are the phase components,                              f                      i                                   are the frequency components,                              A                      i                                   are the amplitude components of the series, and                    j                 is the imaginary unit (                              j                      2                          =        −        1                ).Prony's method is essentially a decomposition of a signal with                     M                 complex exponentials via the following process:
Regularly sample                                                         f              ^                                      (        t        )                 so that the                     n                -th of                     N                 samples may be written as
                              F                      n                          =                                            f              ^                                      (                  Δ                      t                          n        )        =                  ∑                      m            =            1                                M                                                B                                m                                    e                                    λ                              m                                                    Δ                              t                                      n                          ,                n        =        0        ,        …        ,        N        −        1.                If                                                         f              ^                                      (        t        )                 happens to consist of damped sinusoids, then there will be pairs of complex exponentials such that
                                                                                                              B                                                        a                                                                                              =                                                      1                    2                                                                    A                                      i                                                                    e                                                            ϕ                                              i                                                              j                                                  ,                                                                                                          B                                                        b                                                                                              =                                                      1                    2                                                                    A                                      i                                                                    e                                      −                                          ϕ                                              i                                                              j                                                  ,                                                                                      λ                                      a                                                                                              =                                  σ                                      i                                                  +                j                                  ω                                      i                                                  ,                                                                                      λ                                      b                                                                                              =                                  σ                                      i                                                  −                j                                  ω                                      i                                                  ,                                                            where
                                                                                                              B                                                        a                                                                    e                                                            λ                                              a                                                              t                                                  +                                                      B                                                        b                                                                    e                                                            λ                                              b                                                              t                                                                                              =                                                      1                    2                                                                    A                                      i                                                                    e                                                            ϕ                                              i                                                              j                                                                    e                                      (                                          σ                                              i                                                              +                    j                                          ω                                              i                                                              )                    t                                                  +                                                      1                    2                                                                    A                                      i                                                                    e                                      −                                          ϕ                                              i                                                              j                                                                    e                                      (                                          σ                                              i                                                              −                    j                                          ω                                              i                                                              )                    t                                                                                                                                    =                                  A                                      i                                                                    e                                                            σ                                              i                                                              t                                                  cos                                (                                  ω                                      i                                                  t                +                                  ϕ                                      i                                                  )                .                                                            Because the summation of complex exponentials is the homogeneous solution to a linear difference equation, the following difference equation will exist:
                                                        f              ^                                      (                  Δ                      t                          n        )        =        −                  ∑                      m            =            1                                M                                                              f              ^                                      [                  Δ                      t                          (        n        −        m        )        ]                  P                      m                          ,                n        =        M        ,        …        ,        N        −        1.                The key to Prony's Method is that the coefficients in the difference equation are related to the following polynomial:
                              z                      M                          +                  P                      1                                    z                      M            −            1                          +        ⋯        +                  P                      M                          =                  ∏                      m            =            1                                M                                    (          z          −                      e                                          λ                                  m                                                              )                .                These facts lead to the following three steps to Prony's Method:
1) Construct and solve the matrix equation for the                               P                      m                                   values:
                                          [                                                                                F                                          M                                                                                                                    ⋮                                                                                                  F                                          N                      −                      1                                                                                            ]                          =                              [                                                                                F                                          M                      −                      1                                                                                        …                                                                      F                                          0                                                                                                                    ⋮                                                  ⋱                                                  ⋮                                                                                                  F                                          N                      −                      2                                                                                        …                                                                      F                                          N                      −                      M                      −                      1                                                                                            ]                                                [                                                                                P                                          1                                                                                                                    ⋮                                                                                                  P                                          M                                                                                            ]                          .                Note that if                     N        ≠        2        M                , a generalized matrix inverse may be needed to find the values                               P                      m                                  .
2) After finding the                               P                      m                                   values find the roots (numerically if necessary) of the polynomial
                              z                      M                          +                  P                      1                                    z                      M            −            1                          +        ⋯        +                  P                      M                          ,                The                     m                -th root of this polynomial will be equal to                               e                                    λ                              m                                                            .
3) With the                               e                                    λ                              m                                                             values the                               F                      n                                   values are part of a system of linear equations that may be used to solve for the                                           B                                m                                   values:
                                          [                                                                                F                                                                  k                                                  1                                                                                                                                                                  ⋮                                                                                                  F                                                                  k                                                  M                                                                                                                                          ]                          =                              [                                                            (                                      e                                                                  λ                                                  1                                                                                                                          )                                                                  k                                                  1                                                                                                                                      …                                                  (                                      e                                                                  λ                                                  M                                                                                                                          )                                                                  k                                                  1                                                                                                                                                                  ⋮                                                  ⋱                                                  ⋮                                                                              (                                      e                                                                  λ                                                  1                                                                                                                          )                                                                  k                                                  M                                                                                                                                      …                                                  (                                      e                                                                  λ                                                  M                                                                                                                          )                                                                  k                                                  M                                                                                                                                          ]                                                [                                                                                                      B                                                              1                                                                                                                    ⋮                                                                                                                        B                                                              M                                                                                            ]                          ,                where                     M                 unique values                               k                      i                                   are used. It is possible to use a generalized matrix inverse if more than                     M                 samples are used.
Note that solving for                               λ                      m                                   will yield ambiguities, since only                               e                                    λ                              m                                                             was solved for, and                               e                                    λ                              m                                                    =                  e                                    λ                              m                                                  +                        q            2            π            j                                   for an integer                     q                . This leads to the same Nyquist sampling criteria that discrete Fourier transforms are subject to:
                              |          I          m          (                      λ                          m                                )          |                =                  |                      ω                          m                                |                <                              1                          2                              Δ                                  t                                                                    .