Repository: scipy/scipy
Pull request: MAINT: stats: raise error when input to kruskal has >1 dim
#### Reference issue
Closes gh-9981

#### What does this implement/fix?
`scipy.stats.kruskal` only works on 1D arrays (and scalars), but this isn't documented, and no error is raised if this condition is not met. This PR adds a check that the samples have no more than 1 non-singleton dimension

Problem statement:
stats.kruskal : add a warning for an input with 2 or more columns
<!-- 
Thank you for taking the time to report a SciPy issue.
Please describe the issue in detail, and for bug reports
fill in the fields below. You can delete the sections that 
don't apply to your issue. You can view the final output
by clicking the preview button above.
-->

if a dataframe  has two columns, results are weird. For instance let the dataframe x be (first column = indices):
         counts  group
3803673      31      0
3800421      24      0
3800446      21      0
3804785      11      0
3801989      31      0
3801364      15      0
3804710      12      0
3800049      20      0
3801342      11      0
3804642       4      0

```stats.kruskal(x, x)```
gives this result: 
KruskalResult(statistic=205.08158702862053, pvalue=1.6254869278866505e-46)


I expect a warning message that the test is only suited for a series or a dataframe with one column.

Hints from issue discussion:
I can confirm that the behavior is strange. When a 2D array goes in, it's not doing any of the things I might expect:
```
import numpy as np
from scipy.stats import kruskal

np.random.seed(0)
x = np.random.rand(10, 2)

print(kruskal(x, x))  # this is not the same as any of the following
print(kruskal(x.ravel(), x.ravel()))
print(kruskal(*x))
print(kruskal(*x, *x))
print(kruskal(*x.T))
print(kruskal(*x.T, *x.T))
```
```
KruskalResult(statistic=177.4987325456498, pvalue=1.7043908031353597e-40)
KruskalResult(statistic=1.4237566848125315e-14, pvalue=0.9999999047953962)
KruskalResult(statistic=6.799999999999997, pvalue=0.6579333062136292)
KruskalResult(statistic=13.95789473684211, pvalue=0.7861402653468463)
KruskalResult(statistic=0.8228571428571456, pvalue=0.36434612663355204)
KruskalResult(statistic=1.6890225563909933, pvalue=0.6393770842915287)
```

Commit message: MAINT: stats: raise error when input to kruskal has >1 dim

Patch to review:
diff --git a/scipy/stats/stats.py b/scipy/stats/stats.py
index c3cd1f7818bb..46612d33b6cd 100644
--- a/scipy/stats/stats.py
+++ b/scipy/stats/stats.py
@@ -7271,7 +7271,7 @@ def kruskal(*args, nan_policy='propagate'):
     ----------
     sample1, sample2, ... : array_like
        Two or more arrays with the sample measurements can be given as
-       arguments.
+       arguments. Samples may not have more than one non-singleton dimension.
     nan_policy : {'propagate', 'raise', 'omit'}, optional
         Defines how to handle when input contains nan.
         The following options are available (default is 'propagate'):
@@ -7324,6 +7324,7 @@ def kruskal(*args, nan_policy='propagate'):
 
     """
     args = list(map(np.asarray, args))
+
     num_groups = len(args)
     if num_groups < 2:
         raise ValueError("Need at least two groups in stats.kruskal()")
@@ -7331,6 +7332,10 @@ def kruskal(*args, nan_policy='propagate'):
     for arg in args:
         if arg.size == 0:
             return KruskalResult(np.nan, np.nan)
+        elif np.squeeze(arg).ndim > 1:  # 0d (scalar) is ok
+            raise ValueError("Samples may not have more than one "
+                             "non-singleton dimension.")
+
     n = np.asarray(list(map(len, args)))
 
     if nan_policy not in ('propagate', 'raise', 'omit'):
diff --git a/scipy/stats/tests/test_stats.py b/scipy/stats/tests/test_stats.py
index b9b801ba3589..e2fe4398297d 100644
--- a/scipy/stats/tests/test_stats.py
+++ b/scipy/stats/tests/test_stats.py
@@ -5187,6 +5187,19 @@ def test_empty(self):
         z = []
         assert_equal(stats.kruskal(x, y, z), (np.nan, np.nan))
 
+    def test_nd_arrays(self):
+        # N-d arrays are OK...
+        x = [[1]]
+        y = [[2]]  # this is fine
+        h, p = stats.kruskal(x, y)
+        assert_equal(h, 1.0)
+        assert_approx_equal(p, stats.distributions.chi2.sf(h, 1))
+
+        # but only one non-singleton dimension is allowed
+        z = np.random.rand(2, 2)
+        with assert_raises(ValueError, match="Samples may not have more"):
+            stats.kruskal(x, y, z)
+
     def test_kruskal_result_attributes(self):
         x = [1, 3, 5, 7, 9]
         y = [2, 4, 6, 8, 10]


Code review comments:
diff_hunk: @@ -7324,13 +7324,18 @@ def kruskal(*args, nan_policy='propagate'):
 
     """
     args = list(map(np.asarray, args))
