Efficient Nonparametric Density Estimation on the Sphere with Applications in Fluid Mechanics.
The application of nonparametric probability density function estimation for the purpose of data analysis is well established. More recently, such methods have been applied to fluid flow calculations since the density of the fluid plays a crucial role in determining the flow. Furthermore, when the calculations involve directional or axial data, the domain of interest falls on the surface of the sphere. Accurate and fast estimation of probability density functions is crucial for these calculations since the density estimation is performed at each iteration during the computation.  In particular the values fn (X1 ), fn (X2), ... , fn (Xn) of the density estimate at the sampled points Xi are needed to evolve the system. Usual nonparametric estimators make use of kernel functions to construct fn.   We propose a special sequence of weight functions for nonparametric density estimation that is especially suitable for such applications. The resulting method has a computational advantage over kernel methods in certain situations and also parallelizes easily. Conditions for convergence turn out to be similar to those required for kernel-based methods. We also discuss experiments on different distributions and compare the computational efficiency of our method with kernel based estimators.
Introduction
.
esti mati oni s the problem of the esti
mati on of the values of a
li ty
ven samples from the associ
ated
di stri buti on. No
made about the type of the
di stri buti on
from
whi ch the samples are drawn.
Thi si si n contrast to
esti mati on
whi ch the
assumed to come from a
ven
fami ly, and the parameters
are then
esti mated by
ous
stati sti cal methods. Early
contri butors to the theory
of
nonparametri c
esti mati oni nclude
Smi rnov [21], Rosenblatt [16], Parzen [15], and
Chentsov [3].
ve
descri pti ons of
ous approaches to
nonparametri c
on along
wi th a
ve
bi bli ography can be
books by
lverman
[23] and Nadaraya [14]. More recent developments are
presentedi n books by Scott
[18] and Wand and Jones [27]. Results of the
experi mental
compari son of some
dely
used methods
In
addi ti on to data
analysi s,
ani mportant
appli cati on of
nonparametri c
onal
flui d
mechani cs. When the flow
calculati ons are per-
Lagrangi an framework, a set of
space are evolved through
me
usi ng the
ng
ons. In
poi nts that
werei ni ti ally close can move apart,
leadi ng to mesh
di storti on and
cal
di #culti es. Problems
th mesh
di storti on
can be
eli mi nated to a
extent by the use of smoothed
cle
hydrodynami cs
ques [2, 13, 9, 12]. SPH treats the
nts
bei ng tracked as samples comi
ng from an unknown
li ty
di stri buti on. These
calculati ons often
requi re the
computati on of the values of not only the unknown
densi ty,
buti ts
gradi ent as well.
# Received by the editors August 16, 1995; accepted for publication (in revised form) August 3,
1999; published electronically June 13, 2000.
http://www.siam.org/journals/sisc/22-1/29046.html
Department of Computer Science, University of California, Santa Barbara, CA 93106 (omer@cs.
ucsb.edu).
# Department of Mathematics, Indian Institute of Technology, Bombay, India (ashok@math.
iitb.ernet.in).
In contrast to
appli cati ons concerned
wi th the
di splay of the
densi ty,
ci ent to
esti mate the
densi ty on some
gri d,i n these
flui d flow
calculati ons the
requi red at each sample
nt. Another
di #erencei n these two types of
appli cati onsi s that when
deali ng
th data
analysi s,
usually concerned
th
the
opti mal accuracy one can get for a
ven sample
si ze. In
flui d flow
calculati ons,
where
addi ti onal "data" can be
ned
wi thi ncreased
di screti zati on,
usually
more concerned
wi th the
opti mal
vari ati on of the
onal e#ort as a
on
of error.
In some
appli cati ons, for
example,i n
ng
di recti onal data [24], the
samples
li e on the
ci rcle S 1 or along the surface of the
al
case of
di recti onal
al
whi ch the
c about the
center of the
ci rcle or the sphere,
ous methods have been proposed for
nonparametri c
cal
stati sti cs, such as the kernel [15, 1, 28] and the orthogonal
seri es methods
[17, 11]. The kernel method has been
extensi vely
studi ed,
probably the
most popular
appli cati ons such as SPH. In
thi s method, the value of the
densi ty at the
nt
xi s
esti mated as
#, (1)
ni s the
esti mate of the
ven a sample of
are the
posi ti ons of
the samples drawn from a
li ty
di stri buti on
wi th an unknown
on
Ki s a kernel
hi s the
ndow
wi dth, and A
hi s a
normali zati on factor
to make f
ni nto a
li ty
ty. One of the drawbacks of the kernel method
s the
onal
costi nvolved. Even
possi ble to reduce the
the
one-di mensi onal case
usi ng the
expansi on of a
polynomi al kernel and an
ng
strategy [19],
thi s strategy cannot be
ly extended to
hi gher
di mensi ons [5].
ng
methods [5] can be
usedi n any
di mensi on. However,
si nce the
evaluated on a
gri d,
methodi s not
sui table for the
flui d flow
calculati ons
whi ch we
arei nterested, where an
requi red at each sample
nt.
We propose a
cosi ne-based
wei ght
on
nonparametri c
esti mati on,
whi chi s a
al case of the class of
esti mators that form a #
sequence [26, 28].
Thi s
mi lar to the kernel
esti mator but has the ease
of
evaluati on of a
seri es
expansi on. The role of the
ndow
wi dth parameter h of the
kernel
methodi s replaced by a
smoothi ng parameter
our method, and f
ni s now
of the form
(2)
Our
choi ce of c
cularly
sui table for
appli cati onsi n
flui d flow
calculati ons
where the values f n (X 1 at the sampled
di recti ons themselves
are
requi red at each
each
me
stepi n the flow
We show that
th
esti mator the
red n values can be computed
e#ci ently
usi ng only O(m 1+d n)
operati ons for
di recti onal data and O(m d n)
operati ons for
al
di mensi ons, where m need not be large as long
wi thout bound
th n.
Thi si si n contrast to the O(n 2 )
ons
red by the kernel method
for
computati oni n the worst case and an expected
complexi ty of O(h d
th
ng bounded support.
However,i n the
al case of d = 1, the
of the kernel method can be reduced to
li near after
ng step.
GL U AND ASHOK SRINIVASAN
We
deri ve
ons under
whi ch the sequence of
esti mated
ons f n
fashi on converge to the unknown
experi mentally
veri fy the accuracy and the
e#ci ency of our
methodi n
practi cal test cases. Experi
ments and
cal analyses
alsoi ndi cate how m should vary
th n for
opti mal
accuracy.
The
paperi s
organi zed as follows. In
secti on 2 we define the
wei ght
on
ve the
ons for the convergence of the
ntegrated square
error (MISE) when the sample
spacei s S 1 (Theorem 3). The
ons guarantee
that
E(f
as n #. We also present
correspondi ng results for S 2 . In
secti on 3 schemes for
e#ci ent
computati on of these
esti mates on S 1 and S 2 are presented. In
secti ons 4
and 5, we
descri be
experi mental results
th our
esti mator and
wi th the
kernel method for some
di stri buti ons
practi ce. Our
experi mentsi mply
a net
savi ngs on the number of
ons performed over kernel
methodsi n
ons and also
veri fy the formula found for the
opti mal
choi ce of m. The results
show that the kernel method and our
esti mator perform
di #erent
setti ngs, and
thus complement each other. The
conclusi ons are
presentedi n
secti on 6. The
ns
addi ti onal test results.
2. The cosine estimator and the convergence of MISE. In
secti on, we
first
menti on some related work done on
spheri cal data; then we define our
esti mator
and
deri ve
ons
ts convergence for
di recti onal data on the
ci rcle, and
ve
correspondi ng results for
di recti onal and
al data on the sphere and
al data on
the
ci rcle.
The kernel method for
nonparametri c
esti mati on for
di recti onal and
al
di scussedi n [6, 8].
Whi le
deali ng
th
di recti onal data,
Fi sher,
Lewi s, and
Embleton [6] recommend
usi ng the
ng kernel:
exp
For
al data they recommend the kernel
#, (4)
where
normali zes W to a
li ty
on, and C
ni s the
reci procal
of h
usedi n the
defini ti on of kernel
esti mators. x and X i are the
Cartesi an represen-
tati on of
nts P and P i ,
respecti vely, and x T X
ii s
thei nner product of these two
vectors. W n
the role of K(x - X i ) of (1). Hall, Watson, and Cabrera
[8] analyze
esti mators for
di recti onal data
wi th the x - X i
replaced by
. Observe that the term x T X
ii s the
cosi ne of the angle between the
nts
ii s a measure of the
di stance along the surface
of the sphere between
nts P and P i . Inner product plays a
cruci al
rolei n these
esti mators. We
consi der an
esti mator that can be
terms of powers of the
nner product, the power
playi ng the role of the
smoothi ng parameter.
Thi s enables
us to expand the
esti matori n a
seri es and
li tates fast
computati on.
c (x-
c (x-
OOc (x)32
(a) (b)
-x
Fig. 1. The
ns c 32 (x) and c
2.1. The case of S 1 . We first define our
esti mator on S 1 . Assume X
1, 2, .
sequence
ndependently
andi denti cally
di stri buted
(i .i .d.) random
ables
(observati ons) for
di recti onal data on
th
li ty
[-#].
Wei mpose the
addi ti onal
condi ti on that
si nce the
random
ables X j are defined on the
ci rcle S 1 .
As an
esti mator of the
densi ty of
di recti onal data f(x), x #
[-#], we
der
a
nonparametri c
esti mator of the form
ven by (2)
th
cos 2m
on
[-#]. The
normali zati on factor Am
ven below makes c m
(x)i ntegrate to 1 on
[-#]:
-#
cos 2m
2# dx.
Maki ng use of a table
ofi ntegrals such as Gradshteyn and
Ryzhi k [7] and by
usi ng
can be shown that
As examples, the
ons c m (x) and c m
4 are shown on S
Fi gure 1(a) and on
thei nterval
Fi gure 1(b).
We
wi sh to find
su#ci ent
ons under
whi ch the sequence of
esti mators f n
converges to
fi n the MISE sense. In order to do
thi s, we first show the convergence
of the
bi as and then
deri ve the
ons under
whi ch the
ance converges to 0.
We shall then use these results to prove convergence of MISE on S 1 .
Fi rst we show that as m #, the expected value of the
esti mate f n (x) approaches
the actual
uni formly for any
ven n.
Lemma 1. Suppose f # C 2
andl et f n (x) be as given in (2). Then
uniforml y,
independentl y of n.
Proof.
-#
s)f(s)ds, (7)
GL U AND ASHOK SRINIVASAN
as
showni n
lverman [20] and
Whi ttle [29]. By a change of
able
-#
x-#
x-#
usi ng the
peri odi ci ty of c m and f , along
th (7), (8), and
the mean value theorem,
-#
-#
/2)dy, where # x
nt between x and y. Therefore
-#
-#
-#
From (6), the
ntegral evaluates to 1, and
nce yc m
(y)i s an odd
on, the
secondi ntegral evaluates to 0. Let 2M
We then have the
ng
esti mate for the
as:
-#
-#
-#Am
cos 2m (y/2)y 2 dy.
For any # such that 0 < #,
|y|<#Am
cos 2m (y/2)y 2 dy
|y|#Am
cos 2m (y/2)y 2 dy
|y|<#Am cos 2m (y/2)dy
Am cos 2m (#/2)
nce cos y decreases as
|y|i ncreases on
thei nterval under
consi derati on. Furthermore,
bounded above by 1. Therefore,
Am
Am
Am
In order to get a bound, we
wi ll choose # as a
functi on of m. If we take # 0
as m #, then
# as m #, then
thi s term decays
exponenti ally. The second
(10)i s the product of
thi s term and m), and thus the
product approaches 0
si nce the
exponenti al decay
domi nates. In order to get a good
bound on the first term of (10), we
wi sh to choose #
ng the
condi ti on that
# such that the #
2i s as small as
possi ble. We can choose
arbi trari ly small. Thus M #
/mi s an
asymptoti c bound on the
as.
Furthermore, the
ndependent of x; hence, the
uni form.
Lemma2 . Suppose f # C 2
andl et f n (x) be as given in (2). Then
uniforml y as n #, provided m # as n #, and
Proof.
-#
-#
as
showni n
Whi ttle [29]. As a consequence of Lemma 1, the
secondi ntegral approaches
asymptoti cally, and hence, the second term approach 0
nce
bounded. It thus su#ces to show the convergence of the
firsti ntegral to 0.
x#] |f(x)|.
maki ng a change of
able
usi ng
-#
n, where the
expressi on on the
ri ght-hand
a consequence of the
pressi on for
nce m/n 2
#, the
abovei ntegral converges to
nce
mi si ndependent of x, the
uni form. Therefore the
vari ance of
uni formly to 0 under the
ons of the lemma.
Note that the bound on the
bi as for the
cosi ne method
ven by Lemma
1i s of
the form
and the bound for the
ance
ven by Lemma
Therefore, the role played by m for the
cosi ne
methodi s the same as h -2 of the kernel
based methods, where
hi s the
ndow
wi dth of the kernel
esti mator. In other words
the bounds on the
bi as and the
vari ance of the
cosi ne
esti mator
arei n
accordance
wi th the
behavi or of the kernel method
lverman [23].
Such
si mi lari ty of rates of
convergencei s to be expected
si nce the
cosi ne
essenti ally
li ke the kernel
esti mator, though the forms of the
ons
di #er. It
wi ll be
shown later that the
advantage of the
cosi ne
esti mator
li esi ni ts
onal
e#ci ency.
Theorem 3. Suppose f # C 2
[-#] and f n (x) is as given in (2). If m # as
n #, and
-#
E(f
as n #.
Proof.
E(f
as
showni n
Whi ttle [29]. From Lemmas 1 and 2, each of
thei ntegrals approaches 0.
Hence, the MISE converges to 0.
GL U AND ASHOK SRINIVASAN
In fact, the
MISEi s of the form
where the
c bounds on the constants are as
as we shall
explai n later, the exact
asymptoti c constants are not all
mportant
for
practi cal
ons.
Condi ti ons for the convergence of the
esti mates
ts
deri vati ves on the real
li nei nstead of S 1 can be
consi der the case when the
di recti onal
data
li e along the surface of the
2.2. The case of S 2 . Let X
1, 2, .
., n, be a sequence
ables
th values on the surface of the
centered at
the
Suppose that the
li ty
on f(x) of the X j has
bounded second
deri vati ves. We
consi der a
nonparametri c
esti mator of the form
some m to be
determi ned as a
functi on of n. The c m are
thi s case as
follows. If # xX denotes the angle between
nts x and X, then
cos 2m
The
normali zi ng factor
ven below:
Through a
deri vati on along the
li nes of the case of the
ci rcle, the
ng theorem
can be proved for the convergence of the
esti mators.
Theorem 4. Suppose f # C
andl et f n (x) be as given in (14). If m #
as n # and
E(f
Analogous to (13), the form of the
MISEi s found to be
From
expressi on for MISE we see that
asi n the case of S
hi s
the
ndow
wi dth of the kernel
esti mator.
When
deali ng
th
al data, we can
consi der the
ng
al
esti mator for
spheri cal data:
cos 2m (# xX ).
We can also define a
correspondi ng
esti mator on the
ci rcle, where we take the
cosi ne
of the arc length between two
poi nts,i nstead of the
cosi ne of half the arc length
the case of
di recti onal data. The
between the
s the same
for the cases of the
ci rcle and the sphere,
vely.
3. E#cient evaluation of the density estimator. In
secti on, we shall
descri be an
e#ci ent
algori thm for the
computati on of the
esti mates f n (x)
evaluated at a set of n observed
nts X
2, .
on the
ci rcle S 1
case) and on the
We also show
f the value
of f n at some
arbi trary
xi s
desi red, then
ly
accompli shed once
computed. The
e#ci ency of our
methodi s
based on the fact that f n
terms of the
ons c m (x) as
or (15).
Suppose we represent the
posi ti ons of the observed
nts X
2, .
thei r
Cartesi an
coordi nates. We show that for any x, f n (x) can be expressed as a
polynomi al
of total degree
mi n the
coordi nates of x. The
coe#ci ents of
polynomi al can be
determi nedi n turn from the
coordi nates of the X i . Moreover, the
coe#ci ents are the
sums of the
contri buti ons due to each X
ii ndependently.
Fi rst we
consi der the case of
di recti onal data on S 1 . From (2), (5), and the
half-angle formula for
cosi ne we get
Denote the
poi nts on S 1
correspondi ng to the angles x and X i for
1, 2, .
., n, by
Cartesi an
coordi nates and let
# represent the
standardi nner product on R 2 . Then cos(x
.
ng
thi si nto
(17) we get
The
expressi on
(18)i s a
polynomi al of degree
. For a fixed m, we can
compute the
coe#ci ents by
addi ng the
contri buti on of each X i as follows.
Usi ng the
al theorem and (18)
r+s#m
s
x r
where
thei nner
th
Changi ng the order of
summati on
r+s#m
s
where
M(r, s) =n
i2, and
Ami s as
ven by (6). If we use the
expressi on for
onal ease, then (20)
fies to
r+s#m
s
GL U AND ASHOK SRINIVASAN

Table
Co mputatio nal
co mplexityo f the
co sine
estimato r.
Circle
Axial
Directional
for large m. Now we
consi der the number of
ons
requi red for the
evaluati on
of f n (x)
ven the
ons X
2, .
. The powers X r
i1 and X r
i2 for a fixed
1, 2, .
., m can be computed
th O(m)
multi pli cati ons.
ng
thi s for
1, 2, .
res O(mn)
multi pli cati ons. After the
conclusi on of
thi s step, each
of the O(m 2 ) averages
M(r, s) for a
ven r and s can be computed
wi th an
addi ti onal
O(n)
ons.
Si nce there are a total of O(m 2 )
correspondi ng to
the
rs
s
th 0 #
thi s means that the
coe#ci ents of the
polynomi ali n
or (21) can be computed
wi th a total of O(m 2 n)
ons.
Once the
coe#ci ents of f n (x) have been computed, to evaluate f n (x)
th
calculate the powers x r
1 and x r
1, 2, .
., mi n O(m)
ons.
Si nce the
coe#ci ents are already
avai lable, the
remai ni ng
res
only an
multi pli cati ons and
addi ti ons. The results for
di #erent cases
are
summari zedi n Table 1.
Remark. For our MISE convergence
condi ti on for S 1 (Theorem 3) to be
musti ncrease
wi thout bound
th n.
Theoreti cally, we can take m to be as
ng as we
li ke. Then the above
resulti mpli es that the
computati on
of the
densi ty at all of the sample
poi nts can be
accompli shed
usi ng only about O(n)
the
magni tude of m
ves acceptable accuracy for f n (x). The
problem
th
ng m to be too slowly
that the
magni tude of m controls
the error
our convergence proofs.
An
e#ci ent
algori thm for the
evaluati on of f n (#x) for
di recti onal data on S
constructed
larly. When #
2, .
are
observati ons on S 2 drawn from an
unknown
can be shown [4] that
r+s+t#m
where the
th
M(r, s, t) =n
i3 .
Thi s
ti me the
coe#ci ents of the
polynomi ali n (22) can be computed
wi th a total of
O(m 3 n)
ons. After
preprocessi ng, each
evaluati on of f n (x) for
res only an
ons.
ng results can be
deri ved for
al data as
summari zedi n Table 1.
It should also be noted that we needed the
Cartesi an
representati on of the data. If
the data are
cal
coordi nates, then there
wi ll be an
addi ti onal overhead
for
ng the
Cartesi an
representati on. However,
thi s overhead takes only
li near
ti me and so
wi ll be
negli gi ble for
su#ci ently large data.
Furthermore,i t has been
shown [24] that for
ani mportant class of
appli cati ons,
Cartesi an
coordi nates are
preferable to
cal
coordi nates, as the latter
systemi s not
cally stable for
solvi ng the
di #erenti al
equati ons that
se.
In the subsequent parts of
secti on we shall compare the
onal e#-
ci ency of our scheme
wi th that of the kernel method.
3.1. Parallelization. One of the advantages of the
onal strategy de-
bed
abovei s the ease of
paralleli zati on.
Paralleli zati oni s
redi n many
flui d
flow
calculati ons due to the large
si zes of the systems. The kernel
methodi s some-
what
di #cult to
paralleli ze. If we use an
e#ci ent
kerneli mplementati on that performs
kernel
evaluati ons only for those
nts
whi ch are
di stance h of the
ven sample, then an
e#ci enti mplementati on of the
paralleli zati on
res
load
ng and
decomposi ti on so that
poi nts that are close by
remai n on
the same processor, and so that each processor has roughly the same
loadi n terms
of the
onal e#ort. Also, the
communi cati on pattern for the kernel method
i s not very regular. In contrast,
paralleli zati on for the
cosi ne
esti mator can
ly
be
accompli shed by a global
reducti on
operati on, for
whi ch
e#ci ons
are usually
avai lable.
Thi s method
requi res the same
onal e#ort for each
poi nt, and so the
loadi s
ly balanced by
havi ng the same number of
each
processor.
decomposi ti on does not play
ani mportant
si nce the
poi nts can be on any processor.
3.2. Theoretical comparison of the kernel and the cosine estimators.
Now we analyze the
onal
e#ci ency of the kernel and the
cosi ne
esti mati on methods.
Ani mportant measure of the
e#ci ency of the
algori thmsi s not
just the convergence rate of the error
th sample
ze n, but the
of the
onal e#ort C
red as a
functi on of the error E. For the kernel
esti mator,
we can
te the
#, (23)
where
hi s the
smoothi ng parameter,
the
di mensi on, and
ni s the sample
ze. The
onal e#ort
requi red for
nonparametri c
esti mati on can
be expressed as
dependi ng on the
detai ls of the
algori thm used. For
a
ven sample
ves the
mal h as h # n -1/(d+4) . However,
si nce the
equati on for the
onal
also depends on h, we need to
der
the
possi bi li ty that a value of h smaller than
opti mal value may actually result
lower
onal e#ort.
Let us
consi der a
vari ati on of h
th n of the form
From (25), (24), and (23) we
# .
For
mi ni mum error, the exponent of both the terms on the
ri ght should be the same,
otherwi se the error due to the
hi gher term
wi ll
domi nate.
Thi s leads to
whi chi s the same as the value of
mal # for a
ven n. If we let h optn represent
the
opti mal h
mi ni mi zi ng the MISE for a
ven n, and h optC represent the
opti mal
GL U AND ASHOK SRINIVASAN
mi ni mi zi ng the
onal e#ort as a
functi on of the error, then the
expressi on
deri ved above for # does not
necessari lyi mply that h
si nce the
relati on
would
sti ll
sfy the
expressi on for # for some constant k. If k < 1,
then we can choose a
subopti mal value of
hi n order
toi mprove the speed of the
algori thm.
The
opti mal
vari ati on of error
th
onal
usi ng
thi s value of
ven by
4 .
Let us now
consi der the
cosi ne
esti mator. We can
te the
asymptoti c MISE as
follows:
where
Ei s the MISE,
mi s the
smoothi ng parameter,
di s the
di mensi on
the
ci rcle and for the sphere), and
ni s the sample
ze. The
onal
e#ort
requi red for
esti mator can be expressed as
n), where # can be
determi ned from Table 1. The
expressi on for C
abovei s the same as (24)
th
recalli ng that m behaves as h -2 . By an
mi lar to the
previ ous
case we can show that the
opti mal
vari ati on of error
th
onal
ven by
As examples, for the
cosi ne
esti mator on the
ci rcle
th
al data
and
di recti onal data 2), the
onal
complexi ty and the error are related
by
vely.
The
complexi ty of the kernel
esti matori s the same for
al and
di recti onal data.
However, several
di #erent
possi bi li ti es
st
dependi ng on how
e#ci ent
thei mplemen-
tati on of the
esti matori s. If we
der
esti mators of the form
ven by (3) and (4),
then we have
However,i f we
consi der a kernel
th bounded support, and use an
e#ci
mentati on of the
algori thm that computes the kernel only for those
poi nts that have
a nonzero
contri buti on, then the expected value of
for data on the
ci rcle. Note that the worst case
remai ns
(26). For
the case, we can
consi der an
e#ci ent
algori thm
usi ng
polynomi al kernels and
ng [19],
whi ch uses a
li near
me after
ani ni ti al O(n log n)
ng step. In
thi s case
whi ch means that the kernel method has a
better
than the
cosi ne kernel. However there appears to be no
natural
generali zati on of
thi s update strategy to
hi gher
di mensi ons [5].
Results for the
di #erent cases can be
determi nedi n the manner demonstrated
above and are
presentedi n Table 2. We
wi sh to
menti on that the exact
constantsi n Theorem 3 are not
qui te
mportant (compared
wi th the exponent
on E),
nce
asymptoti cally the
slowdowni ncurred by the cache
domi nates the
overall
runni ng
me. We can expect that the
mpler memory access pattern of our
esti mator
wi ll
makei t advantageous over the kernel
methodi n the
asymptoti c case.

Table

theo ptimal
co mputatio nal
e#o rt versus MISE. The numbers in the table represent
#, where the
relatio nship between the
co mputatio nal
es
no t take
into acco unt an initial
so rting step.
Estimator Circle
Cosine, axial data 1.75
Cosine, directional data 2.25
Kernel, worst case 1.25 *
Kernel, expected case 1.25 *
Si nce the worst case
es of the kernel method and the
cosi ne
esti mator
for
di recti onal data on the sphere have the same order, the
ve
e#ci enci es of
the methods can be tested only through
experi ments.
larly,
si nce the worst case
complexi ty of the
cosi ne
esti mator for
al data on the
spherei s the same as the
expected case for an
e#ci enti mplementati on of the kernel
esti mator, we need to
perform
experi ments to test the
ve
meri ts of the two
esti mators.
4. Experimental results. We performed
cal
experi ments for
al and
di recti onal data on the
ci rcle and
spherei n order to test the
e#ecti veness of our
esti mator. We first plot
esti mates for known
di stri buti ons and then demonstrate that
the MISE follows expected trends for
di stri buti ons. We finally compare the
onal
e#ci ency of our
esti mator
wi th that of kernel methods. More
ri cal
results are
presentedi n the
appendi x.
We
consi der the
on
normali zes the
functi on to be a
densi ty on the surface of a sphere and
Si s a known
functi on of U .
The angles # and # are the
azi muth and the
cal
coordi nates.
Thi s
s the
soluti on to a
cular
flui d
mechani cs. In
Fi gure 2(a) we present a
cal
esti mate for the
versi on of the above
on where # was taken
to be In
thi s figure, we take the data to be
di recti onal.
However,
nce
th respect to the center of the
ci rcle, we
can
consi der the data as
al and use the
al
esti mator. We can see from
Fi gure
2(b) that
res a much smaller value of m.
In
Fi gure 3(a) the
MISEi s compared versus m and n for the
one-di mensi onal #
di stri buti on
usi ng the
al
cosi ne
esti mator. We also compare
th one case of the
di recti onal
esti matori n order to show the benefit of
usi ng the
al
esti mator. In
Fi gure 3(b) the
MISEi s compared versus m and n for the
two-di mensi onal
versi on of
the #
di stri buti on on the surface of a sphere
usi ng the
di recti onal
esti mator.
We next present results for
experi ments
compari ng the speed of the
cosi ne and the
kernel
esti mators. We
consi der the
opti mal
vari ati on of the
onal e#ort
th
MISE. In order to get the
opti mal
onal e#ort for a
ven MISE, we allow for
the
possi bi li ty that we may
re
di #erent sample
si zes for the kernel and the
cosi ne
esti mators.
Thi si s
fied
ve
calculati ons one can
ly change
the "sample"
ze by
changi ng the
di screti zati on of the system. We have performed
these
sons only for
spheri cal data. The case of data on the
ci rcle was not
dered because of the
asymptoti c analyses of the
previ ous
secti on
whi ch clearly
cate that the
li near kernel
algori thmi n the
one-di mensi onal case
wi ll outperform
the
cosi ne
esti mator.
However,i n a
paralleli mplementati on, the
ng step for the
li near kernel
algori thm may be slow, and then one may
wi sh to
consi der the
cosi ne
esti mator.
GL U AND ASHOK SRINIVASAN
(a)
(b)
Co sine estimates
theo ne-dimensio nal
caseo f the #
defined
abo ve. The
so lid line represents the true density. (a) The dashed line represents
the
directio nal estimate
(b) The dashed line represent the axial estimate
The
ng kernel was chosen for the
sons:
[1, 2], 0
otherwi se,
(a)
(b)
n=1000
n=2000
Fig. 3. MISE
versusm and n
fo r the # density. (a)
One-dimensio nal
(o n the circle):
exp (US cos 2 (x))/A; the
so lid line
sho ws the results
fo r the axial
co sine
estimato r and the dashed
line
fo r the
directio nal
estimato r (with
Two -dimensio nal
(o n the
exp (US cos 2 (#))/A as defined
abo ve; MISE
fo r the
directio nal
co sine
estimato r.
where
Ai s the
normali zati on constant, and the
rati of the
di stance between two
nts along the surface of the sphere to
hi s
ven as the argument to the kernel
on. The use of
thi s kernel for
compari sons can be
fied
ts popular use
GL U AND ASHOK SRINIVASAN
Time
Fig. 4.
Co mpariso no f time (in
seco nds) versus the MISE
fo r the
co sine and kernel
data sampled
defined
abo ve. The
po ints marked in
represent the kernel estimate. The
po ints marked in x represent the
co sine estimate.
flui d
mechani cs
calculati ons [12 . Furthermore, we cannot expect any other kernel
to
ve a
ficantly better performance, for the
ng reasons:
(i
Iti s well
known that most kernels are equally good [23, Table 3.1
th respect to
"e#ci ency,"
as
Gi ven that the
e#ci encyi s about the same, the only other
consi derati oni s the
onal
e#orti nvolved. Our kernel takes between 6 and 10
ng
nt
operati ons for a nonzero
evaluati on
(i ncludi ng the cost of
computi ng the
square of the
di stance). Any other reasonable kernel would
requi re at least 6
ng
nt
ons. Apart from
thi s, the memory access
ti mes and
zero-evaluati ons
would add the same constant to all kernels.
Fi gures 4 and 5 compare the
onal e#ort
requi red for the
cosi ne
wei ght
on
esti mator and the kernel
esti mator for
es. We
obtai ned the
data
usi ng the
ng procedure. We performed
esti mates for
ous values of n,
m, and h and
obtai ned the MISE and the
ti me for the
calculati ons. For the
cosi ne
and the kernel
esti mates, we separately plotted the data for the
me
requi red for the
calculati ons versus the error. We chose the lower envelope of the data as the curve for
that
cular
esti mator
si nce the values of
m, n for the data on the lower envelope
ve the best
nable speed for a
ven MISE.
In
thei mplementati on of the kernel
esti mator, we
di vi ded the
spherei nto cells
such that the
si des of the cells had length at least 2h
nce the kernel defined above
has
ndow
dth 2h, rather than h). We placed each
samplei n the
ate cell.
When
computi ng the
densi ty for a
cular cell, we need to search over
only a few cells. The expected
complexi ty of
We first
der data
esti mated
usi ng the
al
cosi ne
esti mator and an
al
vari ant of the kernel
esti mator.
Fi gure 4 shows the results for the
two-di mensi onal
di stri buti on.
Thi si s an example of a
hi ghly
nonuni form
di stri buti on. We can see
that the kernel and the
cosi ne
esti mators are about equally fast.
We next
consi der a more
di stri buti on
ven by
where
#i s the
azi muth and
#i s the
on. The results
presentedi n
Fi gure 5(a)
show that the
cosi ne
esti mator outperforms the kernel
esti mator by more than an
order of
magni tude when the
al forms of both the
esti mators are used.
After
thi s we
consi dered the two
es
menti oned above but treated the data
as
di recti onal and
esti mated them
usi ng the
di recti onal
vari ants of the kernel and the
cosi ne
esti mators. For the #
di stri buti on, the
cosi ne
esti mator performed very poorly
onal
e#ci ency, and we do not present results for
thi s case.
Fi gure 5(b)
presents the results for the
di stri buti on
ven by
cos(#)+1/(8#). We can see
that the
cosi ne
esti mator
sti ll outperforms the kernel
esti mator, though only
sli ghtly.
5. Discussion. The
compari son of the
esti mates
wi th the true
densi tyi ndi cate
that the
cosi ne
esti mator produces accurate results for the
di stri buti ons tested. Plots
of MISE versus m and n follow the expected trends. As the sample
zei ncreases, the
error for the
decreases.
Besi des, the
opti mal value of
mi ncreases as the
sample
zei ncreases. It can also be seen that as the number of
poi ntsi ncreases, the
range of m over
whi ch the
esti mate performs well
alsoi ncreases. We can use
thi s to
our advantage by
choosi ng a
subopti mal value of m
whi ch decreases the
onal
e#ort
ficantly
buti ncreases the error only
sli ghtly.
The
experi ments
compari ng the
onal
e#ci enci es show that the
cosi ne
esti mator outperforms the kernel
esti mator for
al data when the
di stri buti oni s
moderately
uni form. If the
di stri buti oni s
hi ghly
nonuni form, then the two
esti ma-
tors
ve comparable performance for
al data. The
cosi ne
esti mator outperforms
the kernel
esti mator
sli ghtly for
di recti onal data when the
di stri buti oni s moderately
uni form. However, the
ng results for the
cosi ne
esti mator are poor for
hi ghly
nonuni form
di recti onal data. In general, when the
datai s not very
nonuni form,
smoother
wei ght
ons are used.
Thi s
ves a low value of m
whi chi mpli es a
fast
evaluati on
usi ng the
cosi ne
esti mator. However,
thi s leads to a
hi gher h for the
kernel
esti mator,
whi chi mpli es that more samples
contri bute to the kernel
evaluati on
of each sample
nt and, hence,
thi s leads to more
onal e#ort. Conversely,
when the
di stri buti oni s
hi ghly
nonuni form,
especi ally for
di recti onal data, the ker-
nel
methodi s to be preferred. The
ri cal test results
presentedi n the
further demonstrate
nt.
We also analyzed our
experi mental data to
esti mate an
opti mum
vari ati on of m
th n.
Usi ng the results of our
experi ments, we can perform a least squares
and
mate m as kn 1/2.5 for
one-di mensi onal
esti mati on
whi chi s the
same as that expected based on the
expressi on for the MISE.
ng
appears to
ve reasonable
esti mates for
densi ty on the surface of a sphere.
Thi s
resulti s also
stent
wi th the
cal
predi cti ons. Here, the
magni tude of
k depends on the
complexi ty of the
functi on. It
es between 1 and 10 for
the
di stri buti ons
dered here.
We also noted the values of m, n, and h,
whi ch gave the
opti mal
onal
e#ort for a
ven MISE, and compared the results for the kernel and the
cosi ne esti
mators. We observed that the values of h were close to the values
whi ch gave the
mi ni mum MISE for the
ven sample
si ze. However, the values of m were
cantly lower than the values
whi ch gave the
mi ni mum MISE for the
ven sample
ze, though the
errori nvolvedi tself was not much
hi gher than the
mi ni mum MISE.
appears that we can choose a
subopti mal
smoothi ng
parameteri n order to
ncrease the
speedi n the case of the
cosi ne
esti mator.
GL U AND ASHOK SRINIVASAN
(b)
(a)
Time
Time
Fig. 5.
Plo to f time (in
seco nds) versus the MISE
fo r the
co sine and the kernel
estimatio no f
data sampled
1/(8#). The
po ints marked in o represent the kernel estimate.
The
po ints marked in x represent the
co sine estimate. (a) Data treated as axial. (b) Data treated
as
directio nal.
Elevation
Density
Fig. 6.
Plo to f density
ns g(#;
di#erent
the
elevatio n. The
so lid line
the dashed line the
dash-do tted line
and the
do tted line
6. Conclu ion . In
thi s paper, we have
descri bed a
wei ght
on
esti mator
for
nonparametri c
esti mati on of
li ty
ons based on
cosi nes, and
we
provi ded
ons under
whi ch the
esti mate
ts
deri vati ves converge to the
actual
ons. We have developed a scheme for the
e#ci ent
computati on of the
densi ty and presented
experi mental results to check the performance of the
esti mator
for
practi cal problems. These results are
cularly relevant to
flui d mechan-
calculati ons
and,i n general, to
ons where the sample
si ze can be controlled,
for example, though refinement of the
di screti zati on. We have also
ven an
empi r-
cal formula for
choosi ng the
wei ght
on exponent parameter of the
esti mator.
Our
experi mental results suggest that the
cosi ne
esti mator outperforms the kernel
esti mator for both
di recti onal and
al data that are moderately
uni form. It
ves
performance comparable to the kernel
esti mator for
hi ghly
nonuni form
al data,
whi le the kernel
methodi s preferable for
hi ghly
nonuni form
di recti onal data. There
potenti al for further
theoreti cal study of our
esti mator.


Appendix

. Further test results. We present more test
resultsi n
secti on
to study the
ve
e#ci enci es of the
cosi ne and the kernel
techni ques, as
the
ed
systemati cally from
bei ng
vely
uni form to
bei ng
sharply peaked on the
For these tests, we chose
ons g(#;
si s a
constant that governs the sharpness of the
#i s the
elevati on, and
normali zes
thi s to a
li ty
on.
Fi gure 6 shows
the
densi ty as a
functi on of the
elevati on alone for
di #erent values of the parameter s.
Thi s
c about the center of the sphere, and thus we can use the
al
esti mators,
on 4. We can
alsoi gnore our knowledge of
and use the general
di recti onal
esti mators.
GL U AND ASHOK SRINIVASAN
GL U AND ASHOK SRINIVASAN
GL U AND ASHOK SRINIVASAN
Time
(a)
(b)
We present the results of the
experi ments as plots of
me versus MISE for the
kernel
esti mator versus the
cosi ne
esti mator, for both
al and
di recti onal data. The
kernel
esti matori s the one
usedi n the
ri cal
on 4. The tests were
performed on an Intel Celeron 300MHz processor
th 64 MB memory. The C code
was
led
wi th the gcc
compi ler at
opti mi zati on level -O3.
We can see from
Fi gures 7, 8, 9, 10, and 11 that when the
on
i s not very sharp, the
cosi ne
esti mator outperforms the kernel
esti mator for both
al and
di recti onal data. As the
densi ty becomes sharper, the kernel method starts
outperformi ng the
cosi ne
esti mator for
di recti onal data, though the
latteri s
sti ll
better for
al data. When the
densi ty becomes extremely sharp, the kernel method
becomes better for both types of data, though for
al data the two methods are
sti ll
comparable to a
extenti n terms of speed. These results follow the
cally
predi cted trends and demonstrate that these two methods complement each other for
di #erent types of data.


Appendix

. We thank the referees for
thei r
detai led comments and
advi ce, espe-
ci ally for
di recti ng our
attenti on to the current
li terature.



--R

On so me glo bal measureso f the deviatio nso f density functio n estimates

Estimatio no f unkno wn pro bability density basedo no bservatio ns

Fast implementatio nso f no nparametric curve estimato rs
Statistical Analysiso f Spherical Data

Kernel density estimatio n with spherical data
TREESPH: A unificatio no f SPH with the hierarchical tree metho d

The estimatio no f pro bability densities and cumulatives by Fo urier series metho ds

thed particle hydro dynamics

On estimatio no f a pro bability density functio n and mo de
Remarkso n so me no n-parametric estimateso f a density functio n
Estimatio no f pro bability density by ano rtho go nal series
Multivariate Density Estimatio n
Fast algo rithms fo r no nparametric curve estimatio n
Kernel density estimatio n using the fast Fo urier transfo rm
On the appro ximatio no f pro bability densitieso f rando m variables
Numerical So lutio no f Partial Di
Density Estimatio n fo r Statistics and Data Analysis
A new co mputatio nal metho d fo r the so lutio no f flo w pro blemso f micro structured fluids.
Pro bability density estimatio n in astro no my
Pro bability density estimatio n using delta sequences
Kernel Smo thing
On the estimatio no f the pro bability density
On the smo ns
--TR

--CTR
Jeff Racine, Parallel distributed kernel estimation, Computational Statistics & Data Analysis, v.40 n.2, p.293-302, 28 August 2002
kernel method;convergence;efficient algorithm;nonparametric estimation;fluid mechanics;probability density