DSP-CIS

Part-III : Optimal & Adaptive Filters

Chapter-9 :
Recursive Least Squares (RLS) Algorithms

Marc Moonen
Dept. E.E./ESAT-STADIUS, KU Leuven
[Link]@[Link]
[Link]/stadius/

Part-III : Optimal & Adaptive Filters

Optimal Filters - Wiener Filters

Chapter-7 General Set-Up / Applications / Wiener Filters
Adaptive Filters

Chapter-8 Least Mean Squares (LMS) Algorithms

Chapter-9 Recursive Least Squares (RLS) Algorithms

Chapter-10 Recursive Least Squares (RLS) Algorithms

(continued)

Chapter-11
Kalman Filters

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 2 / 40

1
 Recap 1/4

Norbert Wiener (1894-1964)

4

Introduction : Optimal and adaptive filters

Prototype optimal filtering set-up :

filter input
Design filter such that for a given
(i.e. ‘statistical info available’)
input signal, filter output signal is filter filter parameters
‘optimally close’ (to be defined)
to a given ‘desired output signal’. filter output
+

error desired signal

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 3 / 40

Recap 2/4
PS: Shorthand notation uk =u[k], yk =y[k], dk =d[k], ek=e[k],
Filter coefficients (‘weights’) are wl (replacing bl of

For adaptive filters wl also have a time index wl[k]

27

1 Will use filtering/ Wiener filters

Optimal

FIR filters (=tapped-delay line filter/‘transversal’ filter)

filter input u[k]

N
! −1
L yk = wl · uk−l u[k] u[k-1] u[k-2] u[k-3]

yk = ∑ wl .uk−l = w .u k = u .w
l=0
T T
k
w0[k] w1[k] w2[k] w3[k]

0
l=0 yk = wT · uk
previous chapters)

e[k]
filter output y[k]
where +

" error# desired signal d[k]

wT = w w w · · · w w

wT = !" w0 0 w01 …2 wL #$N −1

b
a+bw a
" #
uT k = uk uk−1 uk−2 · · · uk−N +1
!
u k = " uk uk−1 … uk−L $
T #

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 4 / 40

2
 Recap 3/4

2 Will use

30

Optimal filtering/ Wiener filters

Quadratic cost function :

