1

Design of Robust Adaptive Beamforming
Algorithms Based on Low-Rank and
Cross-Correlation Techniques

arXiv:1606.01313v1 [] 4 Jun 2016

Hang Ruan and Rodrigo C. de Lamare

Abstract—This work presents cost-effective low-rank techniques for designing robust adaptive beamforming (RAB) algorithms. The proposed algorithms are based on the exploitation
of the cross-correlation between the array observation data and
the output of the beamformer. Firstly, we construct a general
linear equation considered in large dimensions whose solution
yields the steering vector mismatch. Then, we employ the idea of
the full orthogonalization method (FOM), an orthogonal Krylov
subspace based method, to iteratively estimate the steering
vector mismatch in a reduced-dimensional subspace, resulting in
the proposed orthogonal Krylov subspace projection mismatch
estimation (OKSPME) method. We also devise adaptive algorithms based on stochastic gradient (SG) and conjugate gradient
(CG) techniques to update the beamforming weights with low
complexity and avoid any costly matrix inversion. The main
advantages of the proposed low-rank and mismatch estimation
techniques are their cost-effectiveness when dealing with high
dimension subspaces or large sensor arrays. Simulations results
show excellent performance in terms of the output signal-tointerference-plus-noise ratio (SINR) of the beamformer among
all the compared RAB methods.
Index Terms—robust adaptive beamforming, low-rank techniques, low complexity algorithms.

I.

INTRODUCTION

Adaptive beamforming has been one of the most important
research areas in sensor array signal processing. It has also
been recognized that traditional adaptive beamformers are
extremely sensitive to environmental uncertainties or steering
vector mismatches, which may be caused by many different
factors (e.g., imprecise antenna size calibration, signal pointing
errors or local scattering). Furthermore, some radar systems in
advanced applications require antenna arrays with a very large
number of sensor elements in highly dynamic environments,
which leads to the increase of computational complexity
and the decrease of the convergence rate for computing the
parameters of the beamformer.
A. Prior and Related Work
In order to mitigate the effects of uncertainties on adaptive
beamformers, robust adaptive beamforming (RAB) techniques
have been developed. Popular approaches include worst-case
optimization [2], [9], diagonal loading [4], and eigen-subspace
decomposition and projection techniques [6], [8], [11]. However, these RAB approaches have some limitations such as
their ad hoc nature, high probability of subspace swap at
low signal-to-noise ratio (SNR) [3] and high computational

cost due to online optimization or subspace decomposition
techniques.
Furthermore, in the case of large sensor arrays the above
mentioned RAB methods may encounter problems for their
application. This is because in these RAB algorithms, a cubic
or greater computational cost is required to compute the
beamforming parameters. Therefore, dimensionality reduction
(or rank-reduction) methods ([15]-[31], [33]-[36]) have been
employed and developed to reduce the complexity and improve
the convergence rate.
In the recent years, great efforts have been devoted to
the investigation of dimensionality reduction techniques [38],
[39], [40], [41], [42], [43], [44], [45], [59], [47], [48], [49],
[50], [51], [52], [53], [54], [55], [56], [57], [58], [59], [60],
[61], [62], [63], [64], [65], [66], [70], [68], [69], [70], [71],
[72], [73], [74], [75], [76], [77], [78], [79], [80], [81], [82],
[83], [84], [85], [86], [87], [88], [89], [90]. for RAB. The
beamspace approach of [23] projects the data onto a lower
dimension subspace by using a beamspace matrix, whose
columns are determined by linearly independent constrained
optimization problems. A more effective approach (i.e., [16][20],[24],[25]) is based on preprocessing the array observation
data using a Krylov subspace. However, there are different
ways to generate the Krylov subspace and the choice usually
depends on the cost and the performance. The Arnoldi method
[12], [13], [21] and the Lanczos iterations [12], [13], [15]
are typical approaches used to generate orthogonal Krylov
subspaces, whereas [19] also introduces a method to generate
non-orthogonal ones. However, the main challenge in these
techniques is the model order determination. Specifically, the
model order must be properly chosen to ensure robustness to
over-determination of the system model order [17]. Another
effective approach to dimensionality reduction is the joint
iterative optimization (JIO) [26]-[31] techniques, [33]-[36],
which employ a subspace projection matrix and jointly and
iteratively optimize the bases of the subspace and the beamformer weights. The work in [27] has developed a recursive
least squares (RLS) adaptive algorithm based on widelylinear processing using the JIO technique. The study in [29]
has devised efficient stochastic gradient (SG) and RLS RAB
algorithms from a modified JIO (MJIO) scheme.
B. Contributions
In this work, we propose and study novel RAB algorithms
that are based on low-rank and cross-correlation techniques.

2

In the proposed techniques, we exploit the prior knowledge
that the steering vector mismatch of the desired signal is
located within an angular sector which is assumed known.
The proposed algorithms are based on the exploitation of the
cross-correlation between the array observation data and the
output of the beamformer, which avoids costly optimization
procedures. We firstly construct a linear system (considered in
high dimension) involving the mismatched steering vector and
the statistics of the sampled data. Then we employ an iterative
full orthogonalization method (FOM) [12], [13] to compute an
orthogonal Krylov subspace whose model order is determined
by both the minimum sufficient rank [17], which ensures
no information loss when capturing the signal of interest
(SoI) with interferers, and the execute-and-stop criterion of
FOM [12], [13], which automatically avoids overestimating
the number of bases of the computed subspace. The estimated
vector that contains the cross-correlation between the array
observation data and the beamformer output is projected
onto the Krylov subspace, in order to update the steering
vector mismatch, resulting in the proposed orthogonal Krylov
subspace projection mismatch estimation (OKSPME) method.
Furthermore, based on the OKSPME method, we have
also devised adaptive stochastic gradient (SG), conventional
conjugate gradient (CCG) and modified conjugate gradient
(MCG) algorithms derived from the proposed optimization
problems to reduce the cost for computing the beamforming
weights for large sensor arrays [91], [92], [93], [94], [95],
[96], [97], [98], [99], [100], [101], [102], [103], [104], [105],
[106], [107], [108], [109], [110], [111], [112], [114], [115],
[116], [117], [118], [119],[120], [121], [122], [123], [124],
[125], [126], [131], [128], [129], [130], [131], [132], [133],
[134], [135], [136]. , resulting in the proposed OKSPMESG, OKSPME-CCG and OKSPME-MCG RAB algorithms.
We remark that the steering vector is also estimated and
updated using the CG-based recursions to produce an even
more precise estimate. Derivations of the proposed algorithms
are presented and discussed along with an analysis of their
computational complexity.
Moreover, we develop an analysis of the mean squared
error (MSE) between the estimated and the actual steering
vectors for the general approach of using a presumed angular
sector associated with subspace projections. This analysis
mathematically describes how precise the steering vector mismatch can be estimated. Upper and lower bounds are derived
and compared with the approach in [6]. Another analysis on
the computational complexity of the proposed and existing
algorithms is also provided.
In the simulations, we consider local scattering scenarios
(both coherent and incoherent) to model the mismatch effects.
We also study the performance of the proposed algorithms
by testing the output signal-to-interference-plus-noise ratio
(SINR) of the beamformer with respect to training snapshots
and different input SNRs. The number of sensor elements and
interferers is also varied and compared in each scenario to
provide a comprehensive performance study. In summary, the
contributions of this work are:
• The proposed OKSPME RAB method.
• The development of the modified SG and CG type

OKSPME RAB algorithms.
An analysis of the computational complexity and the
MSE performance of the proposed and existing RAB
algorithms.
The remaining sections of this paper are organized as follows: The system model and problem statement are described
in Section II. Section III introduces the proposed OKSPME
method, whereas Section IV introduces the proposed robust
adaptive algorithms. Section V provides the MSE analysis of
the steering vector estimation and the complexity analysis.
Section VI presents and discusses the simulation results.
Section VII gives the conclusion.
•

II. S YSTEM M ODEL AND P ROBLEM S TATEMENT
Let us consider a linear antenna array of M sensors and
K narrowband signals which impinge on the array. The data
received at the ith snapshot can be modeled as
x(i) = A(θ)s(i) + n(i),

(1)

where s(i) ∈ CK×1 are uncorrelated source signals, θ =
[θ1 , · · · , θK ]T ∈ RK is a vector containing the directions
of arrival (DoAs) and [.]T denotes the transpose, A(θ) =
[a(θ1 )+e, · · · , a(θK )] = [a1 , · · · , aK ] ∈ CM×K is the matrix
which contains the steering vector for each DoA and e is the
steering vector mismatch of the desired signal, n(i) ∈ CM×1
is assumed to be complex circular Gaussian noise with zero
mean and variance σn2 . The beamformer output is given by
y(i) = wH x(i),

(2)

where w = [w1 , · · · , wM ]T ∈ CM×1 is the beamformer
weight vector, where (·)H denotes the Hermitian transpose.
The optimum beamformer is computed by maximizing the
SINR and is given by
2

SIN R =

σ12 |wH a1 |
.
wH RI+N w

(3)

where σ12 is the desired signal power, RI+N is the
interference-plus-noise covariance (INC) matrix. The problem
of maximizing the SINR in (3) can be cast as the following
optimization problem:
wH RI+N w

minimize
w

subject to wH a1 = 1,

(4)

which is known as the MVDR beamformer or Capon beamformer [1], [4]. The optimum weight vector is given by
wopt =

−1
RI+N
a1
−1
aH
1 RI+N a1

. Since RI+N is usually unknown in practice, it can be
estimated by the sample covariance matrix (SCM) of the
received data as
i

1X
R̂(i) =
x(k)xH (k).
i

(5)

k=1

Using the SCM for directly computing the weights will lead
to the sample matrix inversion (SMI) beamformer wSMI =

3

R̂−1 a1
.
−1 a
aH
1
1 R̂

However, the SMI beamformer requires a large
number of snapshots to converge and is sensitive to steering
vector mismatches [2], [3]. As previously mentioned, most of
the conventional and existing RAB algorithms are computationally costly especially when encountering arrays with a very
large number of sensors. Therefore, the RAB design problem
we are interested in solving includes the following aspects:
• To design cost-efficient algorithms that are robust against
uncertainties and values of SNRs and interferers in the
presence of uncertainties in the steering vector of a
desired signal.
• The proposed algorithms must preserve their robustness
and low-complexity features for large sensor arrays.

mismatch within the presumed angular sector, say â1 (0) and
set â1 (1) = â1 (0). By adding the snapshot index i, we can
rewrite the array observation data as
x(i) = â1 (i)s1 (i) +

K
X

ak (i)sk (i) + n(i),

(10)

k=2

where â1 (0) and â1 (i) (i = 1, 2, · · · ) designate the initial
guess of the steering vector and its estimate at the ith snapshot,
respectively.
Pre-multiplying the above equation by âH
1 (i) we have
H
âH
1 (i)x(i) = â1 (i)â1 (i)s1 (i) +

K
X

âH
1 (i)ak (i)sk (i) + n(i).

k=2

III. P ROPOSED OKSPME M ETHOD
In this section, the proposed OKSPME method is introduced. This method aims to construct a linear system involving
only known or estimated statistics and then projects an estimated cross-correlation vector between the array observation
data and the beamformer output onto an orthogonal Krylov
subspace, in order to update the steering vector mismatch
with reduced complexity. The SCM of the array observation
data is estimated by (5). The cross-correlation vector between
the array observation data and the beamformer output can
be expressed as d = E[xy ∗ ] (where [.]∗ denotes complex
conjugation) or equivalently as
d = E[(As + n)(As + n)H w].

(6)

H
Assuming that |aH
k w| ≪ |a1 w| for k = 2, · · · , K and all
signals have zero mean, the cross-correlation vector d can be
rewritten as
H
d = E[(As + n)(s∗1 aH
1 w + n w)].

(7)

Note that we also assume that the desired signal is statistically independent from the interferers and the noise, i.e.,
E[sk s∗1 ] = 0 and E[sk ak s∗1 aH
1 w] = 0 for k = 2, · · · , K. With
this assumption the desired signal power is not statistically
affected by the interference and (7) can be rewritten as
H
d = E[σ1 2 aH
1 wa1 + nn w],

(8)

which can be estimated by the sample cross-correlation vector
(SCV) given by
i

d̂(i) =

1X
x(k)y ∗ (k).
i

(9)

k=1

