gh-132893: More accurate CDF computation

Lib/statistics.py

@@ -138,7 +138,7 @@
 from decimal import Decimal
 from itertools import count, groupby, repeat
 from bisect import bisect_left, bisect_right
-from math import hypot, sqrt, fabs, exp, erf, tau, log, fsum, sumprod
+from math import hypot, sqrt, fabs, exp, erfc, tau, log, fsum, sumprod
 from math import isfinite, isinf, pi, cos, sin, tan, cosh, asin, atan, acos
 from functools import reduce
 from operator import itemgetter
@@ -811,8 +811,9 @@ def deco(builder):
 def normal_kernel():
     sqrt2pi = sqrt(2 * pi)
     sqrt2 = sqrt(2)
+    neg_sqrt2 = -sqrt2
     pdf = lambda t: exp(-1/2 * t * t) / sqrt2pi
-    cdf = lambda t: 1/2 * (1.0 + erf(t / sqrt2))
+    cdf = lambda t: 1/2 * erfc(t / neg_sqrt2)
     invcdf = lambda t: _normal_dist_inv_cdf(t, 0.0, 1.0)
     support = None
     return pdf, cdf, invcdf, support
@@ -1257,7 +1258,7 @@ def cdf(self, x):
         "Cumulative distribution function.  P(X <= x)"
         if not self._sigma:
             raise StatisticsError('cdf() not defined when sigma is zero')
-        return 0.5 * (1.0 + erf((x - self._mu) / (self._sigma * _SQRT2)))
+        return 0.5 * erfc((self._mu - x) / (self._sigma * _SQRT2))
 
     def inv_cdf(self, p):
         """Inverse cumulative distribution function.  x : P(X <= x) = p
@@ -1311,7 +1312,7 @@ def overlap(self, other):
         dv = Y_var - X_var
         dm = fabs(Y._mu - X._mu)
         if not dv:
-            return 1.0 - erf(dm / (2.0 * X._sigma * _SQRT2))
+            return erfc(dm / (2.0 * X._sigma * _SQRT2))
         a = X._mu * Y_var - Y._mu * X_var
         b = X._sigma * Y._sigma * sqrt(dm * dm + dv * log(Y_var / X_var))
         x1 = (a + b) / dv

Misc/NEWS.d/next/Library/2025-04-24-21-22-46.gh-issue-132893.KFuxZ2.rst

@@ -0,0 +1 @@
+Improved the accuracy of ``statistics.NormalDist.cdf``.