ARTICLE IN PRESSDTD 5

Neural Networks xx (2005) 1–17
www.elsevier.com/locate/neunet

An incremental network for on-line unsupervisedclassification and topology learning
Shen Furaoa,*, Osamu Hasegawab,caDepartment of Computational Intelligence and Systems Science, Tokyo Institute of Technology
R2-52, 4259 Nagatsuta, Midori-ku, Yokohama 226-8503, Japanb
Imaging Science and Engineering Lab., Tokyo Institute of TechnologycPRESTO, Japan Science and Technology Agency (JST), Japan
Received 12 December 2003; accepted 11 April 2005

Abstract
This paper presents an on-line unsupervised learning mechanism for unlabeled data that are polluted by noise. Using a similarity thresholdbased and a local error-based insertion criterion, the system is able to grow incrementally and to accommodate input patterns of on-line nonstationary data distribution. A definition of a utility parameter, the error-radius, allows this system to learn the number of nodes needed tosolve a task. The use of a new technique for removing nodes in low probability density regions can separate clusters with low-density
overlaps and dynamically eliminate noise in the input data. The design of two-layer neural network enables this system to represent the
topological structure of unsupervised on-line data, report the reasonable number of clusters, and give typical prototype patterns of every
cluster without prior conditions such as a suitable number of nodes or a good initial codebook.
q 2005 Elsevier Ltd. All rights reserved.
Keywords: On-line unsupervised learning; Stationary environment; Non-stationary environment; Clustering; Topology representation

1. Introduction
One objective of unsupervised learning is construction ofdecision boundaries based on unlabeled training data.
Unsupervised classification is also known as data clusteringand is defined as the problem of finding homogeneousgroups of data points in a given multidimensional data set(Jain & Dubes, 1988). Each of these groups is called acluster and defined as a region in which the density ofobjects is locally higher than in other regions.
Clustering algorithms are classifiable into hierarchicalclustering, partitional clustering, fuzzy clustering, nearestneighbor clustering, artificial neural networks for clustering,etc. (Murty, Jain, & Flynn, 1999). Hierarchical clustering
algorithms such as single-link (Sneath & Sokal, 1973),* Corresponding author. Tel.: C81 45 924 5180; fax: C81 45 924 5175.
E-mail addresses: furaoshen@isl.titech.ac.jp (S. Furao), hasegawa@isl.
titech.ac.jp (O. Hasegawa).
URL: http://www.isl.titech.ac.jp/hasegawalab/shen/index.htm.

0893-6080/$ - see front matter q 2005 Elsevier Ltd. All rights reserved.
doi:10.1016/j.neunet.2005.04.006

complete-link (King, 1967), and CURE (Guha, Rastogi, &Shim, 1998) usually find satisfiable clustering results but
suffer from computational overload and the requirement formuch memory space (Murty et al., 1999). Hierarchicalclustering algorithms are therefore unsuitable for large data
sets or on-line data. BIRCH is an extremely efficienthierarchical clustering algorithm (Zhang, Ramakrishnan, &Livny, 1996), but is properly applicable to data setsconsisting only of isotropic clusters: a two-dimensional(2D) circle, or a three-dimensional (3D) sphericity, etc.
Specifically, chain-like and concentric clusters are difficult
to identify using BIRCH (Guha et al., 1998; Oyang, Chen, &Yang, 2001).
Most partitional clustering algorithms run in linear time
and work on large data sets (Lin & Chen, 2002). The
k-means algorithm, a conventionally used partitional
clustering algorithms, suffers from deficiencies such as
dependence on initial starting conditions (Lozano, Pena, &Larranaga, 1999) and a tendency to result in local minima.
Likas, Vlassis, and Verbeek (2003) proposed a global
k-means algorithm, an incremental approach to clusteringthat dynamically adds one cluster center at a time through a

DTD 5ARTICLE IN PRESS

2

S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

deterministic global search consisting of N (the data set size)
execution of the k-means algorithm from suitable initial
positions. Compared to traditional k-means, this algorithm
can obtain equivalent or better results, but it suffers fromhigh computation load. The enhanced LBG algorithm
proposed by Patane and Russo (2001) defines one
parameter-utility of a codeword-to overcome the drawbackof LBG algorithm: the dependence on initial starting
condition. The main difficulties of such methods are how
to determine the number of clusters k in advance and thelimited applicability to data sets consisting only of isotropic
clusters.
Some clustering methods combine features of hierarchical and partitional clustering algorithms, partitioning an
input data set into sub-clusters and then constructing ahierarchical structure based on these sub-clusters (Murty &Krishna, 1981). Representing a sub-cluster as only one
point, however, renders the multilevel algorithm inapplicable to some cases, especially when the dimension of a subcluster is on the same order as the corresponding final
cluster (Lin & Chen, 2002).
Another possible objective of unsupervised learning canbe described as topology learning: given a high-dimensional data distribution, find a topological structure thatclosely reflects the topology of the data distribution. Selforganizing map (SOM) models (Willshaw & von der
Malsburg, 1976; Kohonen, 1982) generate mapping fromhigh-dimensional signal space to lower-dimensional
topological structure. The predetermined structure andsize of Kohonen’s model imply limitations on resulting
mapping (Martinetz & Schulten, 1994). Methods that
identify and repair topological defects are costly(Villmann, Der, Hermann, & Martinetz, 1997). A posterior
choice of class labels for prototypes of the (unsupervised)
SOM causes further problems: class borders are not taken
into account in SOM and several classes may sharecommon prototypes (Villmann, Schleif, & Hammer,
2003). As an alternative, a combination of competitive
Hebbian learning (CHL) (Martinetz, 1993) and neural gas(NG) (Martinetz, Berkovich, & Schulten, 1993) is
effective in constructing topological structure (Martinetz& Schulten, 1994). For each input signal, CHL connects
the two closest centers by an edge, and NG adapts knearest centers whereby k is decreasing from a large initialvalue to a small final value. Problems arise in practicalapplication: it requires an a priori decision about network
size; it must rank all nodes in each adaptation step;furthermore, once adaptation strength has decayed, the useof adaptation parameters ‘freezes’ the network, which
thereby becomes unable to react to subsequent changes insignal distribution. Two nearly identical algorithms areproposed to solve such problems: growing neural gas(GNG) (Fritzke, 1995) and dynamic cell structures(Bruske & Sommer, 1995). Nodes in the network compete
for determining the node with the highest similarity to the
input pattern. Local error measures are gathered during

