RLS algorithm?

classic Classic list List threaded Threaded
3 messages Options
Reply | Threaded
Open this post in threaded view
|

RLS algorithm?

Neal Becker
Anyone have code for RLS (recursive least squares)?  I have one version, but it
seems to be rather unstable.

_______________________________________________
SciPy-User mailing list
[hidden email]
http://mail.scipy.org/mailman/listinfo/scipy-user
Reply | Threaded
Open this post in threaded view
|

Re: RLS algorithm?

Charles R Harris


On Tue, Oct 23, 2012 at 11:20 AM, Neal Becker <[hidden email]> wrote:
Anyone have code for RLS (recursive least squares)?  I have one version, but it
seems to be rather unstable.


A bit more information would be helpful. What are you trying to do and how have you implemented it.

Chuck


_______________________________________________
SciPy-User mailing list
[hidden email]
http://mail.scipy.org/mailman/listinfo/scipy-user
Reply | Threaded
Open this post in threaded view
|

Re: RLS algorithm?

Neal Becker
Charles R Harris wrote:

> On Tue, Oct 23, 2012 at 11:20 AM, Neal Becker <[hidden email]> wrote:
>
>> Anyone have code for RLS (recursive least squares)?  I have one version,
>> but it
>> seems to be rather unstable.
>>
>>
> A bit more information would be helpful. What are you trying to do and how
> have you implemented it.
>
> Chuck

Using Haykin 2002 "Adaptive Filter Theory", pp 443 (table 9.1), I came up with
this:

import numpy as np
from itertools import izip
class rls (object):
    def __init__ (self, w, p, _lambda):
        self.w = w
        self._lambda = _lambda
        self.p = p

    def call1 (self, u, d):
        pi_n = np.dot (u.conj(), self.p)
        kappa = self._lambda + np.dot (pi_n, u)
        k_n = pi_n.conj() / kappa
        z = np.dot (self.w.conj(), u)
        e = d - z
        self.w += k_n * e.conj()
        self.p = 1/self._lambda * (self.p - np.outer (k_n, pi_n))

    def __call__ (self, u, d):
        if hasattr (d, '__len__'):
            for eu, ed in izip (u, d):
                self.call1 (eu, ed)
        else:
            self.call1 (u, d)


_______________________________________________
SciPy-User mailing list
[hidden email]
http://mail.scipy.org/mailman/listinfo/scipy-user