ISR-7
June 196U
XIV.

XIV-1
THE DETERMINATION OF CLUSTERS BY MATRIX ANALYSIS
A. Richard LeSchack

ABSTRACTThe use of clustering techniques related to factor analysis isdiscussed in this section. The application of matrix spectral analysis tothe detection of clusters is based on a matrix model for strong clusteringo
The methods used and the programs developed to implement them on the IBM 7090
computer are described; detailed test results are given; and attempts at
mechanizing the decision procedure are treated.

1.

Introduction
Taxonomy, the science of classification, is an activity fundamental

not only to descriptive biology but to many other areas of scientific andintellectual endeavor as well.

Categorizing individuals - assigning them toone of a tractably small number of subpopulations, the members of which are
sufficiently alike to justify ignoring individual, differences -suppresses

the inessential detail which obscures underlying relationships * The result
is a more economical description of those features of the population which
are truly significant.
A particular case of the classification problem is considered here.
The basic raw data are estimates of the pairwise similarity between membersof the population.

It is desired to identify subsets in which the members are

so similar to one another, and so dissimilar to nonmembers, that it is most

XIV-2

useful to consider them all as a unit. This process has been termed"clumping11 or "clustering."

It is clear that an exact formulation of theclustering problem must be in terms of optimization; for the goal is thebest balance between loss of precision, which is the inevitable consequenceof representing each individual by the typical or average characteristics
of his cluster, and gain in economy, which is the result of replacing many
individuals by a few clusters<>

This general formulation is not pursued

further; we limit our1 selves to those situations of strong clustering inwhich the adherence of individuals to subsets is relatively unambiguous •It is shown that under certain assumptions regarding the nature of thestrong clustering, techniques of matrix spectral analysis will identify
the members of the various clusters. These techniques are extended - with
empirical success, although without rigorous justification - to the generalcase in which the assumption of strong clustering is not so clearly tenable.
A number of techniques for determining clusters have been disussed
in the literature (see Bonner

for a review of many of them):

these include

11"-}
clump theory (Parker-Rhodes and Needham)^ factor analysis (Bonner; Borko
and Bernick)and latent class analysis (Baker).

The present algorithm

bears a close resemblance to factor analysis, but differs from this wellknown statistical technique in that it does not require that the similarity
data be given in terms of the classical product-moment correlation coefficient.
It is based on a different mathematical model (the theory of reducible matrices
of nonnegative elements) and is directed toward a different goal.

XIV-3

A few of the useful shorthand notations of the Iverson programming
language are used i n formulas and flowcharts throughout t h i s section.
reader i s r e f e r r e d to Brooks and Iverson

The

for the general features of theonotation, and to Iverson

for a detailed presentation.

2. A Model for Strong Clustering
The intuitive notion of perfect clustering is best explained byconsidering a population of elements which can be partitioned into nonoverlapping subsets, such that within each subset there is a nonzero linkage
between each pair of elements, but no linkage between a member of the subset
and an element which is not a member.

To generalize to the idea of strong

clustering, assume that we are able to distinguish between strong and weaklinks; the definition is then modified to read:... within each subset there is a strong linkage
between each pair of elements, but at most a weak linkage
between a member of the subset and an element which is nota member.
This part presents the theory of reducible matrices as a model forthe case of perfect clustering, and generalizes the results to the case ofstrong clustering.

The properties of the eigenvectors of such matrices justify

the use of spectral analysis to determine strong clusters.

The exact boundaries

between "strong" and "weak," and the rigorous justification of these techniques in the case of clustering which is not strong, remain open questions.

XIV-U

A square n x n matrix A is reducible (decomposable, uncoupled)
if its rows and columns can be permuted to yield a form partitionable asfollows:

X

;

Y

i

- _T—

0
where 0iZ

is a matrix of all zeros, and X and Z. are square matrices.

More formally, we may state the definition in any of these equivalent
forms:(l)

A is reducible <£=^> there exists a permutation matrix P
such that

T
PAP -

(2)

A is reducible 4-ft there exists a pair of integers (i,j)

such that there exists no chaina• a« • o ci ,1P-L P X P 2
PrJwith all terms nonzero.
(i,j,Pk<L n)(3)

A is reducible <t—> the set {_l,2,...,n} can be partitioned
into mutually exclusive, collectively exhaustive subsets
S and T, such that whenever ifi S and j ^ T, then a. . = 0 e

xiv-5The usual correspondence between a matrix and a directed graph/leads to the conclusion that A is reducible if and only if the corresponding
directed graph is not strongly connected, that is, if there are two nodes
p and q for which no path leads from p to q.
If the matrix A is symmetric, then the definitions take on thisstronger form:

The symmetric matrix A is reducible if and only if there

exists a permutation matrix P such that

r

n-r
r

B = PAPn-rwhere X and Y are square.

Such a matrix is block-diagonal.

that the eigenvectors of B fall into two classes:

It is obvious

those with the last (n-r)

components equal to 0, and those with the first r components equal to 0,
The Jordan canonical formT
B - SASwith(where each £is a normalized eigenvector) reflects the property of

reducibility by taking on the special form
'a.. / 0 in the matrix if and only if an arc connects nodes i and jin the graph.

XIV-6

r

n-r

0 1

-(1)

s=

0

-(2)J n - r

The e i g e n v e c t o r s of A may be expressed a s f o l l o w s :A = PTBP = P T SAS T P = QAQT where Q = PTS
By w r i t i n g -Aj = P -A.P we o b t a i n

