Signal: Least-squares deconvolution. Autoregressive model, impulse responses, yule-walker, etc.

Hi,
I have written some efficient code for learning impulse responses using
least-squares minimization.  First time trying to contribute. Wanted to
know if it seems interesting and if someone wants to review.

*Problem*
I.e. given input x and output y (both N samples long), we wish to learn an
impulse response with M < N samples that minimizes the least squares error
in y.

Given
yp[n] = x[n]h[0]+x[n-1]h[1]+x[n-2]h[2]+... + x[n-M]h[M],

minimize sum( (y-yp)**2)  wrt. h[0],h[1],h[2],... ,

Circular convolution is assumed.

*Applications*
- Learning impulse responses.
- Learn parameters of autoregressive models.

If you use x and y so that x[n] = y[n-1], t*he algorithm will be equal to
learning the parameters of an autoregressive model using the Yule-Walker
equation*.

*Performance*
It will be fast for long time-series, since the A matrix and the b vector
in the normal equation are calculated using FFTs, instead of the regular
dot(X.T,X) and dot(X.T,b).

Also, the toeplitz structure of the normal equation matrix is exploited for
performance. This allows for learning longer impulse responses. And for
using less memory to store the normal equation matrix.

*Notes*

Circular convolution is often not what we want, but can be overcome by
zero-padding.

Alternatively, one wants to avoid that the first M samples of the output to
contributes to the error at all.This can be achieved by first using the
algorithm to learn the impulse response and then fine-tuning it using
scipy.optimize

Can add functionality for this fine-tuning + zero-padding. Should also add
a few alternatives for linear solvers.

Core functionality is only 10 lines of self-written code and uses only
numpy and scipy.linalg.solve_toeplitz.

_______________________________________________
SciPy-Dev mailing list -- [email protected]
To unsubscribe send an email to [email protected]
https://mail.python.org/mailman3/lists/scipy-dev.python.org/
Member address: [email protected]