IEEE Systems, Man and Cybernetics Magazine - January 2021 - 28

Xt k = UYk V T = UU T X k VV T ), the reconstruction square error
2
for the kth image is X k - UU T X k VV T F and the reconstruction mean square error (MSE) is given by
1
MSE = M

	

M

| iX

k

- UU T X k VV T i2F , (7)

k=1

where < $ < F is the Frobenius norm. For MSE in (7), its
square root also appears in the literature, denoted here as
root-mean-square error (RMSE).
In CSA, the iterative process begins by initializing either
the matrix U as some m # m orthogonal matrix and then
retaining the top q columns (eigenvectors) so that the initial
matrix U = U 0 is m # q , or by initializing the matrix V as
an orthogonal n # n matrix and retaining the top p columns (eigenvectors) so the initial matrix V = V0 is n # p.
Using the initial V or U, either scatter matrix S V or S U is
calculated as
1
SV = M

	

M

|X VV
k

T

X Tk (size m # m) (8)

k=1

Algorithm 1. CSA (with U initialized).
Input	
M -centered images X k , k = 1,f, M ; number of
columns p and q of V and U to retain, respectively;
initial orthogonal m # m matrix*, truncated to
size m # q, denoted U 0; error change threshold
f (2 0); maximum iterations t max.
Step	Description
1	
Initialize iteration number t = 1.
2	Compute S U,t = ^1/M h R kM= 1 X kT U t - 1 U Tt - 1 X k .
3	Compute eigenvectors  v i  corresponding to the
p (p 1 n) largest eigenvalues of  S U,t .
Vt = 6v 1 v 2 f v p@; (i th column v i of Vt is the
4	
eigenvector corresponding to the i th largest
eigenvalue of S U,t ).
5	Compute S V,t = ^1/M h R kM= 1 X k Vt V Tt X Tk .
6	Compute eigenvectors u i corresponding to the q
(q 1 m) largest eigenvalues of S V,t .
U t = 6u 1 u 2 f u q@; (i th column u i of U t is the
7	
eigenvector corresponding to the ith largest
eigenvalue of S V,t ).
8	Compute reconstruction RMSE t .†
If t = 1, go to step 9. Else if error change‡
RMSE t -1 - RMSE t # f or maximum iterations
(t = t max), go to step 10. Else continue with step 9.
9	
Increment iteration number: t ! t + 1; go to step 2.
10	
Halt. Output final U = U t ^m # q h,  V = Vt ^n # p h.
The nomenclature and presentations of this algorithm differ notably
in the literature, but the steps above are representative of the core
elements.
*The algorithm with V initialized (denoted V0) is similar, except that steps
5, 6, and 7 would be swapped with steps 2, 3, and 4, respectively; V and
VT subscripts in new step 2 become " t-1, " and U and UT subscripts in new
step 5 become " t. " An orthogonal matrix of proper size (often the identity
matrix, as in [5] and [6]) is specified, then the columns are truncated to
form the initial matrix V0 (or U0).
† RMSE t =

1
M

M

| iX

k

- U t U t X k V t V t iF .
T

T

2

k =1

‡ Theoretically, RMSEt−1  -  RMSEt is monotonic decreasing [6]; however,
|RMSEt−1  -  RMSEt| corrects for the case where a round-off error may
result in an apparent small increase in error change at an iteration.

28	

IEEE SYSTEMS, MAN, & CYBERNETICS MAGAZINE Janu ar y 2021

	

1
SU = M

M

|X

T
k

UU T X k (size n # n). (9)

k=1

For a given V, performing eigendecomposition of S V and
retaining q of the m eigenvectors corresponding to the q largest eigenvalues of S V leads to a new m # q matrix U, where
the eigenvector corresponding to the ith largest eigenvalue
of S V is the ith column of new matrix U. In a similar way,
given U, S U is formed as in (9). Next, eigendecomposition of
S U leads to a new matrix V, where the ith column of V is the
eigenvector of S U corresponding to the ith largest eigenvalue of S U . This process continues, generating a new U and V
for each iteration. This process is summarized as Algorithm 1, for the case where U is initialized.
Distance Measures for Classification
In 2DPCA, B2DPCA, and CSA, a single nearest neighbor
classifier (1-NN) is often used in image recognition applications. The category of the sample image with minimum
distance to the probe image (to be identified) in feature
space, the single nearest neighbor, is assigned to the
probe image. The " distance " d (Yg , Yk) between the feature
matrix Yg of the probe image and the feature matrix Yk of
the kth sample image, often defined for B2DPCA (see [4],
[13]) and CSA [16], is given by
	

d ^Yg, Yk h = Yg - Yk F , (10)

(where once again < $ < F is the Frobenius norm); the value of
k for which d (Yg , Yk) is minimum identifies the nearest
neighbor and corresponding category. We note that in [3], a
slightly different distance is defined for row 2DPCA:
d (Yg , Yk) = R ip= 1 iy ig - y ik i2 , where y ig and y ik are the ith columns of the probe and kth sample image feature matrices
( p columns), respectively. A more general distance measure, the assembled matrix distance, is considered for
B2DPCA in [14].
Some Pros and Cons of PCA
Variants Versus PCA
Previously we noted that the feature matrix for 2DPCA
often has significantly more entries (or coefficients) than
the feature vector in PCA, while B2DPCA and CSA typically have much smaller feature matrices than 2DPCA, and
the number of B2DPCA or CSA coefficients can be comparable to the number of coefficients in the PCA feature vector. In image recognition, where a single nearest neighbor
classifier is often used, computation of distance in feature
space between probe image (to be classified) and each
sample image requires the calculation of M distances since
there are M sample images. Having more feature matrix
coefficients (compared to the PCA feature vector) leads to
increased time determining the best match for the probe
image. This is a disadvantage of 2DPCA.
In PCA, the centered grayscale images with m rows
and n columns of pixels represented as vectors x k of
IEEE Systems, Man and Cybernetics Magazine - January 2021

Table of Contents for the Digital Edition of IEEE Systems, Man and Cybernetics Magazine - January 2021
