[SciPy-User] How to calculate Yulewalk with scipy.optimize.leastsq

[josef.pktd at gmail.com]

On Tue, Jan 17, 2012 at 1:32 PM, Neal Becker <ndbecker2 at gmail.com> wrote:
> klo uo wrote:
>
>> http://www.mathworks.com/help/toolbox/signal/ref/yulewalk.html
>>
>> http://help.scilab.org/docs/5.3.3/en_US/yulewalk.html

Thanks for the references, the paper on the matlab help has some
interesting parts.

I haven't seen anything like this in python, maybe it's work for the
new scikits.signal developers.

The idea of the modified yule walker equations to get the ar
coefficients for an ARMA looked interesting so I coded that part. I
didn't manage to get equation 54 to get the MA coefficients to work.

(Instead I just use a deconvolution to go from autocovariance function
& AR coefficients to MA coefficients. It looks like it's possible to
go from the autocovariance function to recover the ARMA coefficients
with one yule walker, one modified yule walker and one deconvolution.
I just wrote some toy code for this to see if round tripping works. I
still have no idea about filter design. )

------------
import numpy as np
from scipy import linalg

def modified_yule_walker(acov, n_ar, n_ma):
    '''estimate AR part of ARMA from modified Yule-Walker equations

    The autocovariance of an ARMA process for lags greater than n_ma is
    independent of the MA part and can be used to estimate the AR part.
    (If I remember correctly.)

    Parameters
    ----------
    acov : array_like
        autocovariance, needs to have at least n_ar+n_ma terms and start with
        the zero lag, i.e. variance.
    n_ar : int
        length of autoregressive lagpolynomial, including counting the leading
        one.
    n_ma : int
        length of moving-average lagpolynomial, including counting the leading
        one. lagpolynomials with leading (zero lag) value different from one
        have not been tested.

    Returns
    -------
    ar_coefs : ndarray
        the estimated AR coefficients. the autoregressive lag-polynomial is
        [1, -arcoeffs]

    References
    ----------
    Friedlander, B., and B. Porat, 1984, equation (2)

    '''

    p, q = n_ar - 1, n_ma #Friedlander/Porat convention

    r_idx = np.arange(q,q+p+1)  #easier with index than slices
    r2_idx = np.abs(r_idx[0] - 1 - np.arange(p))

    R = linalg.toeplitz(acov[r_idx-1], acov[r2_idx])
    return linalg.lstsq(R, acov[r_idx])[0]
---------

Josef

>>
>>
>> On Mon, Jan 16, 2012 at 11:56 PM, <josef.pktd at gmail.com> wrote:
>>
>>> I have no idea about whether or how it can be used for filter design?
>>>
>>> For estimation of the filter from data, the Yule Walker equations can only
>>> estimate a autoregressive process, not an IIR filter (ARMA process).
>>>
>>>
>
> Many years ago I played around with using YuleWalker as a digital filter design
> technique.  There, I added white noise to the model (adding to the diagonal
> terms of the correlation matrix), and by varying the amount of noise, you could
> vary the filter characteristics.
>
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user