Merging and Splitting Eigenspace Models.
AbstractWe present new deterministic methods that given two eigenspace modelseach representing a set of $n$-dimensional observationswill: 1) merge the models to yield a representation of the union of the sets and 2) split one model from another to represent the difference between the sets. As this is done, we accurately keep track of the mean. Here, we give a theoretical derivation of the methods, empirical results relating to the efficiency and accuracy of the techniques, and three general applications, including the construction of Gaussian mixture models that are dynamically updateable.
Introduction
The contributions of this paper are: (1) a method for
merging eigenspace models; (2) a method for splitting
eigenspace models. These represent an advance in
that previous methods for incremental computation
of eigenspace models considered only the addition (or
subtraction) of a single observation to (or from) an
eigenspace model [1, 2, 3, 4, 5, 6]. Our method also
allows the origin to be updated, unlike most other
methods. Thus, our methods allow large eigenspace
All authors are with the Department of Computer Science,
University of Wales, Cardiff, PO Box 916, Cardiff CF2 3XF,
Wales UK: peter@cs.cf.ac.uk
models to be updated more quickly and more accurately
than when using previous methods. A second
advantage is that very large eigenspace models
may now be reliably constructed using a divide-and-
conquer approach. Limitations of existing techniques
are mentioned in [7].
Eigenspace models have a wide variety of applica-
tions, for example: classification for recognition systems
[8], characterising normal modes of vibration
for dynamic models, such as the heart [9], motion
sequence analysis [10], and the temporal tracking of
signals [4]. Clearly, in at least some of these appli-
cations, merging and splitting of eigenspace models
will be useful, as the data is likely to be presented
or processed in batches. For example, a common
paradigm for machine learning systems is to identify
clusters [11], our methods allow clusters to be
split and merged, and dynamically updated as new
data arrive. As another example, an image database
of university students may require as much as one
quarter of its records to be replaced each year. our
methods permit this without the need to recompute
the eigenspace model ab initio.
This paper is solely concerned with deriving a
theoretical framework for merging and splitting
eigenspaces, and an empirical evaluation of these new
techniques, rather than their particular application in
any area.
An eigenspace model is a statistical description of a
set of N observations in n-dimensional space; such a
model may be regarded as a multi-dimensional Gaussian
distribution. From a geometric point of view,
an eigenspace model can be thought of as a hyperellipsoid
that characterises a set of observations: its
centre is the mean of the observations; its axes point
in directions along which the spread of observations
is maximised, subject to them being orthogonal; the
surface of the hyperellipsoid is a contour that lies at
one standard deviation from the mean. Often, the
hyperellipsoid is almost flat along certain directions,
and thus can be modelled as having lower lower dimension
than the space in which it is embedded.
Eigenspace models are computed using either
eigenvalue decomposition (also called principal component
analysis) or singular-value decomposition.
We wish to distinguish between batch and incremental
computation. In batch computation all observations
are used simultaneously to compute the
eigenspace model. In an incremental computation,
an existing eigenspace model is updated using new
observations.
Previous research in incremental computation of
eigenspace models has only considered adding exactly
one new observation at a time to an eigenspace
model [1, 2, 3, 4, 5, 6]. A common theme of these
methods is that none require the original observations
to be retained. Rather, a description of the
hyperellipsoid is sufficient information for incremental
computation of the new eigenspace model. Each
of previous these approaches allows for a change in
dimensionality of the hyperellipsoid, so that a single
additional axis is added if necessary. Only our
previous work allows for a shift of the centre of the
hyperellipsoid [6]: other methods keep it fixed at the
origin. This proves crucial if the eigenspace model is
to be used for classification, as explained in [6]: a set
of observations whose mean is far from the origin is
clearly not well modelled by a hyperellipsoid centred
at the origin.
When using incremental methods previous observations
need not be kept - thus reducing storage requirements
and making large problems computationally
feasible. Incremental methods must be used if
not all observations are available simultaneously. For
example, a computer may lack the memory resources
required to store all observations. This is true even
if low-dimensional methods are used to compute the
eigenspace [5, 12]. (We will mention low-dimensional
methods later, but they give an advantage when the
number of observations of less than the dimensionality
of the space, N ! n, which is often true when
observations are images.) Even if all observations
are available, it is usually faster to compute a new
eigenspace model by incrementally updating an existing
one rather than using batch computation [3].
This is because the incremental methods typically
compute p eigenvectors, with p - min(n; N ). The
disadvantage of incremental methods is their accuracy
compared to batch methods. When only few incremental
updates are made the inaccuracy is small,
and is probably acceptable for the great majority of
applications [6]. When many thousands of updates
are made, as when eigenspace models are incremented
with a single observation at a time, the inaccuracies
build up, although methods exist to circumvent this
problem [4]. In contrast, our methods allow a whole
new set of observations to be added in a single step,
thus reducing the total number of updates to an existing
model.
Section 2 defines eigenspace models in detail, standard
methods for computing them, and how they
are used for representing and classifying observations.
Section 3 discusses merging of eigenspace models,
while Section 4 treats splitting. Section 5 presents
empirical results, and Section 6 gives our conclusions.
Eigenspace models
In this section, we describe what we mean by
eigenspace models, briefly discuss standard methods
for their batch computation, and how observations
can be represented using them. Firstly, we establish
our notation for the rest of the paper.
Vectors are columns, and denoted by a single un-
derline. Matrices are denoted by a double under-
line. The size of a vector, or matrix, is often im-
portant, and where we wish to emphasise this size,
it is denoted by subscripts. Particular column vectors
within a matrix are denoted by a superscript,
and a superscript on a vector denotes a particular
observation from a set of observations, so we treat
observations as column vectors of a matrix. As an
example, A i
mn is the ith column vector in an (m \Theta n)
matrix. We denote matrices formed by concatenation
using square brackets. Thus [A mn b] is an (m \Theta (n+1)
matrix, with vector b appended to A mn
, as a last column

2.1 Theoretical background
Consider N observations, each a column vector x
We compute an eigenspace model as follows:
The mean of the observations is
and their covariance is
Note that C nn
is real and symmetric.
The axes of the hyperellipsoid, and the spread of
observations over each axis are the eigenvectors and
eigenvalues of the eigenproblem
or, equivalently, the eigenvalue decomposition (EVD)
of C nn is
nn
nn
where the columns of U nn
are eigenvectors, and   nn
is a diagonal matrix of eigenvalues. The eigenvectors
are orthonormal, so that U T
nn
U nn
I nn
, the (n \Theta n)
identity matrix.
The ith eigenvector U i and ith eigenvalue   ii
nn are
associated; the eigenvalue is the length of the eigen-
vector, which is the ith axis of the hyperellipsoid.
Typically, only p - min(n; N) of the eigenvectors
have significant eigenvalues, and hence only p of the
eigenvectors need be retained. This is because
the observations are correlated so that the covariance
matrix is, to a good approximation rank-degenerate:
small eigenvalues are presumed to be negligible. Thus
an eigenspace model often spans a p-dimensional sub-space
of the n-dimensional space in which it is embedded

Different criteria for discarding eigenvectors and
eigenvalues exist, and these suit different applications
and different methods of computation. Three common
methods are: (1) stipulate p as a fixed inte-
ger, and so keep the p largest eigenvectors [5]; (2)
those p eigenvectors whose size is larger than
an absolute threshold [3]; (3) keep the p eigenvectors
such that a specified fraction of energy in the
eigenspectrum (computed as the sum of eigenvalues)
is retained.
Having chosen to discard certain eigenvectors and
eigenvalues, we can recast Equation 4 using block
form matrices and vectors. Without loss of general-
ity, we can permute the eigenvectors and eigenvalues
such that U np are those eigenvectors that are kept,
and   pp their eigenvalues. If nd and
dd
are those discarded. We may rewrite Equation 4
as:
nn
U nd
dd
[U np
U nd
nd (5)
Hence
with error U nd   dd U T
nd , which is small if   dd - 0 nd .
Thus, we define an eigenspace model, \Omega\Gamma as the
mean, a (reduced) set of eigenvectors, their eigen-
values, and the number of observations:
2.2 Low-dimensional computation of
eigenspace models
Low-dimensional batch methods are often used to
compute eigenspace models, and are especially important
when the dimensionality of the observations
is very large compared to their number. Thus, they
may be used to compute eigenspace models that
would otherwise be infeasible. Incremental methods
also use a low dimensional approach.
In principle, computing an eigenspace model requires
that we construct an (n \Theta n) matrix, where n
is the dimension of each observation. In practice, the
model can be computed by using an (N \Theta N) ma-
trix, where N is the number of observations. This
is an advantage in applications like image processing
where, typically, N - n.
We show how this can be done by first considering
the relationship between eigenvalue decomposition
(EVD) and singular value decomposition (SVD). This
leads to a simple derivation for a low-dimensional
batch method for computing the eigenspace model.
The same results were obtained, at greater length,
by [5], see also [12].
Let Y nN be the set of observations shifted to the
mean, so that Y
x. Then a SVD of Y nN
is:
where U nn are the left singular vectors, which are
identical to the eigenvectors previously given; \Sigma nN
is
a matrix with singular values on its leading diagonal,
with   nn
are right singular
vectors. Both U nn
and V NN
are orthonormal matrices

This can now be used to compute eigenspace models
in a low-dimensional way, as follows:
Y nN
is an (N \Theta N) eigenproblem. S NN
is the same as
nn =N , except for the presence of extra trailing zeros
on the main diagonal of   nn
. If we discard the small
singular values, and their singular vectors, following
the above, then remaining eigenvectors vectors are
U np
pp
This result formed the basis of the incremental
technique developed by [5] but they did not allow
for a change in origin, nor does their approach readily
generalise to merging and splitting. Others [3]
observe that a solution based on the matrix product
, as above, is likely to lead to inaccurate
results because of conditioning problems, and they
develop a method for incrementally updating SVD
solutions with a single observation. Although their
SVD method has proven more accurate (see [6]), it
is, again, very difficult to generalise, especially if a
change of origin is to be taken into account.
SVD methods were actually proposed quite early in
the development of incremental eigenproblem analysis
[2]. This early work included a proposal to delete
single observations, but did not extend to merging
and splitting. SVD also formed the basis of a proposal
to incrementally update an eigenspace with several
observations at one step [10]. However, contrary
to our method, a possible change in the dimension of
the solution eigenspace was not considered. Further-
more, none of these methods considered a change in
origin.
Our incremental method is based on the matrix
product C
nN , and specifically its approximation
as in Equation 6. It is a generalisation of our
earlier work [6], which now appears naturally as the
special case of adding a single observation.
2.3 Representing and classifying observation

High-dimensional observations may be approximated
by a low-dimensional vector using an eigenspace
model. Eigenspace models may also be used for clas-
sification. We briefly discuss both ideas, prior to using
them in our results section.
An n-dimensional observation x n is represented using
an eigenspace
N) as a
p-dimensional vector
This shifts the observation to the mean, and then
represents it by components along each eigenvec-
tor. This is also called the Karhunen-Lo'eve transform
[13].
The n-dimensional residue vector is defined by:
and h n is orthogonal to every vector in U np . Thus, h n
is the error in the representation of x n with respect
The likelihood associated with the same observation
is given by:
Clearly, the above definition cannot be used directly
in cases where N - n, as C nn
is then rank
degenerate. In such cases we use an alternative definition
due to Moghaddam and Pentland [8] is appropriate
(a full explanation of which is beyond the
scope of this paper).
Merging Eigenspace models
We now turn our attention to one of the two main
contributions of this paper, merging eigenspace models

We derive a solution to the following problem. Let
and Y nM
be two sets of observations. Let
their eigenspace models
respectively. The problem is
to compute the eigenspace model \Phi (-z; W nr
for Z n(N+M)
Y nM
using
only\Omega and \Psi.
Clearly, the total number of new observations is
The combined mean is:
-z =(N +M)
The combined covariance matrix is:
where C nn
and D nn
are the covariance matrices for
and Y nM
, respectively.
We wish to compute the s eigenvectors and eigen-values
that satisfy:
ns
ss
ns
where some eigenvalues are subsequently discarded
to give r non-negligible eigenvectors and eigenvalues.
The problem to be solved is of size s, and this is
necessarily bounded by
We explain the perhaps surprising additional 1 in
the upper limit later(Section 3.1.1), but briefly, it
is needed to allow for the vector difference between
the means, -
y.
3.1 Method of solution
This problem may be solved in three steps:
1. Construct an orthonormal basis set, \Upsilon ns , that
spans both eigenspace models and -
x\Gamma-y. This basis
differs from the required eigenvectors, W ns ,
by a rotation, R ss , so that:
ns
\Upsilon ns R ss
2. Use \Upsilon ns to derive an intermediate eigenproblem.
The solution of this problem provides the eigen-
values, \Pi ss , needed for the merged eigenmodel.
The eigenvectors, R ss
, comprise the linear transform
that rotates the basis set \Upsilon ns
3. Compute the eigenvectors W ns
, as above, and
discard any eigenvectors and eigenvalues using
the chosen criteria (as discussed above) to yield
and \Pi rr
We now give details of each step.
3.1.1 Construct an orthonormal basis set
To construct an orthonormal basis for the combined
eigenmodels we must chose a set of orthonormal
vectors that span three subspaces: (1) the sub-space
spanned by eigenvectors U np
; (2): the sub-space
spanned by eigenvectors V nq
;(3) the subspace
spanned by (-x \Gamma - y). The last of these is a single vec-
tor. It is necessary because the vector joining the
centre of the two eigenspace models need not belong
to either eigenspace. This accounts for the additional
1 in the upper limit of the bounds of s (Equation 17).
For example, consider the case in which each of the
eigenspaces is a 2D ellipse in a 3D space, and the
ellipses are parallel but separated by a vector perpendicular
to each of them. Clearly, a merged model
should be a 3D ellipse because of the vector between
their origins.
A sufficient spanning set is:
ns nt
nt
is an orthonormal basis set for that component
of the eigenspace of \Psi which is orthogonal to
the eigenspace of \Omega\Gamma and in addition accounts for that
component of (-x \Gamma - y) orthogonal to both eigenspaces;
To construct - nt
we start by computing the
residues of each of the eigenvectors in V nq
with respect
to the eigenspace of \Omega\Gamma
The H nq
are all orthogonal to U np
in the sense that
In general, however, some of
the H nq
are zero vectors, because such vectors represent
the intersection of the two eigenspaces. We also
compute the residue h of -
x with respect to the
eigenspace of \Omega\Gamma using Equation 12.
- nt can now be computed by finding an orthonormal
basis for [H nq ; h], which is sufficient to ensure
that \Upsilon ns is orthonormal. Gramm-Schmidt orthonor-
malisation may be used to do this:
nt
3.1.2 Forming a intermediate eigenproblem
We now form a new eigenproblem by substituting
Equation 19 into Equation and the result together
with Equation 15 into Equation 16 to obtain:
nt
ss
\Pi ss R T
ss
nt
Multiplying both sides on the left by [U nt
on the right by [U np
nt
and using the fact that
[U np
nt
T is a left inverse of [U np
nt
we obtain:
nt
nt ss \Pi ss R T
ss
which is a new eigenproblem whose solution eigenvectors
constitute the R ss
we seek, and whose eigenvalues
provide eigenvalues for the combined eigenspace
model. We do not know the covariance matrices C nn
or D nn , but these can be eliminated as follows:
The first term in Equation 24 is proportional to:
nt
nt
nt C nn U np - T
nt C nn - nt
By Equation 6 U T
U np
. Also, U T
nt
pt by construction, and again, using Equation 6 we
conclude:
[U np
nt
[U np
nt
The second term in Equation 24 is proportional to:
nt
nt
np D nn U np U T
np D nn - nt
nt
nt D nn - nt
We have D nn
, which on substitution
gives:
U np
nt
nt
nt
nt
From Equation 20 we have G pq
.
nt V nq
. We obtain:
nt nt
Now consider the final term in Equation 24:
nt nt
nt
nt
nt
nt
Setting
nt
becomes: "
So, the new eigenproblem to be solved may be approximated
by
G pq
G pq
ss
ss
R T
ss
Each matrix is of size
Thus we have eliminated the need
for the original covariance matrices. Note this also reduces
the size of the central matrix on the left hand
side. This is of crucial computational importance because
it makes the eigenproblem tractable.
3.1.3 Computing the eigenvectors
The matrix \Pi ss
is the eigenvalue matrix we set out
to compute. The eigenvectors R ss
comprise a rotation
for \Upsilon ns
. Hence, we use Equation to compute
the eigenvectors for \Pi ss
. However, not all eigenvectors
and eigenvalues need be kept, and some
of them) may be discarded using a criterion as previously
discussed in Section 2. This discarding of eigen-vectors
and eigenvalues should usually be carried out
each time a pair of eigenspace models is merged.
3.2 Discussion on the form of the so-
lution
We now briefly justify that the solution obtained is of
the correct form by considering several special cases.
First, suppose that both eigenspace models are
null, that is each is specified by (0; 0; 0; 0). Then the
system is clearly degenerate and null eigenvectors and
zero eigenvalues are computed.
If exactly one eigenspace model is null, then the
non-null eigenspace model is computed and returned
by this process. To see this, suppose that \Psi is null.
Then, the second and third matrices on the left-hand
side of Equation 31 both disappear. The first matrix
reduces to   pp
exactly hence the eigen-values
remain unchanged. In this case, the rotation
R ss is the identity matrix, and the eigenvectors are
also unchanged. If
instead\Omega is a null model, then
only the second matrix will remain (as
- nt and V nq will be related by a rotation (or else
identical). The solution to the eigenproblem then
computes the inverse of any such rotation, and the
eigenspace model remain unchanged.
Suppose \Psi has exactly one observation, then it is
specified by (y; 0; 0; 1). Hence the middle term on
the left of Equation 31 disappears, and - nt is the
unit vector in the direction
x. Hence
is a scalar, and the eigenproblem becomes
which is exactly the form obtained when one observation
is explicitly added, as we have proven elsewhere
[6]. This special case has interesting properties
too: if the new observation lies within the subspace
spanned by U np
, then any change in the
eigenvectors and eigenvalues can be explained by rotation
and scaling caused by g p g T
. Furthermore, in
the unlikely event that -
y, then the right matrix
disappears altogether, in which case the eigenvalues
are scaled by N=(N +1), but the eigenvectors are un-
changed. Finally, as
indicating a stable model in
the limit.
If
the\Omega has exactly one observation, then it is specified
by (x; 0; 0; 1). Thus the first matrix on the left
of Equation 31 disappears. The G pq
is a zero matrix,
and - nt = [V nq h], where h is the component of -
which is orthogonal to the eigenspace of \Psi. Hence
the eigenproblem is:
Given that in this case \Gamma tq
, then
has the form of \Delta qq
, but with a row and
column of zeros appended. Also, fl t
Substitution of these terms shows that in this case
too, the solution reduces to the special case of adding
a single new observation: Equation 33 is of the same
form as Equation 32, as can readily shown.
If
the\Omega and \Psi models are identical, then - y.
In this case the third term on the left of Equation 31
disappears. Furthermore, \Gamma tq
is a zero matrix, and
G pq
U np
is the identity matrix, with
Hence, the first and second matrices on the left of
Equation 31 are identical, with
reduce to the matrices of eigenvalues. Hence, adding
two identical eigenmodels yields an identical third.
Finally, notice that for fixed M , as N !1 so the
solution tends to
the\Omega model; for fixed N the reverse
is true; and if M and N tend to 1 simultaneously,
then the final term loses its significance.
3.3 Algorithm
Here, for completeness, we now express the mathematical
results obtained above, for merging models,
in the form of an algorithm for direct computer im-
plementation; see Figure 1
3.4 Complexity
Computing an eigenspace model of size N as a single
batch incurs a computational cost O(N 3 ). Examination
of our merging algorithm shows that it
also requires an intermediate eigenvalue problem to
be solved, as well as other steps; again overall giving
cubic computational requirements. Nevertheless,
let us suppose
with N observations, can be
represented by p eigenvectors, and that \Psi, with M
Function Merge( -
returns (-z; W; \Pi; P )
BEGIN
y
for each column vector of H
discard this column,
if it is of small magnitude.
endfor
of
of \Delta
number of basis vectors in -
construct LHS of Equation 31
eigenvalues of A
eigenvectors of A
discard small eigensolutions, as appropriate
END

Figure

1: Algorithm for merging two eigenspace models

observations, can be represented by q eigenvectors.
Typically p and q are (much) less than N and M ,
respectively.
To compute an overall model with the batch
method requires O((N +M) 3 ) operations. Assuming
that both models to be merged are already known,
our merging method requires at most O((p
the problem to be solved becomes
smaller the greater the amount of overlap between
the eigenspaces
of\Omega and \Psi. (In fact, the number
of operations required is O(s 3 see the end of Section
3.1.1.)
If one, or both, of the models to be merged are unknown
initially, then we incur an extra cost of O(N 3 ),
reduces any advan-
tage. Nevertheless, in one typical scenario, we might
expect\Omega to be known (an existing large database of
N observations), while \Psi is a relatively small batch
of M observations to be added to it. In this case,
the extra penalty, of O(M 3 ), is of little significance
compared to O((N +M) 3 ).
while an exact analysis is complicated and
indeed data dependent, we expect efficiency gains in
both time and memory resources in practice.
Furthermore, if computer memory is limited, sub-division
of the initial set may be unavoidable in order
to manage eigenspace model computation. We have
now provided a tractable solution to this problem.
Splitting eigenspace models
Here we show how to split two eigenspace models.
Given an eigenspace model
remove M) from it to give a third
use \Pi rr
, because
ss
is not available in general. Although splitting is
essentially the opposite of merging, this is not completely
true as it is impossible to regenerate information
which has been discarded when the overall model
is created (whether by batch methods or otherwise).
Thus, if we split one eigenspace model from a larger
one, the eigenvectors of the remnant must still form
some subspace of the larger.
We state the results for splitting without proof.
Clearly, . The new mean is:
As in the case of merging, new eigenvalues and eigen-vectors
are computed via an intermediate eigenprob-
lem. In this case it is:
G rp
rp
r
rr
x).
The eigenvalues we seek are the q non-zero elements
on the diagonal of   rr . Thus we can permute
R rr
and   rr , and write without loss of generality:
R rr
rr
R rt
[R rp
R rt
Hence we need only identify the eigenvectors in R rr
with non-zero eigenvalues, and compute the U np as:
U np
In terms of complexity, splitting must always involve
the solution an eigenproblem of size r. An algorithm
for splitting may readily be derived using a
similar approach to that for merging.
This section describes various experiments that we
carried out to compare the computational efficiency
of a batch method and our new methods for merging
and splitting, and the eigenspace models produced.
We compared models in terms of Euclidean distance
between the means, mean angular deviation
of corresponding eigenvectors, and mean rela-
tiveabsolute difference between corresponding eigen-
values. In doing so, we took care that both models
had the same number of dimensions.
As well as the simple measures above, other
performance measures may be more relevant when
eigenspace models are used for particular ap-
plications, and thus other tests were also per-
formed. Eigenspace models may be used for approximating
high-dimensional observations with a low-dimensional
vector; the error is the size of the residue
vector. The sizes of such residue vectors can readily
be compared for both batch and incremental meth-
ods. Eigenspace models may also be used for classifying
observations, giving the likelihood that an observation
belongs to a cluster. Different eigenspace
models may be compared by relative differences in
likelihoods. We average these differences over all corresponding
observations.
We used a database of 400 face images (each of
10304 pixels) available on-line 1 in the tests reported
here ; similar results were obtained in tests with randomly
generated data. The gray levels in the images
were scaled into the range [0; 1] by division only, but
no other preprocessing was done. We implemented
all functions using commercially available software
(Matlab) on a computer with standard configuration
(Sun Sparc Ultra 10, 300 Hz, 64 Mb RAM).
The results we present used up to 300 images, as
the physical resources of our computer meant that
heavy paging started to occur beyond this limit for
the batch method, although such paging did not affect
the incremental method.
For all tests, the experimental procedure used
was to compute eigenspace models using a batch
method [12], and compare these to models produced
by merging or splitting other models also produced
by the batch method. In each case, the largest of
the three data sets contained 300 images. These were
partitioned into two data sets, each containing a multiple
of 50 images. We included the degenerate cases
when one model contained zero images. Note that we
tested both smaller models merged with larger ones,
and vice-versa.
The number of eigenvectors retained in any model,
including a merged model, was set to be 100 as a
maximum, for ease of comparing results. (Initial
tests using other strategies indicate that the resulting
eigenspace model is little effected.)
5.1 Timing
When measuring CPU time we ran the same code
several times and chose the smallest value, to minimise
the effect of other concurrently running process.
Initially we measured time taken to compute a
model using the batch methods, for data sets of different
sizes. Results are presented in Figure 2 and
show cubic complexity, as predicted.
1 The Olivetti database of faces:
http://www.cam-orl.co.uk/facedatabase.html
Number of input images
time
in
seconds
observed data
cubic fit

Figure

2: Time to compute an eigenspace model with
a batch method versus the number of images, N . The
time is approximated by the cubic: 5:3 \Theta 10 \Gamma4 N
5.1.1 Merging
We then measured the time taken to merge two previously
constructed models. The results are shown
in

Figure

3. This shows that time complexity is approximately
symmetric about the point
half the number of input images. This result may
be surprising because the algorithm given for merging
is not symmetric with respect to its inputs, despite
that fact that the mathematical solution is independent
of order. The approximate symmetry in
time-complexity can be explained by assuming independent
eigenspaces with a fixed upper-bound on
the number of eigenvectors: suppose the numbers of
eigenvectors in the models are N and M . Then complexities
of the main steps are approximately as fol-
lows: computing a new spanning set, - is O(M 3 );
solving an eigenproblem is O(N 3 +M 3 ); rotating the
new eigenvectors is O(N 3 +M 3 ). Thus the time com-
plexity, under the stated conditions, is approximately
O(N 3 +M 3 ), which is symmetric.
Next, the times taken to compute an eigenspace
model from 300 images in total, using the batch
method and our merging method, are compared in
Number of images in first model
time
in
seconds
incremental time
joint

Figure

3: Time to merge two eigenspace models of
images
versus the number of im-
ages, N , in \Omega\Gamma The number of images in \Phi is 300 \Gamma N .
Hence, the total number of different images used to
compute \Phi is constant 300.

Figure

4. The incremental time is the time needed
to compute the eigenspace model to be merged, and
merge it with a pre-computed existing one. The joint
time is the time to compute both smaller eigenmodels
and then merge them. As might be expected, incremental
time falls as the additional number of images
required falls. The joint time is approximately con-
stant, and very similar to the total batch time.
While the incremental method offers no time saving
in the cases above, it does use much less memory.
This could clearly be seen when a model was computed
using 400 images: paging effects set in when
a batch method was used and the time taken rose
to over 800 seconds. The time to produce an equivalent
model by merging two sub-models of size 200,
however, took less than half that.
5.1.2 Splitting
Time complexity for splitting eigenspaces should depend
principally on the size of the large eigenspace
which from which the smaller space is being removed,
and the size of the smaller eigenspace should have
Number of images in first model
time
in
seconds
incremental time
joint
batch time

Figure

4: Time to make a complete eigenspace model
for a database of 300 images. The incremental time
is the addition of the time to construct only the
eigenspace to be added. The joint time is the time to
compute both eigenspace models and merge them.
little effect. This is because the size of the intermediate
eigenproblem to be solved depends on the
size of the larger space, and therefore dominates the
complexity. These expectations are borne out experi-
mentally. We computed a large eigenmodel using 300
images, as before. We then removed smaller models
of sizes between 50 and 250 images inclusive, in
steps of 50 images. At most, 100 eigenvectors were
kept in any model.The average time taken was approximately
constant, and ranged between 9 and 12
seconds, with a mean time of about 11.4 seconds.
These figures are much smaller than those observed
for merging because the large eigenspace contains
only 100 eigenvectors. Thus the matrices involved in
the computation were of size (100 \Theta 100), whereas in
merging the size was at least (150 \Theta 150), and other
computations were involved (such as computing an
orthonormal basis).
5.2 Similarity and performance
The measures used for assessing similarity and performance
of batch and incremental methods were described
above.
5.2.1 Merging
We first compared the means of the models produced
by each method using Euclidean distance. This distance
is greatest when the models to be merged have
the same number of input images (150 in this case),
as fall smoothly to zero when either of the models to
be merged is empty. The value at maximum is typically
very small, and we measured it to be 3:5 \Theta 10 \Gamma14
units of gray level. This compares favourably with
the working precision of Matlab, which is 2:2 \Theta 10 \Gamma16 .
We next compared the directions of the eigenvectors
produced by each method. The error in eigenvector
direction was measured by the mean angular
deviation, as shown in Figure 5. Ignoring the degenerate
cases, when one of the models is empty, we see
that angular deviation has a single minimum when
the eigenspace models were built with about the same
number of images. This may be because when a small
model is added to a large model its information tends
to be swamped.
These results show angular deviation to be very
small on average.
Number of input images
eigenvector
mean
angular
deviation
in
degrees

Figure

5: Angular deviation between eigenvectors
produced by batch and incremental methods versus
the number of images in the first eigenspace model.
The sizes of eigenvalues from both methods were
compared next. In general we observed that the
smaller eigenvalues had larger errors, as might be expected
as they contain relatively little information
and so are more susceptible to noise. In Figure 6
we given a mean absolute difference in eigenvalue.
This rises to a single peak when the number of input
images in both models is the same. Even so, the
maximal value is small, 7 \Theta 10 \Gamma3 units of gray level.
The largest eigenvalue is typically about 100.
Number of input images
eigenvalue
mean
absolute
deviation

Figure

Difference between eigenvalues produced by
batch and incremental methods versus the number of
images in the first eigenspace model.
We now turn to performance measures. The
merged eigenspaces represent the image data with
little loss in accuracy, as measured by the mean difference
in residue error, Figure 7. This performance
measure is typically small, about 10 \Gamma6 units of gray
level per pixel, clearly below any noticeable effect.
Finally we compared differences in likelihood values
(Equation produced by the two methods.
This difference is again small, typically of the order
8 shows; this should be compared
with a mean likelihood over all observations of the
Again the differences in classifications
that would be made by these models would be very
small.
Number of input images
difference
in
mean
residue

Figure

7: Difference in reconstruction errors per pixel
produced by batch and incremental methods versus
the number of images in the first eigenspace model.
Number of input images
mean
class
difference

Figure

8: Difference in likelihoods produced by batch
and incremental methods versus the number of images
in the first eigenspace model.
5.2.2 Splitting
Similar measures for splitting were computed using
exactly those conditions described for testing the timing
of splitting, and for exactly those characteristics
described for merging. In each case a model to be
subtracted was computed by a batch method, and
removed from the overall model by our splitting pro-
cedure. Also, a batch model was made for purposes
of comparison from the residual data set. In all that
follows the phrase "size of the removed eigenspace"
means the number of images used to construct the
eigenspace removed from the eigenspace built from
300 images.
The Euclidean distance between the means of the
models produced by each method grows monotonically
as the size of the removed eigenspace falls,
and never exceeds about 1:5 \Theta 10 \Gamma13 gray-level units.
Splitting is slightly less accurate in this respect than
merging.
The mean angular deviation between corresponding
eigenvector directions rises in similar fashion,
from about 0.6 degrees when the size of the removed
eigenspace is 250, to about 1.1 when the removed
eigenmodel is of size 100. This represented a maximum
in the deviation error, because an error of about
degree was obtained when the removed model is of
size 50. Again, these angular deviations are somewhat
larger than those for merging.
The mean difference in eigenvalues shows the same
general trend. Its maximum is about 0.5 units of
gray level, when the size of the removed eigenspace
is 50. This is a much larger error than in the case
of merging, but is still relatively small compared to
a maximum eigenvalue of about 100. As in the case
of merging, the deviation in eigenvalue grows larger
as the size (importance) of the eigenvalue falls.
Difference in reconstruction error rises as the size
of the removed eigenspace falls. Its size is of the
units of gray level per pixel, which again
is negligible.
The difference in likelihoods is significant, the relative
difference some cases being factors of 10 or more.
After conducting further experiments, we found that
this relative difference is sensitive to the errors introduced
when eigenvectors and eigenvalues are dis-
carded. This is not a surprise, given that likelihood
differences are magnified exponentially. We found
that changing the criteria for discarding eigenvectors
very much reduced: relative difference in likelihood
of the order 10 \Gamma14 were achieved in some cases. We
conclude that should an application require not only
splitting, but also require classification, then eigen-vectors
and eigenvalues must be discarded with care.
We suggest that keeping eigenvectors which exceed
some significance threshold be kept.
Overall the trend is clear; accuracy and performance
grew worse, against any measure we used, as
the size of the eigenmodel being removed falls.
6 Conclusion
We have shown that merging and splitting eigenspace
models is possible, allowing a batch of new observations
to be processed as a whole. This theoretical
result is novel. Our experimental results show that
the methods are wholly practical: computation times
are feasible, the eigenspaces are very similar, and the
performance characteristics differ little enough to not
matter to applications.
time advantage is obtained over batch methods
whenever one of the eigenspace models exists already.
A typical scenario is the addition of a set of observations
(a new year's in take of student faces, say) to
an existing, large, database. Our merging method is
even more advantageous when both eigenspace models
exist already. A typical scenario is dynamic clustering
for classification, in which two eigenspace models
can be merged, perhaps to create a hierarchy of
eigenspace models.



--R



An eigenspace update algorithm for image analysis.

Kumar.
Natural basis functions and topographic memory for face recognition.
Probabilistic visual learning for object representa- tion


3rd Edition.

On the generalised karhunen-loeve expansion
--TR

--CTR
Luis Carlos Altamirano , Leopoldo Altamirano , Matas Alvarado, Non-uniform sampling for improved appearance-based models, Pattern Recognition Letters, v.24 n.1-3, p.521-535, January
Ko Nishino , Shree K. Nayar , Tony Jebara, Clustered Blockwise PCA for Representing Visual Data, IEEE Transactions on Pattern Analysis and Machine Intelligence, v.27 n.10, p.1675-1679, October 2005
Jieping Ye , Qi Li , Hui Xiong , Haesun Park , Ravi Janardan , Vipin Kumar, IDR/QR: an incremental dimension reduction algorithm via QR decomposition, Proceedings of the tenth ACM SIGKDD international conference on Knowledge discovery and data mining, August 22-25, 2004, Seattle, WA, USA
Hoi Chan , Thomas Kwok, An Autonomic Problem Determination and Remediation Agent for Ambiguous Situations Based on Singular Value Decomposition Technique, Proceedings of the IEEE/WIC/ACM international conference on Intelligent Agent Technology, p.270-275, December 18-22, 2006
Tae-Kyun Kim , Ognjen Arandjelovi , Roberto Cipolla, Boosted manifold principal angles for image set-based recognition, Pattern Recognition, v.40 n.9, p.2475-2484, September, 2007
Jieping Ye , Qi Li , Hui Xiong , Haesun Park , Ravi Janardan , Vipin Kumar, IDR/QR: An Incremental Dimension Reduction Algorithm via QR Decomposition, IEEE Transactions on Knowledge and Data Engineering, v.17 n.9, p.1208-1222, September 2005
Xiang Sean Zhou , Dorin Comaniciu , Alok Gupta, An Information Fusion Framework  for Robust Shape Tracking, IEEE Transactions on Pattern Analysis and Machine Intelligence, v.27 n.1, p.115-129, January 2005
John P. Collomosse , Peter M. Hall, Cubist Style Rendering from Photographs, IEEE Transactions on Visualization and Computer Graphics, v.9 n.4, p.443-453, October
principal component analysis;eigenspace models;model merging;model splitting;gaussian mixture models