R. D. De Groat, et. Al. “Subspace Tracking.”
2000 CRC Press LLC. <http://www.engnetbase.com>.
SubspaceTracking
R.D.DeGroat
TheUniversityofTexasatDallas
E.M.Dowling
TheUniversityofTexasatDallas
D.A.Linebarger
TheUniversityofTexasatDallas
66.1Introduction
66.2Background
EVDvs.SVD
•
ShortMemoryWindowsforTimeVarying
Estimation
•
ClassificationofSubspaceMethods
•
Historical
OverviewofMEPMethods
•
HistoricalOverviewofAdaptive,
Non-MEPMethods
66.3IssuesRelevanttoSubspaceandEigenTrackingMethods
BiasDuetoTimeVaryingNatureofDataModel
•
Control-
lingRoundoffErrorAccumulationandOrthogonalityErrors
•
Forward-BackwardAveraging
•

andsubspacetrackingmethods.Fourbasicstrategieshavebeenpursuedtoreducecomputation:
(1)computingonlyafeweigencomponents,(2)computingasubspacebasisinsteadofindividual
eigencomponents,(3)approximatingtheeigencomponentsorbasis,and(4)recursivelyupdatingthe
eigencomponentsorbasis.Themostefficientmethodsusuallyemployseveralofthesestrategies.
In1990,anextensivesurveyofSVDtrackingmethodswaspublishedbyComonandGolub[7].
Theyclassifiedthevariousalgorithmsaccordingtocomplexityandbasicallytwocategoriesemerge:
O(n
2
r)andO(nr
2
)methods,wherenisthesnapshotvectorsizeandristhenumberofextreme
eigenpairstobetracked.Typically,r<norrn,sotheO(nr
2
)methodsinvolvesignificantlyfewer
computationsthantheO(n
2
r)algorithms.However,since1990,anumberofO(nr)algorithmshave
c

