# Numerically compute PDF given a function

Consider $[0,1]$ with the Lebesgue measure $m$ and $f:[0,1]\to \mathbb{R}$, and $x$ a uniformly distributed random variable in $[0,1]$. Then, $f(x)$ itself define a new random variable.

We can then define the Cumulative density function of $f$ as $$P_f (y) :\,= m(f^{-1} (-\infty, y)) \, ,$$ and its PDF as $$\rho _f (y) :\,= \frac{dP_f (y)}{dy} \, ,$$ for all $y\in {\rm range} (f)$.

A straightforward numerical computation of the probability density function (PDF) of $f$ involves inverse derivatives, or, conversely, $1/f'$. As numerical differentiation is numerically ill-conditioned, when $f'$ becomes very low, this is numerically problematic. I'm looking for one of two things:

1. An out of the box function in Matlab that computes a PDF for a given function $f$.
2. A numerical Algorithm for this task.

Thanks

• Can you try to make your question clearer? I don't follow that $f$ is both a random variable and a function, and I don't get the connection between $f$ and the PDF you want to compute. – Wolfgang Bangerth Jan 24 '17 at 13:50
• @WolfgangBangerth is it clearer now? – Amir Sagiv Jan 24 '17 at 19:52
• Why do you say it is "numerically problematic"? Problematic for what purpose? Division is a well-conditioned operation, so dividing by a small number $f'$ will produce an accurate (though large) result, if you have computed $f'$ accurately. Numerical differentiation, on the other hand, is usually ill-conditioned, so what you're asking for might be impossible. – Kirill Jan 24 '17 at 21:00
• @Kirill you are right. I meant to say that the error in numerical differentiation of small derivatives become large when considering $1/f'$. Specifically, this happens a lot in my applications, and so I rather use an $f'$-free algorithm. – Amir Sagiv Jan 24 '17 at 21:27
• @AmirSagiv Then I think there is no reason to think such an algorithm exists. If your problem is ill-conditioned, a well-conditioned numerical algorithm cannot exist (how could it?). Also, since $f'$ is in the definition, you can't really have an $f'$-free algorithm: you could just as easily compute $f'$ by running your algorithm and recovering $f'$. – Kirill Jan 24 '17 at 21:46