Fixed #1077 Add extra check for window size

stumpy/aamp.py

@@ -407,17 +407,17 @@ def aamp(T_A, m, T_B=None, ignore_trivial=True, p=2.0, k=1):
     if T_B.ndim != 1:  # pragma: no cover
         raise ValueError(f"T_B is {T_B.ndim}-dimensional and must be 1-dimensional. ")
 
-    core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
-    ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
-
     n_A = T_A.shape[0]
     n_B = T_B.shape[0]
     l = n_A - m + 1
 
+    ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
     excl_zone = int(np.ceil(m / config.STUMPY_EXCL_ZONE_DENOM))
-    if ignore_trivial:
+    if ignore_trivial:  # self-join
+        core.check_window_size(m, max_size=min(n_A, n_B), n=n_A)
         diags = np.arange(excl_zone + 1, n_A - m + 1, dtype=np.int64)
-    else:
+    else:  # AB-join
+        core.check_window_size(m, max_size=min(n_A, n_B))
         diags = np.arange(-(n_A - m + 1) + 1, n_B - m + 1, dtype=np.int64)
 
     P, PL, PR, I, IL, IR = _aamp(

stumpy/aamped.py

@@ -386,17 +386,17 @@ def aamped(client, T_A, m, T_B=None, ignore_trivial=True, p=2.0, k=1):
     if T_B.ndim != 1:  # pragma: no cover
         raise ValueError(f"T_B is {T_B.ndim}-dimensional and must be 1-dimensional. ")
 
-    core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
-    ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
-
     n_A = T_A.shape[0]
     n_B = T_B.shape[0]
 
+    ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
     excl_zone = int(np.ceil(m / config.STUMPY_EXCL_ZONE_DENOM))
 
     if ignore_trivial:
+        core.check_window_size(m, max_size=min(n_A, n_B), n=n_A)
         diags = np.arange(excl_zone + 1, n_A - m + 1, dtype=np.int64)
     else:
+        core.check_window_size(m, max_size=min(n_A, n_B))
         diags = np.arange(-(n_A - m + 1) + 1, n_B - m + 1, dtype=np.int64)
 
     _aamped = core._client_to_func(client)

stumpy/core.py

@@ -554,11 +554,12 @@ def get_max_window_size(n):
     return max_m
 
 
-def check_window_size(m, max_size=None):
+def check_window_size(m, max_size=None, n=None):
     """
     Check the window size and ensure that it is greater than or equal to 3 and, if
     `max_size` is provided, ensure that the window size is less than or equal to the
-    `max_size`
+    `max_size`. Furthermore, if `n` is provided, then it checks if there is
+    at least one subsequence with no non-trivial neighbour in a self-join.
 
     Parameters
     ----------
@@ -568,6 +569,10 @@ def check_window_size(m, max_size=None):
     max_size : int, default None
         The maximum window size allowed
 
+    n : int, default None
+        The length of the time series in the case of a self-join.
+        `n` should be set to `None` in the case of an `AB-join`.
+
     Returns
     -------
     None
@@ -589,6 +594,41 @@ def check_window_size(m, max_size=None):
     if max_size is not None and m > max_size:
         raise ValueError(f"The window size must be less than or equal to {max_size}")
 
+    if n is not None:
+        # Raise warning if there is at least one subsequence with no
+        # non-trivial neighbour in the case of a self-join.
+
+        # For any time series `T`, an "eligible nearest neighbor" subsequence for
+        # the central-most subsequence must be located outside the `excl_zone`,
+        # and the central-most subsequence will ALWAYS have the smallest index-wise
+        # distance to its furthest neighbour amongs all other subsequences. Therefore,
+        # we only need to check whether the `excl_zone` eliminates all "neighbors" for
+        # the central-most subsequence in `T`. In fact, we just need to verify whether
+        # the `excl_zone` eliminates the "neighbor" that is furthest away (index-wise)
+        # from the central-most subsequence. If it does not, this implies that all
+        # subsequences in `T` will have at least one "eligible nearest neighbors"
+        # outside their respective `excl_zone`
+
+        excl_zone = int(math.ceil(m / config.STUMPY_EXCL_ZONE_DENOM))
+
+        l = n - m + 1
+
+        # If `l` is odd (`l == 2k+1`):
+        # The central subsequence is at index `k`, with furthest neighbors at `0`
+        # and `2k`, both `k == l // 2` indices away from the central-most subsequence.
+
+        # If `l` is even (`l == 2k`):
+        # The central subsequences are at `k-1` and `k`. The furthest neighbor is
+        # at `2k-1` for `k-1`, and `0` for `k`. In both cases, the subsequence
+        # and its furthest neighbor are `k == l // 2` indices away.
+        max_gap = l // 2
+        if max_gap <= excl_zone:
+            msg = (
+                f"The window size, 'm = {m}', may be too large and could lead to "
+                + "meaningless results. Consider reducing 'm' where necessary"
+            )
+            warnings.warn(msg)
+
 
 @njit(fastmath=config.STUMPY_FASTMATH_TRUE)
 def _sliding_dot_product(Q, T):

stumpy/gpu_aamp.py

@@ -536,8 +536,13 @@ def gpu_aamp(T_A, m, T_B=None, ignore_trivial=True, device_id=0, p=2.0, k=1):
             "For multidimensional STUMP use `stumpy.mstump` or `stumpy.mstumped`"
         )
 
-    core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
     ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
+    if ignore_trivial:  # self-join
+        core.check_window_size(
+            m, max_size=min(T_A.shape[0], T_B.shape[0]), n=T_A.shape[0]
+        )
+    else:  # AB-join
+        core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
 
     n = T_B.shape[0]
     w = T_A.shape[0] - m + 1

stumpy/gpu_stump.py

@@ -666,8 +666,13 @@ def gpu_stump(
             "For multidimensional STUMP use `stumpy.mstump` or `stumpy.mstumped`"
         )
 
-    core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
     ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
+    if ignore_trivial:  # self-join
+        core.check_window_size(
+            m, max_size=min(T_A.shape[0], T_B.shape[0]), n=T_A.shape[0]
+        )
+    else:  # AB-join
+        core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
 
     n = T_B.shape[0]
     w = T_A.shape[0] - m + 1

stumpy/maamp.py

@@ -140,7 +140,7 @@ def maamp_subspace(
         returned.
     """
     T = core._preprocess(T)
-    core.check_window_size(m, max_size=T.shape[-1])
+    core.check_window_size(m, max_size=T.shape[-1], n=T.shape[-1])
 
     subseqs, _ = core.preprocess_non_normalized(T[:, subseq_idx : subseq_idx + m], m)
     neighbors, _ = core.preprocess_non_normalized(T[:, nn_idx : nn_idx + m], m)
@@ -269,7 +269,7 @@ def maamp_mdl(
         A list of numpy.ndarrays that contains the `k`th-dimensional subspaces
     """
     T = core._preprocess(T)
-    core.check_window_size(m, max_size=T.shape[-1])
+    core.check_window_size(m, max_size=T.shape[-1], n=T.shape[-1])
 
     if discretize_func is None:
         T_isfinite = np.isfinite(T)
@@ -441,7 +441,7 @@ def maamp_multi_distance_profile(query_idx, T, m, include=None, discords=False,
         err = f"T is {T.ndim}-dimensional and must be at least 1-dimensional"
         raise ValueError(f"{err}")
 
-    core.check_window_size(m, max_size=T.shape[1])
+    core.check_window_size(m, max_size=T.shape[1], n=T.shape[1])
 
     if include is not None:  # pragma: no cover
         include = core._preprocess_include(include)
@@ -933,7 +933,7 @@ def maamp(T, m, include=None, discords=False, p=2.0):
         err = f"T is {T_A.ndim}-dimensional and must be at least 1-dimensional"
         raise ValueError(f"{err}")
 
-    core.check_window_size(m, max_size=min(T_A.shape[1], T_B.shape[1]))
+    core.check_window_size(m, max_size=min(T_A.shape[1], T_B.shape[1]), n=T_A.shape[1])
 
     if include is not None:
         include = core._preprocess_include(include)

stumpy/maamped.py

@@ -389,7 +389,7 @@ def maamped(client, T, m, include=None, discords=False, p=2.0):
         err = f"T is {T_A.ndim}-dimensional and must be at least 1-dimensional"
         raise ValueError(f"{err}")
 
-    core.check_window_size(m, max_size=min(T_A.shape[1], T_B.shape[1]))
+    core.check_window_size(m, max_size=min(T_A.shape[1], T_B.shape[1]), n=T_A.shape[1])
 
     if include is not None:
         include = core._preprocess_include(include)

stumpy/mstump.py

@@ -217,7 +217,7 @@ def subspace(
     array([0, 1])
     """
     T = core._preprocess(T)
-    core.check_window_size(m, max_size=T.shape[-1])
+    core.check_window_size(m, max_size=T.shape[-1], n=T.shape[-1])
     T_subseq_isconstant = core.process_isconstant(T, m, T_subseq_isconstant)
 
     if discretize_func is None:
@@ -409,7 +409,7 @@ def mdl(
     (array([ 80.      , 111.509775]), [array([1]), array([0, 1])])
     """
     T = core._preprocess(T)
-    core.check_window_size(m, max_size=T.shape[-1])
+    core.check_window_size(m, max_size=T.shape[-1], n=T.shape[-1])
     T_subseq_isconstant = core.process_isconstant(T, m, T_subseq_isconstant)
 
     if discretize_func is None:
@@ -1228,7 +1228,9 @@ def mstump(
         err = f"T is {T_A.ndim}-dimensional and must be at least 1-dimensional"
         raise ValueError(f"{err}")
 
-    core.check_window_size(m, max_size=min(T_A.shape[1], T_B.shape[1]))
+    # mstump currently supports self-join. Therefore, the argument `n` should be
+    # passed to the function `core.check_window_size`.
+    core.check_window_size(m, max_size=min(T_A.shape[1], T_B.shape[1]), n=T_A.shape[1])
 
     if include is not None:
         include = core._preprocess_include(include)

stumpy/mstumped.py

@@ -505,7 +505,9 @@ def mstumped(
         err = f"T is {T_A.ndim}-dimensional and must be at least 1-dimensional"
         raise ValueError(f"{err}")
 
-    core.check_window_size(m, max_size=min(T_A.shape[1], T_B.shape[1]))
+    # mstump currently supports self-join. Therefore, the argument `n` should be
+    # passed to the function `core.check_window_size`.
+    core.check_window_size(m, max_size=min(T_A.shape[1], T_B.shape[1]), n=T_A.shape[1])
 
     if include is not None:
         include = core._preprocess_include(include)

stumpy/scraamp.py

@@ -646,10 +646,15 @@ def __init__(
                 "For multidimensional STUMP use `stumpy.mstump` or `stumpy.mstumped`"
             )
 
-        core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
         self._ignore_trivial = core.check_ignore_trivial(
             self._T_A, self._T_B, self._ignore_trivial
         )
+        if self._ignore_trivial:  # self-join
+            core.check_window_size(
+                m, max_size=min(T_A.shape[0], T_B.shape[0]), n=T_A.shape[0]
+            )
+        else:  # AB-join
+            core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
 
         self._n_A = self._T_A.shape[0]
         self._n_B = self._T_B.shape[0]

stumpy/scrump.py

@@ -905,10 +905,15 @@ def __init__(
                 "For multidimensional STUMP use `stumpy.mstump` or `stumpy.mstumped`"
             )
 
-        core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
         self._ignore_trivial = core.check_ignore_trivial(
             self._T_A, self._T_B, self._ignore_trivial
         )
+        if self._ignore_trivial:
+            core.check_window_size(
+                m, max_size=min(T_A.shape[0], T_B.shape[0]), n=T_A.shape[0]
+            )
+        else:
+            core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
 
         self._n_A = self._T_A.shape[0]
         self._n_B = self._T_B.shape[0]

stumpy/stamp.py

@@ -208,13 +208,14 @@ def stamp(
     if T_B.ndim != 1:  # pragma: no cover
         raise ValueError(f"T_B is {T_B.ndim}-dimensional and must be 1-dimensional. ")
 
-    core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
-
     subseq_T_A = core.rolling_window(T_A, m)
     excl_zone = int(np.ceil(m / 2))
 
     # Add exclusionary zone
     if ignore_trivial:
+        core.check_window_size(
+            m, max_size=min(T_A.shape[0], T_B.shape[0]), n=T_A.shape[0]
+        )
         out = [
             _mass_PI(
                 subseq,
@@ -229,6 +230,7 @@ def stamp(
             for i, subseq in enumerate(subseq_T_A)
         ]
     else:
+        core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
         out = [
             _mass_PI(
                 subseq,

stumpy/stomp.py

@@ -81,8 +81,13 @@ def _stomp(T_A, m, T_B=None, ignore_trivial=True):
     if T_B.ndim != 1:  # pragma: no cover
         raise ValueError(f"T_B is {T_B.ndim}-dimensional and must be 1-dimensional. ")
 
-    core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
     ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
+    if ignore_trivial:  # self-join
+        core.check_window_size(
+            m, max_size=min(T_A.shape[0], T_B.shape[0]), n=T_A.shape[0]
+        )
+    else:  # AB-join
+        core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
 
     n = T_A.shape[0]
     l = n - m + 1

stumpy/stump.py

@@ -711,18 +711,17 @@ def stump(
             "For multidimensional STUMP use `stumpy.mstump` or `stumpy.mstumped`"
         )
 
-    core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
-    ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
-
     n_A = T_A.shape[0]
     n_B = T_B.shape[0]
     l = n_A - m + 1
 
+    ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
     excl_zone = int(np.ceil(m / config.STUMPY_EXCL_ZONE_DENOM))
-
-    if ignore_trivial:
+    if ignore_trivial:  # self-join
+        core.check_window_size(m, max_size=min(n_A, n_B), n=n_A)
         diags = np.arange(excl_zone + 1, n_A - m + 1, dtype=np.int64)
-    else:
+    else:  # AB-join
+        core.check_window_size(m, max_size=min(n_A, n_B))
         diags = np.arange(-(n_A - m + 1) + 1, n_B - m + 1, dtype=np.int64)
 
     P, PL, PR, I, IL, IR = _stump(

stumpy/stumped.py

@@ -618,17 +618,17 @@ def stumped(
             "For multidimensional STUMP use `stumpy.mstump` or `stumpy.mstumped`"
         )
 
-    core.check_window_size(m, max_size=min(T_A.shape[0], T_B.shape[0]))
-    ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
-
     n_A = T_A.shape[0]
     n_B = T_B.shape[0]
 
+    ignore_trivial = core.check_ignore_trivial(T_A, T_B, ignore_trivial)
     excl_zone = int(np.ceil(m / config.STUMPY_EXCL_ZONE_DENOM))
 
     if ignore_trivial:
+        core.check_window_size(m, max_size=min(n_A, n_B), n=n_A)
         diags = np.arange(excl_zone + 1, n_A - m + 1, dtype=np.int64)
     else:
+        core.check_window_size(m, max_size=min(n_A, n_B))
         diags = np.arange(-(n_A - m + 1) + 1, n_B - m + 1, dtype=np.int64)
 
     _stumped = core._client_to_func(client)

tests/test_core.py

@@ -192,6 +192,21 @@ def test_check_max_window_size():
             core.check_window_size(m, max_size=3)
 
 
+def test_check_window_size_excl_zone():
+    # To ensure warning is raised if there is at least one subsequence
+    # that has no non-trivial neighbor
+    T = np.random.rand(10)
+    m = 7
+
+    # For `len(T) == 10` and `m == 7`, the `excl_zone` is ceil(m / 4) = 2.
+    # In this case, there are `10 - 7 + 1 = 4` subsequences of length 7,
+    # starting at indices 0, 1, 2, and 3. For a subsequence that starts at
+    # index 1, there are no non-trivial neighbors. So, a warning should be
+    # raised.
+    with pytest.warns(UserWarning):
+        core.check_window_size(m, max_size=len(T), n=len(T))
+
+
 @pytest.mark.parametrize("Q, T", test_data)
 def test_njit_sliding_dot_product(Q, T):
     ref_mp = naive_rolling_window_dot_product(Q, T)