Robust space transformations for distance-based operations.
For many KDD operations, such as nearest neighbor search, distance-based clustering, and outlier detection, there is an underlying &kgr;-D data space in which each tuple/object is represented as a point in the space. In the presence of differing scales, variability, correlation, and/or outliers, we may get unintuitive results if an inappropriate space is used.The fundamental question that this paper addresses is: "What then is an appropriate space?" We propose using a robust space transformation called the Donoho-Stahel estimator. In the first half of the paper, we show the key properties of the estimator. Of particular importance to KDD applications involving databases is the stability property, which says that in spite of frequent updates, the estimator does not: (a) change much, (b) lose its usefulness, or (c) require re-computation. In the second half, we focus on the computation of the estimator for high-dimensional databases. We develop randomized algorithms and evaluate how well they perform empirically. The novel algorithm we develop called the Hybrid-random algorithm is, in most cases, at least an order of magnitude faster than the Fixed-angle and Subsampling algorithms.
INTRODUCTION
For many KDD operations, such as nearest neighbor search,
distance-based clustering, and outlier detection, there is an
underlying k-D data space in which each tuple/object is
represented as a point in the space. Often times, the tuple
represented simply as the point
in the k-D space. More formally, the transformation
from the tuple t to the point p t is the identity
matrix. We begin by arguing that the identity transformation
is not appropriate for many distance-based operations,
Permission to make digital or hard copies of all or part of this work for
personal or classroom use is granted without fee provided that copies are
not made or distributed for profit or commercial advantage and that copies
bear this notice and the full citation on the first page. To copy otherwise, to
republish, to post on servers or to redistribute to lists, requires prior specific
permission and/or a fee.
SIGKDD 2001 San Francisco, California USA
particularly in the presence of variability, correlation, out-
liers, and/or diering scales. Consider a dataset with the
following attributes:
systolic blood pressure (typical range: 100-160 mm of
mercury, with mean =120)
body temperature degrees Celsius, with a very
small standard deviation (e.g., 1-2 degrees for sick pa-
age (range: 20-50 years of age in this example)
Note that dierent attributes have dierent scales and units
(e.g., mm of Hg vs. degree Celsius), and dierent variability
(e.g., high variability for blood pressure vs. low variability
for body temperature). Also, attributes may be correlated
(e.g., age and blood pressure), and there may be outliers.
Example Operation 1 (nearest neighbor search).
Consider a nearest neighbor search using the Euclidean distance
function in the original data space, i.e., the identity
transformation. The results are likely to be dominated by
blood pressure readings, because their variability is much
higher than that of the other attributes. Consider the query
point (blood pressure
Using Euclidean distance, the point (120, 40, 35) is
nearer to the query point than (130, 37, 35) is. But, in terms
of similarity/dissimilarity, this nding is not very meaningful
because, intuitively, a body temperature of 40 degrees is
far away from a body temperature of 37 degrees; in fact, a
person with a body temperature of 40 degrees needs medical
attention immediately!
A simple x to the above problem is to somehow weight
the various attributes. One common approach is to apply
a \normalization" transformation, such as normalizing each
attribute into the range [0,1]. This is usually not a satisfactory
solution because a single outlier (e.g., blood pressure =
could cause virtually all other values to be contained in
a small subrange, again making the nearest neighbor search
produce less meaningful results.
Another common x is to apply a \standardization" trans-
formation, such as subtracting the mean from each attribute
and then dividing by its standard deviation. While this
transformation is superior to the normalization transforma-
tion, outliers may still be too in
uential in skewing the mean
and the standard deviation. Equally importantly, this transformation
does not take into account possible correlation
between attributes. For example, older people tend to have
higher blood pressure than younger people. This means that
we could be \double counting" when determining distances.
Example Operation 2 (Data Mining Operations).
Data clustering is one of the most studied operations in data
mining. As an input to a clustering algorithm, a distance
function is specied. Although some algorithms can deal
with non-metric distance functions (e.g., CLARANS [20]),
most algorithms require metric ones. Among those, a sub-class
of algorithms that has received a lot of attention recently
is the class of density-based algorithms (e.g., DBSCAN
[11] and DENCLUE [14]). The density of a region
is computed based on the number of points contained in
a xed size neighborhood. Thus, density calculation can
be viewed as a xed radius search. Hence, all the concerns
raised above for nearest neighbor search apply just the same
for density-based clustering.
Outlier detection is another important operation in data
mining, particularly for surveillance applications. Many outlier
detection algorithms are distance- or density-based [16,
21, 7]. Again, the issues of diering scale, variability, corre-
lation, and outliers could seriously aect the eectiveness of
those algorithms. At rst glance, the statement that outliers
could impact the eectiveness of outlier detection algorithms
may seem odd. But if attention is not paid to outliers, it is
possible that the outliers may aect the quantities used to
scale the data, eectively masking (hiding) themselves [3].
Contributions of this Paper The fundamental question
addressed in this paper is: \What is an appropriate
space in the presence of diering scale, variability, correla-
tion, and outliers?" So far, we have seen that the spaces
associated with the identity, normalization, and standardization
transformations are inadequate. In this paper, we
focus on robust space transformations, or robust estimators,
so that for distance computation, all points in the space are
treated \fairly". Specically:
Among many robust space estimators that have been
studied in statistics, we propose using the Donoho-
Stahel estimator (DSE). In Section 3, we show two
important properties of the DSE. The rst is the Euclidean
property. It says that while inappropriate in
the original space, the Euclidean distance function becomes
reasonable in the DSE transformed space.
The second, and arguably the more important, property
is the stability property. It says that the transformed
space is robust against updates. That is, in
spite of frequent updates, the transformed space does
not lose its usefulness and requires no re-computation.
Stability is a particularly meaningful property for KDD
applications. If an amount of eort x was spent to set
up an index in the transformed space, we certainly
would not like to spend another amount x after every
single update to the database. In Section 3, we
give experimental results showing that the DSE transformed
space is so stable that it can easily withstand
adding many more tuples to the database (e.g., 50%
of the database size).
Having shown its key properties, in the second half of
this paper, we focus on the computation of the DSE for
high-dimensional (e.g., 10 attributes) databases. The
original DSE algorithm was dened independently by
both Donoho and Stahel [25]; we refer to it as the
Fixed-angle algorithm. In Section 4, we show that the
original algorithm does not scale well with dimension-
ality. Stahel also proposed a version of the algorithm
which uses subsampling (i.e., taking samples of sam-
ples) [25]. However, the number of subsamples to be
used in order to obtain good results is not well known.
We follow the work of Rousseeuw on least median of
squares [22], and come up with a heuristic that seems
to work well, as shown in Section 6. For comparison
purposes, we have implemented this algorithm, applied
some heuristics (e.g., number of subsamples), and evaluated
eectiveness and e-ciency.
Last but not least, in Section 5, we develop a new
algorithm, which we refer to as the Hybrid-random
algorithm, for computing the DSE. Our experimental
results show that the Hybrid-random algorithm
is at least an order of magnitude more e-cient than
the Fixed-angle and Subsampling algorithms. Fur-
thermore, to support the broader claim that the DSE
transformation should be used for KDD operations,
the Hybrid-random algorithm can run very e-ciently
(e.g., compute the estimator for 100,000 5-D tuples in
tens of seconds of total time).
Related Work Space transformations have been studied
in the database and KDD literature. However, they are from
the class of distance-preserving transformations (e.g., [12]),
where the objective is to reduce the dimensionality of the
space. As far as space transformations go, our focus is not so
much on preserving distances, but on providing robustness
and stability.
Principal component analysis (PCA) is useful for data re-
duction, and is well-studied in the statistics literature [17,
15, 10]. The idea is to nd linear combinations of the at-
tributes, while either maximizing or minimizing the variabil-
ity. Unfortunately, PCA is not robust since a few outliers
can radically aect the results. Outliers can also be masked
(hidden by other points). Moreover, PCA lacks the stability
requirements that we desire (cf: Section 3). SVD is
not robust either; it, too, may fail to detect outliers due to
masking.
Many clustering algorithms have been proposed in recent
years, and most are distance-based or density-based [11, 26,
1, 14]. The results presented in this paper will improve
the eectiveness of all these algorithms in producing more
meaningful clusters.
Outlier detection has received considerable attention in
recent years. Designed for large high-dimensional datasets,
the notion of DB-outliers introduced in [16] is distance-
based. A variation of this notion is considered in [21]. The
notion of outliers studied in [7] is density-based. Again, all
of these notions and detection algorithms will benet from
the results presented in this paper.
Developing eective multi-dimensional indexing structures
is the subject of numerous studies [13, 4, 6]. However, this
paper is not about indexing structures. Instead, we focus
on determining an appropriate space within which an index
is to be created.
In [24], nearest neighbor search based on quadratic form
distance functions is considered, that is, distances computed
using some matrix A. That study assumes prior knowledge
of A. For some applications, A may be data-independent
and well-known. For example, for a distance between two
color histograms, each entry in A represents the degree of
perceptual similarity between two colors [12]. However, for
most applications, it is far from clear what a suitable A
could be. The focus of this paper is to propose a meaningful
way of picking A in a data-dependent fashion.
2. BACKGROUND: DONOHO-STAHEL
If two similar attributes are being compared, and those attributes
are independent and have the same scale and vari-
ability, then all points within distance D of a point P lie
within the circle of radius D centered at P . In the presence
of diering scales, variability, and correlation, all points
within distance D of a point P lie within an ellipse. If there
is no correlation, then the major and minor axes of the ellipse
lie on the standard coordinate axes; but, if there is
correlation, then the ellipse is rotated through some angle
. (See Fig. 4 later.) This generalizes to higher dimen-
sions. In 3-D, the ellipsoid resembles a football, with the
covariance determining the football's size, and the correlation
determining its orientation.
An estimator A, also called a scatter matrix, is a k  k
square matrix, where k is the dimensionality of the original
data space. An estimator is related to an ellipsoid as follows.
Suppose x and y are k-dimensional column vectors. The
Euclidean distance between x and y can be expressed as
where T denotes the transpose operator. A quadratic form
distance function can be expressed as dA(x;
y). For x 6= 0, and x T Ax > 0, A is called
a positive denite matrix, and x T Ax yields an ellipsoid.
Donoho-Stahel Estimator and Fixed-angle Algorithm
The DSE is a robust multivariate estimator of location and
scatter. Essentially, it is an \outlyingness-weighted" mean
and covariance, which downweights any point that is many
robust standard deviations away from the sample in some
univariate projection [19, 22]. This estimator also possesses
desirable statistical properties such as a-ne equivariance,
which not all transformations (e.g., principal component
analysis) possess.
The DSE is our estimator of choice, although there are
many estimators to choose from [22]. For example, we chose
the DSE over the Minimum Volume Ellipsoid estimator because
the DSE is easier to compute, scales better, and has
much less bias (especially for dimensions > 2) [18]. In our
extended work, we consider other robust estimators, but the
one that seems to perform the best (based on numerous simulations
and analyses) is the DSE. In the interest of space,
we only deal with the DSE in this paper. Although the application
we focus on here is outlier detection, we add that
the DSE is a general transformation that is useful for many
applications, such as those described in Section 1.
Fig. 1 gives a skeleton of the initial algorithm proposed by
Stahel [25] for computing the estimator in 2-D. Let us step
through the algorithm to understand how the estimator is
dened. The input is a dataset containing N 2-D points
of the form y In step 1, we iterate through
the unit circle, to consider a large number of possible an-
gles/directions , on which to project. We iterate through
degrees rather than 360 degrees since the 180-360 degree
range is redundant. Hereafter, we call this algorithm
The Fixed-angle Donoho Stahel Algorithm
1. For to  (i.e., 0   < ) using some small
increment (e.g., 1 degree), do:
(a) For
is the unit vector)
(b) Compute
(c) Compute (The MAD is
dened to be 1:4826  (medianjx i () m()j)).
(d) For
2. For
3. Compute the robust multivariate centre ^
and the weighting
function w(t) is dened as follows:
2:5
4. Compute the robust covariance matrix ^
5. Return the Donoho-Stahel estimator of location and
scatter:

Figure

1: The DSE Fixed-angle Algorithm for 2-D
the Fixed-angle algorithm.
For each , each point is projected onto the line corresponding
to rotating the x-axis by , giving the value x i ().
Mathematically, this is given by the dot product between y i
and u, which is the unit vector We call u
the projection vector.
In step 1(b), we compute m(), which is the median of
all the x i () values. MAD is an acronym for median absolute
deviation from the median. It is a better estimator of
scatter than the standard deviation in the presence of out-
liers. Finally, step 1(d) yields d i (), which measures how
outlying the projection of y i is with respect to . Note that
d i () is analogous to classical standardization, where each
value x i () is standardized to x i ()
, with  and  being
the mean and the standard deviation of x i (), respectively.
By replacing the mean with the median, and the standard
deviation with the MAD, d i () is more robust against the inuence
of outliers than the value obtained by classical standardization

Robustness is achieved by rst identifying outlying points,
and then downweighting their in
uence. Step 1 computes for
each point and each angle , the degree of outlyingness of
the point with respect to . As a measure of how outlying
each point is over all possible angles, step 2 computes, for
each point, the maximum degree of outlyingness over all
possible 's. In step 3, if this maximum degree for a point
is too high (our threshold is 2.5), the in
uence of this point
is weakened by a decreasing weight function. Finally, with
all points weighted accordingly, the location center ^ R and
the covariance matrix ^
R are computed.
3. KEY PROPERTIES OF THE DSE
In this section, we examine whether the estimator is useful
for distance-based operations in KDD applications. In
Section 6, we provide experimental results showing the difference
the estimator can make. But rst, in this section,
we conduct a more detailed examination of the properties
of the estimator. We show that the estimator possesses the
Euclidean property and the stability property, both of which
are essential for database applications.
Euclidean Property In this section, we show that once
the DSE transformation has been applied, the Euclidean
distance function becomes readily applicable. This is what
we call the Euclidean property.
Lemma 1. The Donoho-Stahel estimator of scatter, ^
R ,
is a positive denite matrix.
The proof is omitted for brevity. According to standard
matrix algebra [2], the key implication of the above lemma is
that the matrix ^
R can be decomposed into: ^
where  is a diagonal matrix whose entries are the eigen-
values, and Q is the matrix containing the eigenvectors of
R . This decomposition is critical to the following lemma.
It says that the quadratic form distance wrt ^
R between
two vectors x and y is the same as the Euclidean distance
between the transformed vectors in the transformed space.
Lemma 2. Let x; y be two vectors in the original space.
Suppose they are transformed into the space described by ^
R ,
Then, the quadratic form distance wrt ^
R is equal to the
Euclidean distance between xR and yR .
Proof:
R
The proof is rather standard, but we include it to provide
a context for these comments:
For each vector x in the original space (or tuple in the
relation), each vector is transformed only once, i.e.,
Future operations do not require any
extra transformations. For example, for indexing, all
tuples are transformed once and can be stored in an
indexing structure. When a query point z is given, z
is similarly transformed to zR . From that point on,
the Euclidean distance function can be used for the
transformed vectors (e.g., xR and zR ).
Furthermore, many existing distance-based structures
are the most e-cient or eective when dealing with
Euclidean-based calculations. Examples include R-trees
and variants for indexing [13, 4], and the outlier
detection algorithm studied in [16].
The key message here is that space transformation wrt ^
R
is by itself not expensive to compute, and can bring further
e-ciency/eectiveness to subsequent processing.
Stability Property The second property we analyze
here for the DSE concerns stability. A transformation is
stable if the transformed space does not lose its usefulness|
even in the presence of frequent updates. This is an important
issue for database applications. If an amount of
eort x was spent in setting up an index in the transformed
space, we certainly would not like to spend another amount
x after every single update to the index. In statistics, there
is the notion of a breakdown point of an estimator, which
quanties the proportion of the dataset that can be contaminated
without causing the estimator to become \arbitrarily
absurd" [27]. But we do not pursue this formal approach
regarding breakdown points; instead, we resort to experimental
evaluation.
In our experiments, we used a real dataset D and computed
the
R (D). We then inserted or deleted tuples
from D, thereby changing D to Dnew . To measure stabil-
ity, we compared matrix ^
R (D) with ^
R(Dnew ). In the
numerical computation domain, there are a few heuristics
for measuring the dierence between matrices; but there is
no universally agreed-upon metric [9]. To make our comparison
more intuitive, we instead picked a distance-based
operation|outlier detection|and compared the results. Section
6 gives the details of our experiments, but in brief, we
proceeded as follows: (a) We used the old estimator ^
to transform the space for Dnew and then found all the
outliers in Dnew ; and (b) We used the updated estimator
(Dnew ) to transform the space for Dnew , and then found
all the outliers in Dnew .
To measure the dierence between the two sets of detected
outliers, we use standard precision and recall [23], and we
dene: (i) the answer set as the set of outliers found by
a given algorithm, and (ii) the target set as the \o-cial"
set of outliers that are found using a su-ciently exhaustive
search (i.e., using the Fixed-angle algorithm with a relatively
small angular increment). Precision is the percentage of the
answer set that is actually found in the target set. Recall
is the percentage of the target set that is in the answer set.
Ideally, we want 100% precision and 100% recall.
Fig. 2 shows the results when there were 25%, 50%, 75%
and 100% new tuples added to D, and when 25%, 50% and
75% of the tuples in D were deleted from D. The new tuples
were randomly chosen and followed the distribution of
the tuples originally in D. The deleted tuples were randomly
chosen from D. The second and third columns show
the number of outliers found, with the third column giving
the \real" answer, and the second column giving the
\approximated" answer using the old transformation. The
fourth and fth columns show the precision and recall. They
clearly show that the DSE transformation is stable. Even
a 50% change in the database does not invalidate the old
transformation, and re-computation appears unnecessary.
For the results shown in Fig. 2, the newly added tuples
followed the same distribution as the tuples originally in D.
For the results shown in Fig. 3, we tried a more drastic sce-
nario: the newly added tuples, called junk tuples, followed
a totally dierent distribution. This is re
ected by the relatively
higher numbers in the second and third columns of
Fig. 3. Nevertheless, despite the presence of tuples from two
distributions, the precision and recall gures are still close
to 100%. This again shows the stability of the DSE.
% Change in # of Outliers in
Dnew , using: Precision Recall
(Dnew
25% Inserts 17 15 88.2% 100%
50% Inserts 17
75% Inserts
100% Inserts 37 29 78.4% 100%
25% Deletes 13 13 100% 100%
50% Deletes
75% Deletes 15 19 100% 78.9%

Figure

2: Precision and Recall: Same Distribution
of Outliers in
Inserted when: Dnew , using: Precision Recall
(Dnew
25% 53 52 94.3% 96.2%
37.5% 74 70 91.9% 97.1%
50% 95 90 92.6% 97.8%
62.5% 108 100 90.7% 98.0%

Figure

3: Precision and Recall: Drastically Dierent
Distribution
4. K-D SUBSAMPLING ALGORITHM
In the previous section, we showed that the DSE possesses
the desirable Euclidean and stability properties for KDD
applications. The remaining question is whether the associated
cost is considerable. Let us consider how to compute
R e-ciently, for k > 2 dimensions.
Complexity of the Fixed-angle Algorithm Recall
that Fig. 1 gives the 2-D Fixed-angle algorithm proposed
by Donoho and Stahel. The extension of this algorithm to
3-D and beyond is straightforward. Instead of using a unit
circle, we use a unit sphere in 3-D. Thus, there are two
angles|1 and 2|through which to iterate. Similarly, in
k-D, we deal with a unit hypersphere, and there are k 1
angles through which to iterate:
To understand the performance of the Fixed-angle algo-
rithm, we conduct a complexity analysis. In step 1 of Fig. 1,
each angle  requires nding the median of N values, where
N is the size of the dataset. Finding the median takes O(N)
time, which is the time that a selection algorithm can partition
an array to nd the median entry. (Note that sorting is
not needed.) Thus, in 2-D, if there are a increments to iterate
through, the complexity of the rst step is O(aN ). For
k-D, there are k 1 angles to iterate through. If there are
a increments for each of these angles, the total complexity
of the rst step is O(a k 1 kN ).
In step 2, in the k-D case, there are a k 1 projection
vectors to evaluate. Thus, the complexity of this step is
nds a robust center, which can be done
in O(kN) time. Step 4 sets up the k  k robust covariance
matrix, which takes O(k 2 N) time. Hence, the total complexity
of the Fixed-angle algorithm is O(a k 1 kN ). Su-ce
it to say that running this basic k-D algorithm is impractical
for larger values of a and k.
Intuition behind the Subsampling Algorithm in 2-D
The rst two steps of the Fixed-angle algorithm compute,
for each point y i , the degree of \outlyingness" d i . The value
of d i is obtained by taking the maximum value of d i () over
all 's, where d i () measures how outlying the projection
of y i is wrt . In the Fixed-angle algorithm, there is an
exhaustive enumeration of 's. For high dimensions, this
approach is infeasible.
Let us see if there is a better way to determine \good" projection
vectors. Consider points A, B, and C in Fig. 4(a),
which shows a 2-D scenario involving correlated attributes.
Fig. 4(b) shows the projection of points onto a line orthogonal
to the major axis of the ellipse. (Not all points are
projected in the gure.) Note that B's projection appears
to belong to the bulk of the points projected down from
the ellipse; it does not appear to be outlying at all in this
projection. Also, although A is outlying on the projection
Algorithm Subsampling
For is the number of iterations
chosen, do:
(a) Select k 1 random points from the dataset. Together
with the origin, these points form a hyper-plane
through the origin, and a subspace V .
i. Compute a basis for the orthogonal complement
of V .
ii. Choose a unit vector u i from the orthogonal
complement to use as vector u in the Fixed-
angle algorithm.
(b) For
(c) (continue with step 1(b) and beyond of the Fixed-
angle algorithm shown in Fig. 1, where i takes the
role of )

Figure

5: DSE Subsampling Algorithm for k-D
line, C is not. In Fig. 4(c), A and C are not outlying, but B
clearly is. Fig. 4(d) shows yet another projection.
As can be seen, the projection vectors which are chosen
greatly in
uence the values for d i in the Donoho-Stahel al-
gorithm. In applying subsampling, our goal is to use lines
orthogonal to the axes of an ellipse (or ellipsoid, in k-D),
to improve our odds of obtaining a good projection line.
While there may be better projection vectors in which to
identify outliers, these are good choices. There is an increased
chance of detecting outliers using these orthogonal
lines because many outliers are likely to stand out after the
orthogonal projection (see Fig. 4). Non-outlying points, especially
those within the ellipsoids, are unlikely to stand
out because they project to a common or relatively short
interval on the line.
If we knew what the axes of the ellipse were, then there
would be no need to do subsampling. However, since: (a)
we do not know the parameters of the ellipsoid, and (b) in
general, there will be too many points and too many dimensions
involved in calculating the parameters of the ellipsoid,
we use the following approach called subsampling. In 2-D,
the idea is to rst pick a random point P from the set of
N input points. Then compute a line orthogonal to the line
joining P and the origin. Note that, with reasonable proba-
bility, we are likely to pick a point P in the ellipse, and the
resulting orthogonal line may approximate one of the axes
of the ellipse. This is the essence of subsampling.
More Details of the Subsampling Algorithm in k-D
In k-D, we rst nd a random sample of k 1 points. Together
with the origin, they form a subspace V . Next, we
need to nd a subspace that is orthogonal to V , which is
called the orthogonal complement of V [2]. From this point
on, everything else proceeds as in the Fixed-angle algorithm.
Fig. 5 outlines the Subsampling algorithm for k-D DSE computation

One key detail of Fig. 5 deserves elaboration: how to compute
m. To determine m, we begin by analyzing the probability
of getting a \bad" subsample. For each subsample,
are randomly chosen. A subsample is likely to
be good if all k 1 points are within the ellipsoid. Let  (a
user-chosen parameter) be the fraction of points outside the
ellipsoid. Typically,  varies between 0.01 to 0.5; the bigger
the value, the more conservative or demanding the user is
on the quality of the subsamples.
A
(a)
A
(b)
A
(c)
A
(d)

Figure

4: Bivariate Plots Showing the Eect of Dierent Projection Lines. (a) Data points only. (b)
Projection onto a line orthogonal to the major axis of the ellipse. (c) Projection onto a line orthogonal to
the minor axis. (d) Projection onto another line.
Let us say that m can be the smallest number of subsamples
such that there is at least a 95% probability that we
get at least one good subsample out of the m subsamples.
Given , the probability of getting a \good" subsample is
the probability of picking all k 1 random points within the
ellipsoid, which is (1 Conversely, the probability of
getting a bad subsample is 1 (1 . Thus, the probability
of all m subsamples being bad is (1 (1
Hence, we can determine a base value of m by solving the
following inequality for m: 1 (1 (1 0:95. For
example, In Section 6,
we show how the quality of the estimator varies with m.
Complexity of the Subsampling Algorithm In k-D,
we can determine a basis for the orthogonal complement of
a hyperplane through the origin and through k 1 non-zero
points in O(k 3 ) time, using Gauss-Jordan elimination [2, 9].
Using this basis, we simply pick any unit vector u as our projection
vector, and then continue with the basic Fixed-angle
algorithm. Recall from Section 4 that the basic algorithm
runs in O(a k 1 kN) time for step 1 and O(k 2 N) time for the
remaining steps. For the Subsampling algorithm, however,
we perform a total of m iterations, where each iteration consists
of k 1 randomly selected points, and thus step 1 of
the Subsampling algorithm runs in O(mk 3 ) time. Thus, following
the analysis in Section 4, the entire algorithm runs
in O(mk 3
5. K-D RANDOMIZED ALGORITHMS
The Subsampling algorithm is more scalable with respect
to k than the Fixed-angle algorithm is, but the mk 3 complexity
factor is still costly when the number of subsamples
m is large (i.e., for a high quality estimator). Thus, in
this section, we explore how the k-D DSE estimator can be
computed more e-ciently. First, we implement a simple alternative
to the Fixed-angle algorithm, called Pure-random.
After evaluating its strengths and weaknesses, we develop a
new algorithm called Hybrid-random, which combines part
of the Pure-random algorithm with part of the Subsampling
algorithm. In Section 6, we provide experimental results,
showing eectiveness and e-ciency.
Pure-random Algorithm Recall from Fig. 1 that in
the Fixed-angle algorithm, the high complexity is due to the
a k 1 factor, where a k 1 denotes the number of projection
unit vectors examined. However, for any given projection
unit vector, the complexity of step 1 reduces drastically to
O(kN ). It is certainly possible for an algorithm to do well
if it randomly selects r projections to examine, and if some
of those projections happen to be \good" or in
uential pro-
jections. A skeleton of this algorithm called Pure-random
Algorithm Pure-random
For where r is the number of projection
vectors chosen, do:
(a) Select a k-D projection unit vector u i randomly
(i.e., pick k 1 random angles)
(b) For
(c) (continue with step 1(b) and beyond of the Fixed-
angle algorithm, where i takes the role of )

Figure

is presented in Fig. 6. Following the analysis shown in Section
4, it is easy to see that the complexity of the Pure-
random algorithm is O(rkN
randomization is also used in the Subsampling algorithm.
But there, each random \draw" is a subspace V formed by
1 points from the dataset, from which the orthogonal
complement of V is computed. In the Pure-random case,
however, each random draw is a projection vector. In order
for the Pure-random algorithm to produce results comparable
to that of the Subsampling algorithm, it is very likely
that r  m.
Hybrid-random Algorithm Conceptually, the Pure-
random algorithm probes the k-D space blindly. This is
the reason why the value of r may need to be high for acceptable
quality. The question is whether random draws
of projection vectors can be done more intelligently. More
specically, are there areas of the k-D space over which the
randomization can skip, or equivalently, are there areas on
which the randomization should focus?
In a new algorithm that we develop called Hybrid-random,
we rst apply the Subsampling algorithm for a very small
number of subsamples. Consider the orthogonal complement
of V that passes through the origin. Imagine rotating
this line through a small angle anchored at the origin, thus
creating a cone. This rotation yields a \patch" on the surface
of a k-D unit hypersphere. From the Fixed-angle al-
gorithm, we know that projection vectors too close to each
other do not give markedly dierent results. So, in the second
phase of the Hybrid-random algorithm, we will restrict
the random draws of projection vectors to stay clear of previously
examined cones/patches.
Using the Euclidean inner product and the Law of Cosines,
a collision between two vectors a and b occurs if dist 2 (a;
is the radius of a
patch on the surface of the k-D unit hypersphere. To determine
-, we used the following heuristic. We say that vectors
a and b are too close to each other if cos   0:95, where  is
the angle between the vectors. Thus, (2-)
hence, as an upper bound,
we use
0:1
Two observations are in order.
First, patches that are too large are counterproductive because
many promising projection vectors may be excluded.
Second, although increasing the number of patches improves
accuracy, favourable results can be obtained with relatively
few patches (e.g., 100), as will be shown in Section 6.
Fig. 7 gives a skeleton of the Hybrid-random algorithm.
Steps 1 to 3 use the Subsampling algorithm to nd some
initial projection vectors (including the eigenvectors of the
scatter matrix) and keep them in S. In each iteration of
step 4, a new random projection vector is generated in such
a way that it stays clear of existing projection vectors.
Algorithm Hybrid-random
1. Run the Subsampling algorithm for a small number m
of iterations (e.g.,
2. Compute the k eigenvectors of the resulting scatter
matrix. This gives us an approximation for the axes
of the ellipsoid.
3. Initialize the set S of previously examined projection
vectors to consist of the m projection vectors from step
1 and the k eigenvectors from step 2.
4. For where r is the number of extra random
patches desired, do:
(a) From S, randomly select 2 unique vectors a and
b that are at least 2- radians apart.
(b) Compute a new vector u i that is a linear combination
of a and b. In particular, u
)b,
where
is randomly chosen between [-, 1-].
(c) If u i is within - radians from an existing vector
in S, then redo the previous step with a new
If these two vectors are still too close during the
second attempt, then go back to step 4(a).
(d) Normalize u i so that it is a unit vector, and add
it to S.
(f) (continue with step 1(b) and beyond of the Fixed-
angle algorithm, where i takes the role of )

Figure

7: DSE Hybrid-random Algorithm for k-D
Recall from our earlier discussion that the complexity
of the Subsampling algorithm is O(m1k 3
m1 is the number of subsamples taken. As for the Pure-
random algorithm, the complexity is O(r1kN ), where r1 is
the number of random projections probed. It is easy to see
that the Hybrid-random algorithm requires a complexity of
We expect that m2  m1 , and r2  r1 .
Experimental results follow.
6. EXPERIMENTAL EVALUATION
Experimental Setup To evaluate the Donoho-Stahel
transformation, we picked the distance-based outlier detection
operation described in [16]. As explained in Section 3,
we use precision and recall [23], to compare the results.
Our base dataset is an 855-record dataset consisting of
1995-96 National Hockey League (NHL) player performance
statistics. These publicly available statistics can be down-loaded
from sites such as the Professional Hockey Server
at http://maxwell.uhh.hawaii.edu/hockey/. Since this real-life
dataset is quite small, we created a number of synthetic
datasets mirroring the distribution of statistics within the
NHL dataset. Specically, we determined the distribution of
each attribute in the original dataset by using a 10-partition
histogram. Then, we generated datasets containing up to
100,000 tuples|whose distribution mirrored that of the base
dataset. As an optional preprocessing step, we applied the
Box and Cox transformation to normality [8] to nd appropriate
parameters p and D for the distance-based outliers
implementation. Unless otherwise stated, we used a 5-D
case of 100,000 tuples as our default, where the attributes
are goals, assists, penalty minutes, shots on goal, and games
played.
Our tests were run on Sun Microsystems Ultra-1 proces-
sor, running SunOS 5.7, and having 256 MB of main mem-
ory. Of the four DSE algorithms presented, only the Fixed-
angle algorithm is deterministic. The other three involve
randomization, so we used the median results of several runs.
Precision was almost always 100%, but recall often varied.
Usefulness of Donoho-Stahel Transformation In
the introduction, we motivated the usefulness of the Donoho-
Stahel transformation by arguing that the identity transformation
(i.e., raw data), as well as the normalization and
standardization transformations, may not give good results.
In the experiment reported below, we show a more concrete
situation based on outlier detection. Based on the 1995-96
NHL statistics, we conducted an experiment using the two
attributes: penalty-minutes and goals-scored. We note
that the range for penalty-minutes was [0,335], and the
range for goals-scored was [0,69].
Fig. 8 compares the top outliers found using the identity,
standardization, and Donoho-Stahel transformations. Also
shown are the actual penalty-minutes and goals-scored
by the identied players. With the identity transformation
(i.e., no transformation), players with the highest penalty-
minutes dominate. With classical standardization, the dominance
shifts to the players with the highest goals-scored
(with Matthew Barnaby appearing on both lists). How-
ever, in both cases, the identied outliers are \trivial", in
the sense that they are merely extreme points for some at-
tribute. Barnaby, May, and Simon were all in the top-5
for penalty-minutes; Lemieux and Jagr were the top-2 for
goals-scored.
With the Donoho-Stahel transformation, the identied
outliers are a lot more interesting and surprising. Donald
Transform- Top Outliers Penalty-mins. Goals-scored
ation Found (raw data) (raw data)
Matthew Barnaby 335 15
Chris Simon 250
Matthew Barnaby 335 15
Standard- Jaromir Jagr 96 62
ization Mario Lemieux 54 69
Matthew Barnaby 335 15
Donoho- Donald Brashear 223 0
Stahel Jan Caloun 0 8
Joe Mullen 0 8

Figure

8: Identied Outliers: Usefulness of Donoho-
Stahel Transformation
Brashear was not even in the top-15 as far as penalty-
minutes goes, and his goals-scored performance was unim-
pressive, that is, penalty-minutes = 223 and goals-scored
Yet, he has a unique combination. This is because to
amass a high number of penalty minutes, a player needs to
play a lot, and if he plays a lot, he is likely to score at least
some goals. (Incidentally, 0 goals is an extreme univariate
point; however, well over 100 players share this value.)
Similar comments apply to Jan Caloun and Joe Mullen;
both had 0 penalty-minutes but 8 goals-scored. While
their raw gures look unimpressive, the players were exceptional
in their own ways. 1 The point is, without an appropriate
space transformation, these outliers would likely be
missed.
Internal Parameters of the Algorithms Every algorithm
presented here has key internal parameters. In the
Fixed-angle case, it is the parameter a, the number of angles
tested per dimension. For the randomization algorithms,
there are m, the number of subsamples, and r, the number
of random projection vectors. Let us now examine how the
choices of these parameters aect the quality of the estimator
computed. Precision and recall will be used to evaluate
quality. However, for the results presented below, precision
was always at 100%. Thus, we only report the recall values.
The four graphs in Fig. 9 each contrast: (i) CPU times,
(ii) recall values, and (iii) number of iterations (or patches
used) for one of the four algorithms. The left hand y-axis
denes CPU times (in minutes for the top two graphs, and in
seconds for the bottom two graphs). The right hand y-axis,
in conjunction with the recall curve (see each gure's legend)
denes recall values. Note, however, that the recall range
varies from one graph to another. Fig. 9(a) measures CPU
time in minutes, and shows that the Fixed-angle algorithm
can take a long time to nish, especially as the number of
random angles a tested increases. The horizontal axis is in
tens of thousands of iterations. Recall that a small decrease
in the angle increment for each dimension can cause a very
large number of additional iterations to occur. For many of
our datasets, it was necessary to use increments as small as
degrees (e.g., 75 hours of CPU time, for 100,000 tuples
in 5-D), before determining the number of outliers present.
We omit these very long runs from our graphs, to allow us
to more clearly contrast CPU times and recall values.
Compared to the Fixed-angle algorithm, the Pure-random
algorithm achieves a given level of recall more quickly, al-
though, as Fig. 9(b) shows, it can still take a long time to
achieve high levels of recall.
Recall that, for the Subsampling algorithm, a key issue
was how many subsamples to use. Based on the heuristic
presented in Section 4, the base value of m was determined
to be 47, and multiples of 47 subsamples were used. From
the recall curve in Fig. 9(c), it is clear that below 47 sub-
samples, the recall value is poor. But even with 3 *
141 subsamples, the recall value becomes rather acceptable.
This is the strength of the Subsampling algorithm, which
1 Actually, we did not even hear of Jan Caloun before our experiment.
During 1995-96, Caloun played a total of 11 games, and scored 8
goals|almost a goal per game, which is a rarity in the NHL. A search
of the World Wide Web reveals that Caloun played a grand total of
13 games in the NHL|11 games in 1995-96, and 2 games in 1996-
97|before disappearing from the NHL scene. We also learned that
he scored on his rst four NHL shots to tie an NHL record.
can give acceptable results in a short time. But, the recall
curve has a diminishing rate of return, and it may take
a very long time for Subsampling to reach a high level of
recall, as conrmed in Fig. 10.
Since the Hybrid-random algorithm uses the Subsampling
algorithm in its rst phase (with d 47
expected that the Hybrid-random algorithm behaves about
as well as the Subsampling algorithm, at the beginning, for
mediocre levels of recall, such as 70-75% (cf: Fig. 10). But,
as shown in Fig. 9(d), if the Hybrid-random algorithm is allowed
to execute longer, it steadily and quickly improves the
quality of its computation. Thus, in terms of CPU time, we
start with the Subsampling curve, but quickly switch to the
Pure-random curve to reap the benets of a fast algorithm
and pruned randomization.
Achieving a Given Rate of Recall The above experiment
shows how each algorithm trades o e-ciency with
quality. Having picked a reasonable set of parameter values
for each algorithm, let us now compare the algorithms head-
to-head. Specically, for xed recall rates, we compare the
time taken for each algorithm to deliver that recall rate. Because
the run time of the Fixed-angle algorithm is typically
several orders of magnitude above the others (for comparable
quality), we omit the Fixed-angle algorithm results from
now on.
Fig. 10 compares the Hybrid-random algorithm with both
the Pure-random and Subsampling algorithms, for higher
rates of recall. In general, the Subsampling algorithm is
very eective for quick, consistent results. However, to improve
further on the quality, it can take a very long time. In
contrast, when the Hybrid-random algorithm is allowed to
run just a bit longer, it can deliver steady improvement on
quality. As a case in point, to achieve about 90% recall in
the current example, it takes the Subsampling algorithm almost
14 hours to achieve the same level of recall produced by
the Hybrid-random algorithm in about two minutes. Never-
theless, we must give the Subsampling algorithm credit for
giving the Hybrid-random algorithm an excellent base from
which to start its computation.
In Fig. 10, the Pure-random algorithm signicantly out-performs
the Subsampling algorithm, but this is not always
the case. We expect the recall rate for Pure-random to be
volatile, and there are cases where the Pure-random algorithm
returns substantially dierent outliers for large numbers
of iterations. The Hybrid-random algorithm tends to
be more focused and consistent.
Scalability in Dimensionality and Dataset Size Fig.
11(a) shows scalability of dimensionality for the Subsampling
and Hybrid-random algorithms. We used moderate
levels of recall (e.g., 75%) and 60,000 tuples for this anal-
ysis. High levels of recall would favor the Hybrid-random
algorithm. The results shown here are for 282 iterations for
the Subsampling algorithm, and 90 Patches for the Hybrid-
random algorithm. Our experience has shown that these
numbers of iterations and patches are satisfactory, assuming
we are satised with conservative levels of recall. Fig. 11(a)
shows that both algorithms scale well, and this conrms our
complexity analysis of Section 4.
Fig. 11(b) shows how the Subsampling and Hybrid-random
algorithms scale with dataset size, in 5-D, for conservative
levels of recall. Again, both algorithms seem to scale well,
and again the Hybrid-random algorithm outperforms the
Run Time and Recall for the Fixed-angle Algorithm: 5-D, 100,000 Tuples
Number of Angles Tested (Iterations)
CPU
Time
in
Minutes
CPU Time
Recall
Recall
Run Time and Recall for the Pure-random Algorithm: 5-D, 100,000 Tuples
Number of Random Angles (Iterations)
CPU
Time
in
Minutes
CPU Time
Recall
Recall
5002060100Run Time and Recall for the Subsampling Algorithm: 5-D, 100,000 Tuples
Number of Subsamples (Iterations)
CPU
Time
in
Seconds
CPU Time
Recall
Recall
700100Run Time and Recall for the Hybrid-random Algorithm: 5-D, 100,000 Tuples; delta=0.0800
Number of Patches
CPU
Time
in
Seconds
CPU Time
Recall
Recall

Figure

9: Plots of Run Time and Recall: (a) Top left: Fixed-angle. (b) Top right: Pure-random. (c) Bottom
left: Subsampling. (d) Bottom right: Hybrid-random.
9050150250350Run Times to Achieve a Given Level of Recall, for 3 Algorithms, 5-D, & 100,000 Tuples
Percent Recall
CPU
Time
in
Seconds
Pure-random
Subsampling
Hybrid-random

Figure

10: Run Time vs. Recall for Subsampling,
Pure-random, and Hybrid-random Algorithms
Subsampling algorithm. High levels of recall would favor
the Hybrid-random algorithm, even more so than shown.
7.


AND CONCLUSION
The results returned by many types of distance-based
KDD operations/queries tend to be less meaningful when
when no attention is paid to scale, variability, correlation,
and outliers in the underlying data. In this paper, we presented
the case for robust space transformations to support
operations such as nearest neighbor search, distance-based
clustering, and outlier detection. An appropriate space is
one that: (a) preserves the Euclidean property, so that ecient
Euclidean distance operations can be performed without
sacricing quality and meaningfulness of results, and (b)
is stable in the presence of a non-trivial number of updates.
We saw that distance operations which ordinarily would be
inappropriate when operating on the raw data (and even on
normalized or standardized data), are actually appropriate
in the transformed space. Thus, the end user sees results
which tend to be more intuitive or meaningful for a given
application. We presented a data mining case study on the
detection of outliers to support these claims.
After considering issues such as eectiveness (as measured
by precision and recall, especially the latter) and efciency
(as measured by scalability both in dimensionality
and dataset size), we believe that the Hybrid-random algorithm
that we have developed in this paper is an excellent
choice among the Donoho-Stahel algorithms. In tens of
seconds of CPU time, a robust estimator can be computed
which not only accounts for scale, variability, correlation,
and outliers, but is also able to withstand a signicant number
of database updates (e.g., 50% of the tuples) without
losing eectiveness or requiring re-computation. For many
cases involving high levels of recall, the randomized algo-
rithms, and in particular, the Hybrid-random algorithm can
be at least an order of magnitude faster (and sometimes
several orders of magnitude faster) than the alternatives. In
conclusion, we believe that our results have shown that robust
estimation has a place in the KDD community, and can
nd value in many KDD applications.
8.



--R

Automatic Subspace Clustering of High Dimensional Data for Data Mining Applications.
Elementary Linear Algebra: Applications Version.
Outliers in Statistical Data.



LOF: Identifying Density-Based Local Outliers
Box and D.
Numerical Analysis.
A fast algorithm for robust principal components based on projection pursuit.
A Density-based Algorithm for Discovering Clusters in Large Spatial Databases with Noise

R-trees: a dynamic index structure for spatial searching.


Algorithms for Mining Distance-Based Outliers in Large Datasets

Bias robust estimation of scale.
The behaviour of the Stahel-Donoho robust multivariate estimator


Robust Regression and Outlier Detection.
Introduction to Modern Information Retrieval.

Breakdown of Covariance Estimators.
STING: A statistical information grid approach to spatial data mining.
High breakdown point estimates of regression by means of the minimization of an e-cient scale
--TR
Robust regression and outlier detection
The R*-tree: an efficient and robust access method for points and rectangles
Efficient and effective querying by image content
Distance-based indexing for high-dimensional metric spaces
Automatic subspace clustering of high dimensional data for data mining applications
Density-based indexing for approximate nearest-neighbor queries
Efficient algorithms for mining outliers from large data sets
Introduction to Modern Information Retrieval
R-trees
Algorithms for Mining Distance-Based Outliers in Large Datasets
Efficient and Effective Clustering Methods for Spatial Data Mining
Efficient User-Adaptable Similarity Search in Large Multimedia Databases

--CTR
Peng Sun , Robert M. Freund, Computation of Minimum-Volume Covering Ellipsoids, Operations Research, v.52 n.5, p.690-706, Sep. - Oct. 2004
S. Cateni , V. Colla , M. Vannucci, A fuzzy logic-based method for outliers detection, Proceedings of the 25th conference on Proceedings of the 25th IASTED International Multi-Conference: artificial intelligence and applications, p.561-566, February 12-14, 2007, Innsbruck, Austria
Leejay Wu , Christos Faloutsos, Making every bit count: fast nonlinear axis scaling, Proceedings of the eighth ACM SIGKDD international conference on Knowledge discovery and data mining, July 23-26, 2002, Edmonton, Alberta, Canada
Jaideep Vaidya , Chris Clifton, Privacy-preserving k-means clustering over vertically partitioned data, Proceedings of the ninth ACM SIGKDD international conference on Knowledge discovery and data mining, August 24-27, 2003, Washington, D.C.
data mining;robust statistics;space transformations;outliers;robust estimators;distance-based operations