two array combinations fix

Author: ryanhammonds

neurodsp/rhythm/phase.py

@@ -63,15 +63,15 @@ def pairwise_phase_consistency(pha0, pha1=None, return_pairs=True, progress=None
 
     else:
 
-        n_combs = int((len(pha0) * (len(pha0) + 1)) / 2)
+        # Include all combinations
+        n_combs = int(len(pha0) ** 2)
 
-        # Include self-combinations
-        iterable = enumerate(combinations_with_replacement(np.arange(len(pha0)), 2))
+        iterable = enumerate((row, col) for row in range(len(pha0)) for col in range(len(pha1)))
 
     # Initialize variables
     if return_pairs:
         cumulative = None
-        distances = np.zeros(n_combs)
+        distances = np.ones((len(pha0), len(pha0)))
     else:
         cumulative = 0
         distances = None
@@ -88,7 +88,7 @@ def pairwise_phase_consistency(pha0, pha1=None, return_pairs=True, progress=None
     # Compute distance indices
     for idx, pair in iterable:
 
-        phi0= pha0[pair[0]]
+        phi0 = pha0[pair[0]]
 
         if pha1 is None:
             phi1 = pha0[pair[1]]
@@ -108,8 +108,11 @@ def pairwise_phase_consistency(pha0, pha1=None, return_pairs=True, progress=None
         # Pairwise circular distance index (PCDI)
         distance = (np.pi - 2 * abs_dist) / np.pi
 
-        if isinstance(distances, np.ndarray):
-            distances[idx] = distance
+        if isinstance(distances, np.ndarray) and pha1 is None:
+            distances[pair[0], pair[1]] = distance
+            distances[pair[1], pair[0]] = distance
+        elif isinstance(distances, np.ndarray) and pha1 is not None:
+            distances[pair[0], pair[1]] = distance
         else:
             cumulative += distance
 

neurodsp/tests/rhythm/test_phase.py

@@ -32,7 +32,7 @@ def test_pairwise_phase_consistency(tsig_sine, return_pairs, phase_shift):
         dist_avg, dists = dist_avg[0], dist_avg[1]
 
         assert isinstance(dists, np.ndarray)
-        assert len(dists) == (len(peaks) * (len(peaks) + 1)) / 2
+        assert len(dists[0]) * len(dists[1]) == len(peaks) ** 2
         assert np.mean(dists) == dist_avg
 
     # Expected consistency
@@ -50,7 +50,7 @@ def test_pairwise_phase_consistency(tsig_sine, return_pairs, phase_shift):
     dist_avg, dists = pairwise_phase_consistency(pha0[peaks], return_pairs=True)
 
     assert dist_avg == 1
-    assert len(dists) == (len(peaks) * (len(peaks) - 1)) / 2
+    assert len(dists[0]) == len(dists[1]) == len(peaks)
 
     # Cases where arrays are invalid sizes
     try: