diff --git a/stumpy/motifs.py b/stumpy/motifs.py index 0c8246dcf..4585dea15 100644 --- a/stumpy/motifs.py +++ b/stumpy/motifs.py @@ -121,7 +121,7 @@ def _motifs( T, M_T=M_T, Σ_T=Σ_T, - max_matches=None, + max_matches=max_matches, max_distance=max_distance, atol=atol, query_idx=candidate_idx, @@ -492,7 +492,12 @@ def match( D = np.empty((d, n - m + 1)) for i in range(d): D[i, :] = core.mass( - Q[i], T[i], M_T[i], Σ_T[i], T_subseq_isconstant=T_subseq_isconstant[i] + Q[i], + T[i], + M_T[i], + Σ_T[i], + T_subseq_isconstant=T_subseq_isconstant[i], + query_idx=query_idx, ) D = np.mean(D, axis=0) diff --git a/tests/test_motifs.py b/tests/test_motifs.py index ca3766d1b..6913c51f2 100644 --- a/tests/test_motifs.py +++ b/tests/test_motifs.py @@ -6,6 +6,68 @@ from stumpy import core, match, motifs +def naive_motifs(T, m, max_motifs, max_matches): + # To avoid complexity, this naive function is written + # such that each array in the ouput has shape + # (max_motif, max_matches). + + # To this end, the following items are considered: + # 1. `max_distance` and `cutoff` are both hardcoded and + # set to np.inf + # 2. If the number of subsequence, i.e. `len(T)-m+1`, is + # not less than `m * max_motifs * max_matches`, then the + # output definitely has the shape (max_motif, max_matches). + + l = len(T) - m + 1 + excl_zone = int(np.ceil(m / 4)) + + output_shape = (max_motifs, max_matches) + motif_distances = np.full(output_shape, np.NINF, dtype=np.float64) + motif_indices = np.full(output_shape, -1, dtype=np.int64) + + D = naive.distance_matrix(T, T, m) + for i in range(D.shape[0]): + naive.apply_exclusion_zone(D[i], i, excl_zone, np.inf) + + P = np.min(D, axis=1) + for i in range(max_motifs): + distances = [] + indices = [] + + idx = np.argmin(P) + + # self match + distances.append(0) + indices.append(idx) + naive.apply_exclusion_zone(P, idx, excl_zone, np.inf) + + # Explore distance profile D[idx] till `max_matches` are found. + naive.apply_exclusion_zone(D[idx], idx, excl_zone, np.inf) + for _ in range(l): + if len(distances) >= max_matches: + break + + nn = np.argmin(D[idx]) + distances.append(D[idx, nn]) + indices.append(nn) + + # Update D[idx] to avoid finding matches that are trivial to + # each other. + naive.apply_exclusion_zone(D[idx], nn, excl_zone, np.inf) + + # Update P after the discovery of each match so that the + # match cannot be selected as the motif next time. + naive.apply_exclusion_zone(P, nn, excl_zone, np.inf) + + # Note that a discovered match cannot be selected as motif but + # it can still be selected again as a match for another motif. + + motif_distances[i] = distances + motif_indices[i] = indices + + return motif_distances, motif_indices + + def naive_match(Q, T, excl_zone, max_distance, max_matches=None): m = Q.shape[0] D = naive.distance_profile(Q, T, m) @@ -45,7 +107,7 @@ def test_motifs_one_motif(): ) npt.assert_array_equal(left_indices, right_indices) - npt.assert_almost_equal(left_profile_values, right_distance_values, decimal=4) + npt.assert_almost_equal(left_profile_values, right_distance_values) def test_motifs_two_motifs(): @@ -102,7 +164,7 @@ def test_motifs_two_motifs(): # We ignore indices because of sorting ambiguities for equal distances. # As long as the distances are correct, the indices will be too. - npt.assert_almost_equal(left_profile_values, right_distance_values, decimal=6) + npt.assert_almost_equal(left_profile_values, right_distance_values) # Reset seed np.random.seed(None) @@ -112,8 +174,9 @@ def test_motifs_max_matches(): # This test covers the following: # A time series contains motif A at four locations and motif B at two. - # If `max_motifs=2` the result should contain only the top two matches of motif A - # and the top two matches of motif B as two separate motifs. + # If `max_moitf=2` and `max_matches=3`, the result should contain + # (at most) two sets of motifs and each motif set should contain + # (at most) the top three matches T = np.array( [ 0.0, # motif A @@ -137,11 +200,12 @@ def test_motifs_max_matches(): 2.3, 2.0, # motif A 3.0, - 2.02, + 3.0, ] ) m = 3 - max_motifs = 3 + max_motifs = 2 + max_matches = 3 left_indices = [[0, 7], [4, 11]] left_profile_values = [ @@ -160,14 +224,79 @@ def test_motifs_max_matches(): T, mp[:, 0], max_motifs=max_motifs, - max_distance=0.1, + max_matches=max_matches, + max_distance=0.05, cutoff=np.inf, - max_matches=2, ) # We ignore indices because of sorting ambiguities for equal distances. # As long as the distances are correct, the indices will be too. - npt.assert_almost_equal(left_profile_values, right_distance_values, decimal=4) + npt.assert_almost_equal(left_profile_values, right_distance_values) + + +def test_motifs_max_matches_max_distances_inf(): + # This test covers the following: + + # A time series contains motif A at two locations and motif B at two. + # If `max_moitf=2` and `max_matches=2`, the result should contain + # (at most) two sets of motifs and each motif set should contain + # (at most) two matches. + T = np.array( + [ + 0.0, # motif A + 1.0, + 0.0, + 2.3, + -1.0, # motif B + -1.0, + -2.0, + 0.0, # motif A + 1.0, + 0.0, + -2.0, + -1.0, # motif B + -1.03, + -2.0, + -0.5, + 2.0, + 3.0, + 2.04, + 2.3, + 2.0, + 3.0, + 3.0, + ] + ) + m = 3 + max_motifs = 2 + max_matches = 2 + max_distance = np.inf + + left_indices = [[0, 7], [4, 11]] + left_profile_values = [ + [0.0, 0.0], + [ + 0.0, + naive.distance( + core.z_norm(T[left_indices[1][0] : left_indices[1][0] + m]), + core.z_norm(T[left_indices[1][1] : left_indices[1][1] + m]), + ), + ], + ] + + # set `row_wise` to True so that we can compare the indices of motifs as well + mp = naive.stump(T, m, row_wise=True) + right_distance_values, right_indices = motifs( + T, + mp[:, 0], + max_motifs=max_motifs, + max_distance=max_distance, + cutoff=np.inf, + max_matches=max_matches, + ) + + npt.assert_almost_equal(left_indices, right_indices) + npt.assert_almost_equal(left_profile_values, right_distance_values) def test_naive_match_exclusion_zone(): @@ -301,3 +430,32 @@ def test_match_mean_stddev_isconstant(Q, T): ) npt.assert_almost_equal(left, right) + + +def test_motifs(): + T = np.random.rand(64) + m = 3 + + max_motifs = 3 + max_matches = 4 + max_distance = np.inf + cutoff = np.inf + + # naive + # `max_distance` and `cutoff` are hard-coded, and set to np.inf. + ref_distances, ref_indices = naive_motifs(T, m, max_motifs, max_matches) + + # performant + mp = naive.stump(T, m, row_wise=True) + comp_distance, comp_indices = motifs( + T, + mp[:, 0].astype(np.float64), + min_neighbors=1, + max_distance=max_distance, + cutoff=cutoff, + max_matches=max_matches, + max_motifs=max_motifs, + ) + + npt.assert_almost_equal(ref_indices, comp_indices) + npt.assert_almost_equal(ref_distances, comp_distance)