# Levinson Method

i have just one question about Levinson, my question is: can we apply this method to determine the coefficients for FIR filter? if yes, please help me ! this is an example:

v = 0.4;
w = sqrt(v)*randn(15000,1);
x = filter(1,a,w);
[r,lg] = xcorr(x,'biased');
r(lg<0) = [];
ar = levinson(r,numel(a)-1) %% this is Gives me good results, but me i search to do this :

v = 0.4;
w = sqrt(v)*randn(15000,1);
x = filter(a,1,w);
[r,lg] = xcorr(x,'biased');
r(lg<0) = [];

ar = levinson(r,numel(a)-1);


The difference between the two codes is in the filter.

• Are you trying to estimate the parameters of an FIR filter given the output of the convolution? What can be assumed on the input signal, is it known, is it white noise or is it unknown (Can be anything)? – Royi Jul 7 '17 at 17:38
• I came across this video a few days ago. I think it helps to clarify what the "levinson" function in Matlab is doing. – Michael_RW Jul 16 '17 at 14:03

## 2 Answers

If you have access to Statistical Digital Signal Processing and Modeling by Monson H. Hayes (Wiley, 1996), see Section 4.7.3 - Moving Average Models for a more detailed description of characterizing moving average signals.

In this section two methodologies are covered, spectral factorization and Durbin's method. Both of these methods have their own innate limitations, but can produce relatively good models for a given moving average process.

The limitation of the spectral factorization method is that it produces a model that might not be unique. Durbin's method, on the other hand, first estimates an all-pole model for the moving average signal and then computes a corresponding set of coefficients that define the zeros for the moving average model. The order of the all-pole model and the order of the all-zero model are adjustable, but it is suggested that the order of the all-pole model be at least four times the order of the all-zero model. In addition to setting these orders, scaling of the computed all-zero coefficients is required.

Revisiting the example from the Matlab "levinson" script-file posted in your original question, consider a moving average process generated using the coefficients,

b=[1, 0.1, -0.8]; a=1;


The frequency response of this system is, Using the Durbin method with an all-zero order of 2 and an all-pole order of 12 without gain scaling, the estimated all-zero coefficients are,

b=[257.9, 15.9, -188.57]; a=1


and the frequency response is, There is good agreement between the original model response and the estimated model response, but the level of gain is not the same. Additionally, the estimated model response has minor deviations in its magnitude and phase response at higher normalized frequency values.

Using the Durbin method with the same orders, but with gain scaling, the all-zero coefficients are,

b=[1.0, 0.056, -0.74]; a=1


and the frequency response is, Characterizing the moving average signal using Durbin's method with gain scaling produces a very good model.

Finally, the frequency response of the model generated using the Matlab "armax" function is, Unlike the frequency response of the model produced using the Durbin method, the "armax" model's responses do not show the high normalized frequency deviations and are almost identical to the frequency response of the original model.

I hope this helps.

It has been about one month since you asked your question, so it is likely you've addressed it and are comfortable with the answer. Regardless, my reply is that the implementation of the Levinson-Durbin algorithm in Matlab (i.e. the "levinson" command) is coded to estimate only an autoregressive process. The comment section of the "levinson" script-file explicitly states it solves the Hermitian Toeplitz system of equations (see "levinson" script-file for a depiction of the set of equations), which are also known as the Yule-Walker AR equations, using the Levinson-Durbin recursion.

Having a look through the time-series model functions in Matlab, it appears there are methods of estimating a moving-average model. For example, the "armax" function can be used.

Using the example from the Matlab "levinson" script-file you originally posted, consider the following Matlab example,

close all; clear; clc;

b=[1, 0.1, -0.8]; a=1;

white_noise = sqrt(1.0)*randn(2^14, 1);
ma_process = filter(b, a, white_noise);

numerator_order = 2;  denominator_order = 0;
ma_model_estimate = armax(ma_process, [denominator_order numerator_order])


The estimated model is,

ma_model_estimate =
Discrete-time MA model:  y(t) = C(z)e(t)
C(z) = 1 + 0.1004 z^-1 - 0.8075 z^-2
`

The estimated model parameters are in good agreement with the values in the "b" vector.

There are a lot of subtleties in the area of time-series modelling and I am not well versed in the field, but I hope this helps.

• Thanks for both answers. I am unsure about OP, but it helped me a lot. Thanks again! – Failed Scientist Jul 13 '17 at 12:34