# [R] followup: lookup function for density(...) objects

Mon May 14 21:13:42 CEST 2001

```Thanks to Ross Ihaka and Bob Wheeler for responding to my earlier question.  I
looked into the Johnson system functions in SuppDists package.  For now, I want
to stick with the density(...) estimator, and so still need the variate lookup
function.

As per Ross' suggestion, I just did a numerical integration on the density
object and used approxfun/splinefun to "lookup" the variate value for specified
cumulative probability.  Here's the function.  Comments are welcome.

--------
#
# return the variate(s) for specified cumulative
# probability(ies) for a density(...) object
#
#
RPDensityLookup <- function(density.obj, p=0.5, interp="l")
{
if(missing(density.obj))
{
cat("\nUsage: RPDensityLookup(density.obj, p=0.5, interp=\"l\")")
cat("\n  density.obj: an object returned from the density(...) function")
cat("\n            p: cumulative probability value(s) for which variate")
cat("\n               values are needed")
cat("\n       interp: \"l\" for linear interpolation between data points")
cat("\n               \"s\" for spline fit using data points\n\n")
}

n <- length(density.obj\$x)
dx <- density.obj\$x[2:n] - density.obj\$x[1:(n-1)]
x <- rep(NA, length(p))

# midpoints
cumx <- (density.obj\$x[2:n] + density.obj\$x[1:(n-1)]) / 2
# numerical integration
cumy <- cumsum((density.obj\$y[2:n] + density.obj\$y[1:(n-1)]) / 2 * dx)

if(interp == "l")
{
Fx <- approxfun(cumy, cumx)  # linear interpolation function
x <- Fx(p)                   # the variates corresponding to p's
}
else if(interp == "s")
{
Fx <- splinefun(cumy, cumx)  # spline fit
x <- Fx(p)                   # the variates corresponding to p's
}
else
{
cat(paste("\nError: unrecognized interp method (\"", interp, "\").\n\n",
sep=""))
}

x
}
--------

Example:

> x _ rnorm(1000)
> d.x _ density(x)
> p _ seq(0.001,0.999,by=0.001)
> x.l _ RPDensityLookup(d.x, p, "l")
> x.s _ RPDensityLookup(d.x, p, "s")
> hist(x,probability=T,ylim=c(0,1)); box()
> lines(d.x,col="blue")
> lines(x.l,p,col="green")
> lines(x.s,p,col="red")
> RPDensityLookup(d.x, 0.5, "l")
 -0.001036816
> RPDensityLookup(d.x, 0.5, "s")
 -0.001036673
> RPDensityLookup(d.x, 0.95, "l")
 1.751115
> RPDensityLookup(d.x, 0.95, "s")
 1.751107

Rajiv.
--------
Postdoctoral Research Associate, Hydrology Group
Pacific Northwest National Laboratory, P.O. Box 999, MSIN K9-33
Richland, WA 99352
Voice: (509) 375-2096  Fax: (509) 372-6089  Email: rajiv.prasad at pnl.gov

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._

```