This page provides stand-alone Haskell code for computing the quantile function, or inverse CDF, for the standard normal probability distribution. The accuracy is modest, but the code is brief and completely self-contained; you can simply copy and paste it into your project.
If you’re willing to take an external dependence, you can get higher accuracy from math-functions and use the formulas here to implement the inverse of the normal CDF in terms of that module’s function invErf.
The code below is based on Formula 26.2.23 from Abramowitz and Stegun.
The code is simply as follows:
rational_approx :: Double -> Double
rational_approx t = t - ((0.010328*t + 0.802853)*t + 2.515517) /
(((0.001308*t + 0.189269)*t + 1.432788)*t + 1.0);
phi_inverse :: Double -> Double
phi_inverse p =
if p < 0.5
then - rational_approx( sqrt (-2.0*log(p)) )
else rational_approx( sqrt (-2.0*log(1.0 - p)) )
The test code is longer than the code itself.
test_phi_inverse :: Bool
test_phi_inverse = maximum error < epsilon
where
error = [ abs(phi_inverse x - y) | (x, y) <- zip xs ys ]
xs = [0.0000001,
0.00001,
0.001,
0.05,
0.15,
0.25,
0.35,
0.45,
0.55,
0.65,
0.75,
0.85,
0.95,
0.999,
0.99999,
0.9999999]
ys = [-5.199337582187471,
-4.264890793922602,
-3.090232306167813,
-1.6448536269514729,
-1.0364333894937896,
-0.6744897501960817,
-0.38532046640756773,
-0.12566134685507402,
0.12566134685507402,
0.38532046640756773,
0.6744897501960817,
1.0364333894937896,
1.6448536269514729,
3.090232306167813,
4.264890793922602,
5.199337582187471]
epsilon = 4.5e-4 -- Accuracy advertised in A&S 26.2.23