核密度估计的峰值

Question

我需要尽可能精确地找到核密度估计的峰值（连续随机变量的模态值）。我能找到的近似值：

x<-rlnorm(100) 
d<-density(x) 
plot(d) 
i<-which.max(d$y) 
d$y[i] 
d$x[i] 

但是，当计算d$y精确的功能是已知的。我怎样才能找到模式的确切值？

Answer 1

如果我理解你的问题，我想你只是想要一个更精细的分离x和y。为此，您可以在density函数中更改n的值（默认值为n=512）。

例如，比较

set.seed(1) 
x = rlnorm(100) 
d = density(x) 
i = which.max(d$y) 
d$y[i]; d$x[i] 
0.4526; 0.722 

有：

d = density(x, n=1e6) 
i = which.max(d$y) 
d$y[i]; d$x[i] 
0.4525; 0.7228 

Answer 2

下面是处理模式的两个功能。 dmode函数查找具有最高峰值的模式（主宰模式），n.mode识别模式的数量。

dmode <- function(x) { 
     den <- density(x, kernel=c("gaussian")) 
     (den$x[den$y==max(den$y)]) 
    } 

    n.modes <- function(x) { 
     den <- density(x, kernel=c("gaussian")) 
     den.s <- smooth.spline(den$x, den$y, all.knots=TRUE, spar=0.8) 
     s.0 <- predict(den.s, den.s$x, deriv=0) 
     s.1 <- predict(den.s, den.s$x, deriv=1) 
     s.derv <- data.frame(s0=s.0$y, s1=s.1$y) 
     nmodes <- length(rle(den.sign <- sign(s.derv$s1))$values)/2 
     if ((nmodes > 10) == TRUE) { nmodes <- 10 } 
      if (is.na(nmodes) == TRUE) { nmodes <- 0 } 
     (nmodes) 
    } 

# Example 
x <- runif(1000,0,100) 
    plot(density(x)) 
    abline(v=dmode(x)) 

Answer 3

我想你需要两个步骤来存档你需要什么：

1）找到KDE峰值

2的X轴值）得到峰的desnity值

因此，使用hdrcde包的解决方案是这样的（如果你使用一个包不记）：

require(hdrcde) 

x<-rlnorm(100) 
d<-density(x) 

# calcualte KDE with help of the hdrcde package 
hdrResult<-hdr(den=d,prob=0) 

# define the linear interpolation function for the density estimation 
dd<-approxfun(d$x,d$y) 
# get the density value of the KDE peak 
vDens<-dd(hdrResult[['mode']]) 

编辑：您还可以使用

hdrResult[['falpha']] 

，如果它是你不够精确！