A. Desired Signal Power Estimation
In this subsection, we describe an iterative method for the
desired signal power (σ12 ) estimation based on our prior work
in [32], which can be accomplished by directly using the
desired signal steering vector. Alternatively, a designer can
employ a maximum likelihood (ML) or a minimum variance (MV) estimator for computing the desired signal power.
However, the approach described here has lower complexity
than the ML and the MV estimators. In the adopted method,
we need to choose an initial guess for the steering vector

(11)
Here we assume that each of the interferers is orthogonal or
approximately orthogonal to the desired signal. Specifically,
the steering vector of each of the interferers is orthogonal (âH
1 (i)ak (i) = 0, k = 2, 3, · · · , K), or approximately
H
orthogonal (âH
1 (i)ak (i) ≪ â1 (i)â1 (i), k = 2, 3, · · · , K)
to the desired signal steering vector (i.e., â1 (i)), so that
âH
1 (i)ak (i) (k = 2, 3, · · · , K) approaches zero and the term
K
P
âH
1 (i)ak (i)sk (i) in (11) can be neglected, resulting in
k=2

H
H
âH
1 (i)x(i) = â1 (i)â1 (i)s1 (i) + â1 (i)n(i).

(12)

2
Taking the expectation of |âH
1 (i)x(i)| , we obtain
2
H
H
∗
E[|âH
1 (i)x(i)| ] = E[(â1 (i)â1 (i)s1 (i) + â1 (i)n(i))
H
(âH
1 (i)â1 (i)s1 (i) + â1 (i)n(i))].

(13)

Assuming that the noise is statistically independent from
the desired signal, then we have
2
H
2
2
E[|âH
1 (i)x(i)| ] = |â1 (i)â1 (i)| E[|s1 (i)| ]
H
+ âH
1 (i)E[n(i)n (i)]â1 (i),

(14)

where E[n(i)nH (i)] represents the noise covariance matrix
Rn (i) that can be replaced by σn2 IM , where the noise variance
σn2 can be easily estimated by a specific estimation method.
A possible approach is to use a Maximum Likelihood (ML)
based method as in [14]. Replacing the desired signal power
E[|s1 (i)|2 ] and the noise variance σn2 by their estimates σ̂12 (i)
and σ̂n2 (i), respectively, we obtain
σ̂12 (i) =

2
H
2
|âH
1 (i)x(i)| − |â1 (i)â1 (i)|σ̂n (i)
.
H
2
|â1 (i)â1 (i)|

(15)

The expression in (15) has a low complexity (O(M )) and
can be directly implemented if the desired signal steering
vector and the noise level are accurately estimated.
B. Orthogonal Krylov Subspace Approach for Steering Vector
Mismatch Estimation
An orthogonal Krylov subspace strategy is proposed in order
to estimate the mismatch with reduced cost and deal with
situations in which the model order is time-varying. Our idea
is based on constructing a linear system, which considers the
steering vector mismatch as the solution, and solving it by

4

using an iterative Krylov subspace projection method. To this
end, consider a general high-dimensional linear system model
given by
Ba1 = b,
(16)
where B ∈ CM×M and b ∈ CM×1 . Then we need to express
B and b only using available information (known statistics or
estimated parameters), so that we can solve the linear system
with the Krylov subspace of order m (m ≪ M ) described by
Km = span{b, Bb, B2 b, · · · , Bm b}.

(17)

TABLE I
A RNOLDI - MODIFIED G RAM -S CHMIDT ALGORITHM
For j = 1, 2, · · · do:
Compute uj = R̂tj
For l = 1, 2, · · · , j, do:
hl,j =< uj , tl >
uj = uj − hl,j tl
End do.
Compute hj,j+1 = kuj k.
If hj,j+1 = 0 or j ≥ K + 1,
set m = j;
break;
u
Else compute tj+1 = h j .
j,j+1
End do.

Taking the complex conjugate of (12), we have
∗
H
xH (i)â1 (i) = âH
1 (i)â1 (i)s1 (i) + n (i)â1 (i).

(18)

Pre-multiplying both sides of (18) by the terms of (10), then
adding an extra term δIâ1 (i) (where δ is a small positive
number defined by the user) and simplifying the terms, we
obtain

T̂(i) = [t1 (i), t2 (i), · · · , tm (i)],

∗
(x(i)xH (i) + δI)â1 (i) = â1 (i)âH
1 (i)â1 (i)s1 (i)s1 (i)

+ n(i)nH (i)â1 (i) + δâ1 (i). (19)
Replacing x(i)xH (i) + δI by R̂(i), s1 (i)s∗1 (i) by σ̂12 (i) and
n(i)nH (i) by σ̂n2 (i)IM , we obtain
R̂(i)â1 (i) =

2
â1 (i)âH
1 (i)â1 (i)σ̂1 (i)

+

{z

|

(σ̂n2 (i)

b̂(i)

+ δ)â1 (i),
}

(20)
in which by further defining the expression on the right-hand
side as b̂(i), we can rewrite (20) as
R̂(i)â1 (i) = b̂(i).

(21)

As can be seen (21) shares the same form as the linear system of equations in (16) and b̂(i) can be expressed in terms of
â1 (i), σ̂12 (i) and σ̂n2 (i) whereas R̂(i) can be estimated by (5).
In the following step, we employ the Arnoldi-modified GramSchmidt algorithm from the FOM method [12], [13] associated
with the minimum sufficient rank criterion discussed in [17] to
compute an orthogonal Krylov subspace. We define a residue
vector to represent the estimation error in the ith snapshot as
r̂(i) = b̂(i) − R̂(i)â1 (i),
and let
t1 (i) =

r̂(i)
.
kr̂(i)k

the number of signal sources), so that no more subspace
components are necessary for capturing the SoI from all
the existing signal sources.
Now by inserting the snapshot index, we have

(22)

(24)

and the Krylov subspace projection matrix is computed by
P̂(i) = T̂(i)T̂H (i).

(25)

It should be emphasized that the Krylov subspace matrix
T̂(i) obtained here is constructed by starting with the residue
vector r̂(i). In other words, T̂(i) is constructed with the
estimation error of the steering vector. In order to extract the
estimation error information and use it to update the steering
vector mismatch, we can project the SCV d̂(i) in (9) onto
P̂(i) and add the estimation error to the current estimate of
â1 (i) as
P̂(i)d̂(i)
.
(26)
â1 (i + 1) = â1 (i) +
kP̂(i)d̂(i)k
C. INC Matrix and Beamformer Weight Vector Computation
Since we have estimated both the desired signal power
σ̂12 (i) and the mismatched steering vector in the previous
subsections, the INC matrix can be obtained by subtracting
the desired signal covariance matrix out from the SCM as
R̂I+N (i) = R̂(i) − σ̂12 (i)â1 (i)âH
1 (i).

(27)

The beamformer weight vector is computed by
(23)

