Recursive least squares filter

From Wikipedia, the free encyclopedia

Recursive least squares (RLS) algorithm is used in adaptive filters to find the filter coefficients that relate to producing the recursively least squares of the error signal (difference between the desired and the actual signal)

1 Discussion
2 Recursive Algorithm
3 RLS algorithm summary
4 Example: Plant Identification
5 References
6 See also

[edit] Discussion

The idea behind RLS filters is to minimize a weighted least squares error function. To stay within the adaptive filter terminology, the weighted least squares error function is the cost function defined as

$C(\mathbf{w}_{n})=\sum_{i=0}^{n}\lambda^{n-i}|e(i)|^{2}=\sum_{i=0}^{n}\lambda^{n-i}e(i)\,e^{*}(i)$

where $0<\lambda\le 1$ is an exponential weighting factor which effectively limits the number of input samples based on which the cost function is minimized. The error signal $e (n)$ is defined in the block diagram section of the adaptive filter page. The cost function is minimized by taking the partial derivatives for all entries $k$ of the coefficient vector $\mathbf{w}_{n}$ and setting the results to zero

$\frac{\partial C(\mathbf{w}_{n})}{\partial w^{*}_{n}(k)}=\sum_{i=0}^{n}\lambda^{n-i}e(i)\,\frac{\partial e^{*}(i))}{\partial w^{*}_{n}(k)}=\sum_{i=0}^{n}\lambda^{n-i}e(i)\,x^{*}(i-k)=0$

Next, replace $e (n)$ with the definition of the error signal

$\sum_{i=0}^{n}\lambda^{n-i}\left[d(i)-\sum_{l=0}^{p}w_{n}(l)x(i-l)\right]x^{*}(i-k)= 0$

Rearranging the equation yields

$\sum_{l=0}^{p}w_{n}(l)\left[\sum_{i=0}^{n}\lambda^{n-i}\,x(i-l)x^{*}(i-k)\right]= \sum_{i=0}^{n}\lambda^{n-i}d(i)x^{*}(i-k)$

This form can be expressed in terms of matrices

$\mathbf{R}_{x}(n)\,\mathbf{w}_{n}=\mathbf{r}_{dx}(n)$

where $\mathbf{R}_{x}(n)$ is the weightend autocorrelation matrix for $x (n)$ and $\mathbf{r}_{dx}(n)$ is the cross-correlation between $d (n)$ and $x (n)$ . Based on this expression we find the coefficients which minimize the cost function as

$\mathbf{w}_{n}=\mathbf{R}_{x}^{-1}(n)\,\mathbf{r}_{dx}(n)$

This is the main result of the discussion.

[edit] Recursive Algorithm

The discussion resulted in a single equation to determine a coefficient vector which minimizes the cost function. In this section we want to derive a recursive solution of the form

$\mathbf{w}_{n}=\mathbf{w}_{n-1}+\Delta\mathbf{w}_{n-1}$

where $\Delta\mathbf{w}_{n-1}$ is a correction factor at time $n - 1$ . We start the derivation of the recursive algorithm by expressing the cross correlation $\mathbf{r}_{dx}(n)$ in terms of $\mathbf{r}_{dx}(n-1)$

$\mathbf{r}_{dx}(n)$	$=\sum_{i=0}^{n}\lambda^{n-i}d(i)\mathbf{x}^{*}(i)$
	$=\sum_{i=0}^{n-1}\lambda^{n-i}d(i)\mathbf{x}^{}(i)+\lambda^{0}d(n)\mathbf{x}^{}(n)$
	$=\lambda\mathbf{r}_{dx}(n-1)+d(n)\mathbf{x}^{*}(n)$

where $\mathbf{x}(i)$ is the $p + 1$ dimensional data vector

$\mathbf{x}(i)=[x(i), x(i-1), ..., x(i-p) ]^{T}$

Similarly we express $\mathbf{R}_{x}(n)$ in terms of $\mathbf{R}_{x}(n-1)$ by

$\mathbf{R}_{x}(n)$	$=\sum_{i=0}^{n}\lambda^{n-i}\mathbf{x}^{*}(i)\mathbf{x}^{T}(i)$
	$=\lambda\mathbf{R}_{x}(n-1)+\mathbf{x}^{*}(n)\mathbf{x}^{T}(n)$

