9
$\begingroup$

See this Wikipedia page:

http://en.wikipedia.org/wiki/Binomial_proportion_confidence_interval#Agresti-Coull_Interval

To get the Agresti-Coull Interval, one needs to calculate a percentile of the normal distribution, called $z$. How do I calculate tha percentile? Is there a ready-made function that does this in Wolfram Mathematica and/or Python/NumPy/SciPy?

$\endgroup$
  • 1
    $\begingroup$ The integral expression in the "normal cdf I got exactly from Wiki" is unfortunately off by a factor of $1/\sqrt{\pi}$. There is no known exact formula for the normal cdf or its inverse using a finite number of terms involving standard functions ($\exp, \log, \sin \cos$ etc) but both the normal cdf and its inverse have been studied a lot and approximate formulas for both are programmed into many calculator, spreadsheets, not to mention statistical packages. I am not familiar with R but I would be astounded if it did not have what you are looking for built in already. $\endgroup$ – Dilip Sarwate Feb 14 '12 at 20:52
  • $\begingroup$ @DilipSarwate, it's fixed! I am doing this using inverse tranformation, also "not allowed" to use too much built in. It's for the sake of learning I suppose. $\endgroup$ – user1061210 Feb 14 '12 at 21:26
  • 1
    $\begingroup$ @Dilip: Not only is there no known exact formula, better yet, it is known that no such formula can exist! $\endgroup$ – cardinal Feb 14 '12 at 21:30
  • 1
    $\begingroup$ The Box-Muller method generates samples from a joint distribution of independent standard normal random variables. So histograms of the values generated will resemble standard normal distributions. But the Box-Muller method is not a method for computing values of $\Phi(x)$ except incidentally as in "I generated $10^4$ standard normal samples of which $8401$ has value $1$ or less, and so $\Phi(1) \approx 0.8401$, and $\Phi^{-1}(0.8401) \approx 1$. $\endgroup$ – Dilip Sarwate Feb 14 '12 at 22:29
  • 1
    $\begingroup$ I just chose $8401$ as an example of the kinds of numbers you might expect. $\Phi(1) = 0.8413\ldots$ and so if you generate $10^4$ samples of a standard normal distribution, you should expect close to $8413$ of the $10000$ samples to have value $\leq 1$. You are implementing the Box-Muller method correctly, but are not understanding the results that you are getting and are not relating them to the cdf etc. $\endgroup$ – Dilip Sarwate Feb 15 '12 at 0:27
3
$\begingroup$

For Mathematica $VersionNumber > 5 you can use

Quantile[NormalDistribution[μ, σ], 100 q]

for the q-th percentile.

Otherwise, you have to load the appropriate Statistics package first.

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ (I have version 7.) I have no problem loading the Statistics package. But what's the function in there called? Because I get the impression that this Quantile line will do the calculation manually instead of using a formula. $\endgroup$ – Ram Rachum Oct 9 '10 at 14:13
  • $\begingroup$ Evaluate it with symbolic parameters (i.e. don't assign values to mu, sigma, and q); you should get an expression involving the inverse error function. $\endgroup$ – J. M. is not a statistician Oct 9 '10 at 14:24
16
$\begingroup$

John Cook's page, Distributions in Scipy, is a good reference for this type of stuff:

In [15]: import scipy.stats

In [16]: scipy.stats.norm.ppf(0.975)
Out[16]: 1.959963984540054
| cite | improve this answer | |
$\endgroup$
4
$\begingroup$

Well, you didn't ask about R, but in R you do it using ?qnorm

(It's actually the quantile, not the percentile, or so I believe)

> qnorm(.5)
[1] 0
> qnorm(.95)
[1] 1.644854
| cite | improve this answer | |
$\endgroup$
  • 1
    $\begingroup$ Quantile vs. percentile (it's merely a matter of terminology), j.mp/dsYz9z. $\endgroup$ – chl Oct 9 '10 at 14:22
  • 1
    $\begingroup$ While we are in, in R Wald-adjusted CIs (e.g. Agresti-Coull) are available in the PropCIs package. Wilson's method is the default in Hmisc::binconf (as suggested by Agresti and Coull). $\endgroup$ – chl Oct 9 '10 at 14:36
3
$\begingroup$

In Python, you can use the stats module from the scipy package (look for cdf(), as in the following example).

(It seems the transcendantal package also includes usual cumulative distributions).

| cite | improve this answer | |
$\endgroup$
0
$\begingroup$

You can use the inverse erf function, which is available in MatLab and Mathematica, for instance.

For the normal CDF, starting from

$$y=\Phi\left(x\right)=\frac{1}{2}\left[1+\text{erf}\left(\frac{x}{\sqrt{2}}\right)\right]$$

We get

$$x=\sqrt{2}\ \text{erf}^{-1}\left(2y-1\right)$$

For the log-normal CDF, starting from

$$y=F_{x}(x;\mu,\sigma)=\frac{1}{2}\text{erfc}\left(\frac{-\log x-\mu}{\sigma\sqrt{2}}\right)$$

We get

$$-\log \left(x\right)=\mu+\sigma\sqrt{2}\ \text{erfc}^{-1}\left(2y\right)$$

| cite | improve this answer | |
$\endgroup$
  • 2
    $\begingroup$ isn't this more of a comment than an answer? $\endgroup$ – Macro Feb 15 '12 at 22:18
  • $\begingroup$ My idea was that if you have inverses for the erf and erfc functions, then the problem is solved. MatLab, for instance, has such preprogrammed functions. $\endgroup$ – Jean-Victor Côté Mar 5 '12 at 19:31
  • $\begingroup$ @Jean-VictorCôté Please, develop your ideas in your reply. Otherwise, it merely looks like a comment as suggested above. $\endgroup$ – chl Mar 5 '12 at 22:31
  • $\begingroup$ The lognormal calculation doesn't look right. After all, its inverse CDF should be identical to the inverse CDF for the normal apart for the use of $\log(x)$ instead of $x$. $\endgroup$ – whuber Mar 6 '12 at 19:15

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.