minimum mean-square error (MMSE) criterion
= minimize 2
JJM (w) = Ε {ek22}}==ΕE {|d
SE (w) = E {e
MSE
k
{ 2
} {
( dkk −− yykk)|2} == ΕE {|d( dkk−−wuTTk w
uk)|2} }
E {x} is ‘expected value’ (mean) of x

Leads to ‘Wiener Filter’ solution

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 5 / 40

Recap 4/4

Introduction : Optimal and adaptive filters

Prototype adaptive filtering set-up :

filter input
Basic operation involves 2 processes :
1. filtering process adaptive filter parameters
2. adaptation process filter

adjusting filter parameters to filter output

(time-varying) environment +
adaptation is steered by error signal
error desired signal

• Depending on the application, either the filter parameters, the filter output
or the error signal is of interest

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 6 / 40

3
 Overview

Recursive Least Squares (RLS) Basics

• Least Squares (LS) Estimation
• Recursive Least Squares

Square-Root RLS Algorithms

• QRD-Based LS Estimation
• QRD-RLS
• Residual Extraction

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 7 / 40

Least Squares & RLS Estimation 5

1. Least Squares (LS) Estimation

Quadratic cost function

MMSE : (see Lecture 8)
2
MSE (w)
JMSE {ek2k}}==ΕE{|d
(w)==ΕE{e
2
{
( dkk −− yk |)2} == E{|d
2
} {(
Ε dk k−−wuTk uwk |2} T
)}
Least-squares(LS) criterion :
if statistical info is not available, may use an alternative ‘data-based’ criterion...

Lk
! Lk
! k
2 "L 2
JJLS
LS (w) =
(w) = e 2 2=
∑ |d − y |2=
∑ dlk −−uwTlTw
ekl = ( kdl − kyl ) = k=1 |d uk |2 ∑( )
k=1
l=1 k=1
l=1 l=1

Interpretation? : see below

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 8 / 40

4
 Least Squares & RLS Estimation 6

1. Least Squares (LS) Estimation

filter input sequence : u1, u2, u3, . . . , uL k

corresponding desired response sequence is : d1, d2, d3, . . . , ddLk
! $ ! $  ! TT

$! $
ee11 dd11 # uu
11  & # ww
# & # &  00 &
# e & # d & 
 e22   d2  # uuT2 T2  &# ww
. 1 1&

# . & ==# . 2 &−−  # .   ·& # . &
# .! & # .! &  # .!  & # . ! &
# eeL #" ddL # TTT & #ww
'"" k k *&
% kk &% " uukuL ' ()NL *&%
()
#$ ' ()
"#$ '"%* ()
#
k
%*$% ""#$ w
error
errorsignal
signale e dd UU w
+2 2 2 2
k

(w)==∑ eL
2
cost function :JJLSLS (w) l =e e=2 "e"
k=1 k
= d=−"d
2
Uw − 2Uw"22
l=1

→ linear least squares problem : minw "d − Uw"22

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 9 / 40

Least Squares & RLS Estimation 7

1. Least Squares (LS) Estimation

Lk
!
JJLS (w) =
LS (w) =
ee2kl2==!e!
∑ 22 T e = !d − Uw!22
e 2 ==eeT· .e = d −Uw 2
2 2
k=1
l=1

minimum obtained by setting gradient = 0 :

∂JLS (w) ∂
0=[ ]w=wLS = [ (dT d + wT U T Uw − 2wT U T d)]w=wL
∂w ∂w
T U w − 2 U T d]
= [2 "U#$ % " #$ % w=w LS
Xuu Xdu

Xuu · wLS = Xdu → wLS = X−1

uu Xdu ‘normal equations’.
This is the
‘Normal equations’ ‘Least Squares Solution’
(L+1 equations in L+1 unknowns)

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 10 / 40

10

5
 Least Squares & RLS Estimation 9

1. Least Squares (LS) Estimation

Note : correspondences with Wiener filter theory ?

♣ estimate X̄uu and X̄du by time-averaging (ergodicity!)

L
1 !1 k T 1 TT 11 1
estimate{X̄uuℵ
estimate {
}= } k=1 k ∑
uk. · uku=. uU · U = .Xuu T
uuL = l L l = L U U = .ℵuu
l=1 k k
L k
1 !1 1 T 11 1
estimate } =du } = uk.∑
estimate{X̄du{ℵ · dk u U T d = .ℵdu
U l ·=d = .Xdu
=l . d
L L L
k=1 k l=1 k k
leads to same optimal filter :

estimate{wW F } = ( L1 Xuu)−1 · ( L1 Xdu) = X−1

uu · Xdu = wLS

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 11 / 40

11

Least Squares & RLS Estimation 10

1. Least Squares (LS) Estimation

Note : correspondences with Wiener filter theory ? (continued)

♣ Furthermore (for ergodic processes!) :

L k
1 !1 1 1
lim k→∞ u.k∑
ℵX̄uuuu = lim uTl =Xuulim k→∞ .ℵuu
· uTku=l . lim
L→∞ L L
k=1k l=1 k
L→∞

k
1 !1 1
L
lim k→∞ u.k∑
= lim
ℵX̄dudu = dl 1=Xdulim k→∞ .ℵdu
· dku=l . lim
L→∞ L
k=1
k l=1 L→∞ L k

so that

lim
lim k→∞ w
L→∞ wLSLS ==ww F.
WWF

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 12 / 40

12

6
 Least Squares (LS) Estimation

In words:
Whenever statistical info (autocorrelation and
crosscorrelation) is missing, this can be estimated
from observed data (assuming ergodicity)
The Wiener filter solution, with true statistical
quantities replaced by estimated quantities, then
turns out to be the same as the LS solution
LS approach in itself optimizes a different (LS)
criterion, without any need for statistical
assumptions (e.g. ergodicity..)

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 13 / 40

13

Overview

Recursive Least Squares (RLS) Basics

• Least Squares (LS) Estimation
• Recursive Least Squares

Square-Root RLS Algorithms

• QRD-Based LS Estimation
• QRD-RLS
• Residual Extraction

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 14 / 40

14

7
 Least Squares & RLS Estimation 11

2. Recursive Least Squares (RLS)

k least squares problem is

For a fixed data segment 1...L,
minw !d − Uw!22
2

Matrices and vectors now

with time index added
 !# d1 $& !# uT1T $& !# w0 $&
d1 # d & # uu1T & # w &
min   2 − # T2  1
&&. # ! &
 2 ## ! && #
d
w[k] u !
d =  . # dU &= # 2  # &
 . ""#k$% " u..Tk &% #" wL &%
 "#$
"%#%$
dL dk uUTLk w[k]
2

−1
T #U
$ . k U k %& .U k d k
T T
w[k]
w
−1
= ℵ T[ k ] −1
.ℵ [ k ] =
LS = [Uuu U] du· U d
Wanted : recursive/adaptive algorithms
L → L +@
Can LS solution 1, time
slidingk be
windows, etc... from solution @ time k-1 ?
computed

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 15 / 40

15

Least Squares & RLS Estimation1 13 -

nk
Ra ates
2.1 Standard RLS d
up

T
Xuu
It is observed that ℵ [k] =ℵ= [k (L) +u [Link]+1
Xuu−1]+ k u(and
T ℵ [k] =ℵ [k −1]+ u k .dk )
uu(L + 1) uu L+1 du du
The matrix inversion lemma states that (check ‘matrix inversion lemma’ in Wikipedia)
1[Xuu(L)]−1u uT T [Xuu(L)]−1
[k]−1 =ℵ 1)][k−1−1]
ℵ [X −1
= [X−uu( (L)]T−1 − 1+u ).k kL+1
L+1 kk with k k =ℵ [k −1]−1 u k
uu uu(L +uu uu
1+ u kℵ [k −1]−1TL+1 u k[Xuu(L)]−1uL+1
uu
With this it is proved that:
Result : 'Kalman
!###gain "### vector'
$ ! 'a ##
priori residual'
#" ### $
w LS [k ] = w LS [k −1] + −1
ℵ [k ]−1 u k . (dk T− u Tk w LS [k −1])
wLS (L + 1) = wLS (L) +% [X###
uu (L
uu+ 1)]
& uL+1 · (d'
### L+1 − uL+1wLS (L))
! "# $ ! 1
"# $
=(
Kalman
1+u T gain −1 u k ). k k
k ℵuu [ k−1]
a priori residual

= standard recursive least squares (RLS) algorithm

O(N22) instead
Remark : O(L operations
of per 3
O(Ltime
) update
Remark : square-root algorithms with better numerical properties
see below

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 16 / 40

16

8
 Next to a mechanism for adding new observations, also need a mechanism
Least Squares & RLS Estimation
for removing old observations. First approach is as follows… 14

e
2.2 Sliding Window RLS

id
sl
Sliding window RLS
!
JLS (w) = L 2
k=L−M+1 ek
M = length of the data window

is
th
wLS (L) = [U (L)T U (L)]−1 [U (L)T d(L)]
" #$ %" #$ %
[Xuu(L)]−1 Xdu (L)

p
with

i
   
dL−M+1 uTL−M+1

k
   
 dL−M+2   uT 
d(L) =   U (L) =  L−M+2 .

S
 ..   . 
   
dL uTL

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 17 / 40

17

Least Squares & RLS Estimation 15

e
2.2 Sliding Window RLS

id
sl
It is observed that Xuu(L + 1) = Xuu(L) − uL−M+1uTL−M+1 +uL+1uTL+1
! "# $
L|L+1
Xuu

• leads to : updating (cfr supra) + similar downdating

is
th
• downdating is not well behaved numerically, hence to be avoided...

• simple alternative : exponential weighting

kip
S
DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 18 / 40

18

9
 Next to a mechanism for adding new observations, also need a mechanism
Least Squares & RLS Estimation
for removing old observations. Simpler approach is as follows… 16

2.3 Exponentially Weighted RLS

Exponentially weighted RLS : Goal is to give a smaller weight to ‘older’ data, i.e.
L!k
JLS (w)== ∑
LS (w) k=1
2(L−k)
λλ2(k−l ) 22
eel k
l=1

0 < λ < 1 is weighting factor or forget factor

1
1−λ is a ‘measure of the memory of the algorithm’
T −1 T 2
Which leads w
to…
LS (L) = [U (L) U (L)] [U (L) d(L)]
" #$ %" #$ %
" % −1"
k−1[Xuu (L)]
$ λ d1 ' $ λ k−1Xudu1T(L) %' " w0 %
$ '
with $ k−2 ' $ k−2 T ' $ w1 '
min $ λ d2 ' − $ λ u 2 ' .$ '
w[k] L−1
$ λ !d1 ' $ ! ' $uT1 ! '
λL−1
$ L−20  ' $
# λ λ d2k  & # λ 0 u T
k λ%
' $ TwL '
L−2
& #u &
d(L) =  $ "#%
"$# % " $$ =#$$
U (L) 2
.d
.  U  .. w[k]

 k  k 2

λ0dL λ0−1
uTL
w[k] = ℵuu [ k ]−1 .ℵdu [ k ] = #$U kTU k %& .U kT d k
DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 19 / 40

19

Least Squares & RLS Estimation 17

2.3 Exponentially Weighted RLS

It is observed that X
ℵuu(L 2
[k]+=1)λ = .ℵλ2X[k .uTk uTL+1
(L) +uukL+1
−1]+ (and ℵ [k] = λ 2 .ℵ [k −1]+ u k .dk )
uu uu uu du du

hence
1 1 1
ℵ [k]−1 = 2 ℵ−1 [k −1] −1 − ( ℵ [k −1]−1 u k
1 [X −1 T 1T[X (L)]−1
2 uu (L)] uL+1 u ).k kk
uu with kk =
uu[Xuu(L +λ 1)]uu = λ12 [Xuu(L)]−1 λ2k
L+1
1 −T λ 1+ 1 uT−1 [X (L)] −1u λ 2 uu
1+ 2 u kℵ [k −1] 2 L+1u kuu L+1
λ uu λ

wwLSLS[k]
(L + = wwLS
1) = (L) +
LS[k [Xuu+ℵ uL+1−1
(L + 1)]−1[k] · (dL+1 uTL+1w
u k−.(d T .
LS (L))
−1] k − u k w LS [k −1])
uu
i.e. exponential weighting hardly changes RLS formulas.. (easy!)

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 20 / 40

20

10
 Recursive Least Squares

Computational Complexity:
Standard RLS algorithm has O(L2) computational
complexity per update
Compare to O(L) for LMS (=cheaper, but slow convergence)

In Chapter-10, will present ‘Fast RLS’ algorithms with

O(L) computational complexity (and without compromising
convergence properties)

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 21 / 40

21

Recursive Least Squares

Numerical Analysis/Stability:
Standard RLS algorithm has been shown to have
unstable quantization error propagation (in low-
precision implementation)
i.e. when an infinite precision version is run next to a finite precision version
(both fed with the same input signals), then after xx iterations the finite precision
version produces results (far) away from the infinite precision results

Will present ‘Square-Root RLS’ algorithm (‘QRD-RLS’) which

is shown to be perfectly stable numerically

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 22 / 40

22

11
 Overview

Recursive Least Squares (RLS) Basics

• Least Squares (LS) Estimation
• Recursive Least Squares

Square-Root RLS Algorithms

• QRD-Based LS Estimation
• QRD-RLS
• Residual Extraction

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 23 / 40

23

Least Squares & RLS Estimation 20

3.1 QRD-based RLS Algorithms

QR Decomposition for LS estimation

least squares problem

minw !d − Uw!22

‘square-root algorithms’ based on QR decomposition (QRD) :

R $ !
U == !"#$
!"#$ QQ · .R !
U! ! !"#$# &= Q
! . R
!
L×N
kx ( L+1) kxk " 0
!%
N×N Q(:,1:L+1) ( L+1) x ( L+1)
L×N
kx ( L+1)
square * rectangular rectangular * square
QT · Q = I, Q is orthogonal R is upper triangular

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 24 / 40

24

12
 Least Squares & RLS Estimation 21

Everything you need to know about QR decomposition

Example :
U Q
Q
!
 "# $ ! "# $ ! R
1 6 10 0.182 0.816 0.174  "# $
2    5.477 14.605 −5.112
 7 −11   0.365 0.408 −0.619  
3 = · 0 4.082 8.981 
8 12   0.547 0 0.716 
0 0 20.668
4 9 −13 0.730 −0.408 −0.270

Remark : QRD ≈ Gram-Schmidt

Remark : U T · U = RT · R
R is Cholesky factor or square-root of U T · U
→ ‘square-root’ algorithms !

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 25 / 40

25

Least Squares & RLS Estimation 22

3.1 QRD-based RLS Algorithms

QRD for LS estimation

if
(**) orthogonal transformation preserves norm

! R $
U Q .# !
Q . R
! = ! &= ! !
= Q · R kxk
Ukx ( L+1) "!0 % Q(:,1:L+1) ( L+1) x ( L+1)

kx ( L+1)

then
2
! " 2!(**) " 2 " z % " R %
QT w· dU −Uw
T
min d =2 =Rmin
z w Q (d −Uw) = min w $ '−$ 'w
2
# * & # 0 & 2

with this triangular

backsubstitution

Rw⋅LSw=LSR−1=· zz ⇒ w LS = R −1 ⋅ z = [Q! TU ]−1 ⋅ Q! T d

−1
This is a numerically better way of computing the LS solution, better than w LS = !"U U #$ ⋅U d
T T

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 26 / 40

26

13
 Square Root RLS Algorithms
PS: This (= QRD + backsubstitution) is also the way Matlab™
solves LS problems ( cfr “w=U\d” or “w=mldivide(U,d)” )

Now back to recursive least squares (RLS) estimation...

This will be based on ‘recursive QRD’, i.e ‘QRD-updating’
DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 27 / 40

27

Least Squares & RLS Estimation 23

3.1 QRD-based RLS Algorithms

QR-updating for RLS estimation

Assume we have computed the QRD at time kk-1
! " !!
" R[kT−1] z[k −1] $ = Q[k
"
#Q̃(k) · U(k) d(k) % = R(k) ⋅"
−1]Tz(k)
# U k−1 d k−1
$
%

w LS [k(k)−1] −1]−1 ⋅ z[k −1]

= R[k −1
The corresponding LS solution is w LS = [R(k)] · z(k)
Our aim is to update the QRD into
! " ! "
Q̃(k + 1)T · Uz[k]
! R[k] (k + 1) #d(k
=+ !1) =T ⋅R(k
Q[k] ! U d+ 1)#
+ 1) z(k
" $ " k k $
−1 z(k + 1)
and then compute wLS(k + 1) = w
[R(k[k] = R[k]· −1
+ 1)]
LS ⋅ z[k]

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 28 / 40

28

14
 Least Squares & RLS Estimation 24

3.1 QRD-based RLS Algorithms

QR-updating for RLS estimation

It is proved that the relevant QRD-updating problem is

"! R[k] z[k]" % "! R[k −1]" z[k −1] %

R(k + 1) z(k + 1) T $ R(k) z(k) '
$ Q[k]
'←=Q(k + 1)T ⋅·$ T T
# 0!0
$ 0 !∗ '
& # uk+1udkk+1 dk '
&

PS: This is based on a QR-factorization as follows:

" R[k −1] % " %
$ ' = Q[k] . $ R[k] '
$# uTk ! $ 0 '
'& (L+2)×(L+2) # &
!#"#$ ! #"#
$
(L+2)×(L+1) (L+2)×(L+1)

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 29 / 40

29

Least Squares & RLS Estimation 25

3.1 QRD-based RLS Algorithms

QR-updating for RLS estimation

"! R[k] z[k] % " " R[k!−1] "
z[k −1] %
$ R(k + 1) z(k + 1) Q[k]
'=← T $
⋅ R(k) z(k) '
Q(k +$ 1)T ·u T T dk '
# 0!0
$ 0 ∗! ' & # kuk+1 dk+1 &
−1
w LS [k]
R(k + 1)= ·R[k]
wLS(k⋅ z[k] z(k + 1) . backsubstitution’
+ 1) ==‘triangular
= square-root (information matrix) RLS = QRD-RLS

Remark . with exponential weighting

! " ! "
" R[k] 1) z(k%+ 1) T "$ λ ⋅ R[kT −1]λ R(k)
R(k + z[k] −1] %
λ z(k)
λ ⋅ z[k
'
$ ' = Q[k]←⋅Q(k + 1)T ·
# 0!00
$ ∗ ' &! $
# uk uTk+1dk dk+1 '&

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 30 / 40

30

15
 Least Squares & RLS Estimation 26

3.1 QRD-based RLS Algorithms

Will now look into the details of how the Q[k] can be constructed
QRD updating
! " " R[k!−1] "
" R(k + 1) z(k
R[k] z[k]+ 1)% T $ T z[k −1] %
R(k) z(k)
$ ← Q(k
' = Q[k] ⋅ + 1) · T '
# 0!0
$ 0 ∗! ' & $ u Tkuk+1 dk+1
# dk ' &
basic tool is Givens rotation
i j
 ↓ ↓ 
Ii−1 0 0 0 0
def
Gi,j,θ =   0 cos θ 0 sin θ 0 ← i
 0 0 Ij−i−1 0 0
 
 0 − sin θ 0 cos θ 0 ← j
0 0 0 0 Im−j

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 31 / 40

31

Least Squares & RLS Estimation 27

3.1 QRD-based RLS Algorithms

QRD updating

Givens rotation applied to a vector x̃ = Gi,j,θ · x :

x̃i = cos θ · xi + sin θ · xj

x̃j = − sin θ · xi + cos θ · xj

x̃l = xl for l "= i, j

x
x̃j = 0 iff tan θ = xj !
i

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 32 / 40

32

16
 Least Squares & RLS Estimation 28

3.1 QRD-based RLS Algorithms

QRD updating

"! R[k] z[k] % " " R[k!−1] "

z[k −1] %
$ R(k + 1) z(k + 1) Q[k]T $ R(k) z(k) '
'=← Q(k +$ 1)T ·u T T
⋅
dk '
# 0!0
$ 0 ∗! ' & # kuk+1 dk+1 &
Q(k[k]
+1) is constructed as a product/sequence of
sequence of Givens transformations
       
3-by-3 example

x x x x x x x x xx x x xxxx
 x x x  x x x  x x x
  x x x
 → → 
 
→


 x x  x x  x x   x x
x x x x x x x x x !

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 33 / 40

33

Square Root RLS Algorithms

A graphical representation (i.e. a ‘realization’) of the

QRD updating process is presented in the next slide

This is also referred to as a ‘signal flow graph’ (SFG)

The SFG in the next slide will be further developed in

later slides, and also used explicitly for the (graphical)
derivation of a ‘fast’ RLS algorithm (Chapter 10)

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 34 / 40

34

17
 29

Least Squares
!"
R[k]&z[k] RLS %" Estimation!
" R[k
$ R(k + 1) z(k + 1) ' = Q[k] ⋅ $
−1]
R(k)
"
z(k) z[k −1] T
%
'
→ Q(k + 1)$T · uT T
# 0!0
$ 0 ∗
! '
& # uk+1
k dk+1 dk '
&
u(1) u(2) u(3) u(4) d
u[k] u[k-1] u[k-2] u[k-3] d[k]
rotation cell
R11 R12 R13 R14 z11 a a’

0
b b’
R22 R23 R24 z21

0
memory cell
R33 R34 z31
(delay)
4-by-4 example

0
R44 z41

0
*

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 35 / 40

35

Overview

Recursive Least Squares (RLS) Basics

• Least Squares (LS) Estimation
• Recursive Least Squares

Square-Root RLS Algorithms

• QRD-Based LS Estimation
• QRD-RLS
• Residual Extraction

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 36 / 40

36

18
 Square Root RLS Algorithms

Residual Extraction

So far, the “star“ (*) in the update equation

(also appearing at the bottom of the SFG)
has not been considered/defined/used

It will turn out that this “star” can be used to compute

least squares residuals
(without explicitly computing the least squares filter vector!)

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 37 / 40

37

Least Squares & RLS Estimation 30

3.1 QRD-based RLS Algorithms

Residual extraction

!! R[k] z[k] $ " ! R[k! −1] z[k "−1] $

# R(k + 1) z(k + 1) Q[k]
&== T #
⋅ T R(k) z(k) &
Q(k +#1) · u T T dk &
" 0!0
# 0 ε ! &
% " uk+1 dk+1
k %

From this it is proved that the ‘a posteriori residual’ is

#
L+1
ddk+1
T
k −−uukk+1
Tw w [k]
LSLS (k = ∏
+ ε1)⋅= ! · i=1Ncos(θ )
i=1 cos(θii )

and the ‘a priori residual’ is

ddk+1 T w
u Tkuw ! ε
k −− LS [k
k+1 LS−1] = #N L+1
(k) =
i=1 cos(θ i) θ ) Hence ε is geometric mean of
∏ cos(
i=1 i
a posteriori & a priori residual

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 38 / 40

38

19
 32

Least
Residual Squares
extraction (v1.1)&: RLS Estimation
u[k]
u(1) u[k-1]
u(2) u[k-2]
u(3) u[k-3]
u(4) d[k]
d 1
0 rotation cell
R11 R12 R13 R14 z11 a a’

0 0 b b’
R22 R23 R24 z21

0 0
memory cell
R33 R34 z31
(delay)

0 0
R44 z41

0
ε cos

L+1 !N
ddkk+1
−u T
−kuw
T ⋅∏
w[k](k=+ε1) =!· cos( θ i )i)
cos(θ
LS
k+1 LS i=1 output1
i=1

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 39 / 40

39

Least Squares & RLS Estimation

33

far-end signal
Example
Acoustic Echo Cancellation
u[k] u[k-1] u[k-2] u[k-3]

near-end signal
+ residual echo
d

1
0

0 0

0 0

0 0

DSP-CIS 2022-2023 / Chapter-9: RLS Algorithms 40 / 40

40

20