It would be better to do the (a1*K + a2*K^2 + a3*K^3 + a4*K^4 +
a5*K^5) portion as follows, so that the power operator is not called.
( ( ( (a5 * K + a4) * K + a3) * K + a2) * K + a1) * K
Not only will this be faster (possibly much faster) but it should
also be more likely. This is just the application of Horner's rule
to the polynomial.
Cheers,
Malcolm Smith
On Mar 10, 2006, at 2:46 AM, Paul Rehill wrote:
OK. My numerical techniques are a bit rusty but the following
function converges to at least 4 decimal places for all except the
extremes for the inverse of the standard normal cumulative
distribution. And the extremes are well extreme, and I'm using the
approximation to the Normal CDF anyway. Any comments to improve it
would be appreciated.
Function InvStdCumNormDist (p As Double) As Double
Const a1 = 0.31938153
Const a2 = -0.356563782
Const a3 = 1.781477937
Const a4 = -1.821255978
Const a5 = 1.330274429
Dim L , K, w, dLow,dHigh, d1 As double
If p<=0 Then p=0.0001
If p>=1 Then p=0.9999
If p<0.5 Then
d1=1-p
Else
d1=p
End If
dLow=0.0
dHigh=5.0
Do
L=(dLow+dHigh)/2
K=1.0 / (1.0 + 0.2316419 * L)
w=1.0 - 1.0 / Sqrt(2.0 * Pi) * Exp(-0.5*L^2) * (a1*K + a2*K^2 +
a3*K^3 + a4*K^4 + a5*K^5)
If d1-w>0.0 Then
dLow=L
Else
dHigh=L
End If
Loop Until Abs(d1-w)<0.000001
If p<0.5 Then
Return -L
Else
Return L
End If
End Function
Regards
Paul
On 3/10/06, Paul Rehill <[EMAIL PROTECTED]> wrote:
Karen
One more question. Do you have a good function for the inverse of
the standard normal cumulative distribution?
_______________________________________________
Unsubscribe or switch delivery mode:
<http://www.realsoftware.com/support/listmanager/>
Search the archives of this list here:
<http://support.realsoftware.com/listarchives/lists.html>
_______________________________________________
Unsubscribe or switch delivery mode:
<http://www.realsoftware.com/support/listmanager/>
Search the archives of this list here:
<http://support.realsoftware.com/listarchives/lists.html>