Then the Krylov subspace bases can be computed using
the modified Arnoldi-modified Gram-Schmidt algorithm as in
Table I (the snapshot index i is omitted here for simplicity).
In Table I, <, > denotes the inner product and the parameters hl,j (l, j = 1, 2, · · · , m) are real-valued coefficients,
the model order is determined once if one of the following
situations is satisfied:
• The execute-and-stop criterion of the original Arnoldimodified Gram-Schmidt algorithm is satisfied (i.e.,
hj,j+1 = 0).
• The minimum sufficient rank for addressing the SoI and
interferers is achieved (i.e., j ≥ K + 1, where K is

ŵ(i) =

R̂−1
I+N (i)â1 (i)
−1
âH
1 (i)R̂I+N (i)â1 (i)

,

(28)

which has a computationally costly matrix inversion R̂−1
I+N (i).
The proposed OKSPME method is summarized in Table II. In
the next section, we will introduce adaptive algorithms to avoid
matrix inversions and reduce the complexity.
IV. P ROPOSED A DAPTIVE A LGORITHMS
This section presents adaptive strategies based on the OKSPME robust beamforming method, resulting in the proposed OKSPME-SG, OKSPME-CCG and OKSPME-MCG algorithms, which are especially suitable for dynamic scenarios.
In the proposed adaptive algorithms, we estimate the desired

5

TABLE II
P ROPOSED OKSPME METHOD

By substituting λL in (30) again, the weight update equation
for OKSPME-SG is obtained as

Initialization:
ŵ(1) = 1;
w(i + 1) = (I − µσ̂12 (i)â1 (i)âH
1 (i))w(i)
Choose an initial guess â1 (0) within the sector and set â1 (1) = â1 (0);
âH (i)x(i)â1 (i)
For each snapshot i = 1, 2, · · · :
− µ(σ̂12 (i)â1 (i) + y ∗ (i)(x(i) − 1 H
)).
i
P
â1 (i)â1 (i)
R̂(i) = 1i
x(k)xH (k)
(32)
k=1

i
P

The adaptive SG recursion circumvents a matrix inversion
when
computing the weights using (28), which is unavoidable
k=1
Step 1. Compute the desired signal power
in
OKSPME.
Therefore, the computational complexity is
2
|âH (i)x(i)|2 −|âH
3
1 (i)â1 (i)|σ̂n (i)
reduced
from
O(M
) in OKSPME to O(M 2 ) in OKSPMEσ̂12 (i) = 1
2
|âH
1 (i)â1 (i)|
SG. It is also important that the step size µ should satisfy
Step 2. Determine the Krylov subspace
2
2
0 < µ < σ̂21(i) to guarantee that I − µσ̂12 (i)â1 (i)âH
b̂(i) = â1 (i)âH
(i)â
(i)σ̂
(i)
+
σ̂
(i)â
(i)
1
1
1 (i)
1
1
n
1
r̂(i) = b̂(i) − R̂(i)â1 (i)
is always a positive-definite matrix so that (32) is ensured
r̂(i)
t1 (i) = kr̂(i)k
converging to a solution. To implement OKSPME-SG, Step
Apply the algorithm in Table I to determine m and t1 (i),· · · ,tm (i)
1, Step 2 and Step 3 from Table II and (32) are required.
1
i

d̂(i) =

x(k)y ∗ (k)

T̂(i) = [t1 (i), t2 (i), · · · , tm (i)]
Step 3. Update the steering vector
P̂(i) = T̂(i)T̂H (i)
d̂(i)
â1 (i + 1) = â1 (i) + kP̂(i)
P̂(i)d̂(i)k
â1 (i + 1) = â1 (i + 1)/kâ1 (i + 1)k
Step 4. Compute the weight vector
R̂I+N (i) = R̂(i) − σ̂12 (i)â1 (i)âH
1 (i)
ŵ(i) =

B. OKSPME-CCG Adaptive Algorithm
In this subsection, the proposed OKSPME-CCG algorithm
is introduced. In CG-based approaches, we usually employ a
forgetting factor (e.g. λ) to estimate the second-order statistics
of the data or the SCM [1], [10], which can be expressed by

R̂−1
(i)â1 (i)
I+N

R̂(i) = λR̂(i − 1) + x(i)xH (i),

−1
âH
1 (i)R̂I+N (i)â1 (i)

End snapshot

whereas the SCV d̂(i) can be estimated with the same forgetting factor as described by

signal power and its steering vector with the same recursions
as in OKSPME, whereas the estimation procedure of the
beamforming weights is different. In particular, we start from a
reformulated optimization problem and use SG and CG-based
adaptive recursions to derive the weight update equations,
which reduce the complexity by an order of magnitude as
compared to that of OKSPME.
A. OKSPME-SG Adaptive Algorithm
We resort to an SG adaptive strategy and consider the
following optimization problem:
minimize
w(i)

subject to

(29)
wH (i)â1 (i) = 1,

∂L
,
w(i + 1) = w(i) − µ
∂w(i)

(30)

where L = wH (i)(x(i)xH (i) − σ̂12 (i)â1 (i)âH
1 (i))w(i) +
λL (wH (i)â1 (i) − 1) and µ is the step size.
By substituting L into the SG equation (30) and letting
wH (i + 1)â1 (i + 1) = 1, λL is obtained as
H
2(σ̂12 (i)âH
1 (i)â1 (i) − y(i)x (i)â1 (i))
.
âH
1 (i)â1 (i)

d̂(i) = λd̂(i − 1) + x(i)y ∗ (i).

(34)

The proposed optimization problem that leads to the
OKSPME-CCG algorithm is described by
minimizeJ = vH (i)(R̂(i) − R̂1 (i))v(i) − âH
1 (i)v(i), (35)
â1 (i),v(i)

where v(i) is the CG-based weight vector. In OKSPME-CCG,
we require N iterations for each snapshot. In the nth iteration,
â1,n (i) and vn (i) are updated as follows
â1,n (i) = â1,n−1 (i) + αâ1 ,n (i)pâ1 ,n (i),

(36)

vn (i) = vn−1 (i) + αv,n (i)pv,n (i),

(37)

where pâ1 ,n (i) and pv,n (i) are direction vectors updated by

wH (i)(R̂(i) − R̂1 (i))w(i)

where R̂(i) can be written as x(i)xH (i) and R̂1 (i) represents
the desired signal covariance matrix and can be written as
σ̂12 (i)â1 (i)âH
1 (i).
Then we can express the SG recursion as

λL =

(33)

(31)

pâ1 ,n+1 (i) = gâ1 ,n (i) + βâ1 ,n (i)pâ1 ,n (i),

(38)

pv,n+1 (i) = gv,n (i) + βv,n (i)pv,n (i),

(39)

where gâ1 ,n (i) and gv,n (i) are the negative gradients of the
cost function in terms of â1 (i) and v(i), respectively, which
are expressed as
gâ1 ,n (i) = −

gv,n (i) = −

∂J
= σ̂12 (i)vn (i)vnH (i)â1,n (i) + vn (i),
∂â1,n (i)
(40)
∂J
∂vn (i)

= gv,n−1 (i) − αv,n (i)(R̂(i) − σ̂12 (i)x(i)xH (i))pv,n (i).
(41)
The scaling parameters αâ1 ,n (i), αv,n (i) can be obtained by
substituting (36) and (37) into (35) and minimizing the cost

6

function with respect to αâ1 ,n (i) and αv,n (i), respectively.
The solutions are given by
αâ1 ,n (i) = −

αv,n (i) =

gâH1 ,n−1 (i)pâ1 ,n (i)
H
σ̂12 (i)pH
â1 ,n (i)vn (i)vn (i)pâ1 ,n (i)

,

(42)

H
gv,n−1
(i)pv,n (i)

pH
v,n (i)(R̂(i)

−

σ̂12 (i)â1,n (i)âH
1,n (i))pv,n (i)

βv,n (i) =

gâH1 ,n (i)gâ1 ,n (i)
gâH1 ,n−1 (i)gâ1 ,n−1 (i)

,

H
gv,n
(i)gv,n (i)
.
H
gv,n−1 (i)gv,n−1 (i)

(44)

(45)

vN (i)
,
H
â1,N (i)vN (i)

(46)

The computational cost of OKSPME-CCG algorithm is
O(N M 2 ), which is higher than the cost required in OKSPMESG due to the inner iterations at every snapshot. The proposed
OKSPME-CCG is summarized in Table III.
C. OKSPME-MCG Adaptive Algorithm
In OKSPME-MCG, we let only one iteration be performed
per snapshot, which further reduces the complexity compared
to OKSPME-CCG. Here we denote the CG-based weights
and steering vector updated by snapshots rather than inner
iterations as
â1 (i) = â1 (i − 1) + αâ1 (i)pâ1 (i),

(47)

v(i) = v(i − 1) + αv (i)pv (i).

(48)

As can be seen, the subscripts of all the quantities for inner
iterations are eliminated. Then, we employ the degenerated
scheme to ensure αâ1 (i) and αv (i) satisfy the convergence
bound [10] given by
H
0 ≤ pH
â1 (i)gâ1 (i) ≤ 0.5pâ1 (i)gâ1 (i − 1),

0≤

pH
v (i)gv (i)

≤

0.5pH
v (i)gv (i

− 1).

(49)
(50)

Instead of updating the negative gradient vectors gâ1 (i) and
gv (i) in iterations, now we utilize the forgetting factor to reexpress them in one snapshot as
gâ1 (i) =(1 − λ)v(i) + λgâ1 (i − 1)
+ σ̂12 (i)αâ1 (i)v(i)vH (i)pâ1 (i)

(51)

− x(i)xH (i)â1 (i),
gv (i) =(1 − λ)â1 (i) + λgv (i − 1)
− αv (i)(R̂(i) − σ̂12 (i)â1 (i)âH
1 (i))pv (i)
H

− x(i)x (i)v(i − 1).

k=1

d̂(i) =

1
i

i
P

x(k)y ∗ (k)

k=1

Step 1 from Table II
Step 2 from Table II
Step 3 from Table II
Steering Vector and Weight Vector Estimations
â1,0 (i) = â1 (i)
gâ1 ,0 (i) = σ̂12 (i)v0 (i)v0H (i)â1,0 (i) + v0 (i)
gv,0 (i) = â1,0 (i) − R̂(i)v0 (i)
pâ1 ,0 (i) = gâ1 ,0 (i); pv,0 (i) = gv,0 (i)
For each iteration index n = 1, 2, · · · , N :
gH

After â1,n (i) and vn (i) are updated for N iterations, the
beamforming weight vector w(i) can be computed by
w(i) =

Initialization:
ŵ(1) = v0 (1) = 1; λ;
Choose an initial guess â1 (0) within the sector and set â1 (1) = â1 (0);
For each snapshot i = 1, 2, · · · :
i
P
R̂(i) = 1i
x(k)xH (k)

.

(43)
The parameters βâ1 ,n (i) and βv,n (i) should be chosen to
provide conjugacy for direction vectors [10], which results in
βâ1 ,n (i) =

TABLE III
P ROPOSED OKSPME-CCG ALGORITHM

(52)

(i)pâ ,n (i)
1
H
n (i)vn (i)pâ1 ,n (i)
â1 ,n
H
gv,n−1
(i)pv,n (i)
2
H
pH
(i)(
R̂(i)−σ̂
v,n
1 (i)â1,n (i)â1,n (i))pv,n (i)

αâ1 ,n (i) = − σ̂ 2 (i)pHâ1 ,n−1
(i)v
1

αv,n (i) =

â1,n (i) = â1,n−1 (i) + αâ1 ,n (i)pâ1 ,n (i)
vn (i) = vn−1 (i) + αv,n (i)pv,n (i)
gâ1 ,n (i) = σ̂12 (i)vn (i)vnH (i)â1,n (i) + vn (i)
gv,n (i) = gv,n−1 (i) − αv,n (i)(R̂(i) − σ̂12 (i)x(i)xH (i))pv,n (i)
βâ1 ,n (i) =
βv,n (i) =

H
gâ
(i)gâ ,n (i)
1
1 ,n
H
gâ ,n−1 (i)gâ ,n−1 (i)
1
1
H
gv,n (i)gv,n (i)
H
gv,n−1
(i)gv,n−1 (i)

pâ1 ,n+1 (i) = gâ1 ,n (i) + βâ1 ,n (i)pâ1 ,n (i)
pv,n+1 (i) = gv,n (i) + βv,n (i)pv,n (i)
End iteration
v0 (i + 1) = vN (i)
N (i)
w(i) = âH v(i)v
(i)
1,N

N

End snapshot

H
Pre-multiplying (51) and (52) by pH
â1 (i) and pv (i), respectively, and taking expectations we obtain
H
H
E[pH
â1 (i)gâ1 (i)] = E[pâ1 (i)(v(i) − x(i)x (i)â1 )(i)]
H
+ λE[pH
â1 (i)gâ1 (i − 1)] − λE[pâ1 (i)v(i)]
2
H
+ E[αâ1 (i)pH
â1 (i)σ̂1 (i)v(i)v (i)pâ1 (i)],

(53)

H
H
E[pH
v (i)gv (i)] = λE[pv (i)gv (i − 1)] − λE[pv (i)â1 (i)]
2
H
− E[αv (i)pH
v (i)(R̂(i) − σ̂1 (i)â1 (i)â1 (i))pv (i)],

(54)

where in (54) we have E[R̂(i)v(i − 1)] = E[â1 (i)]. After
substituting (54) in (50) we obtain the bounds for αv (i) as
follows
H
(λ − 0.5)E[pH
v (i)gv (i − 1)] − λE[pv (i)â1 (i)]

≤E[αv (i)]
2
H
E[pH
v (i)(R̂(i) − σ̂1 (i)â1 (i)â1 (i))pv (i)]
H
λE[pH
v (i)gv (i − 1)] − λE[pv (i)â1 (i)]
≤
. (55)
2
H
E[pH
v (i)(R̂(i) − σ̂1 (i)â1 (i)â1 (i))pv (i)]

Then we can introduce a constant parameter ηv ∈ [0, 0.5]

7

TABLE IV
P ROPOSED OKSPME-MCG ALGORITHM

to restrict αv (i) within the bounds in (55) as
αv (i) =
H
H
λ(pH
v (i)gv (i − 1) − pv (i)â1 (i)) − ηv pv (i)gv (i − 1)
2
H
pH
v (i)(R̂(i) − σ̂1 (i)â1 (i)â1 (i))pv (i)

.

(56)

Similarly, we can also obtain the bounds for αâ1 (i).
For simplicity let us define E[pH
â1 (i)gâ1 (i − 1)] = A,
H
H
E[pH
(i)v(i)]
=
B,
E[p
(i)x(i)x
(i)â1 (i)] = C and
â1
â1
2
H
H
E[pâ1 (i)σ̂1 (i)v(i)v (i)pâ1 (i)] = D. Substituting (53) into
(49) gives
λ(B − A) − B + C
≤E[αâ1 (i)]
D
λ(B − A) − B + C + 0.5A
≤
, (57)
D
in which we can introduce another constant parameter ηâ1 ∈
[0, 0.5] to restrict αâ1 (i) within the bounds in (57) as
E[αâ1 (i)] =

λ(B − A) − B + C + ηâ1 A
,
D

(58)

or
H
H
αâ1 (i) = [λ(pH
â1 (i)v(i) − pâ1 (i)gâ1 (i − 1)) − pâ1 (i)v(i)
H
H
+ pH
â1 (i)x(i)x (i)â1 (i) + ηâ1 pâ1 (i)gâ1 (i − 1)]
2
H
/[σ̂1 (i)pâ1 (i)v(i)vH (i)pâ1 (i)]. (59)

Then we can update the direction vectors pâ1 (i) and pv (i)
by
(60)
pâ1 (i + 1) = gâ1 (i) + βâ1 (i)pâ1 (i),
pv (i + 1) = gv (i) + βv (i)pv (i),

(61)

where βâ1 (i) and βv (i) are updated by
H

βâ1 (i) =
βv (i) =

[gâ1 (i) − gâ1 (i − 1)] gâ1 (i)
,
gâH1 (i − 1)gâ1 (i − 1)

(62)

[gv (i) − gv (i − 1)]H gv (i)
.
gvH (i − 1)gv (i − 1)

(63)

Finally we can update the beamforming weights by
w(i) =

v(i)
,
âH
1 (i)v(i)

(64)

The MCG approach employs the forgetting factor λ and
constant η for estimating α(i), which means its performance
may depend on a suitable choice of these parameters. The
proposed OKSPME-MCG algorithm requires a complexity of
O(M 2 ). However, the cost is usually much lower compared
to CCG approach for the elimination of inner recursions and it
presents a similar performance in most studied scenarios. From
an implementation point of view, the choice of using the CCG
and MCG algorithms is based on the stationarity of the system:
the CCG algorithm is more suitable for scenarios in which
the system is stationary and we can compute the beamformer
with K iterations while the MCG algorithm is suggested for
non-stationary scenarios as we only run one iteration per
snapshot and can track variations in the environment. Table
IV summarizes the OKSPME-MCG algorithm.

Initialization:
ŵ(1) = v(0) = 1; λ; ηv = ηâ1 ;
Choose an initial guess â1 (0) within the sector and set â1 (1) = â1 (0);
gv (0) = pv (1) = â1 (1); gâ1 (0) = pâ1 (1) = v(0);
For each snapshot i = 1, 2, · · · :
i
P
x(k)xH (k)
R̂(i) = 1i
k=1

d̂(i) =

1
i

i
P

x(k)y ∗ (k)

k=1

Step 1 from Table II
Step 2 from Table II
Step 3 from Table II
Steering Vector and Weight Vector Estimations
H
H
αâ1 (i) = [λ(pH
â1 (i)v(i) − pâ1 (i)gâ1 (i − 1)) − pâ1 (i)v(i)
H
H
+pH
(i)x(i)x
(i)â
(i)
+
η
p
(i)g
(i
−
1)]
1
â1 â1
â1
â1
H
/[σ̂12 (i)pH
â1 (i)v(i)v (i)pâ1 (i)]
λ(pH (i)g (i−1)−pH (i)â (i))−η pH (i)g (i−1)

v
v v
v
1
v
v
αv (i) =
2
H
pH
v (i)(R̂(i)−σ̂1 (i)â1 (i)â1 (i))pv (i)
â1 (i) = â1 (i − 1) + αâ1 (i)pâ1 (i)
v(i) = v(i − 1) + αv (i)pv (i)
gâ1 (i) = (1 − λ)v(i) + λgâ1 (i − 1)
+σ̂12 (i)αâ1 (i)v(i)vH (i)pâ1 (i) − x(i)xH (i)â1 (i)
gv (i) = (1 − λ)â1 (i) + λgv (i − 1) − αv (i)(R̂(i)
H
−σ̂12 (i)â1 (i)âH
1 (i))pv (i) − x(i)x (i)v(i − 1)

βâ1 (i) =

[gâ (i)−gâ (i−1)]H gâ (i)
1
1
1
H (i−1)g
gâ
â (i−1)
1

1
H

v (i−1)] gv (i)
βv (i) = [gv (i)−g
H (i−1)g (i−1)
gv
v
pâ1 (i + 1) = gâ1 (i) + βâ1 (i)pâ1 (i)
pv (i + 1) = gv (i) + βv (i)pv (i)
w(i) = âHv(i)
1 (i)v(i)
End snapshot

V. A NALYSIS
In this section, we present an analysis of the following
aspects of the proposed and existing algorithms:
• An analysis of the MSE between the estimated and actual
steering vectors for the general approach that employs a
presumed angular sector.
• MSE analysis results of the proposed OKSPME method
and the SQP method in [6] and their relationships and
differences.
• A complexity analysis for the proposed and existing
algorithms.
A. MSE analysis
Firstly, we present the MSE analysis of the general approach
that employs a presumed angular sector. Since we have the
steering vector estimate â1 (i) in the ith snapshot, by denoting
the true steering vector as a1 , we can express the MSE of the
estimate â1 (i) as
MSE{â1 (i)} = tr(E[(â1 (i) − a1 )(â1 (i) − a1 )H ])
= E[(â1 (i) − a1 )H (â1 (i) − a1 )].

(65)

In the approach that employs an angular sector, we usually
choose an initial guess (i.e., â1 (0)) from the presumed sector.
Let us express the accumulated estimation error as
ê(i) = â1 (i) − â1 (0),

(66)

8

then (65) can be rewritten as
H

MSE{â1 (i)} = E[(â1 (0) + ê(i) − a1 ) (â1 (0) + ê(i) − a1 )].
(67)
The initial guess â1 (0) can be described as the true steering
vector plus a guess error vector (i.e., ǫ):
â1 (0) = a1 + ǫ.

(68)

Taking expectation of both sides of the above, we have
E[â1 (0)] = a1 + E[ǫ].

(69)

Substituting (68) into (67), taking into account that the accumulated estimation error is uncorrelated with the initial guess
error and simplifying the expression, we obtain

equal probability, in other words, the choice of τ is uniformly
distributed within [0, θ]. If the sample size of the selected ǫ
is large enough, we can approximately describe its probability
density function (pdf) as a continuous function given by
1
.
(76)
θ
Meanwhile, we are also able to calculate the chord length kǫk
from a simple geometric criterion as
τ
(77)
kǫk = 2ka1 k sin .
2
Then the expectation of kǫk can be computed by
f(τ ) =

E[kǫk] =
MSE{â1 (i)}
= E[ǫH ǫ] + E[ǫH ]E[ê(i)] + E[êH (i)]E[ǫ] + E[êH (i)ê(i)].
(70)
Furthermore, it should be emphasized that both ǫ and ê(i) are
in vector forms, which means that their second-order statistics
can be re-expressed in terms of their first-order statistics of
their Euclidean norms. Then we can re-express (70) as

2

8ka1 k sin2 θ4
.
(79)
θ
At this point, we can also compute the variance of kǫk by
using (79) as
E[kǫk] =

Var[kǫk] =

(72)

E[kê(i)k2 ] = Var[kê(i)k] + E 2 [kê(i)k].

(73)

At this stage, we can employ Popoviciu’s inequality [37] to
obtain the upper bounds for the variances of the norms of the
random vectors ǫ and ê(i), which are given by
2

(supkǫk − infkǫk)
,
4

(74)

(supkê(i)k − infkê(i)k)2
.
(75)
4
However, the last term in (71) is not analytical when conducting a norm analysis. Actually, E[ǫ] depends on how the
presumed sector is chosen: if the sector is chosen in an
unbiased manner (i.e., the true steering vector lies in the centre
of the sector), then we have E[ǫ] = 0 by symmetry criterion, in
which case we can omit the last terms of (71). For convenience
of carrying out the norm analysis as the next step, we focus
on the unbiased case only, so that the MSE only depends on
the expectation, the infimum and the supremum of kǫk and
kê(i)k. In Fig. 1, we utilize Euclidean geometry to illustrate
the relationships among the norms of the errors and the norm
of the steering vector, which is a fixed parameter due to the
re-normalization procedure after it is estimated each time.
According to Fig. 1, we can use θ (i.e., half of the angular
sector, assumed less than π/4) and ka1 k to obtain E[kǫk] by
the following (any angular parameter appeared in the equations
should be measured in radians rather than degrees): kǫk is
equivalent to the chord length which corresponds to the arc
of a variable τ , which can be any value from 0 to θ with
Var[kê(i)k] ≤

Zθ

(kǫk − E[kǫk])2 f(τ )dτ,

(80)

0

from which after simplification we obtain

E[kǫk ] = Var[kǫk] + E [kǫk],

Var[kǫk] ≤

(78)

from which after simplification we obtain

= E[kǫk2 ] + E[kê(i)k2 ] + 2E[ǫH ]E[ê(i)]. (71)

2

kǫkf(τ )dτ,

0

MSE{â1 (i)}
Since both kǫk and kê(i)k are scalars we have

Zθ

sin θ 32 sin4 4θ
).
(81)
−
θ
θ2
In addition, it is clear that we have infkǫk = 0 and supkǫk =
2ka1 k sin θ2 , which can be substituted in (74) and result in
Var[kǫk] = 2ka1 k2 (1 −

θ
(82)
Var[kǫk] ≤ ka1 k2 sin2 .
2
We can see that the right-hand side of (81) satisfies the
inequality in (82). After substituting (79) and (81) in (72),
we obtain
sin θ
E[kǫk2 ] = 2ka1 k2 (1 −
).
(83)
θ

a1

â1(i)

ǫ

â1(0)

ê(i)

||a1||
θi
θ

Fig. 1.

Euclidean norm interpretation of the MSE

9

...

ê(i) − ê(i − 1)

and the norm of the accumulated estimation error is no greater
than the sum of the norms of all the piecewise estimation
errors, which is given by the inequality

ê(1)
â1(1)

â1(i − 1)

kê(i)k ≤

kê(i)kmax = ka1 k tan θi .
Fig. 2.

... â1(1)
â1(i − 1)

ê(1)

θi

Fig. 3.

lim kê(i)k = θi ka1 k,

â1(0)

â1(i)

i→∞

||a1||

Regarding the computation of the norm of the accumulated
estimation error kê(i)k, we need to emphasize that even
though the steering vector is always re-normalized each time
after it is updated, the piecewise estimation error in each
snapshot does not directly update the steering vector to its normalized version, which means it is inappropriate to calculate
the estimation error by geometric methods directly from Fig. 1
because the accumulated estimation error partially comes from
the unnormalized steering vectors. However, we can obtain
the infimum and supremum values for kê(i)k if we assume
the update scheme is unidirectional (i.e., the steering vector is
updated from â1 (0) to â1 (i) in a single direction within the
sector), with the unnormalized steering vectors considered.
We firstly look at the SQP method scenario in [6]. The
steering vector update scheme is shown in Fig. 2. It is
necessary to emphasize that now we focus on the angular
sector range of θi (i.e., the angle difference between the
initially guessed steering vector and its estimate in the ith
snapshot) rather than θ. In [6], an online-optimization program
was used to iteratively solve for the piecewise estimation error
in every snapshot, which was always orthogonal to the current
steering vector estimate. Let us consider that at each time
instant the steering vector is updated, its direction changes
by θi,k , where i is the snapshot index and k (1 ≤ k ≤ i) is
the index for the kth update. Since the total direction change
in a snapshot is θi , then we have
i
X

k=1

θi,k ,

(87)

which is actually the infimum of kê(i)k and cannot be
achieved since the number of snapshots or iterations are always
limited in practical situations. By combining (86) and (87),
kê(i)k is bounded by

update scheme of the OKSPME method

θi =

(86)

On the other hand, we notice that the piecewise estimation
error vector can never enter into the angular sector, but at most
move along with the arc if the number of iterations is large
enough. In this case, we can approximately and geometrically
illustrate the arc length corresponding with θi as the lower
bound by taking the limit i → ∞, i.e.,

update scheme of the SQP method

ê(i) − ê(i − 1)

(85)

If we assume θi is less than π/2, then the right-hand side of
(85) achieves its maximum value when θi,k = tan θi , which
is also the supremum of kê(i)k and equals

||a1||

θi

ka1 k tan θi,k .

k=1

â1(0)

â1(i)

i
X

(84)

infkê(i)k = θi ka1 k < kê(i)k ≤ ka1 k tan θi = supkê(i)k.
(88)
Different from the SQP method, the proposed OKSPME
method utilizes the Krylov subspace and the cross-correlation
vector projection approach to extract the error information then
use it to update the steering vector. From (9) we have
i

i

k=1
i

k=1

1X
1X
x(k)y ∗ (k) =
x(k)(wH (k)x(k))∗
d̂(i) =
i
i
i

1X
1X
=
x(k)xH (k)w(k) =
R̂(k)w(k).
i
i
k=1

(89)

k=1

Note that an initialization for vector d̂ or matrix R̂ should be
considered to ensure R̂ is full-rank and invertible, which can
be done by either setting d̂(0) = δIw(0) or R̂(0) = δI. We
also know that
w(k) =

R̂−1 (k)â1 (k)
H
â1 (k)R̂−1 (k)â1 (k)

=

R̂−1 (k)â1 (k)
.
σ̂12 (k)

(90)

Pre-multiplying (90) by R̂(k) on both sides we obtain
R̂(k)w(k) =

â1 (k)
,
σ̂12 (k)

(91)

which is then substituted in (89) and results in
i

d̂(i) =

1 X â1 (k)
,
i
σ̂12 (k)

(92)

k=1

where σ̂12 (k) is a scalar, which means the SCV contains the
direction of the desired signal steering vector. Projecting d̂(i)
onto the Krylov subspace represented by P̂(i) is therefore
similar to projecting â1 (i). In our method, the estimation
of d̂(i) is separate from the update of â1 (i), which means

10

the steering vector estimation error used for the updates is
obtained from d̂(i), so that in the kth (1≤k < i) snapshot, the
error does not have to be orthogonal to â1 (k), but should
be orthogonal to another potentially better estimate â1 (j)
(1≤k < j≤i), resulting in a situation where the error is located
inside the sector (see Fig. 3). There are two benefits in the case
which the error is inside the sector: faster convergence rate
and smaller estimation error. We can obtain the infimum and
supremum values in a similar way. By applying the inequality
that the norm of the accumulated estimation error is no greater
than the sum of the norms of all the piecewise estimation
errors, we have
kê(i)k ≤

i
X

ka1 k sin θi,k ,

(93)

k=1

where the parameters θi,k (k = 1, 2, · · · , i) satisfy the constraint in (84). However, the right-hand side of (93) achieves
its maximum value when all these parameters are equal (i.e.,
θi,1 = θi,2 = · · · = θi,i = θii ) and it is given by
θi
,
(94)
i
The right-hand side of (94) can be treated as a function of i
which is an increasing function on i = 1, 2, · · · , ∞. Therefore,
we can take the limit of it to obtain the upper bound of
kê(i)kmax , and so as to kê(i)k. In fact, when i → ∞, the
piecewise estimation error moves along the arc corresponding
with θi , resulting in the upper bound obtained is the same as
the lower bound of the SQP method case, which is given by the
right-hand side expression of (87) and defines the supremum
of kê(i)k in this case. Since we have already assumed that θi is
less than π/2 so that ê(i) must be inside of the angular sector
but its Euclidean norm cannot be smaller than the orthogonal
distance between â1 (0) to â1 (i), so this orthogonal distance
can define the lower bound of kê(i)k, which is actually the
infimum and calculated by ka1 k sin θi . Then, in the OKSPME
method, kê(i)k is bounded by
kê(i)kmax = ika1 k sin

infkê(i)k = ka1 k sin θi ≤ kê(i)k < θi ka1 k = supkê(i)k.
(95)
By taking expectations of both (88) and (95), we obtain
E[θi ]ka1 k < {E[kê(i)k]}SQP ≤ ka1 k tan(E[θi ]),

(96)

ka1 k sin(E[θi ]) ≤ {E[kê(i)k]}OKSP ME < E[θi ]ka1 k. (97)
On the other side, by substituting (88) and (95) in (75), we
obtain
0 ≤ {Var[kê(i)k]}SQP ≤

ka1 k2 (tan θi − θi )2
,
4

(98)

ka1 k2 (θi − sin θi )2
. (99)
4
Substituting (96), (98) and (97), (99) in (73), respectively, we
have
0 ≤ {Var[kê(i)k]}OKSP ME ≤

2

E 2 [θi ]ka1 k2 < {E[kê(i)k ]}SQP
ka1 k2 (tan θi − θi )2
+ ka1 k2 tan2 (E[θi ]), (100)
≤
4

ka1 k2 sin2 (E[θi ]) ≤ {E[kê(i)k2 ]}OKSP ME
ka1 k2 (θi − sin θi )2
+ E 2 [θi ]ka1 k2 . (101)
<
4
However, E[θi ] also has its lower and upper bounds. Since
our analysis focuses on the unbiased case only as mentioned,
the true steering vector is located in the center of the angular
sector and the estimate â1 (i) is always closer to the center
than â1 (0). Let us assume that even if the estimate â1 (i)
always happens to be very close to either edge of the sector,
no matter how â1 (0) is chosen within the sector, θi will vary
from 0 to 2θ with equal probability, or equivalently, uniformly
distributed within [0, 2θ), in which case we can obtain the
upper bound for E[θi ] by taking the average between 0 to
2θ, which is obtained as θ. On the other hand, if we assume
that the estimate â1 (i) always happens to be exactly at the
center of the sector, resulting in that θi can only vary from
0 to θ, or uniformly distributed within [0, θ] in which case
E[θi ] = θ/2, resulting in the lower bound of E[θi ] is θ/2.
Therefore, the upper and lower bounds for MSE{â1 (i)} can be
further obtained by substituting E[θi ]max → θ, [θi ]max → 2θ
and E[θi ]min = θ/2, [θi ]min = 0 into the upper and lower
bounds of (100) and (101) respectively, resulting in
θ2
2
ka1 k2 < {E[kê(i)k ]}SQP
4
ka1 k2 (tan 2θ − 2θ)2
+ ka1 k2 tan2 θ, (102)
<
4

ka1 k2 sin2

θ
≤ {E[kê(i)k2 ]}OKSP ME
2
ka1 k2 (2θ − sin 2θ)2
<
+ θ2 ka1 k2 . (103)
4

Finally, by combining the expectation of the mean-squared
initial guess error E[kǫk2 ] in (83) with (102) and (103), we
obtain the bounds for the MSE of the steering vector estimate
MSE{â1 (i)} as
2 sin θ θ2
+ )ka1 k2 < {MSE{â1 (i)}}SQP
θ
4
2 sin θ (tan 2θ − 2θ)2
+
+ tan2 θ)ka1 k2 , (104)
< (2 −
θ
4

(2 −

(2 −

2 sin θ
θ
+ sin2 )ka1 k2 ≤ {MSE{â1 (i)}}OKSP ME
θ
2
2 sin θ (2θ − sin 2θ)2
< (2 −
+
+ θ2 )ka1 k2 . (105)
θ
4

From (104) and (105), we can see that the MSEs now only
depend on two parameters: the norm of the true steering vector
and the angular sector spread. The lower and upper bounds of
the proposed OKSPME method are lower than those of the
SQP method. As mentioned before, it is important that the
presumed angular sector spread 2θ must be less than π/2 (i.e.,
90◦ ) to ensure the previous assumption of θi < π/2 is always
valid.

11

TABLE V
C OMPLEXITY C OMPARISON

OKSPME
OKSPME-SG
OKSPME-CCG
OKSPME-MCG

Flops
4M 3 + 3M 2 + 20M
2M 3 + 11M 2
O(M 3.5 )
3
2M + 4M 2 + 5M
2nM 2 + 7nM
M 3 + (4m + 11)M 2
+(3m2 + 5m + 20)M
(4m + 7)M 2
+(3m2 + 5m + 33)M
(4m + 8n + 8)M 2
+(3m2 + 5m + 33n + 29)M
(4m + 14)M 2
+(3m2 + 5m + 86)M

7

10

6

10

flops

RAB Algorithms
LOCSME [32]
RCB [4]
SQP [6]
LOCME [7]
LCWC [9]

8

10

5

10

LOCSME [32]
RCB [4]
SQP [6]
LOCME [7]
LCWC [9]
OKSPME
OKSPME−SG
OKSPME−CCG
OKSPME−MCG

4

10

3

10

2

10

0

10

20

30

40

50

60

70

80

90

100

M

B. Complexity Analysis
The computational complexity analysis is discussed in this
subsection. We measure the total number of additions and multiplications (i.e., flops) in terms of the number of sensors M
performed for each snapshot for the proposed algorithms and
the existing ones and list them in Table V. Note that the SQP
method in [6] has a highly-variant computational complexity
in different snapshots, due to the online optimization program
based on random choices of the presumed steering vector.
However, it is usually in O(M 3.5 ). The complexity of the
LCWC algorithm in [9] often requires a much larger n than
that in the proposed LOCSME-CCG algorithm. It is obvious
that all of the proposed algorithms have their complexity
depending on the Krylov subspace model order m, which is
determined from Table I and is no larger than K + 1. For
the convenience of comparison, we eliminate all parameters
except M by setting them to common values (the values of n
in LCWC and OKSPME-CCG is set to 50 and 5 respectively,
m = K + 1 where K = 3) and illustrate their complexity
with M varying from 10 to 100 as shown in Fig. 4. As can
be seen that the proposed OKSPME-SG and OKSPME-MCG
algorithms have lower complexity than the other algorithms.
VI. S IMULATIONS
In this section, we present and discuss the simulation results
of the proposed RAB algorithms by comparing them to some
of the existing RAB techniques. We consider uniform linear
arrays (ULA) of omnidirectional sensors with half wavelength
spacing. To produce all the figures (if unspecified in a few
scenario), 100 repetitions are executed to obtain each point of
the curves and a maximum of i = 300 snapshots are observed.
The desired signal is assumed to arrive at θ1 = 10◦ . The
signal-to-interference ratio (SIR) is fixed at 0dB. As the prior
knowledge, the angular sector in which the desired signal is
assumed to be located is chosen as [θ1 −5◦ , θ1 +5◦ ]. The results
focus on the beamformer output SINR performance versus the
number of snapshots, or a variation of input SNR (−10dB
to 30dB) and both coherent and incoherent local scattering
mismatch [5] scenarios are considered.

Fig. 4.

Complexity Comparison

A. Mismatch due to Coherent Local Scattering
All simulations in this subsection consider coherent local
scattering. With time-invariant coherent local scattering, the
steering vector of the desired signal is modeled as
a1 = p +

4
X

ejϕk b(θk ),

(106)

k=1

where p corresponds to the direct path while b(θk )(k =
1, 2, 3, 4) corresponds to the scattered paths. The angles
θk (k = 1, 2, 3, 4) are randomly and independently drawn in
each simulation run from a uniform generator with mean 10◦
and standard deviation 2◦ . The angles ϕk (k = 1, 2, 3, 4) are
independently and uniformly taken from the interval [0, 2π] in
each simulation run. Both θk and ϕk change from trials while
remaining constant over snapshots.
We firstly compare our proposed methods with some classical RAB methods (i.e., standard diagonal loading method with
a fixed loading factor equal to 10 times the noise variance,
the RCB method in [4] which estimates the loading factor
iteratively, and the method that solves an online quadratic
optimization programming, which refers to the SQP method
[6]). The numbers of sensors and signal sources (including the
desired signal) are set to M = 10 and K = 3, respectively.
For this case only, we set the INR (interferences-to-noise
ratio) to 20dB and illustrate the SINR performance versus
snapshots within 100 snapshots in Fig. 5. The two interferers
are arranged to be in the directions of θ2 = 30◦ and θ3 = 50◦ ,
respectively. The other user-defined parameters, if unspecified,
(e.g. the step size µ and the forgetting factor λ) are manually
optimized to give the best algorithm performance, which is
also applied for the other simulation scenarios.
We then set the number of sensors to M = 12, the number
of signal sources as (including the desired signal) K = 3
and illustrate the SINR versus snapshots and the SINR versus
input SNR performance in Fig. 6 and Fig. 7 respectively.

12

15

15

10

10

5

SINR (dB)

SINR (dB)

5

0
Optimum SINR
fixed diagonal loading
RCB [4]
SQP [6]
OKSPME
OKSPME−SG
OKSPME−CCG
OKSPME−MCG

−5

−10

−15

0

20

40

60

0

−5
Optimum SINR
LOCSME [32]
LCWC [9]
SQP [6]
OKSPME
OKSPME−SG
OKSPME−CCG
OKSPME−MCG

−10

80

100

snapshots
−15

0

50

100

150

200

250

300

snapshots

Fig. 5. Coherent local scattering, SINR versus snapshots, M = 10, K = 3,
INR = 20dB
Fig. 6.

Coherent local scattering, SINR versus snapshots, M = 12, K = 3

TABLE VI
C HANGES OF I NTERFERERS

0 − 150

Number of Interferers
(K − 1)
2

151 − 300

5

Snapshots

DoAs
◦

30

◦

θ2 = 30 , θ3 = 50 .
θ2 = 20◦ , θ3 = 30◦ , θ4 = 40◦ ,
θ5 = 50◦ , θ6 = 60◦ .

25

20

The two interferers are arranged to be in the directions of
θ2 = 30◦ and θ3 = 50◦ , respectively. In either Fig. 6
or Fig. 7, we can see that the proposed OKSPME method
has a very similar or slightly better performance compared
to the LOCSME algorithm of [32] and both of them have
the best performance. Furthermore, the proposed OKSPMECCG and OKSPME-MCG algorithms also achieve very close
performance to OKSPME.
In Fig. 8, we assess the SINR performance versus snapshots
of those selected algorithms in a specific time-varying scenario
which encounters a halfway redistribution of the interferers at
a certain snapshot. In this case, the number of sensors is kept
at M = 12, whereas the details of the interferers are given in
Table VI.
In Figs. 9 and 10, we set the number of signal sources to
K = 3, but increase the number of sensors from M = 12 to
M = 40 and study the SINR versus snapshots and the SINR
versus input SNR performance of the selected and proposed
dimensionality reduction RAB algorithms, respectively. We
set the reduced-dimension as D = 4 for the beamspace
based algorithm [23] in all simulations. This time, it is clear
that the proposed OKSPME, OKSPME-SG, OKSPME-CCG
and OKSPME-MCG algorithms all have a certain level of
performance degradation compared to the scenario where
M = 12. The proposed OKSPME based algorithms achieve
better performances than the beamspace approach.
B. Mismatch due to Incoherent Local Scattering
In this case, the desired signal affected by incoherent local
scattering has a time-varying signature and its steering vector

SINR (dB)

15

10

5

0
Optimum SINR
LOCSME [32]
LCWC [9]
SQP [6]
OKSPME
OKSPME−SG
OKSPME−CCG
OKSPME−MCG

−5

−10

−15

−5

0

5

10

15

20

25

30

SNR (dB)

Fig. 7.

Coherent local scattering, SINR versus SNR, M = 12, K = 3

is modeled by
a1 (i) = s0 (i)p +

4
X

sk (i)b(θk ),

(107)

k=1

where sk (i)(k = 0, 1, 2, 3, 4) are i.i.d zero mean complex
Gaussian random variables independently drawn from a random generator. The angles θk (k = 0, 1, 2, 3, 4) are drawn
independently in each simulation run from a uniform generator
with mean 10◦ and standard deviation 2◦ . At this time, sk (i)
changes both from run to run and from snapshot to snapshot.
In order to show the effects caused by incoherent scattering
only, we set the parameters M = 40 and K = 3, study the
SINR versus SNR performance of the selected algorithms in
Fig. 11 and compare the results with Fig. 10. As a result,
a performance degradation is observed for all the studied
algorithms. This is because the time-varying nature of incoherent scattering results in more dynamic and environmental

13

15

30

25
10

20

15

SINR (dB)

SINR (dB)

5

0

−5

Optimum SINR
LOCSME [32]
LCWC [9]
SQP [6]
OKSPME
OKSPME−SG
OKSPME−CCG
OKSPME−MCG

−10

−15

0

50

100

150

200

250

10

5

0

−5
Optimum SINR
OKSPME
OKSPME−SG
OKSPME−CCG
OKSPME−MCG
Beamspace [23]

−10

−15

−20
−10

300

−5

0

5

snapshots

Fig. 8.

10

15

20

25

30

SNR (dB)

Coherent local scattering, SINR versus snapshots, M = 12

Fig. 10.

15

Coherent local scattering, SINR versus SNR, M = 40, K = 3

30

25
10

20

15

SINR (dB)

SINR (dB)

5

0

10

5

0

−5

−5
Optimum SINR
OKSPME
OKSPME−SG
OKSPME−CCG
OKSPME−MCG
Beamspace [23]

−10

−15

0

50

100

150

200

250

−15

−20
−10

300

snapshots

Fig. 9.

Coherent local scattering, SINR versus snapshots, M = 40, K = 3

uncertainties in the system, which increases the steering vector
mismatch.
VII. C ONCLUSION
We have proposed the OKSPME algorithm based on the
exploitation of cross-correlation mismatch estimation and the
orthogonal Krylov subspace. In addition, low complexity RAB
algorithms, OKSPME-SG, OKSPME-CCG and OKSPMEMCG have been developed to enable the beamforming weights
to be updated recursively without matrix inversions. A detailed
steering vector estimation MSE analysis for the general RAB
design approach that relies on a presumed angular sector as
prior knowledge has been provided. The computational complexity of the proposed and some of the existing algorithms
have been compared and discussed. Simulation results have
shown that the proposed algorithms have robustness against

Optimum SINR
OKSPME
OKSPME−SG
OKSPME−CCG
OKSPME−MCG
Beamspace [23]

−10

−5

0

5

10

15

20

25

30

SNR (dB)

Fig. 11.

Incoherent local scattering, SINR versus SNR, M = 40, K = 3

different choices of user-defined parameters and environmental effects, and achieved excellent output SINR performance
especially in medium-high input SNR values.
R EFERENCES
[1] H. L. Van Trees, Optimum Array Processing, New York: Wiley, 2002.
[2] S. A. Vorobyov, A. B. Gershman and Z. Luo, “Robust Adaptive Beamforming using Worst-Case Performance Optimization: A Solution to the
Signal Mismatch Problem,” IEEE Transactions on Signal Processing,
Vol. 51, No. 4, pp 313-324, Feb 2003.
[3] A. Khabbazibasmenj, S. A. Vorobyov and A. Hassanien, “Robust Adaptive Beamforming Based on Steering Vector Estimation with as Little as
Possible Prior Information,” IEEE Transactions on Signal Processing,
Vol. 60, No. 6, pp 2974-2987, June 2012.
[4] J. Li, P. Stoica and Z. Wang, “On Robust Capon Beamforming and
Diagonal Loading,” IEEE Transactions on Signal Processing, Vol. 57,
No. 7, pp 1702-1715, July 2003.
[5] D. Astely and B. Ottersten, “The effects of Local Scattering on Direction of Arrival Estimation with Music,” IEEE Transactions on Signal
Processing, Vol. 47, No. 12, pp 3220-3234, Dec 1999.

14

[6] A. Hassanien, S. A. Vorobyov and K. M. Wong, “Robust Adaptive
Beamforming Using Sequential Quadratic Programming: An Iterative
Solution to the Mismatch Problem,” IEEE Sig. Proc. Letters., Vol. 15,
pp 733-736, 2008.
[7] L. Landau, R. de Lamare, M. Haardt, “Robust Adaptive Beamforming
Algorithms Using Low-Complexity Mismatch Estimation,” Proc. IEEE
Statistical Signal Processing Workshop, 2011.
[8] J. Zhuang and A. Manikas, “Interference cancellation beamforming
robust to pointing errors,” IET Signal Process., Vol. 7, No. 2, pp. 120127, April 2013.
[9] A. Elnashar, “Efficient implementation of robust adaptive beamforming
based on worst-case performance optimization,” IET Signal Process.,
Vol. 2, No. 4, pp. 381-393, Dec 2008.
[10] L. Wang and R. C. de Lamare, “Constrained adaptive filtering algorithms
based on conjugate gradient techniques for beamforming,” IET Signal
Process., Vol. 4, No. 6, pp. 686-697, Feb 2010.
[11] J. P. Lie, W. Ser and C. M. S. See, “Adaptive Uncertainty Based
Iterative Robust Capon Beamformer Using Steering Vector Mismatch
Estimation,” IEEE Transactions on Signal Processing, Vol. 59, No. 9,
Sep 2011.
[12] A. Filimon, “Krylov Subspace Iteration Methods,” May 2008.
[13] H. A. van der Vorst, “Iterative Krylov Methods for Large Linear
Systems,” Cambridge University Press, 2003.
[14] M. Morelli, L. Sangunietti and U. Mengali, “Channel Estimation for
Adaptive Frequency-Domain Equalization,” IEEE Transactions on Wireless Communications, Vol. 4, No. 5, pp. 53-57, Sep 2005.
[15] Z. Bai, “Krylov subspace techniques for reduced-order modeling of
large-scale dynamical systems,” Applied Numerical Mathematics, Department of Computer Science, University of California, 2002.
[16] R. C. de Lamare, M. Haardt and R. Sampaio-Neto, “Blind Adaptive
Constrained Reduced-Rank Parameter Estimation Based on Constant
Modulus Design for CDMA Interference Suppression,” IEEE Transactions on Signal Processing, Vol. 56, No. 6, June 2008.
[17] H. Ge, I. P. Kirsteins and L. L. Scharf, “Data Dimension Reduction
Using Krylov Subspaces Making Adaptive Beamformers Robust to
Model Order-Determination,” in Int. Conf. Acoustics, Speech, and Signal
Processing (ICASSP), Vol. 4, 2006.
[18] I. P. Kirsteins and H. Ge, “Performance Analysis of Krylov Space
Adaptive Beamformers,” Sensor Array and Multichannel Processing,
Fourth IEEE Workshop, pp. 16-20, July 2006.
[19] S. D. Somasundaram, N. H. Parsons, P. Li and R. C. de Lamare, “DataAdaptive Reduced-Dimension Robust Capon Beamforming,” Int. Conf.
Acoustics, Speech, and Signal Processing (ICASSP), pp. 4159-4163,
May 2013.
[20] S. D. Somasundaram, N. H. Parsons, P. Li and R. C. de Lamare,
“Reduced-Dimension Robust Capon Beamforming Using KrylovSubspace Techniques,” IEEE Transactions on Aerospace and Electronic
Systems, Vol. 51, No. 1, pp. 270-289, Jan 2015.
[21] V. Druskin and V. Simoncini, “Adaptive Rational Krylov Subspaces for
Large-Scale Dylymical Systems,” Systems and Control Letters, Vol. 60,
No. 8, pp. 546-560, 2011.
[22] K. Carlberg and C. Farhat, “An Adaptive POD-Krylov Reduced-Order
Model for Structural Optimization,” 8th World Congress on Structural
and Multidisciplinary Optimization, June 2009.
[23] A. Hassanien and S. A. Vorobyov, “A Robust Adaptive Dimension
Reduction Technique With Application to Array Processing,” IEEE Sig.
Proc. Letters. Vol. 16, No. 1, pp. 22-25, Jan 2009.
[24] M. Yukawa, R. C. de Lamare and I. Yamada, “Robust ReducedRank Adaptive Algorithm Based on Parallel Subgradient Projection and
Krylov Subspace,” IEEE Transactions on Signal Processing, Vol. 57,
No. 12, Dec 2009.
[25] C. Dumard and T. Zemen, “Double Krylov Subspace Approximation for
Low Complexity Iterative Multi-User Decoding and Time-Variant Channel Estimation,” 6th IEEE Workshop on Signal Processing Advances in
Wireless Communications, pp. 328-332, June 2005.
[26] R. C. de Lamare and R. Sampaio-Neto, “Reduced-Rank Adaptive
Filtering Based on Joint Iterative Optimization of Adaptive Filters,”
IEEE Sig. Proc. Letters. Vol. 14, No. 12, Dec 2007.
[27] N. Song, W. U. Alokozai, R. C. de Lamare and M. Haardt, “Adaptive
Widely Linear Reduced-Rank Beamforming Based on Joint Iterative
Optimization,” IEEE Sig. Proc. Letters. Vol. 21, No. 3, Mar 2014.
[28] S. D. Somasundaram, “A Framework for Reduced Dimension Robust
Capon Beamforming,” IEEE Statistical Signal Processing Workshop
(SSP), 2011.
[29] P. Li and R. C. de Lamare, “Low-Complexity Robust Data-Dependant
Dimensionality Reduction Based on Joint Iterative Optimization of Pa-

[30]

[31]

[32]

[33]

[34]

[35]

[36]

[37]
[38]

[39]

[40]

[41]

[42]

[43]

[44]

[45]

[46]

[47]

[48]

[49]

[50]

[51]

rameters,” IEEE 5th International Workshop on Computational Advances
in Multi-Sensor Adaptive Processing (CAMSAP), pp. 49-52, Dec 2013.
R. C. de Lamare and R. Sampaio-Neto, “Adaptive Reduced-Rank
Processing Based on Joint and Iterative Interpolation, Decimation, and
Filtering,” IEEE Transactions on Signal Processing, Vol. 57, No. 7, July
2009.
L. Wang, R. C. de Lamare and M. Haardt, “Direction finding Algorithms
with Joint Iterative Subspace Optimization,” IEEE Transactions on
Aerospace and Electronic Systems, Vol. 50, No. 4, pp. 2541-2553, Oct
2014.
H. Ruan and R. C. de Lamare, “Robust Adaptive Beamforming Using
a Low-Complexity Shrinkage-Based Mismatch Estimation Algorithm,”
IEEE Sig. Proc. Letters., Vol. 21, No. 1, pp. 60-64, Nov 2013.
C. Zhang, R. Tian, X. Tan, “A Leakage-Based Dimensionality Reduction
Beamforming for MIMO Cognitive Radio Networks,” 8th International
Conference on Communications and Networking in China (CHINACOM), 2013.
R. Fa, R. C. de Lamare and D. Zanatta-Filho, “Reduced-Rank STAP
Algorithm for Adaptive Radar Based on Joint Iterative Optimization of
Adaptive Filters,” 42nd Asilomar Conference on Signals, Systems and
Computers, pp. 533-537, Oct 2008.
R. Fa, R. C. de Lamare and P. Clarke, “Reduced-Rank STAP for
MIMO Radar Based on Joint Iterative Optimization of KnowledgeAided Adaptive Filters,” Conference Record of the Forty-Third Asilomar
Conference on Signals, Systems and Computers, pp. 496-500, Nov 2009.
S.D. Somasundaram, “Reduced dimension robust Capon beamforming
for large aperture passive sonar arrays,” IET Radar, Sonar and Navigation, Vol. 5, No. 7, pp. 707-715, 2011.
Popoviciu, “Sur Les Equations Algebriques Ayant Toutes Leurs Racines
Relles,” Mathematica 9, pp. 129-145, 1935.
L. L. Scharf and D. W. Tufts, “Rank reduction for modeling stationary
signals,” IEEE Transactions on Acoustics, Speech and Signal Processing,
vol. ASSP-35, pp. 350-355, March 1987.
A. M. Haimovich and Y. Bar-Ness, “An eigenanalysis interference
canceler,” IEEE Trans. on Signal Processing, vol. 39, pp. 76-84, Jan.
1991.
D. A. Pados and S. N. Batalama ”Joint space-time auxiliary vector filtering for DS/CDMA systems with antenna arrays” IEEE Transactions
on Communications, vol. 47, no. 9, pp. 1406 - 1415, 1999.
J. S. Goldstein, I. S. Reed and L. L. Scharf ”A multistage representation
of the Wiener filter based on orthogonal projections” IEEE Transactions
on Information Theory, vol. 44, no. 7, 1998.
Y. Hua, M. Nikpour and P. Stoica, ”Optimal reduced rank estimation
and filtering,” IEEE Transactions on Signal Processing, pp. 457-469,
Vol. 49, No. 3, March 2001.
M. L. Honig and J. S. Goldstein, “Adaptive reduced-rank interference
suppression based on the multistage Wiener filter,” IEEE Transactions
on Communications, vol. 50, no. 6, June 2002.
E. L. Santos and M. D. Zoltowski, “On Low Rank MVDR Beamforming
using the Conjugate Gradient Algorithm”, Proc. IEEE International
Conference on Acoustics, Speech and Signal Processing, 2004.
Q. Haoli and S.N. Batalama, “Data record-based criteria for the selection
of an auxiliary vector estimator of the MMSE/MVDR filter”, IEEE
Transactions on Communications, vol. 51, no. 10, Oct. 2003, pp. 1700
- 1708.
R. C. de Lamare and R. Sampaio-Neto, “Reduced-Rank Adaptive
Filtering Based on Joint Iterative Optimization of Adaptive Filters”,
IEEE Signal Processing Letters, Vol. 14, no. 12, December 2007.
Z. Xu and M.K. Tsatsanis, “Blind adaptive algorithms for minimum
variance CDMA receivers,” IEEE Trans. Communications, vol. 49, No.
1, January 2001.
R. C. de Lamare and R. Sampaio-Neto, “Low-Complexity Variable
Step-Size Mechanisms for Stochastic Gradient Algorithms in Minimum
Variance CDMA Receivers”, IEEE Trans. Signal Processing, vol. 54,
pp. 2302 - 2317, June 2006.
C. Xu, G. Feng and K. S. Kwak, “A Modified Constrained Constant
Modulus Approach to Blind Adaptive Multiuser Detection,” IEEE
Trans. Communications, vol. 49, No. 9, 2001.
Z. Xu and P. Liu, “Code-Constrained Blind Detection of CDMA Signals
in Multipath Channels,” IEEE Sig. Proc. Letters, vol. 9, No. 12,
December 2002.
R. C. de Lamare and R. Sampaio Neto, ”Blind Adaptive CodeConstrained Constant Modulus Algorithms for CDMA Interference
Suppression in Multipath Channels”, IEEE Communications Letters,
vol 9. no. 4, April, 2005.

15

[52] L. Landau, R. C. de Lamare and M. Haardt, “Robust adaptive beamforming algorithms using the constrained constant modulus criterion,”
IET Signal Processing, vol.8, no.5, pp.447-457, July 2014.
[53] R. C. de Lamare, “Adaptive Reduced-Rank LCMV Beamforming Algorithms Based on Joint Iterative Optimisation of Filters”, Electronics
Letters, vol. 44, no. 9, 2008.
[54] R. C. de Lamare and R. Sampaio-Neto, “Adaptive Reduced-Rank
Processing Based on Joint and Iterative Interpolation, Decimation and
Filtering”, IEEE Transactions on Signal Processing, vol. 57, no. 7, July
2009, pp. 2503 - 2514.
[55] R. C. de Lamare and Raimundo Sampaio-Neto, “Reduced-rank Interference Suppression for DS-CDMA based on Interpolated FIR Filters”,
IEEE Communications Letters, vol. 9, no. 3, March 2005.
[56] R. C. de Lamare and R. Sampaio-Neto, “Adaptive Reduced-Rank MMSE
Filtering with Interpolated FIR Filters and Adaptive Interpolators”, IEEE
Signal Processing Letters, vol. 12, no. 3, March, 2005.
[57] R. C. de Lamare and R. Sampaio-Neto, “Adaptive Interference Suppression for DS-CDMA Systems based on Interpolated FIR Filters with
Adaptive Interpolators in Multipath Channels”, IEEE Trans. Vehicular
Technology, Vol. 56, no. 6, September 2007.
[58] R. C. de Lamare, “Adaptive Reduced-Rank LCMV Beamforming Algorithms Based on Joint Iterative Optimisation of Filters,” Electronics
Letters, 2008.
[59] R. C. de Lamare and R. Sampaio-Neto, “Reduced-rank adaptive filtering
based on joint iterative optimization of adaptive filters”, IEEE Signal
Process. Lett., vol. 14, no. 12, pp. 980-983, Dec. 2007.
[60] R. C. de Lamare, M. Haardt, and R. Sampaio-Neto, “Blind Adaptive
Constrained Reduced-Rank Parameter Estimation based on Constant
Modulus Design for CDMA Interference Suppression”, IEEE Transactions on Signal Processing, June 2008.
[61] M. Yukawa, R. C. de Lamare and R. Sampaio-Neto, “Efficient Acoustic
Echo Cancellation With Reduced-Rank Adaptive Filtering Based on
Selective Decimation and Adaptive Interpolation,” IEEE Transactions
on Audio, Speech, and Language Processing, vol.16, no. 4, pp. 696710, May 2008.
[62] R. C. de Lamare and R. Sampaio-Neto, “Reduced-rank space-time adaptive interference suppression with joint iterative least squares algorithms
for spread-spectrum systems,” IEEE Trans. Vehi. Technol., vol. 59, no.
3, pp. 1217-1228, Mar. 2010.
[63] R. C. de Lamare and R. Sampaio-Neto, “Adaptive reduced-rank equalization algorithms based on alternating optimization design techniques
for MIMO systems,” IEEE Trans. Vehi. Technol., vol. 60, no. 6, pp.
2482-2494, Jul. 2011.
[64] R. C. de Lamare, L. Wang, and R. Fa, “Adaptive reduced-rank LCMV
beamforming algorithms based on joint iterative optimization of filters:
Design and analysis,” Signal Processing, vol. 90, no. 2, pp. 640-652,
Feb. 2010.
[65] R. Fa, R. C. de Lamare, and L. Wang, “Reduced-Rank STAP Schemes
for Airborne Radar Based on Switched Joint Interpolation, Decimation
and Filtering Algorithm,” IEEE Transactions on Signal Processing,
vol.58, no.8, Aug. 2010, pp.4182-4194.
[66] L. Wang and R. C. de Lamare, ”Low-Complexity Adaptive Step Size
Constrained Constant Modulus SG Algorithms for Blind Adaptive
Beamforming”, Signal Processing, vol. 89, no. 12, December 2009, pp.
2503-2513.
[67] L. Wang and R. C. de Lamare, “Adaptive Constrained Constant Modulus
Algorithm Based on Auxiliary Vector Filtering for Beamforming,” IEEE
Transactions on Signal Processing, vol. 58, no. 10, pp. 5408-5413, Oct.
2010.
[68] L. Wang, R. C. de Lamare, M. Yukawa, ”Adaptive Reduced-Rank
Constrained Constant Modulus Algorithms Based on Joint Iterative
Optimization of Filters for Beamforming,” IEEE Transactions on Signal
Processing, vol.58, no.6, June 2010, pp.2983-2997.
[69] L. Wang, R. C. de Lamare and M. Yukawa, “Adaptive reduced-rank
constrained constant modulus algorithms based on joint iterative optimization of filters for beamforming”, IEEE Transactions on Signal
Processing, vol.58, no. 6, pp. 2983-2997, June 2010.
[70] L. Wang and R. C. de Lamare, “Adaptive constrained constant modulus
algorithm based on auxiliary vector filtering for beamforming”, IEEE
Transactions on Signal Processing, vol. 58, no. 10, pp. 5408-5413,
October 2010.
[71] R. Fa and R. C. de Lamare, “Reduced-Rank STAP Algorithms using
Joint Iterative Optimization of Filters,” IEEE Transactions on Aerospace
and Electronic Systems, vol.47, no.3, pp.1668-1684, July 2011.
[72] Z. Yang, R. C. de Lamare and X. Li, “L1-Regularized STAP Algorithms
With a Generalized Sidelobe Canceler Architecture for Airborne Radar,”

[73]

[74]

[75]

[76]

[77]

[78]

[79]

[80]

[81]

[82]

[83]

[84]

[85]

[86]

[87]

[88]

[89]

[90]

[91]
[92]

[93]

IEEE Transactions on Signal Processing, vol.60, no.2, pp.674-686, Feb.
2012.
Z. Yang, R. C. de Lamare and X. Li, “Sparsity-aware spacetime adaptive
processing algorithms with L1-norm regularisation for airborne radar”,
IET signal processing, vol. 6, no. 5, pp. 413-423, 2012.
Neto, F.G.A.; Nascimento, V.H.; Zakharov, Y.V.; de Lamare, R.C.,
”Adaptive re-weighting homotopy for sparse beamforming,” in Signal
Processing Conference (EUSIPCO), 2014 Proceedings of the 22nd
European , vol., no., pp.1287-1291, 1-5 Sept. 2014
Almeida Neto, F.G.; de Lamare, R.C.; Nascimento, V.H.; Zakharov,
Y.V.,“Adaptive reweighting homotopy algorithms applied to beamforming,” IEEE Transactions on Aerospace and Electronic Systems, vol.51,
no.3, pp.1902-1915, July 2015.
L. Wang, R. C. de Lamare and M. Haardt, “Direction finding algorithms
based on joint iterative subspace optimization,” IEEE Transactions on
Aerospace and Electronic Systems, vol.50, no.4, pp.2541-2553, October
2014.
S. D. Somasundaram, N. H. Parsons, P. Li and R. C. de Lamare,
“Reduced-dimension robust capon beamforming using Krylov-subspace
techniques,” IEEE Transactions on Aerospace and Electronic Systems,
vol.51, no.1, pp.270-289, January 2015.
S. Xu and R.C de Lamare, , Distributed conjugate gradient strategies for
distributed estimation over sensor networks, Sensor Signal Processing
for Defense SSPD, September 2012.
S. Xu, R. C. de Lamare, H. V. Poor, “Distributed Estimation Over
Sensor Networks Based on Distributed Conjugate Gradient Strategies”,
IET Signal Processing, 2016 (to appear).
S. Xu, R. C. de Lamare and H. V. Poor, Distributed Compressed
Estimation Based on Compressive Sensing, IEEE Signal Processing
letters, vol. 22, no. 9, September 2014.
S. Xu, R. C. de Lamare and H. V. Poor, “Distributed reduced-rank
estimation based on joint iterative optimization in sensor networks,”
in Proceedings of the 22nd European Signal Processing Conference
(EUSIPCO), pp.2360-2364, 1-5, Sept. 2014
S. Xu, R. C. de Lamare and H. V. Poor, “Adaptive link selection
strategies for distributed estimation in diffusion wireless networks,” in
Proc. IEEE International Conference onAcoustics, Speech and Signal
Processing (ICASSP), , vol., no., pp.5402-5405, 26-31 May 2013.
S. Xu, R. C. de Lamare and H. V. Poor, “Dynamic topology adaptation
for distributed estimation in smart grids,” in Computational Advances
in Multi-Sensor Adaptive Processing (CAMSAP), 2013 IEEE 5th International Workshop on , vol., no., pp.420-423, 15-18 Dec. 2013.
S. Xu, R. C. de Lamare and H. V. Poor, “Adaptive Link Selection
Algorithms for Distributed Estimation”, EURASIP Journal on Advances
in Signal Processing, 2015.
N. Song, R. C. de Lamare, M. Haardt, and M. Wolf, “Adaptive Widely
Linear Reduced-Rank Interference Suppression based on the Multi-Stage
Wiener Filter,” IEEE Transactions on Signal Processing, vol. 60, no. 8,
2012.
N. Song, W. U. Alokozai, R. C. de Lamare and M. Haardt, “Adaptive
Widely Linear Reduced-Rank Beamforming Based on Joint Iterative
Optimization,” IEEE Signal Processing Letters, vol.21, no.3, pp. 265269, March 2014.
R.C. de Lamare, R. Sampaio-Neto and M. Haardt, ”Blind Adaptive
Constrained Constant-Modulus Reduced-Rank Interference Suppression
Algorithms Based on Interpolation and Switched Decimation,” IEEE
Trans. on Signal Processing, vol.59, no.2, pp.681-695, Feb. 2011.
Y. Cai, R. C. de Lamare, “Adaptive Linear Minimum BER ReducedRank Interference Suppression Algorithms Based on Joint and Iterative
Optimization of Filters,” IEEE Communications Letters, vol.17, no.4,
pp.633-636, April 2013.
R. C. de Lamare and R. Sampaio-Neto, “Sparsity-Aware Adaptive
Algorithms Based on Alternating Optimization and Shrinkage,” IEEE
Signal Processing Letters, vol.21, no.2, pp.225,229, Feb. 2014.
L. Qiu, Y. Cai, R. C. de Lamare and M. Zhao, “Reduced-Rank DOA Estimation Algorithms Based on Alternating Low-Rank Decomposition,”
in IEEE Signal Processing Letters, vol. 23, no. 5, pp. 565-569, May
2016.
R. C. de Lamare, “Massive MIMO Systems: Signal Processing Challenges and Future Trends”, Radio Science Bulletin, December 2013.
W. Zhang, H. Ren, C. Pan, M. Chen, R. C. de Lamare, B. Du and J.
Dai, “Large-Scale Antenna Systems With UL/DL Hardware Mismatch:
Achievable Rates Analysis and Calibration”, IEEE Trans. Commun.,
vol.63, no.4, pp. 1216-1229, April 2015.
M. Costa, ”Writing on dirty paper,” IEEE Trans. Inform. Theory, vol.
29, no. 3, pp. 439-441, May 1983.

16

[94] R. C. de Lamare and A. Alcaim, ”Strategies to improve the performance
of very low bit rate speech coders and application to a 1.2 kb/s codec”
IEE Proceedings- Vision, image and signal processing, vol. 152, no. 1,
February, 2005.
[95] P. Clarke and R. C. de Lamare, ”Joint Transmit Diversity Optimization
and Relay Selection for Multi-Relay Cooperative MIMO Systems Using
Discrete Stochastic Algorithms,” IEEE Communications Letters, vol.15,
no.10, pp.1035-1037, October 2011.
[96] P. Clarke and R. C. de Lamare, ”Transmit Diversity and Relay Selection
Algorithms for Multirelay Cooperative MIMO Systems” IEEE Transactions on Vehicular Technology, vol.61, no. 3, pp. 1084-1098, March
2012.
[97] T. Peng, R. C. de Lamare and A. Schmeink, “Adaptive Distributed
Space-Time Coding Based on Adjustable Code Matrices for Cooperative
MIMO Relaying Systems,” IEEE Transactions on Communications, vol.
61, no. 7, pp. 2692-2703, July 2013.
[98] T. Peng and R. C. de Lamare, “Adaptive Buffer-Aided Distributed SpaceTime Coding for Cooperative Wireless Networks,” IEEE Transactions
on Communications, vol. 64, no. 5, pp. 1888-1900, May 2016.
[99] Y. Cai, R. C. de Lamare, and R. Fa, “Switched Interleaving Techniques
with Limited Feedback for Interference Mitigation in DS-CDMA Systems,” IEEE Transactions on Communications, vol.59, no.7, pp.19461956, July 2011.
[100] Y. Cai, R. C. de Lamare, D. Le Ruyet, “Transmit Processing Techniques
Based on Switched Interleaving and Limited Feedback for Interference
Mitigation in Multiantenna MC-CDMA Systems,” IEEE Transactions
on Vehicular Technology, vol.60, no.4, pp.1559-1570, May 2011.
[101] T. Wang, R. C. de Lamare, and P. D. Mitchell, “Low-Complexity
Set-Membership Channel Estimation for Cooperative Wireless Sensor
Networks,” IEEE Transactions on Vehicular Technology, vol.60, no.6,
pp.2594-2607, July 2011.
[102] T. Wang, R. C. de Lamare and A. Schmeink, ”Joint linear receiver
design and power allocation using alternating optimization algorithms
for wireless sensor networks,” IEEE Trans. on Vehi. Tech., vol. 61, pp.
4129-4141, 2012.
[103] R. C. de Lamare, “Joint iterative power allocation and linear interference suppression algorithms for cooperative DS-CDMA networks”, IET
Communications, vol. 6, no. 13 , 2012, pp. 1930-1942.
[104] T. Peng, R. C. de Lamare and A. Schmeink, “Adaptive Distributed
Space-Time Coding Based on Adjustable Code Matrices for Cooperative
MIMO Relaying Systems”, IEEE Transactions on Communications, vol.
61, no. 7, July 2013.
[105] K. Zu, R. C. de Lamare, “Low-Complexity Lattice Reduction-Aided
Regularized Block Diagonalization for MU-MIMO Systems”, IEEE.
Communications Letters, Vol. 16, No. 6, June 2012, pp. 925-928.
[106] K. Zu, R. C. de Lamare, “Low-Complexity Lattice Reduction-Aided
Regularized Block Diagonalization for MU-MIMO Systems”, IEEE.
Communications Letters, Vol. 16, No. 6, June 2012.
[107] K. Zu, R. C. de Lamare and M. Haart, “Generalized design of
low-complexity block diagonalization type precoding algorithms for
multiuser MIMO systems”, IEEE Trans. Communications, 2013.
[108] M. Tomlinson, ”New automatic equaliser employing modulo arithmetic,” Electronic Letters, vol. 7, Mar. 1971.
[109] C. T. Healy and R. C. de Lamare, “Decoder-optimised progressive edge
growth algorithms for the design of LDPC codes with low error floors”,
IEEE Communications Letters, vol. 16, no. 6, June 2012, pp. 889-892.
[110] A. G. D. Uchoa, C. T. Healy, R. C. de Lamare, R. D. Souza, “LDPC
codes based on progressive edge growth techniques for block fading
channels”, Proc. 8th International Symposium on Wireless Communication Systems (ISWCS), 2011, pp. 392-396.
[111] A. G. D. Uchoa, C. T. Healy, R. C. de Lamare, R. D. Souza, “Generalised Quasi-Cyclic LDPC codes based on progressive edge growth
techniques for block fading channels”, Proc. International Symposium
Wireless Communication Systems (ISWCS), 2012, pp. 974-978.
[112] A. G. D. Uchoa, C. T. Healy, R. C. de Lamare, R. D. Souza, “Design
of LDPC Codes Based on Progressive Edge Growth Techniques for
Block Fading Channels”, IEEE Communications Letters, vol. 15, no.
11, November 2011, pp. 1221-1223.
[113] C. T. Healy and R. C. de Lamare, “Design of LDPC Codes Based
on Multipath EMD Strategies for Progressive Edge Growth”, IEEE
Transactions on Communications, 2016.
[114] H. Harashima and H. Miyakawa, ”Matched-transmission technique for
channels with intersymbol interference,” IEEE Trans. Commun., vol. 20,
Aug. 1972.
[115] K. Zu, R. C. de Lamare and M. Haardt, “Multi-branch tomlinsonharashima precoding for single-user MIMO systems,” in Smart Antennas

(WSA), 2012 International ITG Workshop on , vol., no., pp.36-40, 7-8
March 2012.
[116] K. Zu, R. C. de Lamare and M. Haardt, “Multi-Branch TomlinsonHarashima Precoding Design for MU-MIMO Systems: Theory and
Algorithms,” IEEE Transactions on Communications, vol.62, no.3,
pp.939,951, March 2014.
[117] L. Zhang, Y. Cai, R. C. de Lamare and M. Zhao, “Robust Multibranch
TomlinsonHarashima Precoding Design in Amplify-and-Forward MIMO
Relay Systems,” IEEE Transactions on Communications, vol.62, no.10,
pp.3476,3490, Oct. 2014.
[118] B. Hochwald, C. Peel and A. Swindlehurst, ”A vector-perturbation
technique for near capacity multiantenna multiuser communication - Part
II: Perturbation,” IEEE Trans. Commun., vol. 53, no. 3, Mar. 2005.
[119] C. B. Chae, S. Shim and R. W. Heath, ”Block diagonalized vector
perturbation for multiuser MIMO systems,” IEEE Trans. Wireless Commun., vol. 7, no. 11, pp. 4051 - 4057, Nov. 2008.
[120] R. C. de Lamare, R. Sampaio-Neto, “Adaptive MBER decision feedback multiuser receivers in frequency selective fading channels”, IEEE
Communications Letters, vol. 7, no. 2, Feb. 2003, pp. 73 - 75.
[121] A. Rontogiannis, V. Kekatos, and K. Berberidis,” A Square-Root
Adaptive V-BLAST Algorithm for Fast Time-Varying MIMO Channels,”
IEEE Signal Processing Letters, Vol. 13, No. 5, pp. 265-268, May 2006.
[122] R. C. de Lamare, R. Sampaio-Neto, A. Hjorungnes, “Joint iterative
interference cancellation and parameter estimation for CDMA systems”,
IEEE Communications Letters, vol. 11, no. 12, December 2007, pp. 916
- 918.
[123] Y. Cai and R. C. de Lamare, ”Adaptive Space-Time Decision Feedback
Detectors with Multiple Feedback Cancellation”, IEEE Transactions on
Vehicular Technology, vol. 58, no. 8, October 2009, pp. 4129 - 4140.
[124] J. W. Choi, A. C. Singer, J Lee, N. I. Cho, “Improved linear softinput soft-output detection via soft feedback successive interference
cancellation,” IEEE Trans. Commun., vol.58, no.3, pp.986-996, March
2010.
[125] R. C. de Lamare and R. Sampaio-Neto, “Blind adaptive MIMO
receivers for space-time block-coded DS-CDMA systems in multipath
channels using the constant modulus criterion,” IEEE Transactions on
Communications, vol.58, no.1, pp.21-27, January 2010.
[126] R. Fa, R. C. de Lamare, “Multi-Branch Successive Interference Cancellation for MIMO Spatial Multiplexing Systems”, IET Communications,
vol. 5, no. 4, pp. 484 - 494, March 2011.
[127] R.C. de Lamare and R. Sampaio-Neto, “Adaptive reduced-rank equalization algorithms based on alternating optimization design techniques
for MIMO systems,” IEEE Trans. Veh. Technol., vol. 60, no. 6, pp.
2482-2494, July 2011.
[128] P. Li, R. C. de Lamare and R. Fa, “Multiple Feedback Successive
Interference Cancellation Detection for Multiuser MIMO Systems,”
IEEE Transactions on Wireless Communications, vol. 10, no. 8, pp. 2434
- 2439, August 2011.
[129] R.C. de Lamare, R. Sampaio-Neto, “Minimum mean-squared error
iterative successive parallel arbitrated decision feedback detectors for
DS-CDMA systems,” IEEE Trans. Commun., vol. 56, no. 5, May 2008,
pp. 778-789.
[130] R.C. de Lamare, R. Sampaio-Neto, “Minimum mean-squared error
iterative successive parallel arbitrated decision feedback detectors for
DS-CDMA systems,” IEEE Trans. Commun., vol. 56, no. 5, May 2008.
[131] R.C. de Lamare and R. Sampaio-Neto, “Adaptive reduced-rank equalization algorithms based on alternating optimization design techniques
for MIMO systems,” IEEE Trans. Veh. Technol., vol. 60, no. 6, pp.
2482-2494, July 2011.
[132] P. Li, R. C. de Lamare and J. Liu, “Adaptive Decision Feedback Detection with Parallel Interference Cancellation and Constellation Constraints
for Multiuser MIMO systems”, IET Communications, vol.7, 2012, pp.
538-547.
[133] J. Liu, R. C. de Lamare, “Low-Latency Reweighted Belief Propagation
Decoding for LDPC Codes,” IEEE Communications Letters, vol. 16, no.
10, pp. 1660-1663, October 2012.
[134] P. Li and R. C. de Lamare, Distributed Iterative Detection With
Reduced Message Passing for Networked MIMO Cellular Systems,
IEEE Transactions on Vehicular Technology, vol.63, no.6, pp. 29472954, July 2014.
[135] A. G. D. Uchoa, C. T. Healy and R. C. de Lamare, “Iterative
Detection and Decoding Algorithms For MIMO Systems in BlockFading Channels Using LDPC Codes,” IEEE Transactions on Vehicular
Technology, 2015.
[136] R. C. de Lamare, “Adaptive and Iterative Multi-Branch MMSE Decision Feedback Detection Algorithms for Multi-Antenna Systems”, IEEE
Trans. Wireless Commun., vol. 14, no. 10, October 2013.