the learning process to determine where to insert newnodes, and new node is inserted near the node with thehighest accumulated error. The major drawbacks of thesemethods are their permanent increase in the number ofnodes and drift of the centers to capture the input
probability density (Hamker, 2001). Thresholds such as a
maximum number of nodes predetermined by the user, oran insertion criterion depending on overall error or on
quantization error are not appropriate, because appropriate
figures for these criteria are unknown a priori.
For the much more difficult problems of non-stationary
data distributions, on-line learning or life-long learning
tasks, the above-mentioned methods are not suitable. The
fundamental issue for such problems is how a learning
system can adapt to new information without corrupting orforgetting previously learned information—the StabilityPlasticity Dilemma (Carpenter & Grossberg, 1988). Using autility-based removal criterion, GNG-U (Fritzke, 1997)
deletes nodes that are located in regions of low input
probability density. GNG-U uses a network of limited size
to track the distribution in each moment, and the target of
GNG-U is to represent the actual state. Life-long learning(Hamker, 2001) emphasizes learning through the entire
lifespan of a system. For life-long learning, the ‘dead nodes’
removed by GNG-U can be interpreted as a kind of memorythat may be useful again, for example, when the probability
distribution takes on a previously held shape. Therefore,
those ‘dead nodes’ preserve the knowledge of previous
situations for future decisions and play a major role in lifelong learning. The GNG-U serves to follow a non-stationary
input distribution, but the previously learned prototype
patterns are destroyed. For that reason, GNG-U is unsuitable
for life-long learning tasks.
Lim and Harrison (1997) propose a hybrid network thatcombines advantages of Fuzzy ARTMAP and probabilisticneural networks for incremental learning. Hamker (2001)
proposes a life-long learning cell structure (LLCS) that isable to learn the number of nodes needed to solve a current
task and to dynamically adapt the learning rate of each node
separately. Both methods work for supervised on-line
learning or life-long learning tasks, but how to process
unsupervised learning remains controversial.
The goal of the present study is to design an autonomous
learning system for unsupervised classification and topology representation tasks. The objective is to develop a
network that operates autonomously, on-line or life-long,and in a non-stationary environment. The network grows
incrementally, learns the number of nodes needed to solve a
current task, learns the number of clusters, accommodates
input patterns of on-line non-stationary data distribution,and dynamically eliminates noise in input data.
Below, we describe the proposed algorithm, and then use
artificial and real-world data sets to illustrate the learning
process and observe details. A comparison with typicalincremental networks, GNG and GNG-U, elucidates thelearning behavior.

DTD 5ARTICLE IN PRESS
S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

2. Proposed method
For unsupervised on-line learning tasks, we separate
unlabeled non-stationary input data into different classes
without prior knowledge such as how many classes exist.
We also intend to learn input data topologies. Briefly, thetargets of the proposed algorithm are:
† To process on-line or life-long learning non-stationarydata.
† Using no prior conditions such as a suitable number ofnodes or a good initial codebook or knowing how many
classes exist, to conduct unsupervised learning, report asuitable number of classes and represent the topologicalstructure of input probability density.
† To separate classes with low-density overlap and detect
the main structure of clusters that are polluted by noise.
To realize these targets, we emphasize the key aspects oflocal representation, insertion of new nodes, similarity
threshold, adaptive learning rate, and deletion of low
probability density nodes.
2.1. Overview of the proposed method
In this study, we adopt a two-layer neural network
structure to realize our targets. The first layer is used togenerate a topological structure of input pattern. We obtainsome nodes to represent the probability density of an input
pattern when we finish the first-layer learning. For thesecond layer, we use nodes identified in the first layer as the
input data set. We report the number of clusters and give
typical prototype nodes of every cluster when we finish thesecond-layer learning.
For unsupervised classification task, we must determine
if an input sample belongs to previously learned clusters or
to a new cluster. Suppose we say that two samples belong to
the same cluster if the Euclidean distance between them isless than threshold distance T. If T is too large, all samples
will be assigned to one cluster. If T is too small, each sample
will form an isolated, singleton cluster. To obtain ‘natural’clusters, T must be greater than the typical within-cluster
distance and less than the typical between-cluster distance(Duda, Hart & Stork, 2001).
For both layers, we must calculate the thresholddistance T. In the first layer, we set the input signal as a
new node (the first node of a new cluster) when the
distance between the signal and the nearest node (or thesecond nearest node) is greater than a threshold T that ispermanently adapted to the present situation. In the second
layer, we calculate the typical within-cluster distances andtypical between-cluster distances based on those nodes
generated in the first layer, then give a constant thresholddistance Tc according to the within-cluster distances andbetween-cluster distances.

3

To represent the topological structure, in on-line orlife-long learning tasks, growth is an important featurefor decreasing task error and adapting to changing
environments while preserving old prototype patterns.
Therefore, the insertion of new nodes is a very usefulcontribution to plasticity of the Stability-Plasticity
Dilemma without interfering with previous learning(stability). Insertion must be stopped to prohibit apermanent increase in the number of nodes and to
avoid overfitting. For that reason, we must decide whenand how to insert a new node within one cluster andwhen insertion is to be stopped.
For within-cluster insertion, we adopt a scheme used bysome incremental networks (such as GNG (Fritzke, 1995),
GCS (Fritzke, 1994)) to insert a node between node q withmaximum accumulated error and node f, which is among theneighbors of q with maximum accumulated error. Current
incremental networks (GNG, GCS) have no ability to learn
whether further insertion of a node is useful or not. Nodeinsertion leads to catastrophic allocation of new nodes.
Here, we suggest a strategy: when a new node is inserted,
we evaluate insertion by a utility parameter, the errorradius, to judge if insertion is successful. This evaluation
ensures that the insertion of a new node leads to decreasing
error and controls the increment of nodes, eventuallystabilizing the number of nodes.
We adopt the competitive Hebbian rule proposed byMartinetz in topology representing networks (TRN) to build
connections between neural nodes (Martinetz & Schulten,
1994). The competitive Hebbian rule can be described as:
for each input signal, connect the two closest nodes(measured by Euclidean distance) by an edge. It is provedthat each edge of the generated graph belongs to the
Delaunay triangulation corresponding to the given set ofreference vectors, and that the graph is optimally topologypreserving in a very general sense. In on-line or life-longlearning, the nodes change their locations slowly but
permanently. Therefore, nodes that are neighboring at an
early stage might not be neighboring at a more advanced
stage. It thereby becomes necessary to remove connections
that have not been recently refreshed.
In general, overlaps exist among clusters. To detect the
number of clusters precisely, we assume that input data areseparable: the probability density in the centric part ofevery cluster is higher than the density in intermediateparts between clusters; and overlaps between clusters have
low probability density. We separate clusters by removing
those nodes whose position is in a region with very low
probability density. To realize this, Fritzke (1994) designed
an estimation to find the low probability-density region; if
the density is below threshold h, the node is removed.
Here, we propose a novel strategy: if the number of input
signals generated so far is an integer multiple of aparameter, remove those nodes with only one or no
topological neighbor. We infer that, if the node has only
one or no neighbor, during that period, the accumulated

