1
$\begingroup$

I used R to find kernel density estimates of my dataset (for experiment I used 1000 samples generated from a known distribution in this step).

I used code density() to find the kernel density estimate. Now I need to find the cdf of this pdf in order to find quantiles of the distribution. I used these codes for normal:

x=rnorm(1000,8,3)
pdf=density(x)
f=approxfun(pdf$x,pdf$y,yleft=0,yright=0)
FF2=function(x) integrate(f,-Inf,x)$v-0.25
uniroot(FF2,c(2,10))

By the last code I could find the first quantile of the distribution and this is what I want. But my question is about another distribution like gamma. These codes do not run for this distribution as here:

x1=rgamma(1000,shape=2,scale=1)
pdf1=density(x1)
f1=approxfun(pdf1$x,pdf1$y,yleft=0,yright=0)
FF1=function(x) integrate(f1,-Inf,x)$value-0.25

The result of the last code is not logical!

$\endgroup$
  • 3
    $\begingroup$ Which quantiles do you want? I see 0.25 in your code and wonder if you want quartiles. (The "first quantile" makes no sense.) If so, or even if not so, why would quartiles from a kernel density estimate be better than, or in any sense preferable to, quartiles from the raw data? There could be a good reason, but please tell us what it is. You could also smooth the quantile function directly. $\endgroup$ – Nick Cox Jun 15 '18 at 8:21
  • $\begingroup$ To second @NickCox comment, the only way that a density estimate would provide more precision in quantile estimates would be if you had injected prior knowledge into the density estimation process. If density estimation used 3 or less effective degrees of freedom you'd expect it to be more efficient, if the structure imposed on the density happened to be right. $\endgroup$ – Frank Harrell Jun 15 '18 at 11:15
0
$\begingroup$

The true value is

> qgamma(.25,shape=2,scale=1)
[1] 0.9612788

so searching over the support $(2,10)$ does not seem like a good idea. For your normal example, the true quantile is

> qnorm(.25,8,3)
[1] 5.976531

so that the support $(2,10)$ contains the solution, explaining why the first example works but the second does not.

Also, the gamma distribution has support on $(0,\infty)$, so that the integration in FF1 ought to be started at 0.

That is,

x1=rgamma(1000,shape=2,scale=1)
pdf1=density(x1)
f1=approxfun(pdf1$x,pdf1$y,yleft=0,yright=0)
FF1=function(x) integrate(f1,0,x)$value-0.25
uniroot(FF1,c(0,10))

works for me (starting the integration at 0 seems optional, as the distribution has no mass in the negative values).

$\endgroup$
  • $\begingroup$ thanks for your reply, but I am not searching for the real quantiles of the gamma distribution. I want to find the approximated distribution. I mean I would like to fit a kernel density distribution and then find this density's quantiles so that I can compare them with the real ones. that is why I am doing so. $\endgroup$ – fattmim Jun 15 '18 at 6:48
  • $\begingroup$ I think I got that - I printed the "real" quantiles in order to emphasize issues I see in your code, namely that the support you specify in uniroot does not seem well-chosen - see also my edit to my answer. $\endgroup$ – Christoph Hanck Jun 15 '18 at 7:34

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.