Chapter 3 - Extrapolation Theory

 

  1. Introduction

    2. The general goal of extrapolation is to estimate a signal, f(t), from a limited known segment, g(t). If there is a transform that can describe g(t) from f(t), g(t) = T{f(t)}, then its inverse could estimate f(t) from g(t), f(t) = T-1{g(t)}. The known segment is a time-limited portion of the real signal,

    g(t) = w(t) f(t)

    ( 3.1)

    where, w(t) = 1, between times [-T, T]

    w(t) = 0, elsewhere

    If the signal is band-limited between frequencies -s and s , B = [-s ,s ], then

    g(t) = w(t) B(t) f(t).

    ( 3.2)

    where, B(t) = 1, between frequencies [-s ,s ]

    B(t) = 0, elsewhere

    If a signal is continuous, deterministic, and band-limited, it is analytic, and therefore can be represented as a Taylor series expansion [8]. In theory, this expansion proves that since all derivatives can be determined within the known segment, the signal can be determined outside of this interval. Since the derivatives are extremely noise-sensitive, the smallest noise would strongly corrupt the extrapolation, and, thus, it is impractical.

     

  2. Minimum Energy Extrapolation

    3.  

    1. Gerchberg – Super-Resolution through Error Energy Reduction

      2. Gerchberg [9] showed that resolution of a data object could be greatly improved through error energy reduction. He explained that a finite object function has an analytic spectrum. Since the sum of analytic functions is itself analytic, a finite segment of an analytic function could determine the whole function in a unique way. He concluded that if a segmented (band-limited) spectrum is continued, then increased resolution, or as he termed super-resolution, of the time-domain object could be achieved. He proposed an iterative algorithm where the energy outside the known segment is zeroed at each iteration.

      wpe8.jpg (17304 bytes)

      Figure 3.1 Block diagram of computation procedure for Gerchberg Super-Resolution

      Since only error energy existed beyond object’s known limit, then zeroing this energy necessarily had to reduce the error energy. The extent of the true object is required to be known a priori in order for the algorithm to be successful. If the extent is overestimated, then there is less reduction in the error energy. If the extent is underestimated, the algorithm fails.

    2. Papoulis – Minimum Energy Extrapolation

Based on Gerchburg’s findings, Papoulis [10] developed an algorithm for extrapolating band-limited signals. His iterative algorithm consists of the following:

It can be written in equation form as

( 3.3)

The first iteration of the extrapolation method is illustrated in Figure 3.2. The figure displays the full time-domain signal and its spectrum. Below the full signal, the known segment of the full signal is shown. The first iteration shows the known time-domain segment transferred to the frequency domain. This frequency-domain representation is band-limited and then transferred back into the time domain. The known time-domain segment is then replaced in the new time-domain signal. This process is repeated until it converges to one signal.

wpe9.jpg (25729 bytes)