DTD 5ARTICLE IN PRESS

4

S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

error of this node has a very low probability of becoming
maximum and the insertion of new nodes near this node isdifficult: the probability density of the region containing
the node is very low. However, for one-dimensional input
data, the nodes form a chain, the above criterion willremove the boundary nodes repeatedly. In addition, if the
input data contain little noise, it is not good to delete those
nodes having only one neighbor. Thus, we use anotherparameter, the local accumulated number of signals of thecandidate-deleting node, to control the deletion behavior. If
this parameter is greater than an adaptive threshold, i.e. if
the node is the winner for numerous signals, the node will
not be deleted because the node does not lie in a lowdensity area. This strategy works well for removing nodes
in low-density regions without added computation load. In
addition, the use of this technique periodically removes
nodes caused by noise because the probability density ofnoise is very low.
If two nodes can be linked with a series of edges, a pathexists between the two nodes. Martinetz and Schulten(1994) prove some theorems and reach the conclusion thatthe competitive Hebbian rule is suitable for forming a pathpreserving representations of a given manifold. If thenumber of nodes is sufficient for obtaining a dense
distribution, the connectivity structure of the network
corresponds to the induced Delaunay triangulation thatdefines both a perfectly topology-preserving map and apath-preserving representation. In TRN, neural gas (NG) isused to obtain weight vectors Wi2M, iZ1,.,N,M is thegiven manifold. The NG is an efficient vector quantization
procedure. It creates a homogeneous distribution of theweight vectors Wi on M. After learning, the probability
distribution of the weight vectors of the NG isrðWi Þ wPðWi Þa

A constant similarity threshold is used in the second layer
instead of the adaptive similarity threshold used in the firstlayer.
Notations to be used in the algorithm
A node set, used to store nodes
NA number of nodes in A
C connection set (or edge set), used to store
connections (edges) between nodes
NC number of edges in C
Wi n-dimension weight vector of node i
Ei local accumulated error of node i; it is updated whennode i is the nearest node (winner) from the input
signalMi local accumulated number of signals of node i; thenumber is updated when node i is the winner
Ri inherited error-radius of node i; the error-radius ofnode i is defined by the mean of accumulated error,
Ei/Mi. Ri serves as memory for the error-radius ofnode i at the moment of insertion. It is updated ateach insertion, but only for affected nodes
Ci cluster label. This variable is used to judge whichcluster node i belongs to.
Q number of clusters.
Ti similarity threshold. If the distance between an
input pattern and node i is larger than Ti, the input
pattern is a new node.
Ni set of direct topological neighbors of node i.
Li number of topological neighbors of node i.
age(i,j) age of the edge that connects node i and node j.
path given a series of nodes xi2A, iZ1,2,.,n, makes (i,x1),(x1,x2),.,(xnK1,xn),(xn,j)2C. We say that a‘path’ exists between node i and node j.

(1)

With the magnification factor aZ Deff =ðDeff C 2Þ, theintrinsic data dimension is Deff (Villmann, 2000). Therefore,
with the path-preserving representation, the competitive
Hebbian rule allows the determination of which parts of agiven pattern manifold are separated and form different
clusters. In the proposed algorithm, we adopt a scheme likeneural gas with only nearest neighbor learning for theweight vector adaptation; the competitive Hebbian rule isused to determine topology neighbors. Therefore, we canuse the conclusions of TRN to identify clusters, i.e. if two
nodes are linked with one path, we say the two nodes belong
to one cluster.
2.2. Complete algorithm

Algorithm 2.1: basic algorithm(1) Initialize node set A to contain two nodes, c1 and c2:
A Z fc1 ; c2 g(2)with weight vectors chosen randomly from the input
pattern. Initialize connection set C, C3A!A, to theempty set
C ZF

(3)
n

(2) Input new pattern x2R .
(3) Search node set A to determine the winner s1, andsecond-nearest node (second winner) s2 by
s1 Z arg min kxKWc k

(4)

c2A

Using the analysis presented in Section 2.1, we give thecomplete algorithm here. The same algorithm is used totrain both the first layer and the second layer. The difference
between the two layers is that the input data set of thesecond layer is the nodes generated by the first layer.

s2 Z arg min jxKWc jj(5)

c2Anfs1 g

If the respective distances separating x and s1 or s2 aregreater than similarity thresholds Ts1 or Ts2 , the input

ARTICLE IN PRESS

DTD 5

S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

signal is a new node; add the new node to A and go tostep (2) to process the next signal, i.e. if kxKWs1 kO
Ts1 or kxKWs2 kO Ts2 , then AZAgr and WrZx
(4) If a connection between s1 and s2 does not exist
already, create it and add it to connection set C.
C Z Cg ðs1 ; s2 Þ(7)(5) Increase the age of all edges emanating from s1.
ageðs1 ;iÞ Z ageðs1 ;iÞ C 1ðci 2Ns1 Þ(8)(6) Add the Euclidian distance between the input signal
and the winner to local accumulated error Es1
Es1 Z Es1 C kxKWs1 k:(9)(7) Add 1 to the local accumulated number of signals Ms1 :
Ms1 Z Ms1 C 1:

Rr Z a3 ðRq C Rf Þ(18)

† Decrease accumulated error variables of q and f by
fraction b
Eq Z bEq ; Ef Z bEf

(19)

(6)

Set the age of the connection between s1 and s2 to zeroageðs1 ;s2 Þ Z 0

5

(10)

(8) Adapt the weight vectors of the winner and its direct
topological neighbors by fraction 31(t) and 32(t) of thetotal distance to the input signal† Decrease the accumulated number of signal
variables of q and f by fraction g
Mq Z gMq ;Mf Z gMf(20)† Judge whether or not insertion is successful. If theerror-radius is larger than the inherited error-radius
Ri (ci2{q,r,f}), in other words, if insertion cannot
decrease the mean error of this local area, insertion
is not successful; else, update the inherited errorradius. I.e. if Ei/MiORi (ci2{q,r,f}), insertion is
not successful, new node r is removed from set A,
and all parameters are restored; else, RqZEq/Mq,
RfZEf/Mf, and RrZEr/Mr.
† If insertion is successful, insert edges connecting
new node r with nodes q and f, and remove theoriginal edge between q and f.

DWs1 Z 31 ðtÞðxKWs1 Þ(11)

C Z Cg fðr; qÞ; ðr; f Þg(21)

DWi Z 32 ðtÞðxKWi Þ ðci 2Ns1 Þ(12)

C Z Cnfðq; f Þg

