[R] Antw: Re: Density estimate with bounds

[Duncan Murdoch]

On 11/5/2009 8:36 AM, Justine Rochon wrote:
> 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")

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

> 
> 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 
>  
>  
>>>> Duncan Murdoch <murdoch at stats.uwo.ca> 05.11.2009 12:36 >>>
> On 05/11/2009 4:35 AM, Justine Rochon wrote:
>> 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.
> 
> 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