how R implement qnorm()
On Fri, Oct 19, 2012 at 12:21 PM, Sheng Liu <sheng.liu at live.ca> wrote:
Thanks a lot. It's very helpful. I've read through the c code. Just FYI and for my completion of the question, I post some of my thought on it: To me it looks like the algorithm is actually inquiring through an approximation table (the approximations, at least for pnom, is "derived from those in "Rational Chebyshev approximations for the error function" by W. J. Cody, Math. Comp., 1969, 631-637."), it is not using any function such as intergration or inverse erf to compute the value based on an equation as I thought, though lack a bit subtlety but the up side is the speed.
It's not an approximation table, it's a rational Chebyshev approximation, a ratio of polynomials. qnorm also uses a rational polynomial approximation. At some level this *has to be* what is going on: computers don't implement pnorm/qnorm or erfc in hardware, and there is no closed-form expression for them in terms of quantities that are implemented in hardware, so the functions must be some sort of approximate expression using arithmetic and other hardware computations.
From the point of view of R, the only relevant issues are precision,
speed, and portability, and these approximations are portable,
accurate, and fast.
-thomas
Hope this can help others who had similar questions as well. Sheng On Thu, Oct 18, 2012 at 2:32 AM, peter dalgaard <pdalgd at gmail.com> wrote:
On Oct 18, 2012, at 09:55 , Prof Brian Ripley wrote:
R is a bit confusing as it requires inverse error function (X = - sqrt(2)* erf-1 (2*P)), while R doesn't have a build in one. The InvErf function most people use is through qnorm( InvErf=function(x)
I think you are wrong about 'most people': this is the notation used by
a small group of non-statisticians (mainly physicists, I think). Well, he's right in the sense that it is erf/erfc that are commonly documented in collections of special functions (like Abramowitz & Stegun), and also those that appear in common C math libraries. In both cases of course because physicists have dominated in their development. Of course "most people" use Excel these days (which has the inverse normal distribution but AFAICS not the inverse error function). -- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Email: pd.mes at cbs.dk Priv: PDalgd at gmail.com
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[[alternative HTML version deleted]]
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Thomas Lumley Professor of Biostatistics University of Auckland