gsignal (version 0.3-7)

levinson: Durbin-Levinson Recursion

Description

Use the Durbin-Levinson algorithm to compute the coefficients of an autoregressive linear process.

Usage

levinson(acf, p = NROW(acf))

Value

A list containing the following elements:

a

vector or matrix containing (p+1) autoregression coefficients. If x is a matrix, then each row of a corresponds to a column of x. a has p + 1 columns.

e

white noise input variance, returned as a vector. If x is a matrix, then each element of e corresponds to a column of x.

k

Reflection coefficients defining the lattice-filter embodiment of the model returned as vector or a matrix. If x is a matrix, then each column of k corresponds to a column of x. k has p rows.

Arguments

acf

autocorrelation function for lags 0 to p, specified as a vector or matrix. If r is a matrix, the function finds the coefficients for each column of acf and returns them in the rows of a.

p

model order, specified as a positive integer. Default: NROW(acf) - 1.

Author

Paul Kienzle, ,
Peter V. Lanspeary, .
Conversion to R by Geert van Boxtel, .

Details

levinson uses the Durbin-Levinson algorithm to solve: $$toeplitz(acf(1:p)) * x = -acf(2:p+1)$$ The solution c(1, x) is the denominator of an all pole filter approximation to the signal x which generated the autocorrelation function acf.

From ref [2]: Levinson recursion or Levinson–Durbin recursion is a procedure in linear algebra to recursively calculate the solution to an equation involving a Toeplitz matrix. Other methods to process data include Schur decomposition and Cholesky decomposition. In comparison to these, Levinson recursion (particularly split Levinson recursion) tends to be faster computationally, but more sensitive to computational inaccuracies like round-off errors.

References

[1] Steven M. Kay and Stanley Lawrence Marple Jr. (1981). Spectrum analysis – a modern perspective. Proceedings of the IEEE, Vol 69, 1380-1419.
[2] https://en.wikipedia.org/wiki/Levinson_recursion

Examples

## Estimate the coefficients of an autoregressive process given by
## x(n) = 0.1x(n-1) - 0.8x(n-2) - 0.27x(n-3) + w(n).
a <- c(1, 0.1, -0.8, -0.27)
v <- 0.4
w <- sqrt(v) * rnorm(15000)
x <- filter(1, a, w)
xc <- xcorr(x, scale = 'biased')
acf <- xc$R[-which(xc$lags < 0)]
lev <- levinson(acf, length(a) - 1)