Hi David,
Thanks a lot for the reply.
I might not have stated the problem clearly. Let me try again.
Given a set of observations X, I want to find out the estimated density
values for the observations X?
I believe that the values "x" returned from "density" function is not the
observations
that are fed into the function and the returned "y" values are estimated
density values for "x".
Below are in the R manual
x
the n coordinates of the points where the density is estimated.
y
the estimated density values. These will be non-negative, but can be zero.
We can also check this using the code below.
X <- rnorm(100)
density(X)-> den0
den0
X[1:10]
(den0$x)[1:10]
(den0$y)[1:10]
round(dnorm((den0$x)[1:10]), 6)
round(dnorm(X[1:10]), 6)
Thank you.
Hannah
2012/8/8 David L Carlson <dcarl...@tamu.edu>
> The numbers are there, they just aren't listed by the default print method
> for density. When you type the object name, den0, R runs
> print.density(den0) which provides summary statistics. You need to look at
> the manual page (?density) and pay close attention to the section labeled
> "Value" which provides information about what values are returned by the
> function.
>
> str(den0)
> den0$x
> den0$y
> plot(den0$x, den0$y, typ="l")
>
> ----------------------------------------------
> David L Carlson
> Associate Professor of Anthropology
> Texas A&M University
> College Station, TX 77843-4352
>
>
> > -----Original Message-----
> > From: r-help-boun...@r-project.org [mailto:r-help-bounces@r-
> > project.org] On Behalf Of li li
> > Sent: Wednesday, August 08, 2012 9:03 PM
> > To: r-help
> > Subject: [R] Density
> >
> > 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
> >
> > > den0
> >
> >
> > Call:
> >
> > 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)}
> >
> > >
> >
> > > est
> >
> >
> >
> >
> > Thanks in advance.
> >
> > Hannah
> >
> > [[alternative HTML version deleted]]
> >
> > ______________________________________________
> > R-help@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-<http://www.r-project.org/posting->
> > guide.html
> > and provide commented, minimal, self-contained, reproducible code.
>
>
[[alternative HTML version deleted]]
______________________________________________
R-help@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.