+
     num_groups = len(args)
     if num_groups < 2:
         raise ValueError("Need at least two groups in stats.kruskal()")
 
     for arg in args:
         if arg.size == 0:
             return KruskalResult(np.nan, np.nan)
+        elif np.squeeze(arg).ndim > 1:  # 0d (scalar) is ok
path: scipy/stats/stats.py
text: @author:
```suggestion
        elif arg.ndim != 1: 
```

diff_hunk: @@ -5187,6 +5187,19 @@ def test_empty(self):
         z = []
         assert_equal(stats.kruskal(x, y, z), (np.nan, np.nan))
 
+    def test_nd_arrays(self):
+        # N-d arrays are OK...
+        x = [[1]]
+        y = [[2]]  # this is fine
path: scipy/stats/tests/test_stats.py
text: @author:
```suggestion
        x = [1]
        y = [2]
```

diff_hunk: @@ -7324,13 +7324,18 @@ def kruskal(*args, nan_policy='propagate'):
 
     """
     args = list(map(np.asarray, args))
+
     num_groups = len(args)
     if num_groups < 2:
         raise ValueError("Need at least two groups in stats.kruskal()")
 
     for arg in args:
         if arg.size == 0:
             return KruskalResult(np.nan, np.nan)
+        elif np.squeeze(arg).ndim > 1:  # 0d (scalar) is ok
path: scipy/stats/stats.py
text: @user1:
in a few other tests (e.g. in _hypotests.py), an error is raised if `arg.ndim > 1`. I wonder if we should follow the same approach here, I did not understand the statement 'more than one non-singleton dimension' right away. it is a minor point, I am also fine with the current version.

diff_hunk: @@ -7324,13 +7324,18 @@ def kruskal(*args, nan_policy='propagate'):
 
     """
     args = list(map(np.asarray, args))
+
     num_groups = len(args)
     if num_groups < 2:
         raise ValueError("Need at least two groups in stats.kruskal()")
 
     for arg in args:
         if arg.size == 0:
             return KruskalResult(np.nan, np.nan)
+        elif np.squeeze(arg).ndim > 1:  # 0d (scalar) is ok
+            raise ValueError("Samples may not have more than one "
+                             "non-singleton dimension.")
path: scipy/stats/stats.py
text: @author:
```suggestion
            raise ValueError("Samples must be one-dimensional.")
```

diff_hunk: @@ -5187,6 +5187,19 @@ def test_empty(self):
         z = []
         assert_equal(stats.kruskal(x, y, z), (np.nan, np.nan))
 
+    def test_nd_arrays(self):
+        # N-d arrays are OK...
+        x = [[1]]
+        y = [[2]]  # this is fine
+        h, p = stats.kruskal(x, y)
+        assert_equal(h, 1.0)
+        assert_approx_equal(p, stats.distributions.chi2.sf(h, 1))
+
+        # but only one non-singleton dimension is allowed
+        z = np.random.rand(2, 2)
+        with assert_raises(ValueError, match="Samples may not have more"):
path: scipy/stats/tests/test_stats.py
text: @author:
```suggestion
        with assert_raises(ValueError, match="Samples must be one-dimensional."):
```

diff_hunk: @@ -7271,7 +7271,7 @@ def kruskal(*args, nan_policy='propagate'):
     ----------
     sample1, sample2, ... : array_like
        Two or more arrays with the sample measurements can be given as
-       arguments.
+       arguments. Samples may not have more than one non-singleton dimension.
path: scipy/stats/stats.py
text: @author:
```suggestion
       arguments. Samples must be one-dimensional.
```