(22)Here, we call 31(t) the learning rate of the winner, and
32(t) the learning rate of the neighbor.
(9) Remove edges with an age greater than a predefined
threshold agedead, i.e. if (i,j)2C, and age(i,j)O
agedead(ci,j2A), then CZC\{(i,j)}.
(10) If the number of input signals generated so far is an
integer multiple of parameter l, insert a new node andremove nodes in low probability density as follows:† Determine node q with maximum accumulated
error E:q Z arg max Ec

(13)

c2A

† Determine, among the neighbors of q, node f withmaximum accumulated error:
f Z arg max Ec

(14)

c2Nq

† Add new node r to the network and interpolate its
weight vector from q and f:
A Z Ag frg; Wr Z ðWq C Wf Þ=2:0(15)† Interpolate accumulated error Er, accumulated
number of signals Mr, and the inherited errorradius Rr from Eq, Ef, Mq, Mf, and Rq, Rf by:
Er Z a1 ðEq C Ef Þ

(16)

Mr Z a2 ðMq C Mf Þ

(17)

† For all nodes in A, search for nodes having only
one neighbor, then compare the accumulated
number of signals of these nodes with the averageaccumulated number of all nodes. If a node has
only one neighbor and the accumulated number ofsignals is less than an adaptive threshold, remove itfrom the node set, i.e. if LiZ1 (ci2A) andNAP
Mi ! c Mj =NA , then AZA\{i}. Here, c isjZ1determined by the user and 1RcO0. If muchnoise exists in the input data, c will be larger andvice versa.
† For all nodes in A, search for isolated nodes, thendelete them, i.e. if LiZ0 (ci2A), then AZA\{i}.
(11) After a long constant time period LT, report the
number of clusters, output all nodes belonging todifferent clusters. Use the following method to classify
nodes into different classes:† Initialize all nodes as unclassified.
† Loop: Randomly choose one unclassified node ifrom node set A. Mark node i as classified and label
it as class Ci.
† Search A to find all unclassified nodes connected tonode i with a ‘path.’ Mark these nodes as classified
and label them as the same class as node i (Ci).
† If unclassified nodes exist, go to Loop to continue
the classification process until all nodes areclassified.

DTD 5ARTICLE IN PRESS

6

S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

(12) Go to step (2) to continue the on-line unsupervised
learning process.

2.3. Parameter discussion
We determine parameters in Algorithm 2.1 as follows.
2.3.1. Similarity threshold Ti of node i
As discussed in Section 2.1, similarity threshold Ti is avery important variable. In step (3) of Algorithm 2.1, thethreshold is used to judge if the input signal belongs topreviously learned clusters or not.
For the first layer, we have no prior knowledge of input
data and therefore adopt an adaptive threshold scheme. For
every node i, the threshold Ti is adopted independently.
First, we assume that all data come from the same cluster;thus, the initial threshold of every node will be CN. After aperiod of learning, the input pattern is separated into
different small groups; each group comprises one node andits direct topological neighbors. Some of these groups canbe linked to form a big group; the big groups are thenseparated from each other. We designate such big groups asclusters. The similarity threshold must be greater than
within-cluster distances and less than between-cluster
distances. Based on this idea, we calculate similarity
threshold Ti of node i with the following algorithm.
Algorithm 2.2: Calculation of similarity threshold T for the
first layer.
(1) Initialize the similarity threshold of node i to CN whennode i is generated as a new node.
(2) When node i is a winner or second winner, update
similarity threshold Ti by† If the node has direct topological neighbors (LiO0),
Ti is updated as the maximum distance between
node i and all of its neighbors.
Ti Z max kWi KWc k
c2Ni

(23)

† If node i has no neighbors (LiZ0), Ti is updated as
the minimum distance of node i and all other nodes
in A.
Ti Z min kWi KWc k(24)

c2Anfig

For the second layer, the input data set contains theresults of the first layer. After learning of the first layer, we
obtain coarse clustering results and a topological structure.
Using this knowledge, we calculate the within-cluster
distance dw as
1 X
kWi KWj k;
(25)
dw Z
NC ði;jÞ2C
and calculate the between-cluster distance db(Ci,Cj) ofclusters Ci and Cj as

db ðCi ; Cj Þ Z max kWi KWj k(26)

i2Ci ;j2Cj

That is, the within-cluster distance is the mean distanceof all edges, and the between-cluster distance is theminimum distance between two clusters. With dw and db,
we can give a constant threshold Tc for all nodes.
The threshold must be greater than the within-cluster
distance and less than the between-cluster distances.
Influenced by overlap or noise, some between-cluster
distances are less than within-cluster distances, so we use
the following algorithm to calculate the threshold distanceTc for the second layer.
Algorithm 2.3: Calculation of similarity threshold T for thesecond layer.
(1) Set Tc as the minimum between-cluster distance.
Tc Z db ðCi1 ; Cj1 Þ Z

mink;lZ1;.;Q;ksl

db ðCk ; Cl Þ

(27)(2) If Tc is less than within-cluster distance dw, set Tc as thenext minimum between-cluster distance.
Tc Z db ðCi2 ; Cj2 Þ Z

mink;lZ1;.;Q;ksl;ksi1 ;lsj1

db ðCk ; Cl Þ

(28)(3) Go to step (2) to update Tc until Tc is greater than dw.

2.3.2. Adaptive learning rate
In step (8) of Algorithm 2.1, the learning rate 3(t)determines the extent to which the winner and the neighborsof the winner are adapted towards the input signal.
A constant learning rate is adopted by GNG and GNG-U,
i.e. 3(t)Zc0, 1Rc0O0. With this scheme, each reference
vector Wi represents an exponentially decaying average of
those input signals for which the node i has been a winner.
However, the most recent input signal always determines afraction c0 of the current Wi. Even after numerous iterations,the current input signal can cause a considerable change inthe reference vector of the winner.
Ritter et al. (1991) proposed a decaying adaptation
learning rate scheme. This scheme was adopted in the neuralgas (NG) network (Martinetz et al., 1993). The exponential
decay scheme is3ðtÞ Z 3i ð3f =3i Þt=tmax(29)

wherein 3i and 3f are the initial and final values of thelearning rate and tmax is the total number of adaptation steps.
This method is less susceptible to poor initialization. For
many data distributions, it gives a lower mean-square error.
However, we must choose parameter 3i and 3f by ‘trial anderror’ (Martinetz et al., 1993). For on-line or life-longlearning tasks, the total number of adaptation steps tmax is
not available.
In this study, we adopt a scheme like k-means to adapt
the learning rate over time by

DTD 5ARTICLE IN PRESS
S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

1
1
31 ðtÞ Z ; 32 ðtÞ Z
t
100t

(30)

