ENH: Add rng to cluster permutation tests

Author: bruAristimunha

mne/stats/cluster_level.py

@@ -12,6 +12,7 @@
 from ..utils import (
     ProgressBar,
     _check_option,
+    _check_rng_compat,
     _pl,
     _validate_type,
     check_random_state,
@@ -1120,6 +1121,8 @@ def permutation_cluster_test(
     check_disjoint=False,
     buffer_size=1000,
     verbose=None,
+    *,
+    rng=None,
 ):
     """Cluster-level statistical permutation test.
 
@@ -1151,7 +1154,9 @@ def permutation_cluster_test(
     %(stat_fun_clust_f)s
     %(adjacency_clust_n)s
     %(n_jobs)s
-    %(seed)s
+    %(rng)s
+    seed : None | int | instance of ~numpy.random.RandomState
+        Deprecated. Use ``rng`` instead.
     %(max_step_clust)s
     %(exclude_clust)s
     %(step_down_p_clust)s
@@ -1181,6 +1186,7 @@ def permutation_cluster_test(
     .. footbibliography::
     """
     stat_fun, threshold = _check_fun(X, stat_fun, threshold, tail, "between")
+    rng = _check_rng_compat(rng, legacy=seed, legacy_name="seed")
     return _permutation_cluster_test(
         X=X,
         threshold=threshold,
@@ -1189,7 +1195,7 @@ def permutation_cluster_test(
         stat_fun=stat_fun,
         adjacency=adjacency,
         n_jobs=n_jobs,
-        seed=seed,
+        seed=rng,
         max_step=max_step,
         exclude=exclude,
         step_down_p=step_down_p,
@@ -1218,6 +1224,8 @@ def permutation_cluster_1samp_test(
     check_disjoint=False,
     buffer_size=1000,
     verbose=None,
+    *,
+    rng=None,
 ):
     """Non-parametric cluster-level paired t-test.
 
@@ -1238,7 +1246,9 @@ def permutation_cluster_1samp_test(
     %(stat_fun_clust_t)s
     %(adjacency_clust_1)s
     %(n_jobs)s
-    %(seed)s
+    %(rng)s
+    seed : None | int | instance of ~numpy.random.RandomState
+        Deprecated. Use ``rng`` instead.
     %(max_step_clust)s
     %(exclude_clust)s
     %(step_down_p_clust)s
@@ -1291,6 +1301,7 @@ def permutation_cluster_1samp_test(
     .. footbibliography::
     """
     stat_fun, threshold = _check_fun(X, stat_fun, threshold, tail)
+    rng = _check_rng_compat(rng, legacy=seed, legacy_name="seed")
     return _permutation_cluster_test(
         X=[X],
         threshold=threshold,
@@ -1299,7 +1310,7 @@ def permutation_cluster_1samp_test(
         stat_fun=stat_fun,
         adjacency=adjacency,
         n_jobs=n_jobs,
-        seed=seed,
+        seed=rng,
         max_step=max_step,
         exclude=exclude,
         step_down_p=step_down_p,
@@ -1328,6 +1339,8 @@ def spatio_temporal_cluster_1samp_test(
     check_disjoint=False,
     buffer_size=1000,
     verbose=None,
+    *,
+    rng=None,
 ):
     """Non-parametric cluster-level paired t-test for spatio-temporal data.
 
@@ -1351,7 +1364,9 @@ def spatio_temporal_cluster_1samp_test(
     %(stat_fun_clust_t)s
     %(adjacency_clust_st1)s
     %(n_jobs)s
-    %(seed)s
+    %(rng)s
+    seed : None | int | instance of ~numpy.random.RandomState
+        Deprecated. Use ``rng`` instead.
     %(max_step_clust)s
     spatial_exclude : list of int or None
         List of spatial indices to exclude from clustering.
@@ -1397,6 +1412,7 @@ def spatio_temporal_cluster_1samp_test(
         adjacency=adjacency,
         n_jobs=n_jobs,
         seed=seed,
+        rng=rng,
         max_step=max_step,
         exclude=exclude,
         step_down_p=step_down_p,
@@ -1425,6 +1441,8 @@ def spatio_temporal_cluster_test(
     check_disjoint=False,
     buffer_size=1000,
     verbose=None,
+    *,
+    rng=None,
 ):
     """Non-parametric cluster-level test for spatio-temporal data.
 
@@ -1450,7 +1468,9 @@ def spatio_temporal_cluster_test(
     %(stat_fun_clust_f)s
     %(adjacency_clust_stn)s
     %(n_jobs)s
-    %(seed)s
+    %(rng)s
+    seed : None | int | instance of ~numpy.random.RandomState
+        Deprecated. Use ``rng`` instead.
     %(max_step_clust)s
     spatial_exclude : list of int or None
         List of spatial indices to exclude from clustering.
@@ -1496,6 +1516,7 @@ def spatio_temporal_cluster_test(
         adjacency=adjacency,
         n_jobs=n_jobs,
         seed=seed,
+        rng=rng,
         max_step=max_step,
         exclude=exclude,
         step_down_p=step_down_p,

mne/stats/tests/test_cluster_level.py

@@ -56,6 +56,15 @@ def _get_conditions():
     return condition1_1d, condition2_1d, condition1_2d, condition2_2d
 
 
+def test_cluster_rng_transition():
+    """Test the transition from seed to rng."""
+    X = np.arange(24.0).reshape(8, 3)
+    with pytest.warns(FutureWarning, match="seed"):
+        permutation_cluster_1samp_test(X, threshold=0, n_permutations=2, seed=0)
+    with pytest.raises(TypeError, match="Specify only one"):
+        permutation_cluster_1samp_test(X, threshold=0, n_permutations=2, seed=0, rng=0)
+
+
 def test_thresholds(numba_conditional):
     """Test automatic threshold calculations."""
     # within subjects
@@ -69,7 +78,7 @@ def test_thresholds(numba_conditional):
     with catch_logging() as log:
         with pytest.warns(RuntimeWarning, match="threshold is only valid"):
             out = permutation_cluster_1samp_test(
-                X, stat_fun=my_fun, seed=0, verbose=True, out_type="mask"
+                X, stat_fun=my_fun, rng=0, verbose=True, out_type="mask"
             )
     log = log.getvalue()
     assert str(want_thresh)[:6] in log
@@ -86,12 +95,12 @@ def test_thresholds(numba_conditional):
     with catch_logging() as log:
         with pytest.warns(RuntimeWarning, match="threshold is only valid"):
             out = permutation_cluster_test(
-                X, tail=1, stat_fun=my_fun, seed=0, verbose=True, out_type="mask"
+                X, tail=1, stat_fun=my_fun, rng=0, verbose=True, out_type="mask"
             )
     log = log.getvalue()
     assert str(want_thresh)[:6] in log
     assert len(out[1]) == 1  # 1 cluster
-    assert_allclose(out[2], 0.031250, atol=1e-6)
+    assert_allclose(out[2], 0.03515625, atol=1e-6)
     with pytest.warns(RuntimeWarning, match='Ignoring argument "tail"'):
         permutation_cluster_test(X, tail=0, out_type="mask")
 
@@ -103,7 +112,7 @@ def test_thresholds(numba_conditional):
         pytest.warns(RuntimeWarning, match="invalid value"),
     ):  # NumPy
         out = permutation_cluster_1samp_test(
-            X, seed=0, threshold=dict(start=0, step=0.1), out_type="mask"
+            X, rng=0, threshold=dict(start=0, step=0.1), out_type="mask"
         )
     assert (out[2] < 0.05).any()
     assert not (out[2] < 0.05).all()
@@ -112,7 +121,7 @@ def test_thresholds(numba_conditional):
         with np.errstate(invalid="ignore"):
             permutation_cluster_1samp_test(
                 X,
-                seed=0,
+                rng=0,
                 threshold=dict(start=0, step=0.1),
                 buffer_size=None,
                 out_type="mask",
@@ -137,7 +146,7 @@ def test_cache_dir(tmp_path, numba_conditional):
                 buffer_size=None,
                 n_jobs=2,
                 n_permutations=1,
-                seed=0,
+                rng=0,
                 stat_fun=ttest_1samp_no_p,
                 verbose=False,
                 out_type="mask",
@@ -152,7 +161,7 @@ def test_cache_dir(tmp_path, numba_conditional):
                 buffer_size=10,
                 n_jobs=2,
                 n_permutations=1,
-                seed=random_state,
+                rng=random_state,
                 stat_fun=stat_fun,
                 verbose=False,
                 out_type="mask",
@@ -175,7 +184,7 @@ def test_permutation_large_n_samples(numba_conditional):
         tails = (0, 1) if n_samples <= 20 else (0,)
         for tail in tails:
             H0 = permutation_cluster_1samp_test(
-                X[:n_samples], threshold=1e-4, tail=tail, seed=0, out_type="mask"
+                X[:n_samples], threshold=1e-4, tail=tail, rng=0, out_type="mask"
             )[-1]
             assert H0.shape == (1024,)
             assert len(np.unique(H0)) >= 1024 - (H0 == 0).sum()
@@ -221,7 +230,7 @@ def test_cluster_permutation_test(numba_conditional):
             [condition1, condition2],
             n_permutations=100,
             tail=1,
-            seed=1,
+            rng=1,
             buffer_size=None,
             out_type="mask",
         )
@@ -235,7 +244,7 @@ def test_cluster_permutation_test(numba_conditional):
             [condition1, condition2],
             n_permutations=100,
             tail=1,
-            seed=1,
+            rng=1,
             n_jobs=2,
             buffer_size=buffer_size,
             out_type="mask",
@@ -268,7 +277,7 @@ def test_cluster_permutation_t_test(numba_conditional, stat_fun):
             condition1,
             n_permutations=100,
             tail=0,
-            seed=1,
+            rng=1,
             out_type="mask",
             buffer_size=None,
         )
@@ -281,7 +290,7 @@ def test_cluster_permutation_t_test(numba_conditional, stat_fun):
             n_permutations=100,
             tail=1,
             threshold=1.67,
-            seed=1,
+            rng=1,
             stat_fun=stat_fun,
             out_type="mask",
             buffer_size=None,
@@ -292,7 +301,7 @@ def test_cluster_permutation_t_test(numba_conditional, stat_fun):
             n_permutations=100,
             tail=-1,
             threshold=-1.67,
-            seed=1,
+            rng=1,
             stat_fun=stat_fun,
             buffer_size=None,
             out_type="mask",
@@ -314,7 +323,7 @@ def test_cluster_permutation_t_test(numba_conditional, stat_fun):
                 tail=-1,
                 out_type="mask",
                 threshold=-1.67,
-                seed=1,
+                rng=1,
                 n_jobs=2,
                 stat_fun=stat_fun,
                 buffer_size=buffer_size,
@@ -347,7 +356,7 @@ def test_cluster_permutation_with_adjacency(numba_conditional, monkeypatch):
     n_pts = condition1_1d.shape[1]
     # we don't care about p-values in any of these, so do fewer permutations
     args = dict(
-        seed=None,
+        rng=None,
         max_step=1,
         exclude=None,
         out_type="mask",
@@ -610,7 +619,7 @@ def test_permutation_adjacency_equiv(numba_conditional):
             n_jobs=2,
             max_step=max_step,
             stat_fun=stat_fun,
-            seed=0,
+            rng=0,
             out_type="mask",
         )
         # make sure our output datatype is correct
@@ -677,7 +686,7 @@ def test_spatio_temporal_cluster_chain_merge():
         max_step=1,
         n_permutations=20,
         out_type="indices",
-        seed=0,
+        rng=0,
         verbose=False,
     )
     assert len(clusters) == 1
@@ -760,7 +769,7 @@ def test_spatio_temporal_cluster_adjacency(numba_conditional):
         adjacency=adj,
         n_permutations=50,
         tail=1,
-        seed=1,
+        rng=1,
         threshold=threshold,
         buffer_size=None,
     )
@@ -770,7 +779,7 @@ def test_spatio_temporal_cluster_adjacency(numba_conditional):
         [data1_2d, data2_2d],
         n_permutations=50,
         tail=1,
-        seed=1,
+        rng=1,
         threshold=threshold,
         n_jobs=2,
         buffer_size=buffer_size,
@@ -783,7 +792,7 @@ def test_spatio_temporal_cluster_adjacency(numba_conditional):
         [data1_2d, data2_2d],
         n_permutations=50,
         tail=1,
-        seed=1,
+        rng=1,
         threshold=threshold,
         n_jobs=2,
         buffer_size=None,
@@ -866,19 +875,19 @@ def test_permutation_test_H0(numba_conditional):
     data = rng.random((7, 10, 1)) - 0.5
     with pytest.warns(RuntimeWarning, match="No clusters found"):
         t, clust, p, h0 = spatio_temporal_cluster_1samp_test(
-            data, threshold=100, n_permutations=1024, seed=rng
+            data, threshold=100, n_permutations=1024, rng=rng
         )
     assert_equal(len(h0), 0)
 
     for n_permutations in (1024, 65, 64, 63):
         t, clust, p, h0 = spatio_temporal_cluster_1samp_test(
-            data, threshold=0.1, n_permutations=n_permutations, seed=rng
+            data, threshold=0.1, n_permutations=n_permutations, rng=rng
         )
         assert_equal(len(h0), min(n_permutations, 64))
         assert isinstance(clust[0], tuple)  # sets of indices
     for tail, thresh in zip((-1, 0, 1), (-0.1, 0.1, 0.1)):
         t, clust, p, h0 = spatio_temporal_cluster_1samp_test(
-            data, threshold=thresh, seed=rng, tail=tail, out_type="mask"
+            data, threshold=thresh, rng=rng, tail=tail, out_type="mask"
         )
         assert isinstance(clust[0], np.ndarray)  # bool mask
         # same as "128 if tail else 64"