In order to generate the coefficient vector we are interested in the inverse of the deterministic autocorrelation matrix. For that task the Woodbury matrix identity comes in handy. With

$A$	$=\lambda\mathbf{R}_{x}(n-1)$ is $(p + 1)$ -by- $(p + 1)$
$U$	$=\mathbf{x}^{*}(n)$ is $(p + 1)$ -by- $1$
$V$	$=\mathbf{x}^{T}(n)$ is $1$ -by- $(p + 1)$
$C$	$= I = 1$ is the $1$ -by- $1$ identity matrix

The Woodbury matrix identity follows

$\mathbf{R}_{x}^{-1}(n)$	$=$	$\left[\lambda\mathbf{R}_{x}(n-1)+\mathbf{x}^{*}(n)\mathbf{x}^{T}(n)\right]^{-1}$
	$=$	$\lambda^{-1}\mathbf{R}_{x}^{1}(n-1)$
		$-\lambda^{-1}\mathbf{R}_{x}^{-1}(n-1)\mathbf{x}^{*}(n)$
		$\left\{1+\mathbf{x}^{T}(n)\lambda^{-1}\mathbf{R}_{x}^{-1}(n-1)\mathbf{x}^{*}(n)\right\}^{-1} \mathbf{x}^{T}(n)\lambda^{-1}\mathbf{R}_{x}^{-1}(n-1)$

To come in line with the standard literature, we define

$\mathbf{P}(n)$	$=\mathbf{R}_{x}^{-1}(n)$
	$=\lambda^{-1}\mathbf{P}(n-1)-\mathbf{g}(n)\mathbf{x}^{T}(n)\lambda^{-1}\mathbf{P}(n-1)$

where

$\mathbf{g}(n)$	$=\lambda^{-1}\mathbf{P}(n-1)\mathbf{x}^{}(n)\left\{1+\mathbf{x}^{T}(n)\lambda^{-1}\mathbf{P}(n-1)\mathbf{x}^{}(n)\right\}^{-1}$
	$=\mathbf{P}(n-1)\mathbf{x}^{}(n)\left\{\lambda+\mathbf{x}^{T}(n)\mathbf{P}(n-1)\mathbf{x}^{}(n)\right\}^{-1}$

Before we move on, it is necessary to bring $\mathbf{g}(n)$ into another form

$\mathbf{g}(n)\left\{1+\mathbf{x}^{T}(n)\lambda^{-1}\mathbf{P}(n-1)\mathbf{x}^{*}(n)\right\}$	$=\lambda^{-1}\mathbf{P}(n-1)\mathbf{x}^{*}(n)$
$\mathbf{g}(n)+\mathbf{g}(n)\mathbf{x}^{T}(n)\lambda^{-1}\mathbf{P}(n-1)\mathbf{x}^{*}(n)$	$=\lambda^{-1}\mathbf{P}(n-1)\mathbf{x}^{*}(n)$

Subtracting the second term on the left side yields

$\mathbf{g}(n)$	$=\lambda^{-1}\mathbf{P}(n-1)\mathbf{x}^{}(n)-\mathbf{g}(n)\mathbf{x}^{T}(n)\lambda^{-1}\mathbf{P}(n-1)\mathbf{x}^{}(n)$
	$=\lambda^{-1}\left[\mathbf{P}(n-1)-\mathbf{g}(n)\mathbf{x}^{T}(n)\mathbf{P}(n-1)\right]\mathbf{x}^{*}(n)$

With the recursive definition of $\mathbf{P}(n)$ the desired form follows

$\mathbf{g}(n)=\mathbf{P}(n)\mathbf{x}^{*}(n)$

Now we are ready to complete the recursion. As discussed

$\mathbf{w}_{n}$	$=\mathbf{P}(n)\,\mathbf{r}_{dx}(n)$
	$=\lambda\mathbf{P}(n)\,\mathbf{r}_{dx}(n-1)+d(n)\mathbf{P}(n)\,\mathbf{x}^{*}(n)$