A = (PTSP)A!(PTSP)T
which involves only a reordering of the eigenvalues.

It follows thatthe eigenvectors of A fall into two classes: one with zero in r
positions, the other with zero in the remaining (n-r) positions.
Let us new make the further assumption that all elements of A
are nonnegative, and that each of the matrices X and Y is irreducible.
Even more generally, assume that A has been completely reduced to theform

B = PAP1 =

fx

—"
0

0
—
Y

0
—
0

| 0

0

z
•

where each of the diagonal blocks is irreducible.
The eigenvalues of A are clearly those of B; B, however, has
as its set of roots

i.\(x)}U[x(Y)ju [x(z)}u.

XIV-7The classical theorem of Perron and Frobenius, applied to the matrices
X,Y, etc., yields the result that each of the irreducible matrices along
the main diagonal has a positive maximum eigenvalue to which corresponds
an eigenvector, all of whose components are nonnegative.

The other eigen-values are all strictly smaller, and the other eigenvectors do not, in general,
possess the property of nonnegativity.

This fact implies that if the eigen-values of the entire matrix are obtained, one at a time, in descending orderaccording to magnitude, the largest eigenvalue will be the maximum eigenvalueof one of the submatrices. The corresponding eigenvector has positiveelements in the positions corresponding to the members of the index set
identifying the submatrix.
language (see Iverson

Adopting the terminology of the Iverson programmingor Brooks and Iversonfor details), we regard thelogical vector (£ > 0) as a selection vector- for membership in one of theclusters.
The next largest eigenvalue may belong to the same submatrix (inwhich case the corresponding eigenvector either has negative elements, orelse selects the same cluster again), but it is much more likely that thissecond eigenvalue is the maximum root for one of the othersubmatrices; inthis case, its eigenvector selects the elements of a second cluster in thesame manner as the first.

In this way the eigenvectors corresponding to the_.

A logical vector u such that u. = 1 if and only if a* ^ 0.
'In the sense that u selects those elements of the vector (1,2,3*•«•,n)corresponding to its nonzero elements.

XIV-8eigenvalues obtained in descending order of magnitude identify systematically the clusters originally present.
For the binary case, where each a. . = 0 or 1, the submatrix

11...1
1
1... 1

has an eigenvalue r, and an (r-l)-fold root of 0.

If there are p clusters,

each with n. members:ln., + n 0 + .. . + n = n
1 2pthen there are exactly p nonzero eigenvalues:

Vn2,...,np
to which correspond eigenvectors selecting the p clusters.
The fully-reducible matrices just discussed provide a model ofperfect or ideal clustering, where each member of a cluster is linked toother members of the same cluster but to no nonmember of the cluster. To
bridge the gap between this ideal model and those actual cases of strong
clustering in which there are large links and small, but not zero, links,
consider a matrix in which the zero elements of the fully reducible caseare replaced by small quantities £-. The result will be called an
G-reducible matrix.

Its properties are deducible from those of the reduci-

ble matrix by an application of matrix perturbation theory, for a systematic2
development of which the reader is referred to Bellman, pp. 60 ff.

XIV-9Let the original reducible matrix be R.

The matrix P describes

the pattern of perturbations. Thusis the perturbed, or £ -reducible, matrix. Let X. be the eigenvalues of
R; the corresponding normalized eigenvectors (which form a complete
orthonormal set as a consequence of the earlier assumption that R be
symmetric) are x^ 'j finally, let U. be the eigenvalues of P.
are considered:

Two cases

the special case in which R and P commute, and then thegeneral case.
If R and P commute, then they have a common basis; that is, any
eigenvector of one is an eigenvector of the other. From this it follows
immediately that the eigenvectors of the perturbed matrix R/e\ are the sameas those of the original matrix R; to each eigenvalue X. of R corresponds
an eigenvalue

V nof

5(e)The general case is treated by expressing the eigenvalues andeigenvectors of the perturbed matrix as power series in the quantity e*
Retaining only first-order terms, we obtain these results:

eigenvalues of R:

X.+e(x^iy> Rx')neigenvectors of R:

X ^ J C Xa)

XIV-10

where(3)^(1)
C

Mote that these reduce to the results of the special case if the x^ 1 'are also eigenvectors of P.
We conclude that, to first order in 6, the elements of any
eigenvector are changed (if at all) by a term of order e, thus preserving
the structure of the eigenvectors noted above.

3.

Spectral Analysis and Factor Analysis
The key procedure in the computer program for determining clusters

is the extraction of the characteristic roots of a matrix R in decreasingorder of magnitude.

The most appropriate algorithm in this case is theso-called power method, using the Rayleigh quotient at each iteration toestimate the root.

(See Fadeeva^

Chapter 3-)

To describe the process briefly, let us assume that the first p
eigenvalues have been determined. We desire to evaluate X _ . An initial
p+1
vector yis constructed, orthogonal to each of the first p eigenvectors.

This is done in practice by subtracting the vector£(e,x(i))x(i)i-1

XIV-llfrom an initial vector a, each component of which is unity. Here x
is the normalized eigenvector corresponding to \..
is clearly orthogonal to all of the x(i)The resulting y. One now forms the sequences:

(0)

Z
y

(i+1)

(i)
= vRy

8 ^

- * *>ands(1*1) . ( z ( i + 1 ) , y ( i ) )(y (i) ,y (i) )where the usual notation (x,y) for dot product is used. When the sequences has converged to the desired degree of accuracy, we take the final valueof s as X + 1 , and the final vector jr 1 +eigenvector corresponding to X
errors reintroduce into yx,•*•*£, suitably normalized, as the, • In actual computations, numerical rounding' components along one or more of the vectors

5 this may cause the process to converge back to an earlier

eigenvalue (usually X^) instead of X +]_.