Here, time parameter t represents the number of input
signals for which this particular node has been a winner thus
far, i.e. tZMi. This algorithm is known as k-means. The
node is always the exact arithmetic mean of the input signals
it has been a winner for.
This scheme is adopted because we are hopeful ofmaking the position of the node more stable by decreasing
the learning rate when the node becomes a winner for more
and more input patterns. After the network size becomesstable, the network is fine tuned by stochastic approximation (Wasan, 1969). This approximation denotes anumber of adaptation steps with a strength 3(t) decaying
NP
slowly, but not too slowly, i.e.
3ðtÞZN, and
N
tZ1
P
2
3 ðtÞ!N. The harmonic series, Eq. (30), satisfies the7

the rate of decrease becomes: a1Z1/6, a2Z1/4, a3Z1/4,
bZ2/3, and gZ3/4.
The above analysis is based on the supposition thatsignals are uniformly distributed over Voronoi regions. This
supposition might be untenable for some tasks. For thatreason, the above parameter set is not the optimal choice forsuch tasks. Fortunately, the choice of these parameters is notsensitive. In Section 3, we use different parameter sets totest an artificial data set and achieve nearly identical results.
With the same parameter set, the system also works well fordifferent real-world data experiments such as face recognition and vector quantization for different images.

3. Experiment3.1. Artificial data set

tZ1

conditions.
2.3.3. Decreasing rate of accumulated variables E (error)
and M (number of signals)
When insertion between node q and f happens, how do
we decrease accumulated error Eq, Ef and the accumulated
number of signals Mq, Mf? How do we allocate theaccumulated error Er and the accumulated number ofsignals Mr to new node r?
In Fig. 1, the Voronoi regions of node q and f before
insertion are shown at left, and Voronoi regions belonging
to q, f, and new node r after insertion are shown at right.
We assume that signals in these Voronoi regions are
distributed uniformly. Comparing left to right reveals that
one-fourth of the accumulated number of signals of q and f
are reallocated to new node r, and that three-fourths of theaccumulated number of signals remaining for q and f. For
accumulated error, it is reasonable to assume that error
attributable to V1 is double the error caused by V2 for nodeq. Consequently, after insertion, the accumulated error of q
and f is two-thirds of the accumulated error before
insertion. We also assume that the error to r caused by
V1 is equal to the error to q attributable to V2.
Consequently, the reallocated accumulated error for r is
one-sixth of Eq and Ef. Based on the above analysis,Fig. 1. Voronoi regions of nodes (Left, before insertion; Right, after
insertion).

We conducted our experiment on the data set shown in
Fig. 2. An artificial 2-D data set is used to take advantage of
its intuitive manipulation, visualization, and resulting
insight into system behavior. The data set is separated into
five parts: A, B, C, D, and E. Data sets A and B satisfy 2-D
Gaussian distribution. The C and D data sets are a famous
single-link example. Finally, E is sinusoidal and separated
into E1, E2, and E3 to clarify incremental properties. We
also add random noise (10% of useful data) to the data set tosimulate real-world data. As shown in Fig. 2, overlaps exist
among clusters; noise is distributed over the entire data set.
As stated in (Guha et al., 1998), neither BIRCH nor singlelink clustering can correctly partition such a data set. The
CSM algorithm proposed by Lin and Chen (2002) partiallysolves this problem, but it cannot eliminate noise, and noise
forms new clusters.
This experiment compares our proposed method withtypical incremental network GNG (Fritzke, 1995) and
GNG-U (Fritzke, 1997) to demonstrate the advantages ofthe proposed method. In all experiments, we set parameters
as lZ100, agedeadZ100, and cZ1. For GNG and GNG-U,
the maximum number of nodes is predefined as 300;

Fig. 2. 2D artificial data set used for the experiment.

DTD 5ARTICLE IN PRESS

8

S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

Fig. 3. Randomly input samples from a stationary environment. GNG
results: 1 cluster.

Fig. 5. Randomly input samples from a stationary environment. Results ofproposed method: second layer, five clusters.

the parameter sets in (Fritzke, 1995) and (Fritzke, 1997) are
adopted for other parameters.

environment changes and patterns from area B are chosen.
At step 40,001 the environment changes again, etc. Table 1
details specifications of the test environments. The environment changes from I to VII. In each environment, areas used
to generate patterns are marked as ‘1’; other areas are marked
as ‘0.’ For each environment, we add 10% noise to test dataand noise is distributed over the whole data space.
Before we test our proposed algorithm, we use GNG and
GNG-U to conduct this experiment. We show the last results
for GNG and GNG-U in Figs. 6 and 7, and do not reportintermediate results. These results show that GNG cannot
represent the topological structure of non-stationary on-line
data well. GNG-U deletes all old learned patterns and only
represents the structure of new input patterns. Neither
method eliminates noise. In GNG-U results, for example,
nodes beyond area E3 are all attributable to noise distributed
over the whole space (Fig. 7).
Fig. 8 shows the first-layer results of the proposed method.
After learning in one environment, we report intermediate
topological structures. Environments I and II test isotropic
data sets A and B, and an overlap exists between A and B. In
environment I, the data set of the 2-D Gaussian distribution istested. This system represents the Gaussian distribution
structure very well. In environment II from 20,001 to 40,000
steps, probability changes to zero in area A. Remaining
nodes of area A, often called ‘dead nodes,’ play a major rolein on-line or life-long learning. They preserve the knowledge

3.1.1. Experiment in a stationary environment
First, we use Fig. 2 as a stationary data set to compare ouralgorithm with GNG; 100,000 patterns are chosen randomly
from areas A, B, C, D, and E. Topological results of GNG
and the proposed method are shown in Figs. 3–5. For
stationary data, GNG can represent the topologicalstructure, but it is affected by noise and all nodes are linked
to form one cluster (Fig. 3). The first-layer of the proposed
method partially eliminates the influence of noise andseparates original data set to some different clusters (Fig. 4).
The second layer of the proposed method efficiently
represents the topology structure and gives the number ofclusters and typical prototype nodes of every cluster (Fig. 5).
3.1.2. Experiment in non-stationary environment
We simulate on-line learning through the use of thefollowing paradigm: from step 1 to 20,000, patterns arechosen randomly from area A. At step 20,001, theTable 1
Experiment environments for on-line learning
Area

Fig. 4. Randomly input samples from a stationary environment. Results ofproposed method: first layer, three clusters.

A
B
C
D
E1
E2
E3

Environment
I

II

III

IV

V

VI

VII

1
0
0
0
0
0
0

0
1
0
0
0
0
0

1
0
1
0
0
0
0

0
1
0
1
0
0
0

0
0
0
1
1
0
0

0
0
1
0
0
1
0

0
0
0
0
0
0
1

DTD 5ARTICLE IN PRESS
S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

Fig. 6. Sequentially input samples from non-stationary environments. GNG
results.

