# 多次元ベクトルデータの超球面表現と密度推定2

• Rのkde2d関数の多次元版…スケーリングとか、bandwidth.nrd関数からの出力の補正とかの点でちょっと怪しい…
kdeNd<-function (x,h, n = 25)
{
max<-apply(x,2,max)
min<-apply(x,2,min)
nx <- ncol(x)#No. axis
ns <- nrow(x)#No. samples
g<-mapply(seq.int,min,max,length.out=n)

if (missing(h))
h <- apply(x,2,bandwidth.nrd)
#h <- h/4
h <-h/2^nx

a<-list(NULL)

for(i in 1:nx){
a[[i]]<-outer(g[,i],x[,i],"-")/h[i]
}
ret<-list(NULL)
counter<-0
kurai<-rep(1,nx)
N<-n^nx
for(i in 1:N){

for(j in 1:(length(kurai)-1)){
if(kurai[j]==n+1){
kurai[j]<-1
kurai[j+1]<-kurai[j+1]+1
}
}
ret[[counter<-counter+1]]<-c(rep(0,nx+1))
#ret[[counter]][1]<-1/ns
tmp<-c(rep(1,ns))
tmp2<-0
for(j in 1:nx){

for(k in 1:ns){
tmp[k]<-tmp[k]*dnorm(a[[j]][kurai[j],k])

}
}
for(j in 1:ns){
tmp2<-tmp2+tmp[j]
}
for(j in 1:nx){

tmp2<-tmp2/h[j]

ret[[counter]][j+1]<-g[kurai[j],j]
}
ret[[counter]][1]<-tmp2/(ns^(nx-1))
kurai[1]<-kurai[1]+1
}

#z <- matrix(dnorm(ax), n, nx) %*% t(matrix(dnorm(ay), n, nx))/(nx * h[1] * h[2])
return(ret)
}
• ちなみにkde2dは
function (x, y, h, n = 25, lims = c(range(x), range(y)))
{
nx <- length(x)
if (length(y) != nx)
stop("data vectors must be the same length")
if (any(!is.finite(x)) || any(!is.finite(y)))
stop("missing or infinite values in the data are not allowed")
if (any(!is.finite(lims)))
stop("only finite values are allowed in 'lims'")
gx <- seq.int(lims[1], lims[2], length.out = n)
gy <- seq.int(lims[3], lims[4], length.out = n)
if (missing(h))
h <- c(bandwidth.nrd(x), bandwidth.nrd(y))
h <- h/4
ax <- outer(gx, x, "-")/h[1]
ay <- outer(gy, y, "-")/h[2]
z <- matrix(dnorm(ax), n, nx) %*% t(matrix(dnorm(ay), n,
nx))/(nx * h[1] * h[2])
return(list(x = gx, y = gy, z = z))
}
<environment: namespace:MASS>