Hi all,
I have a question regarding the density function which gives the
kernel density estimator.
I want to decide the bandwidth when using gaussian kernel, given a set of
observations. I am not familiar with different methods for bandwidth
determination. Below are the different ways in R on deciding the bandwidth.
Can anyone give an idea on which ones are preferred.
Also, how can I take a look at the source code for the density function?
Thank you very much.
Hannah
x <- rnorm(1000)
> bw.nrd(x)
[1] 0.2688588
> bw.nrd0(x)
[1] 0.2282763
> bw.ucv(x)
[1] 0.2112366
> bw.bcv(x)
[1] 0.2890085
Warning message:
In bw.bcv(x) : minimum occurred at one end of the range
> bw.SJ(x)
[1] 0.2716242
> density(x, give.Rkern=T, kernel="gaussian")
[1] 0.2820948
> density(x, kernel="gaussian")
Call:
density.default(x = x, kernel = "gaussian")
Data: x (1000 obs.); Bandwidth 'bw' = 0.2283
x y
Min. :-3.974672 Min. :0.0000199
1st Qu.:-1.987712 1st Qu.:0.0076405
Median :-0.000752 Median :0.0529498
Mean :-0.000752 Mean :0.1256971
3rd Qu.: 1.986208 3rd Qu.:0.2552411
Max. : 3.973168 Max. :0.3883532
>
[[alternative HTML version deleted]]
I can't help you decide which bandwidth method to use, but here's how
you view the density source code...
methods("density")
density.default
On Wed, Jul 25, 2012 at 5:56 PM, li li <hannah.hlx at gmail.com>
wrote:>
> Hi all,
> I have a question regarding the density function which gives the
> kernel density estimator.
> I want to decide the bandwidth when using gaussian kernel, given a set
> of
> observations. I am not familiar with different methods for bandwidth
> determination. Below are the different ways in R on deciding the
> bandwidth.
> Can anyone give an idea on which ones are preferred.
> Also, how can I take a look at the source code for the density function?
> Thank you very much.
> Hannah
>
>
> x <- rnorm(1000)
>
> > bw.nrd(x)
>
> [1] 0.2688588
>
> > bw.nrd0(x)
>
> [1] 0.2282763
>
> > bw.ucv(x)
>
> [1] 0.2112366
>
> > bw.bcv(x)
>
> [1] 0.2890085
>
> Warning message:
>
> In bw.bcv(x) : minimum occurred at one end of the range
>
> > bw.SJ(x)
>
> [1] 0.2716242
>
> > density(x, give.Rkern=T, kernel="gaussian")
>
> [1] 0.2820948
>
> > density(x, kernel="gaussian")
>
>
> Call:
>
> density.default(x = x, kernel = "gaussian")
>
>
> Data: x (1000 obs.); Bandwidth 'bw' = 0.2283
>
>
> x y
>
> Min. :-3.974672 Min. :0.0000199
>
> 1st Qu.:-1.987712 1st Qu.:0.0076405
>
> Median :-0.000752 Median :0.0529498
>
> Mean :-0.000752 Mean :0.1256971
>
> 3rd Qu.: 1.986208 3rd Qu.:0.2552411
>
> Max. : 3.973168 Max. :0.3883532
>
> >
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
Dear all, Given a set of observations X, I want to evaluate the kernel density estimator at these observed values. If I do the following, I do not get the those estimated values directly. Can anyone familiar with this give an idea on how to find out the estimated density values for X?> X <- rnorm(100)> density(X)-> den0> den0Call: density.default(x = X) Data: X (100 obs.); Bandwidth 'bw' = 0.354 x y Min. :-3.2254 Min. :0.0002658 1st Qu.:-1.6988 1st Qu.:0.0359114 Median :-0.1721 Median :0.1438772 Mean :-0.1721 Mean :0.1635887 3rd Qu.: 1.3545 3rd Qu.:0.2866889 Max. : 2.8812 Max. :0.3776935 I did write the code for the kernel density estimator myself. So once I find a proper bandwidth, I can use the following function. However, it would be nicer if there is a more direct way.> fhat <- function(x, X){+ h <- density(X, bw="SJ")$bw + n <- length(X) + 1/(n*h)*sum(dnorm((x-X)/h))}>> est <- numeric(length(X))> for (i in 1:length(X)){est[i] <- fhat(x=X[i], X=X)}>> estThanks in advance. Hannah [[alternative HTML version deleted]]