Documentation

# levinson

Levinson-Durbin recursion

## Syntax

```a = levinson(r) a = levinson(r,n) [a,e] = levinson(r,n) [a,e,k] = levinson(r,n) ```

## Description

The Levinson-Durbin recursion is an algorithm for finding an all-pole IIR filter with a prescribed deterministic autocorrelation sequence. It has applications in filter design, coding, and spectral estimation. The filter that `levinson` produces is minimum phase.

`a = levinson(r)` finds the coefficients of a `length(r)-1` order autoregressive linear process which has `r` as its autocorrelation sequence. `r` is a real or complex deterministic autocorrelation sequence. If `r` is a matrix, `levinson` finds the coefficients for each column of `r` and returns them in the rows of `a`. `n=length(r)-1` is the default order of the denominator polynomial A(z); that is, `a = [1 a(2) ... a(n+1)]`. The filter coefficients are ordered in descending powers of z–1.

`$H\left(z\right)=\frac{1}{A\left(z\right)}=\frac{1}{1+a\left(2\right){z}^{-1}+\cdots +a\left(n+1\right){z}^{-n}}$`

`a = levinson(r,n)` returns the coefficients for an autoregressive model of order n.

`[a,e] = levinson(r,n)` returns the prediction error, e, of order n.

`[a,e,k] = levinson(r,n)` returns the reflection coefficients `k` as a column vector of length `n`.

### Note

`k` is computed internally while computing the `a` coefficients, so returning `k` simultaneously is more efficient than converting `a` to `k` with `tf2latc`.

## Examples

collapse all

Estimate the coefficients of an autoregressive process given by

`$x\left(n\right)=0.1\phantom{\rule{0.16666666666666666em}{0ex}}x\left(n-1\right)-0.8\phantom{\rule{0.16666666666666666em}{0ex}}x\left(n-2\right)+w\left(n\right).$`

`a = [1 0.1 -0.8];`

Generate a realization of the process by filtering white noise of variance 0.4.

```rng('default') v = 0.4; w = sqrt(v)*randn(15000,1); x = filter(1,a,w);```

Estimate the correlation function. Discard the correlation values at negative lags. Use the Levinson-Durbin recursion to estimate the model coefficients.

```[r,lg] = xcorr(x,'biased'); r(lg<0) = []; ar = levinson(r,numel(a)-1)```
```ar = 1×3 1.0000 0.0957 -0.8026 ```

## Algorithms

`levinson` solves the symmetric Toeplitz system of linear equations

where `r = [`r(1) ... r(n + 1)`]` is the input autocorrelation vector, and r(i)* denotes the complex conjugate of r(i). The input `r` is typically a vector of autocorrelation coefficients where lag 0 is the first element, r(1).

### Note

If `r` is not a valid autocorrelation sequence, the `levinson` function might return `NaN`s even if the solution exists.

The algorithm requires O(n2) flops and is thus much more efficient than the MATLAB® backslash command for large `n`. However, the `levinson` function uses `\` for low orders to provide the fastest possible execution.

## References

 Ljung, Lennart. System Identification: Theory for the User. 2nd Ed. Upper Saddle River, NJ: Prentice Hall, 1999.