BUG: fix stats.kendalltau to handle larger lists. Closes #1445.

rgommers · rgommers · commit ce14ddb13d99 · 2011-06-12T20:55:55.000+02:00
diff --git a/scipy/stats/stats.py b/scipy/stats/stats.py
@@ -2624,7 +2624,7 @@ def kendalltau(x, y, initial_lexsort=True):
 
     x = np.asarray(x).ravel()
     y = np.asarray(y).ravel()
-    n = len(x)
+    n = np.int64(len(x))
     temp = range(n) # support structure used by mergesort
     # this closure recursively sorts sections of perm[] by comparing
     # elements of y[perm[]] using temp[] as support
@@ -2707,8 +2707,9 @@ def mergesort(offs, length):
     if tot == u and tot == v:
         return 1    # Special case for all ties in both ranks
 
-    tau = ((tot - (v + u - t)) - 2.0 * exchanges) / \
-                    np.sqrt((tot - u) * (tot - v))
+    # Prevent overflow; equal to np.sqrt((tot - u) * (tot - v))
+    denom = np.exp(0.5 * (np.log(tot - u) + np.log(tot - v)))
+    tau = ((tot - (v + u - t)) - 2.0 * exchanges) / denom
 
     # what follows reproduces the ending of Gary Strangman's original
     # stats.kendalltau() in SciPy
diff --git a/scipy/stats/tests/test_stats.py b/scipy/stats/tests/test_stats.py
@@ -585,8 +585,8 @@ def test_kendalltau():
     # with some ties
     x1 = [12, 2, 1, 12, 2]
     x2 = [1, 4, 7, 1, 0]
-    res = (-0.47140452079103173, 0.24821309157521476)
-    expected = stats.kendalltau(x1, x2)
+    expected = (-0.47140452079103173, 0.24821309157521476)
+    res = stats.kendalltau(x1, x2)
     assert_approx_equal(res[0], expected[0])
     assert_approx_equal(res[1], expected[1])
 
