| Library: | pp |
| See also: | rqua |
| Quantlet: | ppftind | |
| Description: | computes the Friedman Tukey index (integral about the squared density of the projected data) via kernel density estimation |
| Usage: | ind = ppftind (x, h) | |
| Input: | ||
| x | n x p matrix containing the user-data | |
| h | 1 x p matrix of bandwidths | |
| Output: | ||
| ind | scalar, the Friedman Tukey index | |
; loads the library pp
library("pp")
; initialize random generator
randomize(0)
; generate a dataset with mean(x)=0 and cov(x)=I_2
x = normal(100,2)
; compute the index with Scott's rule
h = 2.62 .*100^(-1/6)
ppftind(x, h~h)
Content of object ind [1,] 0.086871
proc() = plotindex1(x)
n = 50
phi = grid(0, 2*pi/n, n)
ind = NaN*matrix(n)
h = 2.62.*rows(x)^(-1/5)
i = 0
while(i<n)
i = i+1
xp = x * #(cos(phi[i]), sin(phi[i]))
ind[i] = ppftind(xp, h)
line(phi~ind)
endo
endp
; load the library pp and plot
library("pp")
library("plot")
; read the banknote data
x = read("bank2")
; select the fourth and sixth column and
; sphere the data
x = transform(x[,4|6], grc.prep.sphere)
; generate the index function for 1D projection
plotindex1(x)
We see the plot of the index function in dependence of the projection vector. The projection vector is parametrized by the angle phi. Since the projection vector v and -v generate the same projection, we see the index function doubled. Note the maxima at phi=1,76 and phi=2,8 which represent the "most" interesting projections.
library("pp")
library("plot")
; read the banknote data
x = read("bank2")
; select the fourth and sixth column
; and sphere the data
x = transform(x[,4|6], grc.prep.sphere)
; generate a display
d = createdisplay(2,1)
; choose bandwidth related to Scott's rule
h = 2.62.*rows(x)^(-1/5)
; create a "white" data set
w = setmask((-4~0)|(4~0), "white")
; generate first projection(global maximum)
px = x* #(cos(1.76), sin(1.76))
show(d, 1, 1, grdot(px), w)
ind = ppftind(px, h)
setgopt(d, 1, 1, "title", string("FT-Index = %.4f", ind))
; generate second projection(local maximum)
px = x* #(cos(2.80), sin(2.80))
show(d, 2, 1, grdot(px), w)
ind = ppftind(px, h)
setgopt(d, 2, 1, "title", string("FT-Index = %.4f", ind))
Shows the plot of the two projections at the maxima of the index curve in the example before. Note that the Friedman Tukey index is attracted by the extreme values rather than by the separation into cluster.