Such parasitic behavior is inhi-

bited by occasional reorthogonalization of yto remove any such components.

At the same time, it is desirable to renormalize the y
to prevent possible computer overflow in case X

from time to time,(for example) is quite large.

The basic process, for the case where a fixed number M of roots is to
be extracted, is outlined in Flowcharts 1 and 2.

XIV-12INPUTINITIALIZATIONLOGICAL CHECKS

Logical Checks
Fail<JPRINT OUT. STOP0

H, LOOP FOR \/ THE VARIOUS h\ EIGENVALUES/(

DONE?)/ MAIN\^ ITERATION/LOOPI FAILURE \(TO'< CON VERGE/

SPECIAL FAILURE
PROCEDURE

1RECORD DETAILS
STOP.
)

c

©

RE-DEFINE Y BY
ORTHOGONALIZING WITHRESPECT TO FIRST ( i - l )COLUMNS OF XNORMALIZE THE
RESULTING Y

•0

Not
MAIN
Converged
CALCULATION(NEXT FIGURE)
Converged H/8

Determination of Dominant Eigenvalues
and Corresponding Eigenvectors(Program Organization)
Flowchart 1

XIV-13
Legend to Flowchart 1

ParametersM

number of eigenvalues to be determinedLmaximum number of iterations

MM

how often to reorthogonalize and normalize

Variables

3.

i:
k:

control index; the i t h eigenvalue i s being determinedi t e r a t i o n counters:

previous value of the Rayleigh quotient

ArraysY:current vector iterateX:

the ith normalized eigenvector, as soon as it has been
determined, is stored as the ith column of the matrix X,

The close relation between this algorithm and the principal factor
solution, or method of principal axes, in factor analysis will be discussed
briefly.

In essence, the present algorithm, if applied to a matrix ofproduct-moment correlations, will yield the principal factor solutions to
the factor analysis problem.
Factor analysis, developed originally as a mathematical model for
certain psychological theories of human abilities and behavior, has becomea standard tool of multivariate statistical analysis. Although no attempt
can be made here to discuss the development of methods or their applications f'For this, refer to Harmon.'NORMALIZE YN.
STORE RESULT AS
COLUMN i OF
THE ARRAY X

RECORDRESULTS®~
Determination of Dominant Eigenvalues
and Corresponding Eigenvalues(Main Calculation)
Flowchart 2

xiv-i5
Legend to Flowchart 2

1.

Parameters
TOL:

2.

convergence criterion for successive iterates

Variabless:

previous value of the Rayleigh quotientr:current value of the Rayleigh quotient3• Arrays
A

the matrix to be analyzed

Y

the previous vector iterate

YN

the current vector iterateXthe eigenvector just determined, after normalization,is stored as the ith column of the matrix X.

a brief outline of the fundamental factor analysis problem may make clear
the relation between factor analysis and spectral analysis.
Consider the following:

given the symmetric, nonnegative definite

matrix E(nx n), find an (n x m) matrix A, with m as small as possible,

for whichR« = AAis sufficiently close to R.

T

If we express R in Jordan canonical form as
R - SAST

XIV-16

where S is orthogonal and A. is a diagonal matrix displaying the eigenvaluesof R, ordering the elements of A and the corresponding columns of S (which are
the normalized eigenvectors of R) so that\^ ^2 ^ •' * — \i ^ °then we may write
R = TTT

where
T = SA1/2

andEach matrix AA m ' is formed, clearly, by taking the first m columns of T.
If the matrix R is of rank p (less than n), then at m = p the productreproduces R exactly, and the process terminates. In the usual applicationsof factor analysis, the given matrix R differs slightly, because of sampling
errors, from an ideal matrix R of rank p, less than n (in fact, usuallyp^^.n).

Both p and R are initially unknown. When successive factors cease

XIV-17
to improve the representation, it may be assumed that p has been reached,
and the current Rf is taken as the best estimate of R.
Finally, it will be seen that the method of spectral analysis
described here may be considered as "using factor analysis to determine
clusters" if, and only if, two conditions are satisfied:(1)

the pairwise similarities are expressed by the productmoment correlation, and(2) the factors turn out to have the properties assumed forthe eigenvectors of ^almost separable'1 matrices.

k*

Organization of Programs for Finding Eigenvalues and Eigenvectors
The fundamental matrix and vector manipulations which form theinnermost loops of the program were coded (for the IBM 7090 computer) in
the FAP language, to ensure greater speed and accuracy.

The remaining

subprograms are coded in FORTRAN.- NORM, 0RTH0G, RAQUO, and SMPY are theentries to the FAP-coded vector-matrix package (DOT), based on a double
precision accumulation of dot products.

Subprograms

NORM:
0RTH0G:

normalizes the input vector to length one.
orthogonalizes the input vector with respect to any
desired number of stored normalized vectors.

RAQUO:

computes the Rayleigh quotient of two input vectors,,

SMPY:

does matrix-vector multiplication.

Main:

the calling program which uses this routine to extract
factors, determine clusters, etc.

BOSS:

controls the basic logic of the successive determinationof eigenvalues and vectors.

CALC

performs one iteration of the "main calculation'1

