2013-04-27 31 views
8

カーネル密度推定のピーク(連続ランダム変数のモーダル値)を可能な限り正確に見つける必要があります。 d$y正確な機能が知られている計算カーネル密度推定のピーク

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

しかしとき:私は近似値を見つけることができます。どのようにしてモードの正確な値を見つけることができますか?

答えて

5

ご質問いただければ、xyというより細かい離散をご希望されていると思います。これを行うには、関数の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 
+0

ありがとうございます! ;)それはかなり正確に動作するようです – 16per9

10

モードに対処するための2つの機能があります。 dmode関数は最高ピーク(優勢モード)を持つモードを検出し、n.modesはモードの数を示します。

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)) 
0

私はあなたが必要なものをアーカイブするには、2つのステップが必要だと思う:KDEピーク

2のX軸値を検索)

1)をピーク

の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']] 

を使用することができます!

関連する問題