Skip to content

Density estimate with bounds

6 messages · Justine Rochon, Duncan Murdoch, Giovanni Petris

#
Dear R users,

I would like to show the estimated density of a (0, 1) uniformly distributed
random variable. The density curve, however, goes beyond 0 and 1 because of the
kernel smoothing. 

Example:

x = runif(10000)
plot(density(x))

Is there a way to estimate the density curve strictly within (0, 1) and still
use some sort of smoothing?

Any help would be greatly appreciated.

Best regards,

Justine Rochon



________________________
Justine Rochon
- Biostatistician -
Center for Clinical Studies
University Hospital Regensburg 
Franz-Josef-Strau?-Allee 11
D-93053 Regensburg
Phone: ++49-(0)941-944-5626
Fax: ++49-(0)941-944-5632
Email: justine.rochon at klinik.uni-regensburg.de
#
On 05/11/2009 4:35 AM, Justine Rochon wrote:
One way is to extend the data by reflection on each end.  That is,

x <- runif(10000)
ex_x <- c(-x, x, 2-x)
den <- density(ex_x)
plot(den$x, 3*den$y, xlim=c(0,1), type="l")

You need the rescaling to 3*den$y because you've tripled the range.

Duncan Murdoch
#
Hi Duncan,

Thank you for your e-mail.

It works for the uniform distribution, but I have trouble with the exponential
distribution: 

x <- rexp(10000)
ex_x <- c(-x, x)
den <- density(ex_x)
plot(den$x, 2*den$y, xlim=c(0,5), type="l")

Best regards,

Justine






________________________
Justine Rochon
- Biostatistician -
Center for Clinical Studies
University Hospital Regensburg 
Franz-Josef-Strau?-Allee 11
D-93053 Regensburg
Phone: ++49-(0)941-944-5626
Fax: ++49-(0)941-944-5632
Email: justine.rochon at klinik.uni-regensburg.de
On 05/11/2009 4:35 AM, Justine Rochon wrote:
distributed
the
still
One way is to extend the data by reflection on each end.  That is,

x <- runif(10000)
ex_x <- c(-x, x, 2-x)
den <- density(ex_x)
plot(den$x, 3*den$y, xlim=c(0,1), type="l")

You need the rescaling to 3*den$y because you've tripled the range.

Duncan Murdoch
#
On 11/5/2009 8:36 AM, Justine Rochon wrote:
Just don't plot the out-of-range values.  For example,

keep <- den$x >= 0
plot(den$x[keep], 2*den$y[keep], type="l")

It doesn't do a good job of estimating the density right near zero; for 
that, you'd need to pre-transform to get it flat, then transform back.
For example, if you knew it was like an exponential near zero, you could 
do the following:

u <- 1-exp(-x)  # transform to uniform
u <- c(-u, u)   # reflect
x <- -log(1-u)  # transform back

If you get the scale wrong (i.e. it isn't Exp(1), but it is like some 
other exponential near zero), this should still be okay if you're close 
or you have lots of data.  If you are way off (e.g. some other shape of 
gamma distribution) it won't work well at all.

Duncan Murdoch
#
I think the cleanest solution would be to use kernels whose support is
(included in) (0,1). For example you could use Beta kernels instead of
Normal kernels. Implementing this approach requires a little work,
though -- but not much if you just want a density estimate.

Best,
Giovanni

ps: please do not post the same question multiple times.

  
    
#
What is the problem that you see?
(I think the trouble is that Laplace density is not smooth at x=0,
while a kernel density estimate returns a smooth density)

It would help to know what you are trying to achieve with this
exercise... 

Giovanni