PARAM

inputs control parameters for the run

CHECK

makes logical and consistency checks on parameters andinput values

LOOP:

controls the iterations to find one eigenvalue

FAIL:

attempts corrective action in case iterations fail toconverge

MATRIX:

computes or stores the matrix to be analyzed

Main
Obtain required no.
of e.values
and vectors

BOSS^ PARAM

Input various parameters

3

Logical checks, etc.

CHECK

Input matrix to be factored

to MainFind next eigenvalue
and vector

If first failure, re-adjust
tolerances and try again

Logical Organization Eigenvalue-vector Program
Flowchart 3

XIV-19

j

LOOP

|

Orthogonalize?

Yes^ ORTHOG

Orthogonalize with respect
to previously found e.vectors

No
Another*| iteration?

Yes

->

CALC

One iteration of the power
methodNo
Converged
yet?

Yes

J Return to BOSS

No
Logical Organization
Eigenvalue-Vector Program (inner Loop)
Flowchart k

5.

Program Test
The program was tested by using it to determine, to six-figure

accuracy, the first several eigenvalues and eigenvectors of three matrices
which have been used as examples in the literature of factor analysis. Eachof the matrices was factored once with unities on the main diagonal, oncewith estimates (or exact values) of the communalities, giving six tests inall.

Results agreed in every case with published values. The following

summarizes the tests*

XIV-20
Size

Test Matrix
Harmon, Table 8*5, p . 1^2,

Number of F a c t o r s&

5 factors(unities)

5 factors(communalities)

3 factors(unities)

2 factors(communalities)

3 factors(unities)c l a s s i c a l example i n v o l v i n g i n t e r r e l a t i o n s h i p s among r e s u l t s of ab a t t e r y of p s y c h o l o g i c a l t e s t s

13 X 13forchildrenHaimon, Table £ . 6 , p . 9 1 , a c o n t r i v e d
6X6"textbook example"
Harmon, Table 9 . 1 , p . l 6 i | , an a c t u a lexample i n v o l v i n g c o r r e l a t i o n s among

8X8

v a r i o u s p h y s i c a l measurements

The

3 factors (communalities)

following statistics, for the first of the cases listed above, aretypical of the performance of the program.

Eigen-value Number

Eigenvalue

Number <Df I t e r ;i t i o n s

15.066688721.801387

27

3

h
5

6.

12

0.857069
0.71991U

2k
66

Spectral Analysis to Determine Clusters
Two examples were chosen to test the process of determining

clusters by eigenvalue-eigenvector analysis. The first is based on a

The references are to Harmon.

XIV-21classical example in the literature of factor analysis, involving ananalysis of results of psychological testing.

This is the 1 3 ^ 13 matrix

used in the first two tests described in the previous part. A number ofauthors, studying this matrix and the underlying data by a variety ofmethods, have determined, on the basis of certain patterns of strong correlations, that the 13 tests may be divided into three groups, with themembers of each group strongly interrelated.

It was expected that an

examination of the eigenvectors belonging to the dominant eigenvalues would
reveal this known "clustering" directly.
The second test, not contrived, but actual, involved a matrix ofword-word correlations generated by the document analysis program underdevelopment by the Information Retrieval Group at the Computation Laboratoryof Harvard University.

Certain tendencies - fairly weak, to be sure -toward clustering of the words had been identified using the simple methods
Qof the type described by Lesk in a previous report from this Laboratory <,It was hoped that the present techniques would be capable of detecting
clusters, if any, in this case.
A.

The 13-variable Case (Psychological Tests)
Two runs were made, one using unities, the other using known

estimates of the communalities, on the main diagonal of the matrix. As
expected from theory, the eigenvectors in the two cases were not significantly different.

The values quoted here are for unities on the main

diagonal, since this is felt to be a more realistic testeXIV-22
Fir s t e i g e n v a l u e :5.067

Corresponding e i g e n v e c t o r :

1

0.271

2

0.172

3

0.201

h
5

0.232
0.3UU

6

0.337

7

0.3UO

8

0.333

9

0.331

10

0.205

11

0.2U8

12

0.209

13

0.295

Examination r e v e a l s f i v e elements of approximately the same magnitude( 0 . 3 3 or 0•3U) which a r e l a r g e r than any o t h e r e l e m e n t s .
vector u' = (x^ ' > 0 . 3 3 ) e x t r a c t s t h i su(l)/(l,2,3,...,13)

The l o g i c a lcluster:= (5,6,7,8,9).

We now e x t r a c t t h e second r o o t .
Second e i g e n v a l u e :Corresponding e i g e n v e c t o r :

1.801
1

0.075

2

-O.Oll;3-0.070k
5-0.083-0.176

6-0.2717-0.25kXIV-238-0.110910-0.300

11

O.363

12

0.517

13

0.331*

0.1*1*9It will be noted, first of all, that all the elements of x(2)

corresponding

to the cluster (5,6,7,8,9) are negative and relatively large in magnitude;
this fact further sequesters these five elements and confirms their grouping
as a cluster. There are four relatively large and positive components of
(2)
x

(2)
. The logical vector (x> 0.33) selects this second cluster:(10,11,12,13).

A third root is now extracted:
1 Mh5

Third eigenvalue:
Corresponding eigenvector:

1

0.1*12

2

0.387

3
1*

0.1*915-0.211*

6-0.161*7
8-0.211*

9
10

0.357

-0.052-0.213-0.31111
12-0.02513

0.136-0.169