Figure 3.2 Diagram of Papoulis extrapolation method

      1. Sabri and Steenaart - Extrapolation Matrix

        2. The Papoulis-Gerchberg algorithm has a slow convergence and requires a large number of iterations. The iterative procedure is replaced with a single operator using the Sabri-Steenaart [11] extrapolation matrix. Since the first iteration could be thought of in the following way,

        {f1} = {f0} + {f0}{B}{X}

        ( 3.4)

        where {f0} = {g},

        {X} is a diagonal matrix zeroing out the known time interval, and

        {B} = {DFT}-1{BL}{DFT},

        where {DFT} is a discrete Fourier Transform matrix and

        {BL} is a matrix of weights: band-limiter.

        Let {H} = {B}{X}. The second iteration would therefore be:

        {f2} = {f0} + {f1}{H} = {f0} + {f0}{H} + {f0}{H}2.

        ( 3.5)

        This follows that

        ( 3.6)

        where {En} is the extrapolation matrix.

        A simple example is now shown with a four sample input vector,, and a bandpass filter, BL = diag([0 1 1 0]). The first iteration is computed as follows:

        ( 3.7 )

        This example converges at iteration 15 where. The extrapolation matrix for this example accelerates the iteration process and obtains the same result.

        ( 3.8)

        The extrapolation matrix can be further written,

        ( 3.9)

        where {I} is the identity matrix.

        As n approaches infinity, the extrapolation matrix becomes

        {E¥ } = ({I} - {H})-1.

        ( 3.10)

        The extrapolation matrix exists if and only if ({I} - {H}) is nonsingular. In practical situations, the {H} matrix needs to be finite (although usually large). Therefore, En does exist, but as n approaches infinity, E¥ becomes ill-conditioned, and therefore unsolvable in a practical sense.

      2. Jain-Ranganath - Minimum Norm Least Squares Extrapolation

        3. Due to the fact that discrete signals are not analytic signals, discrete extrapolation algorithms do not always converge to the original signal. Jain and Ranganath [8] showed that continuous extrapolation algorithms could be converted to discrete algorithms using the minimum norm least squares (MNLS) solution.

        Let A be an arbitrary matrix. If Af = g, then f + = A+g. The MNLS solution of f + is defined as

        ( 3.11)

        The matrix A+ is called the generalized inverse and is written

        ( 3.12)

        Jain and Ranganath showed that since the infinite extrapolation matrix does not exist in practical situations, the correct extrapolation matrix is a pseudo-inverse matrix. Let A = wBT; band-limited and time-limited matrix. When this is substituted into equation ( 3.12), the pseudo-inverse extrapolation matrix is obtained.

        E+ = B wT (w BTB w T) –1

        ( 3.13)

        Where B is the band-limiting operator and

        w is the time-limiting operator.

        Since band-limiting a function twice is the same as band-limiting it once, B is idempotent, BTB = B. Therefore, equation ( 3.13) can be simplified.

        E+ = B wT (w B w T) –1

        ( 3.14)

        This matrix, equation ( 3.14), is the MNLS solution to the extrapolation problem.

      3. Other Minimum Energy Extrapolation Methods

An alternative method of looking at the extrapolation problem was developed by Cadzow [12]. His algorithm claimed not to have the truncation errors found in the Gerchberg-Papoulis method. He proposed a method of breaking up the extrapolated signal into two parts, the known segment and the unknown segment.

f = Bg + (I – Bw) f

( 3.15)

Therefore, the iterative solution is

fn = Bg + (I – Bw) fn-1

( 3.16)

Chamaz and Xu [13] improved the Papoulis-Gerchberg algorithm by speeding up its convergence. They proposed a steepest decent method to be included in the Papoulis -Gerchberg algorithm. A constant is multiplied by the signal outside of the known region at each iteration. This constant, An, is chosen to minimize the energy of the error at each iteration.

( 3.17)

where An is defined as

( 3.18)

with

 

    1. Extrapolation using Time Bandwidth Dimensions

      2.  

      1. Kolba and Parks – BLTC Estimation

        2. Minimum energy extrapolation methods are effective for a limited region beyond the known data. This effective region is TL = TM/2B where TM is the length of the known segment and B is the bandwidth. This region is only true when B £ 0.5, normalized bandwidth. In response to this, Kolba and Parks [14] developed an algorithm that considered time bandwidth dimensions instead of minimum energy criteria. Their method, entitled Band-Limited Time-Concentrated (BLTC) estimation, improves the extrapolation in the effective region by exploiting both the band-limiting and time-concentrating characteristic of the signal. In some instances the effective region increases when using this algorithm.

        The algorithm consists of determining the orthonormal basis functions that could construct the inner product of the signal. These basis functions are derived from a special family of band-limited sequences known as Discrete Prolate Spheroidal (DPS) sequences. The DPS eignenfunctions, y i(t), are defined

        ( 3.19)

        The eigenfunctions are orthogonal in the intervals (-¥ ,¥ ) and (-T,T).

        Once the basis functions are determined, they are then recombined to form the extrapolation.

        ( 3.20)

        with

        y = DPS basis functions

        L = Length of Time bandwidth dimension

         

      2. Dharanipragada and Arun Algorithm

Dharanipragada and Arun [15] improved the Kolba and Parks algorithm by including reduced dimensions and known energy concentration. Instead of using a dimension equal to the size of the observation, this algorithm uses the smallest time dimension for which the worst-case error is within a tolerance. Orthonormal eigenvectors are formed based on the known passbands and concentration regions. A set of Discrete Prolate Spheroidal (DPS) basis functions is then formed using the eigenvectors and the reduced time-bandwidth dimensions. These basis functions are then recombined to form the extrapolation.