1999byCRCPressLLC
been developed. This article will primarily focus on recursive subspace and eigen updating methods
developed since 1990, especially, the O(nr
2
) and O(nr) algorithms.
66.2 Background
66.2.1 EVD vs. SVD
Let X =[x
1
|x
2

is a diagonal matrix whose nonzero entries are positive. It is easy to see that the left
singular vectors of X are the eigenvectors of XX
H
= USS
T
U
H
, and the right singular vectors of X
are the eigenvectors of X
H
X = VS
T
SV
H
. This is so because XX
H
and X
H
X are positive definite
Hermitian matrices (which have orthogonal eigenvectors and real, positive eigenvalues). Also note
that the nonzero singular values of X are the positive square roots of the nonzero eigenvalues of XX
H
and X
H
X. Mathematically, the eigen information contained in the SVD of X or the EVD of XX
H
(or X
H
X) is equivalent, but the dynamic range of the eigenvalues is twice that of the corresponding
singular values. With finite precision arithmetic, the greater dynamic range can result in a loss of

slowly enough within the effective observation window that the process is approximately stationary
and some averaging is desirable. In any event, at time k, a length N moving rectangular data window
results in a rank two modification of the correlation matrix estimate, i.e.,
R
(rect)
k
= R
(rect)
k−1
+
1
N
(x
k
x
H
k
− x
k−N
x
H
k−N
)
(66.1)
where x
k
is the new snapshot vector and x
k−N
is the oldest vector which is being removed from the
estimate. The corresponding data matrix is given by X

An exponentially faded data window produces a rank one modification in
R
(f ade)
k
= αR
(f ade)
k−1
+ (1 − α)x
k
x
H
k
(66.2)
where α is the fading factor with 0 ≤ α ≤ 1. In this case, the data matrix is growing in size, but the
older data is de-emphasized with a diagonal weighting matrix,
X
(f ade)
k
=[x
k
|x
k−1
|...|x
1
] sqrt(diag(1,α,α
2
, ..., α
k−1
))and R
(f ade)

(66.3)
Based on a Fourier analysis of linearly varying frequencies, equal frequency lags occur when [14]
N ≈
(1 + α)
(1 − α)
or equivalent ly α ≈
N − 1
N + 1
.
(66.4)
Either one of these relationships could be used as a rule of thumb for relating the effective observation
window of the two most popular short memory windowing schemes.
66.2.3 Classification of Subspace Methods
Eigenstructure estimation can be classified as (1) block or (2) recursive. Block methods simply
compute an EVD, SVD, or related decomposition based on a block of data. Recursive methods
update the previously computed eigen information using new data as it arrives. Wefocus on recursive
subspace updating methods in this article.
Most subspace tracking algorithms can also be broadly categorized as (1) modified eigen problem
(MEP) methods or (2) adaptive (or non-MEP) methods. With short memory windowing, MEP
methods are adaptive in the sense that they can track time varying eigen information. However,
when we use the word adaptive, we mean that exact eigen information is not computed at each
update, but rather, an adaptive method tends to move towards an EVD (or some aspect of an EVD) at
each update. For example, gradient-based, perturbation-based, and neural network-based methods
are classified as adaptive because on average they move towards an EVD at each update. On the other
hand, rank one, rank k, and sphericalized EVD and SVD updates are, by definition, MEP methods
because exact eigen information associated with an explicit matrix is computed at each update. Both
MEP and adaptive methods are supposed to track the eigen information of the instantaneous, time
varying correlation matrix.
c


In [42], Xu and Kailath present a Lanczos based subspace tracking method with an associated
detection scheme to track the number of sources. A reference list for systolic implementations of
SVD based subspace trackers is contained in [12].
66.2.5 Historical Overview of Adaptive, Non-MEP Methods
Owsley pioneered orthogonal iteration and stochastic-based subspace trackers in [32]. Yang and
KavehextendedOwsley’sworkin[44] by devising a family of constrained gradient-based algorithms.
A highly parallel algorithm, denoted the inflation method, is introduced for the estimation of the
noise subspace. The computational complexity of this family of gradient-based methods varies from
(approximately) n
2
r to
7
2
nr for the adaptation equation. However, since the eigenvectors are only ap-
proximatelyorthogonal, an additional nr
2
flops may be needed if Gram Schmidt orthogonalization is
used. It maybe that a partial orthogonalization scheme (see Section 66.3.2 Controlling Roundoff Error
Accumulation and Orthogonality Errors) can be combined with Yang and Kaveh’s methods to improve
orthogonality enough to eliminate the O(nr
2
) Gram Schmidt computation. Karhunen [22] also ex-
tended Owsley’s work by developing a stochastic approximation method for subspace computation.
Bin Yang [43] used recursive least squares (RLS) methods with a projection approximation approach
to develop the projection approximation subspace tracker (PAST) which tracks an arbitrary basis for
the signal subspace, and PASTd which uses deflation to track the individual eigencomponents. A
multi-vector eigen tracker based on the conjugate gradient method is developed in [18]. Previous
conjugate gradient-based methods tracked a single eigenvector only. Orthogonal iteration, lossless
adaptive filter, and perturbation-based subspace trackers appear in [40][36], and [5] respectively.
A family of non-EVD subspace trackers is given in [16]. An adaptive subspace method that uses a

fashion. However, recursive updating algorithms cannot tolerateeven a linear buildup oferrorif large
(possibly unbounded) numbers of updates are to be performed. For real time processing, periodic
reinitialization is undesirable. Most of the subspace tracking algorithms involve the product of at
least k orthogonal matrices by the time the kth update is computed. According to Parlett [33], the
error propagated by a product of orthogonal matrices is bounded as
|U
k
U
H
k
− I|
E
≤ (k + 1)n
1.5

(66.5)
where the n × n matrix U
k
= U
k−1
Q
k
= Q
k
Q
k−1
...Q
1
is a product of k matrices that are each
orthogonal to working accuracy,  is machine precision, and|.|