of previous situations for future decisions. In environment
III, 40,001 to 60,000 steps, the reappearance of area A does
not raise error and knowledge is preserved completely.
Therefore, nearly no insertion happens and most nodes
remain at their original positions. Environments III and IV
test a difficult situation (data sets such as area C and area D).
The system works well, removing noise between areas andseparating C and D. Environments V, VI, and VII test acomplicated artificial shape (area E). E is separated into threeparts and data come to the system sequentially, notrandomly. We find nodes of area E increasing following
the change of environment from V to VI and VII, but allnodes are linked to form the same class.
Fig. 8 is used as the input data set of the second layer anda constant similarity threshold based on Fig. 8 is calculated.
Fig. 9 is the result of the second layer. It reports the numberof clusters and gives prototype nodes of every cluster.
Figs. 4 and 8 show that the first layer of the proposed
method processes stationary and non-stationary environments well. It detects the main structure from original data,
which is polluted by noise. It controls the number of nodes
needed for the current task. However, some nodes generated

Fig. 7. Sequentially input samples from non-stationary environments.
GNG-U results.

9

by noise remain in results of the first layer, and the nodes are
too many to represent the topology structure efficiently.
Figs. 5 and 9 show that the second layer not only removes
those nodes caused by noise; it also removes someredundant nodes and rearranges the position of nodes torepresent the topology structure efficiently. We say that,using the second layer, the system can ‘understand’ original
raw input data very well.
Fig. 10 shows how the number of nodes changes during
on-line learning (first-layer). The number of nodes increases
when an input signal comes from a new area (see Table 1,environment changes from I to VII). In the sameenvironment, after some learning steps, the increase innodes stops and the number of nodes converges to a constantbecause further insertion cannot engender decreased error.
Noise leads the system to frequently insert and delete nodes.
For that reason, a small fluctuation exists in the number ofnodes in Fig. 10.
Finally, we perform some experiments to test thesensitivity of parameters. Four different parameter sets for{a1,a2,a3,b,g} are used to test the artificial data set (Fig. 2):{1/4, 1/4, 1/4, 2/3, 3/4}, {1/4, 1/4, 1/4, 2/3, 1/2}, {1/6, 1/4,
1/4, 1/2, 1/2} and {1/6, 1/4, 1/4, 1/4, 3/4}. The experiments
are done in stationary and non-stationary environments. In
both environments, for all four parameter sets, the system
can detect the five clusters from the original data set andgive prototype nodes; the last results are nearly the same as
those shown in Fig. 5 for a stationary environment and thoseshown in Fig. 9 for a non-stationary environment.
Consequently, the experiment results indicate that the
choice of parameter set {a1,a2,a3,b,g} is not influential.
In the following, we only use the original parameter set{1/6, 1/4, 1/4, 2/3, 3/4} to do different real-world
experiments and to check whether or not the same parameter
set is suitable to different tasks.
3.2. Real-world data set
Two real-world data sets are used to test the proposed
method. For clustering, 10 classes are taken from ATT_FACE image database to train the proposed network system.
For topology learning, we do vector quantization (VQ) for
Lena (512!512!256) and Boat (512!512!256) to test if
the proposed method can represent the topological structure
well. In both experiments, we compare the proposed methodwith typical incremental networks GNG and GNG-U.
3.2.1. Test performance of clustering: application of facerecognition
For application of clustering with the proposed method,
we use some facial images as the input signal of thefirst layer. The input facial images are taken from the
ATT_FACE database (http://www.uk.research.att.com/).
The database comprises 40 distinct subjects with 10
different images of every subject. For some subjects, theimages were taken at different times with various lighting,DTD 5ARTICLE IN PRESS

10

S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

Fig. 8. Sequentially input samples from non-stationary environments. Results of the proposed method: first layer.

DTD 5ARTICLE IN PRESS
S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

11

facial expressions (open/closed eyes, smiling/not smiling)
and facial details (glasses/no glasses). The size of each
image is 92!112 pixels, with 256 grey levels per pixel. We
selected 10 subjects (Fig. 11(a)) from the database to test
proposed network system. Fig. 11(b) portrays the different
images of the first subject in Fig. 11(a). Feature vectors ofsuch images are taken as follows: first, the original image
with size 92!112 is re-sampled to 23!28 image using thenearest neighbor interpolation method. Then Gaussiansmoothing is used to smooth the 23!28 image with
Gaussian widthZ4, sZ2 to obtain the 23!28 dimensional
input feature vectors (Fig. 12).
The test consists of two processes: learning andevaluation. During learning, we choose vectors randomly(stationary) or sequentially (more difficult task, nonstationary on-line) from an original image vector set for
the first layer, then input the results of the first layer to thesecond layer, and report number of clusters and prototype
vectors of every cluster as the results of the second layer. In
both layers, we set parameters lZ25 and agedeadZ25. In
the first layer, c is set to 1.0 because we want to delete nodes
lying in overlap area; for the second layer, c is set as 0.05 to

focus on clustering and to avoid deleting too many nodes.
For both stationary and non-stationary environments, we do
10,000 training iterations for the first layer and 2000 training
iterations for the second layer.
In the evaluation process, nearest neighbor is used toclassify vectors in the original image vector set to theclusters learned through the learning process. The correct
recognition ratio is reported as the performance index oflearning process.
For a stationary environment, i.e. vectors are randomly
chosen from the original data set, the learning process reportsthat 10 clusters exist in the input patterns and gives prototype
vectors of every cluster. The evaluation process reports thatthe correct recognition ratio is 90%. In addition, GNG is tested
with parameter lZ25 and amaxZ25, maximum number ofnodes is predefined as 50, other parameters are the same as in
Fritzke (1995), and 10,000 samples are used to train GNG. The
GNG results show four clusters in the original data set.
Consequently, GNG is unsuitable for this task.
For a non-stationary on-line environment, i.e. 10,000
vectors are sequentially taken from original data set (oneclass by one class, 1000 samples for each class), the learning
process reports 10 clusters and gives prototype vectors of
every cluster (Fig. 13). Comparison of Fig. 13 with Fig. 12
reveals that the proposed algorithm reports the number ofclusters correctly. It gives the prototype vectors of every
cluster reasonably. The evaluation process reports that thecorrect recognition ratio is 86%. Furthermore, GNG-U istested with parameter lZ25 and amaxZ25. The maximum
number of nodes is predefined as 25 and other parameters
are the same as in (Fritzke, 1997); 10,000 samples are input
sequentially to GNG-U (one class by one class, 1000
samples for each class). The GNG-U results show twoclusters in the original data set. That result indicates that
GNG-U is unsuitable for this task.
The experimental results show that the proposed method
works well for unsupervised classification (learning thenumber of clusters and reporting prototype vectors of every
cluster). Some other typical incremental networks (such as
GNG, GNG-U) are unsuitable for this task.

Fig. 10. Number of nodes during on-line learning (Environment I–
Environment VII).

3.2.2. Test performance of topology learning: application ofvector quantization
One goal of competitive learning systems is theminimization of distortion error. Vector quantization (VQ)
(Linde et al., 1980) is a typical application in which error
minimization is important. In vector quantization, data aretransmitted over limited bandwidth channels by transmitting, for each data vector, only the index of the nearest
reference vector. The set of reference vectors, the so-called
codebook, is assumed to be known both to the sender andreceiver. Therefore, the receiver can use the transmitted
indexes to reconstruct an original vector with thecorresponding reference vector. Information loss occurs
because of the difference of the current data vector andFig. 9. Sequentially input samples from non-stationary environments.
Results of the proposed method: second layer, five clusters.

DTD 5ARTICLE IN PRESS

12

S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

Fig. 11. Facial image (a) 10 subjects, (b) 10 images of subject one.

Fig. 12. Feature vector (a) vector of Fig. 11(a), (b) vector of Fig. 11(b).

DTD 5ARTICLE IN PRESS
S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

13

the nearest reference vector. The main problem in VQ is
how to find a codebook such that the information loss issmall.
This experiment uses vector quantization to code a wellknown image: Lena (512!512!256). The proposed
algorithm is used to determine the codebook. Different
from traditional methods, the codebook size is notpredefined. Instead, it is determined by the algorithm itself.
Because GNG or GNG-U is unsuitable for sequential input
vectors, this experiment only compares the proposed
method with GNG in a stationary environment, i.e.
randomly input vectors from the original vector data set.
For a non-stationary on-line environment, we only reportresults of the proposed method.
In image compression applications, the compression
ratio and the peak signal to noise ratio (PSNR) are often
used to evaluate the compression algorithm performance.
Here, we use bit per pixel (bpp) to measure the compression
ratio, and PSNR is defined as
PSNR Z 10 log10 N
P

2552

(31)

ðf ðiÞKgðiÞÞ =N
2

iZ1

where f(i) is the grey level of the original image, g(i) is thegrey level of the reconstructed image, and N is the total
number of pixels in the image. For VQ, with the same bpp,
higher PSNR means less information loss. Consequently,
the corresponding codebook is better.
In the experiment, first, original image Lena (Fig. 14) is
separated into non-overlapping image blocks of 4!4. These
blocks are input patterns (16 dimensional vectors) of thefirst layer. In the first layer, parameter lZ50, agedeadZ50,
and cZ0 (for that we assume that all blocks are usefulFig. 13. Sequentially input samples from a non-stationary environment.
Results of proposed method: 10 clusters.
Fig. 14. Original image of Lena (512!512!256), 8 bpp.

DTD 5ARTICLE IN PRESS

14

S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17
Table 2Comparison of the proposed method and GNG for Lena in stationary
environment

The first layer
GNG
The second layer
GNG

Number of nodes

bpp

PSNR

130
130
52
52

0.45
0.45
0.34
0.34

30.79
29.98
29.29
28.61

and no block can be deleted as noise). From the original
input vector data set, 200,000 samples are taken randomly.
After we finish the learning of the first layer, we set thenodes generated in the first layer as the codebook (first-layer
codebook) and find the index of the nearest reference vector
for all original vectors. Then we use such indexes andreference vectors in the first-layer codebook to reconstruct

the Lena image. Fig. 15 is the reconstructed image of thefirst-layer codebook.
Then, we set the size of the first-layer codebook as thepredefined maximum number of nodes of GNG and use
GNG to learn the codebook. During learning, parameters
used by GNG are lZ50, amaxZ50. Other parameters are
identical to those in Fritzke (1995), and 200,000 samples areprocessed. Fig. 16 depicts the reconstructed image of the
GNG codebook.
For the second layer, the first-layer codebook is set as the
input data set of the second layer. With the same parameter
set as the first layer, i.e. lZ50, agedeadZ50, and cZ0, 10,000 samples are taken randomly from the first-layer
codebook to train the second layer. We obtain the secondlayer codebook and use this second-layer codebook to
encode and reconstruct the Lena image. With the size of
second-layer codebook, GNG is also used to obtain acodebook. We use this GNG codebook to encode andreconstruct the Lena image, and calculate the PSNR.
Table 2 summarizes a comparison of the proposed method
and GNG in a stationary environment. For both the firstlayer and the second layer, with the same compression ratio,

Fig. 16. Randomly input samples from a stationary environment.
Reconstructed image of Fig. 14 GNG codebook with 130 nodes,
0.45 bpp, PSNRZ29.98 dB.

Fig. 17. Sequentially input samples from a non-stationary environment.
Reconstructed image of Fig. 14 the first-layer codebook with 499 nodes,
0.56 bpp, PSNRZ32.91 dB.

Fig. 15. Randomly input samples from stationary environment. Reconstructed image of Fig. 14 the first-layer codebook with 130 nodes, 0.45 bpp,
PSNRZ30.79 dB.

DTD 5ARTICLE IN PRESS
S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

15

Table 3Comparison of the proposed method and GNG for Boat in stationary
environment

First layer
GNG
The second layer
GNG

Fig. 18. Sequentially input samples from a non-stationary environment.
Reconstructed image of Fig. 14 the second-layer codebook with 64 nodes,
0.375 bpp, PSNRZ29.66 dB.

the proposed method obtains higher PSNR than GNG: theproposed method provides a better codebook than GNG.
From experiments in Sections 3.1.2 and 3.2.1, we knowthat GNG and GNG-U are unsuitable for an on-line nonstationary environment. Therefore, we only test the proposed
method in an on-line non-stationary environment. The
original vector data set is separated into five subsets. These
subsets are input sequentially to the first layer, with 40,000
samples taken from each subset. The reconstructed results ofFig. 19. Original image of Boat (512!512!256), 8 bpp.

Number of nodes

bpp

PSNR

225
225
62
62

0.487
0.487
0.375
0.375

30.05
29.45
28.13
27.59

the first layer (Fig. 17) and the second layer (Fig. 18) showthat the proposed method is useful to learn the codebook ofLena well in an on-line non-stationary environment.
With the same parameter set as the Lena image, we also
compare the proposed method with GNG for anotherwell-known image-Boat (512!512!256) (Fig. 19)-in astationary environment. Table 3 shows the results: theproposed method obtains a better codebook than GNG for
the Boat image. In the on-line non-stationary environment,
Fig. 20 is the reconstructed image of the first-layer
codebook, and Fig. 21 is the reconstructed image of thesecond layer codebook. Both Lena and Boat experiments
show that, for topology learning, the proposed method
works well for stationary and non-stationary environments.

4. Conclusion
In this study, we proposed a new on-line learning method
for unsupervised classification and topology representation.
Using a similarity threshold-based and a locally accumulated error-based insertion criterion, the system can growFig. 20. Sequentially input samples from a non-stationary environment.
Reconstructed image of Fig. 19 the first-layer codebook with 557 nodes,
0.575 bpp, PSNRZ31.46 dB.

DTD 5ARTICLE IN PRESS

16

S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

Fig. 21. Sequentially input samples from a non-stationary environment.
Reconstructed image of Fig. 19 the second-layer codebook with 81 nodes,
0.4 bpp, PSNRZ28.54 dB.

incrementally and accommodate input patterns of on-line
non-stationary data distribution. A novel on-line criterion
for node removal in low probability-density regions enables
this system to separate clusters with low-density overlap andto dynamically eliminate noise from input data. The utility
parameter ‘error-radius’ is used to judge if the insertion is
successful and to control the increase of nodes. Adoption ofa two-layer neural network makes it possible for this system
to ‘understand’ the original data set well. In summary, thealgorithm can cope with difficulties in online or life-longunsupervised learning such as overlap, never-seen inputs,
temporarily non-appearing patterns, and noise.
Other problems remain unsolved. For example, when ahigh-density overlap pertains between clusters, it isextremely difficult for the proposed algorithm to separate
clusters from each other. Three parameters must be
determined by the user: l, agedead, and c. The difficulty of
automatically determining such parameters is based on thefact that, for different tasks, the optimal choice of suchparameters is different: it is difficult to give a standard ofsuch parameters for every task. Although these parameters
are not so sensitive, we remain hopeful that some methods
are useful to automatically deduce optimal choice of suchparameters for the task. Such problems will be addressed insubsequent studies.

References
Bruske, J., & Sommer, G. (1995). Dynamic cell structure learns perfectly
topology preserving map. Neural Computation, 7, 845–865.

Carpenter, G. A., & Grossberg, S. (1988). The ART of adaptive patternrecognition by a self-organizing neural network. IEEE Computer, 21,
77–88.
Duda, R. O., Hart, P. E., & Stork, D. G. (2001). Pattern classification (2nded.). New York: Wiley–Interscience.
Fritzke, B. (1994). Growing cell structures—a self-organizing network
for unsupervised and supervised learning. Neural Networks, 7,1441–1460.
Fritzke, B. (1995). A growing neural gas network learns topologies.
Advances in neural information processing systems (NIPS) pp.
625–632.
Fritzke, B. (1997). A self-organizing network that can follow
non-stationary distributions Proceedings of ICANN-97 pp.
613–618.
Guha, S., Rastogi, R., & Shim, K. (1998). CURE: An efficient clustering
algorithm for large databases. Proceeding of ACM SIGMOD
Conference on Management of Data , 73–84.
Hamker, F. H. (2001). Life-long learning cell structures-continuously
learning without catastrophic interference. Neural Networks, 14,
551–573.
Jain, A. K., & Dubes, R. C. (1988). Algorithms for clustering data.
Englewood Cliffs: Prentice Hall.
King, B. (1967). Step-wise clustering procedures. Journal of American
Statistical Association, 69, 86–101.
Kohonen, T. (1982). Self-organized formation of topologically correctfeature maps. Biological Cybernetics, 43, 59–69.
Likas, A., Vlassis, N., & Verbeek, J. J. (2003). The global k-meansclustering algorithm. Pattern Recognition, 36, 451–461.
Lim, C. P., & Harrison, R. F. (1997). A incremental adaptive network foron-line supervised learning and probability estimation. Neural Networks, 10, 925–939.
Lin, C. R., & Chen, M. S. (2002). A robust and efficient clustering
algorithm based on cohesion self-merging. Proceedings of ACM
SIGKDD, Edmonton, Alberta, Canada .
Linde, Y., Buzo, A., & Gray, R. M. (1980). An algorithm for vector quantizer
design. IEEE Transactions on Communication, COM-28, 84–95.
Lozano, J. A., Pena, J. M., & Larranaga, P. (1999). An empiricalcomparison of four initialization methods for the k-means algorithm.
Pattern Recognition Letters, 20, 1027–1040.
Martinetz, T. M. (1993). Competitive Hebbian learning rule forms perfectly
topology preserving maps. ICANN , 427–434.
Martinetz, T. M., Berkovich, S. G., & Schulten, K. J. (1993). Neuralgas network for vector quantization and its application to timeseries prediction. IEEE Transactions on Neural Networks, 4(4),556–558.
Martinetz, T., & Schulten, K. (1994). Topology representing networks.
Neural Networks, 7(3), 507–552.
Murty, M. N., Jain, A. K., & Flynn, P. J. (1999). Data clustering: A review.
ACM Computing Survey, 31(3), 232–264.
Murty, N. M., & Krishna, G. (1981). A hybrid clustering procedure forconcentric and chain-like clusters. International Journal of Computerand Information Sciences, 10(6), 341–397.
Oyang, Y.-J., Chen, C.-Y., & Yang, T.-W. (2001). A study on thehierarchical data clustering algorithm based on gravity theory.
Proceedings of the Fifth European Conference on Principles andPractice of Knowledge Discovery in Databases , 350–361.
Patane, G., & Russo, M. (2001). The enhanced LBG algorithm. Neural
Networks, 14, 1219–1237.
Ritter, H. J., Martinetz, T. M., & Schulten, K. J. (1991). Neuronale Netze.
Munchen: Addison-Wesley.
Sneath, H. A., & Sokal, R. R. (1973). Numerical taxonomy. London, UK:
Freeman.
Villmann, T. (2000). Controlling strategies for the magnification
factor in the neural gas network. Neural Network World, 10,
739–750.
Villmann, T., Der, R., Herrmann, M., & Martinetz, T. (1997). Topology
preservation in self-organizing feature maps: Exact definition andDTD 5ARTICLE IN PRESS
S. Furao, O. Hasegawa / Neural Networks xx (2005) 1–17

measurement. IEEE Transactions on Neural Networks, 8(2),226–256.
Villmann, T., Schleif, F.-M., & Hammer, B. (2003). Supervised neural gas
and relevance learning in learning vector quantization. Proceedings ofthe Workshop on Self-Organizing Maps (WSOM), Japan .
Wasan, M. T. (1969). Stochastic approximation. Cambridge: Cambridge
University Press.
Willshaw, D. J., & von der Malsburg, C. (1976). How patterned neural17

connections can be set up by self-organization. Proceedings of theRoyal Society of London B, 194, 431–445.
Zhang, T., Ramakrishnan, R., & Livny, M. (1996). BIRCH: An
efficient data clustering method for a very large database.
Proceedings of ACM SIGMOD Conference on Management of
Data , 103–114.