The algorithm can be described in the following steps:

    1. Precompute (off-line) the time-bandwidth dimension M as the largest integer for which

      where
    2. Solve Aa = g for a . A is a matrix comprised of the first M DPS basis functions.
    3. Construct the extrapolation:

 

    1. Energy Concentration Extrapolation

      2.  

      1. Potter and Arun – Energy Concentration Extrapolation

        2. Since minimum energy extrapolation methods assume that the known samples are taken at places of maximum energy, Potter and Arun [16] developed an extrapolation algorithm that includes weights based on energy distribution. Their method is based on optimizing the Hilbert space of the function. The extrapolation sequence, f , is the sequence that minimizes the following Hilbert space inner product.

        ( 3.21)

        with

        Q = Matrix of weights

        E = Operator selecting the indexes of high energy concentration

        F = Spectrum operator

        I = Identity matrix

        The generalized solution to the extrapolation problem is

        ( 3.22)

        where u and v are the solutions to the following linear equations.

        ( 3.23)

        with

        w = Time operator

        X = (I + F)-1B

        B = Band-limiter

         

      2. Cabrera and Parks - Minimum Weighted Norm Extrapolation

The above extrapolation methods are restricted to narrow band signals, B £ 0.5 using normalized bandwidth. Cabrera and Parks [1] introduced an improved algorithm using minimum weighted norm extrapolation. This method uses a priori spectral information to improve the extrapolation and does not constrict the bandwidth.

Instead of using a frequency band-limiter, B, as a weight, this method uses the power spectrum of the known signal as the main criteria for the weight, Q. The frequency weighted norm can be written as

( 3.24)

This is equivalently translated into the discrete weighted inner product:

( 3.25)

In order to minimize the inner product, the power density spectrum is used as the weight, Q(k). The inverse DFT of the weight is the circular autocorrelation.

Revisiting equation ( 3.14) and replacing Q for B, the minimum norm solution of the extrapolation problem is

E+ = Qw* (w Q w*)-1

Where Q is the weight and w is the time limiter.

( 3.26)

Let G = wQw*, and let a=G-1g(n). So,

F(n) = E+g(n) = Qw* (w Q w*)-1 g(n)

( 3.27)

F(n) = Qw* G-1 g(n) = Qw* a

( 3.28)

Since the power density spectrum is being used as the weight, the circular autocorrelation function can be used as the time-domain equivalent. The time-limited circular autocorrelation for this algorithm is defined as follows,

f (n) = h*(-n) * h(n)

where h(n) is a zero-padded version of g(n) to the length of the extrapolation.

( 3.29)

The minimum weighted norm extrapolation therefore becomes

where L is the length of extrapolation.

( 3.30)

By using regularization, a faster convergence can be obtained due to numerical stability. Regularization offers a trade-off between numerical stability and consistency to valid data. The regularization parameter is usually chosen to be small in order to increase the numerical stability without diverging too far from the valid data. The regularized solution of the weight is

G = (wQw* + r I)

Where I is the identity matrix.

( 3.31)

The steps involved in this extrapolation method are as follows.

    1. Define h(n) as g(n) zero-padded to the length of the extrapolation.
    2. Compute the circular autocorrelation:

      3. f (n) = h*(-n) * h(n)
    3. Compute the G matrix:

      4. Gr,c = f ((r-c)mod L)

      where L is the length of f (n).
    4. Compute vector of constants, a:

      5. a = (G + r I)-1 g(n).
    5. Determine the extrapolation:

.

 
    1. Discussion
      This chapter examined several extrapolation algorithms. Although in theory a Taylor series expansion can be used to extrapolate a known segment, it is a not practical solution. Minimum energy extrapolation provides a practical solution to the extrapolation problem by minimizing the error energy. While an analytic signal extrapolation converges to the original signal, a discrete signal extrapolation does not converge to the original signal. Discrete signal extrapolations must use methods such as the MNLS solution to determine the best convergence. Extrapolations can be improved by exploiting both the band-limiting and time-concentrating characteristics of the signal.

      Minimum energy extrapolation methods assume that the known data is taken at positions of maximum energy and that the signal is a narrow band signal. Energy concentration extrapolation methods weigh the signal based on its energy distribution thus eliminating the assumption of where the known data is taken. Minimum weighted norm extrapolation uses the power density spectrum to weigh the signal and does not restrict the input signal to narrow band signals. Since the known segment of an audio signal is not assumed to be taken from positions of maximum energy and is not a narrow band signal, minimum weighted norm extrapolation is chosen as the extrapolation method in the proposed algorithm.

 

Continue

Back to Table of Contents