XIV-21*

The criterion u V ) = (x'3) ^ 0.33) selects the third cluster:u(3)/(l,2,3,...A3) - (1*2,3,h)
The fourth eigenvalue (\, =0.857) is somewhat smaller inmagnitude; the corresponding eigenvector has only one relatively large
component, the second.

This tends to set off element "2" in a classby itself, agreeing with the observation that "2" is more weakly boundin cluster (1,2,3,1;) than the other elements.

The further eigenvaluesof the matrix, and their corresponding eigenvectors, provide no usefulindications.
The clustering

as summarized in Fig. 1 agrees exactly with that

determined previously and reported in the literature.

Links Among the 13 Elements, Showing Clustering
Figure 1

XIV-25

B.

The 2 5 - v a r i a b l e Case (Word C o r r e l a t i o n s )An o r i g i n a l u n e d i t e d t e x t , a s e c t i o n of a j o u r n a l a r t i c l e on thes u b j e c t " e x t e n s i o n s of ALGOL," c o n s i s t e d of 37 s e n t e n c e s , t h e t o t a l l e n g t hbeing somewhat under 1000 words.

P r o c e s s i n g of t h e t e x t

identified

o c c u r r e n c e s of 51 d i s t i n c t t e r m s , each p o t e n t i a l l y " c o n t e n t - r i c h , " whicha r e found i n the a u t o m a t i c t h e s a u r u s of computer-science t e c h n i c a l terms
b u i l t i n t o t h e p r o c e s s i n g system.

\
WordNo\A word-sentence i n c i d e n c e m a t r i x :Sentence No.
1

2

3

...

37

"T~

i r
2
3

I
j

where R. is the number of occurrences of word i in sentence j, describes
the distribution of these words among the sentences in the text.

From thislist of 5>1 terms all were deleted which occurred only once in the text;
after this, any surviving term was deleted if it did not have at least two
co-occurrences with other undeleted terms. At the conclusion of thisselection process, 25 terms remained.

These are listed according to their

thesaurus numbers in Table 1, but the actual terms are not listed to avoid
prejudicing the evaluation of the results.

XIV-26•••—"———•———————

Index
NumberThesaurus
Number

1

16

2

26

3

Thesaurus
NumberIndex
Number

1

IE

'

103
117

27

15
16

u

29

17

13U

5

18

137

19

7
8

35
U9
53
57

21

1U3
1U7
156

9

60

22

178

10

68

23

181

11

77

2U

191

12

80

25

208

13

102

6

121

20

Terms Included in Clustering Study
TABLE 1

The next step is the preparation of a (2£x 2$) term-term
association matrix by row-wise correlation, using the cosine measuredescribed by Salton.

Thus the matrix element
S^ = cos(R1,R;5)

measures the p u t a t i v e s i m i l a r i t y between term i and term j , based on t h e i rtendency to co-occur within sentences of the t e x t .

There were only 116

nonzero terms i n t h i s symmetric (2$ X 25) matrix as shown i n Fig. 2.

XIV-27

Index
Numbers

Similarity
Coefficient

1
1
1
1
2
2
2
2
2
2
3
3
3

6
13
20
22
6
12
13
1U
20
23
7
18
2ii

1*

5

1*
1*
h
U
1*
5
5
5
5
5
5
5
6
6
7

8
9
18
21
22
8
10
li
18
20
21
22
12
23
13

.5000
.2^00
.2673
.5000
.5000

Index
Numbers
7
7
8
8
8
10
10
10
10
10
11
12
12
12
12
13
13
13
13
lli
lli
1)4
15
16
16
18
18
19
20

.S9kk
.5000
.2887

.531*5
.7071
.7071
.5000
.5000
.2860
.17U1
.2132
.2132

.2|26U
.2132
.1826
.1826
.2236
.1*1*72.1195.1*1*72
.2236
.6708
.7071
.1768

SimilarityCoefficient

18
2J4
11
19
21
13
1U

17
19
20
19
13

lU
20
23
lli
19
20
22
17
19
20
16

2k
2$
21

2k
20
22

Nonzero Elements of the (25 X 25) Similarity Matrix S
TABLE 2.3536 "I
.3536
.8165
.1+082
J4O82.201*1.7071.5771;.1*082
.2182
.5000
.2236
.1291
.2390
.91*87.1*330
.2500
.8018
.2500
.8165
.2887
.1*629
.7559
.2673
.2673
.5000
.5000
.2673
.2673

XIV-28

oo
•J
CM

ON

.;1

!

«

H !

*

p

H
CO

H

6

u

0C)

<D
-P
1

u

«• 1©

CD

«

• 4

H

•>

in
CM

*

H
e

41

©

/

« «

M
•H
£H

/

<t

/

a

+3
c\1
r't

<w*

UN
CM

x:

<

/

xrs

+/
/

X
u\

CM

0>
-P
CM
•rl

0)

CM

CO
(D

CO
•H

•H
-P
•H

X!
•P

til

+ ••

•H
10

«

•rl

0)
•H

CO

nal entr

«

+ /

O

be
Cb
•H
-H
CH

•H
4-3
crj

o

o
o
hi

CD

cq

A A A &

• 4+ftXOJ

CM C\J CM n j ^ U \ v O MO o - c o O O H CM CO ro_jrt ~ | U \ l*-ct

ON orH rH rH rH rH rH rH H rH rH H H CMooo o o §vO vT) I— U \ i r \ ON r o r - O 00 l>- O CM <*> O- rH -~t O- CA F - M O OO rH H CO