The second step follows from the recursive definition of $\mathbf{r}_{dx}(n )$ . Next we incorporate the recursive definition of $\mathbf{P}(n)$ together with the alternate form of $\mathbf{g}(n)$ and get

$\mathbf{w}_{n}$	$=\lambda\left[\lambda^{-1}\mathbf{P}(n-1)-\mathbf{g}(n)\mathbf{x}^{T}(n)\lambda^{-1}\mathbf{P}(n-1)\right]\mathbf{r}_{dx}(n-1)+d(n)\mathbf{g}(n)$
	$=\mathbf{P}(n-1)\mathbf{r}_{dx}(n-1)-\mathbf{g}(n)\mathbf{x}^{T}(n)\mathbf{P}(n-1)\mathbf{r}_{dx}(n-1)+d(n)\mathbf{g}(n)$
	$=\mathbf{P}(n-1)\mathbf{r}_{dx}(n-1)+\mathbf{g}(n)\left[d(n)-\mathbf{x}^{T}(n)\mathbf{P}(n-1)\mathbf{r}_{dx}(n-1)\right]$

With $\mathbf{w}_{n-1}=\mathbf{P}(n-1)\mathbf{r}_{dx}(n-1)$ we arrive at the update equation

$\mathbf{w}_{n}$	$=\mathbf{w}_{n-1}+\mathbf{g}(n)\left[d(n)-\mathbf{x}^{T}(n)\mathbf{w}_{n-1}\right]$
	$=\mathbf{w}_{n-1}+\mathbf{g}(n)e(n)$

where

$e(n)=d(n)-\mathbf{x}^{T}(n)\mathbf{w}_{n-1}$

That means we found the correction factor

$\Delta\mathbf{w}_{n-1}=\mathbf{g}(n)e(n)$

[edit] RLS algorithm summary

The RLS algorithm for a $p$ th order algorithm can be summarized as

Parameters:	$p =$ filter order
	$λ =$ forgetting factor
	$δ =$ value to initialize $\mathbf{P}(0)$
Initialisation:	$\mathbf{w}_{n}=0$
	$\mathbf{P}(0)=\delta^{-1}I$ where $I$ is the $(p + 1)$ -by- $(p + 1)$ identity matrix
Computation:	For $n = 0,1,2,...$
	$\mathbf{x}(n) = \left[ \begin{matrix} x(n)\\ x(n-1)\\ \vdots\\ x(n-p) \end{matrix} \right]$
	$e(n) = d(n)-\mathbf{w}_{n}^{T}\mathbf{x}(n)$
	$\mathbf{g}(n)=\mathbf{P}(n-1)\mathbf{x}^{}(n)\left\{\lambda+\mathbf{x}^{T}(n)\mathbf{P}(n-1)\mathbf{x}^{}(n)\right\}^{-1}$
	$\mathbf{P}(n)=\lambda^{-1}\mathbf{P}(n-1)-\mathbf{g}(n)\mathbf{x}^{T}(n)\lambda^{-1}\mathbf{P}(n-1)$
	$\mathbf{w}_{n} = \mathbf{w}_{n-1}+\,e(n)\mathbf{g}(n)$

[edit] Example: Plant Identification

The goal for a plant identification structure is to match the properties of an unknown system (plant) with an adaptive filter. The following figure shows the block diagram of a plant identification system. The Matlab source code is RLSPlantIdent.

In general any adaptive filter can be used for plant identification, however in this example we use an LMS structure. This example allows us to discuss the merits of the step size factor $μ$ . The following figure shows the error signal for $λ = 0.95$ .

[edit] References

Monson H. Hayes Statistical Digital Signal Processing and Modeling, Wiley, 1996, ISBN 0-471-59431-8
Simon Haykin Adaptive Filter Theory, Prentice Hall, 2002, ISBN 0-13-048434-2

Categories: Digital signal processing | Filter theory

Recursive least squares filter

From Wikipedia, the free encyclopedia

Contents

[edit] Discussion

[edit] Recursive Algorithm

[edit] RLS algorithm summary

[edit] Example: Plant Identification

[edit] References

[edit] See also

Views

Navigation

Search

In other languages