. i

o
c

£H

l A ^

r

1
r

H

H

Non zero

C

|

OC

•H
TT,

XIV-29
The nonzero s i m i l a r i t i e s range i n magnitude from 0.1195 to 0.91*87,
with ten e n t r i e s equal to or greater thanVo" ( t h a t is,-£(R ,|N) ± h5°) •The matrix S i s now used as input to the eigenvalue-eigenvector
analysis program.

The dominant eigenvalue i s X_ = 3.659 and the corresponding

eigenvector:(1)0.139151(9(2)

0.1^2920

(10

(3)
(1*)
(5)
(6)
(7)
(8)

0.011*859
0.022021
0.051981
0.313962
0.029859
0.01*2727

(11
(12
(13
(H*
(15
(16

0.00178)4
O.176368
0.01*2998
O.U27966
0.315^35
O.2780I4.6
0.000390
0.001321;

(17)
(18)
(19)
(20)
(21)
(22)
(23)
(2U)
(25)

0.12 381*7
0.021*163
0.13U669
0.327092
0.023593
0.09U985
0.3561*21
0.011611
O.OOOI38

Although i t would be d i f f i c u l t to j u s t i f y a c le a r - c ut d i s t i n c t i o n between"large terms" and "small terms," i t was f e l t t h a t a reasonable thresholdwould bet± = 0.278
since there is, after this point, a relatively large gap in the magnitudesof the elements. Then the selection vectoru(l) = U ( l ) > t l € )selects the cluster

CI = (2,6,12,13,11*,20,23).

XIV-30
The second eigenvalue is X- = 2.8U3> with corresponding vector:

(1)
(2)
(3)

(U)
(5)
(6)
(7)
(8)

O.OOU183-0.1696U9
0.172282
0.152761
0.2iiOl|29-0.251M3
0.160U53
0.239106

(9)
(10)
(11)
(12)
(13)

(HO
(15)
(16)

0.017721
0.308153
0.20893U-0.261*186
0.165931
0.298267
0.01121*6
0.02725k

(17)
(18)
(19)
(20)
(21)
(22)
(23)

(210
(25)

0.228352
0.25U1080.27l67ii
0.1697)47
0.2161UU
0.095270-0.297750
0.151109
0.003977

It will be noted first that elements 2,6,12 and 23 are negative,
confirming the adherence of these four items to cluster CI and expressing
their remoteness or disjointness from cluster C2. The setting of a thresholdfor membership in C2 is, once again, somewhat arbitrary. If we taket 2 = 0.2
then

H (2) -(£ (2) ^t 2S )s e l e c t s the cluster
C2 = (5,8,10,11,11;,17,18,19,21).
I t i s of i n t e r e s t that c l u s t e r s CI and C2 are not d i s j o i n t ; element lU i sa member of both subsets.
Had we wished to impose s t r i c t e r c r i t e r i a , in order to obtain smaller
and more t i g h t l y bound c l u s t e r s , we would have chosen\= 0.325 and t 2 = 0.250

XIV-31
to obtain, r e s p e c t i v e l y ,CI = (2,12,20,23) and £2 = (lO,lU,l8,19)
which are now disjoint.
The analysis was halted at this point; a third root was notextracted.
These results have been compared with those obtained by theclustering routines currently incorporated in the document analysis
system.

The first three clusters of interest (these are actually No. 1,

No. 3> and No. 5 of the computer output) include the following items (inaddition to a few others, not listed, which had not been retained in our
set of 25).

Kl = (8,10,11,Ik,17,19)
K2 = (1,2,6,12,13,20,22,23)

K3 - (3,7,18,210
There can be no doubt of the identification of CI with K2; in fact, thesimilarity between the selection vector u

for cluster CI and the corre-

sponding vector v^ ' associated with K2, iss - cos(u^,v( 2 )) = 0.802
using the cosine measure of similarity.

To answer the question of whether

this agreement is to be regarded as significant, we compute (from the truncatedhypergeometric distribution) the probability of a chance agreement this highor higher:

XIV-32

p = 0.001
Now C2 may be identified with Kl. In this case the similarity measure iss = 0.816$
with an even higher significance level (p = 0.0008).
The set K3 has only a single item (number 18) in common with theunion of CI and C2; this points to a third cluster lying in the complementof CI U C2.
One final test was undertaken to examine the validity of theclusters selected by the eigenvector analysis. This was a hand simulationof a simple clustering technique based on methods discussed by Needham
and by Parker-Rhodes and Needham.

The first step is sorting the nonzero

elements of the matrix J3, as shown in Table 3 •Only the final results will be quoted here. Incorporating links
of strength 0,5 and greater, we build the following clusters:(2,6,12,23)(13,20)(10,1U,17)
The correspondence with clusters CI and C2 is seen immediately.
It may now be concluded that the clusters identified by eigenvector
analysis correspond closely to those detected by the other approaches.
The attempt to seek a final confirmation by determining whether the actual
terms associated with the thesaurus entry numbers are conceptually related

IndexNumbers

5
12

h

7
5
5
10

U
h
1*
10
5
5
12
12
1
13
13
1
16
16
19
20

h

2
111
7
7
8

20
H(
8
13
8
10
13
9
18
22
20
11
22
13
20
13
19
22
20

2l*
25
20
22
5
Hi
19
18
2l*
19

Similarity
Coefficient

Index
Numbers

8
10
1*
13
5
5
lli
1
1
2
2
3
3
11
18
18
2
10
6
2
3
6
10
15
13
8
111
2
12

•1195
.1291
.171*1
.1768
.1826
.1826
.2010.
.2132
.2132
.2132
.2182
.2236
.2236
.2236
.2390
.2500
.2500
.2500
.2673
.2673
.2673
.2673
.2673
.2860
.2887
.2887
.3536
.3536
.1*082

Similarity
Coefficient

21
19
21
111
18
21
20
6
22
6
13
18
2)4
19
21

2l|
20
17
12
23
7
23
111
16
20
11
17
12
23

Nonzero Elements of S in Ascending Order
TABLE 3

,1|082.1*082.ll261l.1*330.Iili72.10*72.1*629
.5000
.5000
.5000
.5000
.5000
.5000
.5000
.5000
.5000

.53hS
.577a
.6708
.7071
.7071
.7071
.7071

.7559
.8018
.8165
.8165
.891*1*.91*87

XIV-3U

reveals no particular relation among the items in any of the clusters.
This is, however, no criticism of the methods used to infer the presenceof clusters from the correlation matrix, but merely an indication thatthe raw correlation data themselves (the matrix R from which S was derived)
probably express little more than chance co-occurrence of words insentences.
The results of the clustering attempts are displayed in Figs.
2 and 3> each of which may be regarded as a portion of a labeled symmetricgraph showing the pairwise correlations among the 25 items. The graphis drawn to represent the strength of a link by the relative nearness ofthe vertices and/or the thickness of the line joining them, according
to these conventions:s > 0.7071=

O . ^ s ^ 0.7071=

s^.0.5
An examination of these graphs suggests a few more remarks.
Cluster CI may seem somewhat unsatisfactory because it includes item 1)4,
but breaks the strong links which join "ll|!l with ".10" and "17".

Had thethreshold t, been lowered only slightly (to t- = 0.176), then item 10
would have been included along with "II4.". On the other hand, had t beenraised slightly to 0.3, then item Ik would have dropped out, leaving amore compact cluster(2,6,12,13,20,23)
which would have been even better correlated with cluster K2.

XIV-35

PQ

o
CO

U
CO

rH

o
CO

ti PL
•N

eg

cd

«

e -^
0)
-P
HI

Pi
O
-P
1A O
OJ 0)

>
0)

C

^

v

CO

CJ
o
•H

-P
0)

o
o

£

XIV-36

fcl

pq
xJ
G
3

o
Pn
(0

^
<1>
P>
to

3

3
M

+>

•H

a

8H
ti

•H

^ -d

to

0

e
(1) ,o
p
M

uC\J\
CD

•c

O
10
<D

P
(0

,£3 TJ
-P 0

bo
£

5

Q)

o X

1

CO

a

o

•H
-P
a}

H

0)

fc
fc
o
o

c<\

XIV-37

Similar considerations apply to cluster C2. A slight lowering ofthreshold t~ would have picked up items 3, 7, and 2k (which are strongly
linked with "18").

7. Evaluation and Conclusions
An evaluation of the clustering technique presented here must
depend on these three questions:(1) Are clusters found which are intuitively reasonable,
and which match clusters found by other methods?
(2) Is the algorithm rapid enough to be of any practicalsignificance?
(3) Can the decision procedure implicit in the selectionof the ,flarge elements" of the normalized eigenvectors
be mechanized fully?
The first question is answered affirmatively, it is felt, by thetheoretical analysis of Part 2 and the results presented in Parts 5 and 60
The tests have not been extensive enough to permit a direct evaluationof the running time of the clustering programs. The actual running times
were so short, in fact, that timing was not feasible„

A theoretical estimateof the efficiency of the algorithm is presented here.
It is reasonable to assume that the bulk of the time will be
expended in the actual determination of eigenvalues and eigenvectorsc

The

computing time for this process depends critically on the rate of convergenceof the basic iterative method.

This, in turn, is a function of the distri-bution of the eigenvalues — in other words, of the structure of the data

XIV-38

matrix under analysis. Assume that the average number of iterations is
M.

In the course of each iteration, (n+2) dot products are computed, eachof which involves n multiplications and additions.

For the first eigen-value we expect, then, an operation count proportional toMn(n +2);

for succeeding eigenvalues, a somewhat greater count (because of orthogonalization)•

Assuming this to be absorbed in an adjusted factor M f , werequire an operation count proportional tokM'n(n + 2)
to extract k eigenvalues and eigenvectors. Both k and M ! depend strongly
on the structure of the matrix: k, the number of eigenvectors required toapproximate the structure adequately, being equal to the approximate rankof the matrix; M 1 , on the other hand, being related to the separation
between eigenvalues, but not necessarily to its size*

Therefore, all things

being equal, the computing time should increase asn(n + 2)p

or approximately as n , with increase in n.
This estimate compares favorably with running times for procedures
involving complete factor analysis^- or latent class analysis, which are
3 For the processing of large matrices,
expected to grow at least as fast as n..
however, it appears that if truly practical methods are ever found, their
computing time should increase no faster than n.

XIV-39
Using the eigenvectors of the similarity matrix to determine
membership in clusters reduces what is basically a two-dimensional problem
to a one-dimensional one. Each of the components of a normalized eigenvector may be regarded as a measure of the adherence of the corresponding
element to the cluster associated with that vector.

Thus a scan over thetwo dimensions of the original matrix (to find the large elements whichindicate that certain members of the index set belong together) is replaced
by a one-dimensional scan to determine cluster membership. This does not,however, eliminate the element of human intervention or personal choice inthe process of deciding which elements are large enough to warrant inclusionin the cluster; in other words, in deciding what the cut-off point should
be between the r,largeM elements and the nsmall" elements of the eigenvector.
An automatic decision procedure which has been proposed is describedin the following paragraphs, and the results of a few tests, by no means
extensive or definitive, are now presented.
The two steps in the procedure are the following:(1) permute the normalized eigenvector x into a vector ^, where

Zl > l2 > 13 - • • • - *nby a sorting process, and(2) find that value of k which maximizeskIAF(k) = i £0<p<1

where only nonnegative y_. are e l i g i b l e for i n c l u s i o n .

XIV-liO
Since F(k) first increases, then decreases, with increasing k,
the maximum is determined sequentially by evaluating F(l), F(2), etc.
until F no longer increases. The parameter p, which may be considereda shaping factor, influences the cut-off between "large*1 and "small"elements.

At the extremes, p = 0 will cause inclusion of all elements;

p = 1 will include only the first (or the first q, if £= y_ = ... = y ).

The method works well, as expected, for contrived cases in which
the distinction between large and small elements is clear-cut andunambiguous. Here, however, much simpler procedures could be programmed.
A more critical test has been based on attempts to reproduce the resultsof the 13 x 13 case described in Part 6.
The maximum attained by F(k), for each of the five values of p
used in the test, is -underlined in Table U.

The corresponding cluster

would contain all elements above the horizontal line. Thus, forp = 0.8

the same cluster is selected as before:(5,6,7,8,9).
When the procedure was applied to the second normalized eigenvectorof the 13 X 13 test matrix, it was found that the value p = 0.8 did notselect precisely the same cluster as our previous methods, although aslightly larger value of p did. It is not known whether one value of theparameter p can be found which will make the correct selection for all threeeigenvectors of this matrix. Further tests will be necessary to settle

XIV-1*1

8

2

3

1*

5

6.

7

O.3W1
0.31*00.11830.1183

O.11830.1183

O.II83 125
70.23390.1651*

0.11*73

0.131*3

0.1170

3

6

0.337

0.3U750.20060.1671

0.11*1*3

0.1158

1*

8

0.333

0.1*581*0.22920.18190.15120.111*6

5

9

0.331

0.56800.251*00.191*3

0.15670.1136

6

13

0.295

0.6550O.267I1

0.1981*0.15620.1092

7

1

0.271

0.728U

0.27530.19910.15360.101*1811

0.21*80.79000.27930.19750.11*91*

0.0988

9

1*

0.232

0.81*370.28120.191*6

0.11*55

0.0937

10

12

0.209

0.887k

0.28060.19120.11*06

0.0887

11

10

0.205

0.92950.28020.18750.13650.081*5

12

3

0.201

0.96990.28000.18500.13280.0808

13

2

0.172

0.99990.27720.18080.1281;0.0769

1
I

1

j

1
Leger id
Colurrin 1 :c u r r e n t index i

2:o r i g i n a l index

3:o r i g i n a l component of normsi l i z e d e i g e n v e c t o r

1*:

F(k), with p = 0

5":

P = 0.5

6:

p = 0.67

7:

p = 0.88:

p = 1.0
Normalized and Sorted Vecto:ITS
TABLE h

XIV-U2
this question, and to explore the more general problems of the properchoice of p and the efficacy of the algorithm in general.
The difficulty here in identifying the "large11 components of avector is felt to be, in large part, a reflection of the more basic
difficulty in deciding cluster membership in those cases where theclustering cannot be regarded as strong.

This is a problem inherent

in the data and in the criteria adopted for cluster selection.

REFERENCES
1.

Baker, F. B., "Information Retrieval Based on Latent Class Analysis,"
JACM, Vol. 9, No. k (October 1962).

2.

Bellman, R., Introduction to Matrix Analysis, McGraw-Hill, New York(I960).

3.

Bonner, R. E., "On Some Clustering Techniques," IBM Journal of
Research and Development, Vol. 8, No. 1 (January 1961|) .

U.

Borko, H., and Bernick, M. D., "Automatic Document Classification,"
JACM, Vol. 10, No. 2 (April 1963).

5.

Brooks, F. P., Jr., and Iverson, K. E., Automatic Data Processing,Wiley, New York (1963).

6.

Faddeeva, V. N., Computational Methods of Linear Algebra (Translated
by Benster, C. D.), Dover Publications, New York (1959).

7.

Harmon, H. H., Modern Factor Analysis, University of Chicago Press(I960).

8.

Iverson, K. E., A Programming Language, Wiley, New York (1962).

9.

Lesk, M., "Attempts to Cluster Documents with Citation Data,"
Information Storage and Retrieval, Report ISR-3, Sec. VI, The
Computation Laboratory of Harvard University (1963) -XIV-U3
10.

Needhain, R. M., "A Method for Using Computers in Information Classification," Proc. IFIP Congress 1962, Vol. VIII, No. 3 (1962)
pp. 121-123.

11.

Parker-Rhodes, A. F., and Needham, R. M., "The Theory of Clumps: A
New Concept of Classification and Selection," Report ML 126, Cambridge
Language Research Unit, Cambridge, England (i960).

12.

Salton, G. A., "Associative Document Retrieval Techniques Using
Bibliographic Information," JACM, Vol. 10, No. k (1963), pp. hkP-UZl *