Negative Values in Scrump Matrix Profile

[zeydabadi]

Hi,
What am I doing wrong here?

ch2_motif = stumpy.stump(T_A, m, T_B, ignore_trivial=False)
print(np.min(ch2_motif))
print(np.max(ch2_motif))

-1
194.886

Is a negative value expected in the matrix profile?

[seanlaw]

So, when you look at the documentation for stumpy.stump, it explains that output that gets returned (i.e., ch2_motif) is:

The first column consists of the matrix profile, the second column consists of the matrix profile indices, 
the third column consists of the left matrix profile indices, and the fourth column consists of the right 
matrix profile indices.

Thus, the first thing to understand is that ch2_motif contains a lot of information packed into different columns. It looks like what you are looking to do is to retrieve the max/min from the matrix profile value column:

print(np.min(ch2_motif[:, 0]))
print(np.max(ch2_motif[:, 0]))

The reason why you are seeing -1 in print(np.min(ch2_motif)) is because there are matrix profile indices (in columns 2, 3, or 4 of ch2_motif) that don't have a valid nearest neighbor. For "invalid" or unavailable nearest neighbors, the matrix profile index is automatically set to -1, which indicates to you that it is invalid.

[seanlaw]

First suggestion is to try and use triple back-ticks when writing out code blocks instead of a single back-tick. :)

Second, I'm not sure why you are doing approx_P = (approx_P - np.min(approx_P))/(np.max(approx_P)-np.min(approx_P)). This doesn't seem necessary and I don't know what affect this will have on the sorting.

Note that all subsequences are z-normalized before they are compared and so while the raw data may not look similar, they may be similar when z-normalized. If you don't want z-normalization then you may want to consider setting normalize=False.

Since you have the indices, you can verify for yourself that the subsequences at those two locations produces the same z-normalized distance as scrump does.

[zeydabadi]

Hi Sean,
Thanks for your response. I modified my code but still getting the same results as before:

def plot_motifs(x, m):
    print("\nm: ", m)
    seed = np.random.randint(100000)
    np.random.seed(seed)
    approx = stumpy.scrump(x, m, percentage=0.1, pre_scrump=True, s=None)
    approx.update()
    approx_P = approx.P_ # updated matrix profile
    
    plt.rcParams['figure.figsize'] = (24, 8)
    fig, ax = plt.subplots()
   
    ax.plot(x, label="True signal")
    sort_index = np.argsort(approx_P)

    for idx in sort_index[0:2]:
        print("idx: ", idx)
        rect = Rectangle((idx, -300), m, 600, facecolor='lightgrey')
        ax.add_patch(rect)

    plt.grid("both")
    plt.xlabel("Time steps")
    plt.ylabel("Amplitude")
    plt.legend()
    plt.show()

for m in [10000]:
    plot_motifs(bp_df["channel1"][600000:960000].values, m)

Please note that I have not normalized anything myself. So, I still have the same question that the two gray rectangles do not seem very similar. What should I do here to improve the outcome?

[zeydabadi]

I thought this might clarify the matter. For example, visually, it seems that the segments in red rectangles are more similar than what stumpy has found:

[seanlaw]

@zeydabadi There are several conflating factors here. Firstly, you are computing an approximate matrix profile and so you are not guaranteed that the minima corresponds to the global motif (we can only hope that it is "good enough"). Secondly, I recommend plotting the matrix profile. Lastly, when I said:

Note that all subsequences are z-normalized before they are compared and so while the raw data may not look similar, they may be similar when z-normalized. If you don't want z-normalization then you may want to consider setting normalize=False.

I was explaining that scrump automatically z-normalizes each subsequence by default before computing pairwise distances. Thus, looking strictly at your raw time series without subsequence-level normalization can be misleading (there is no way to visualize all of the sliding window z-normalized subsequences in one plot). Thus, you shouldn't be comparing the two raw subsequences because they may "look" pretty close in z-normalized space since, in fact, they may be close after z-normalization.

Unlike stump, I think the two smallest matrix profile values computed by scrump may not both correspond to the same motif (since not all distances are computed). So, to be safe, you'd want to do:

idx = np.argsort(approx_P)[0]
print("idx: ", idx)
rect = Rectangle((idx, -300), m, 600, facecolor='lightgrey')
ax.add_patch(rect)

nn_idx = approx.I_[idx]
print("nearest neighbor idx: ", nn_idx)
rect = Rectangle((nn_idx, -300), m, 600, facecolor='lightgrey')
ax.add_patch(rect)

Since you have the indices, you can verify for yourself that the subsequences at those two locations produces the same z-normalized distance as scrump does.

Here, you should confirm that the distance between the two subsequences is the same as the minima in your matrix profile:

import numpy as np

def z_norm(a, axis=0, threshold=1e-7):
    std = np.std(a, axis, keepdims=True)
    std[np.less(std, threshold, where=~np.isnan(std))] = 1.0

    return (a - np.mean(a, axis, keepdims=True)) / std


def distance(a, b, axis=0):
    return np.linalg.norm(a - b, axis=axis)

D = distance(
        z_norm(x[:, idx : idx + m]),
        z_norm(x[:, nn_idx : nn_idx + m]),
    )

if D == np.min(approx.P_):
    print("This is expected")
else:
    print("This is unexpected")

Either way, you'd want an overlaid plot of z_norm(x[:, idx : idx + m]) and z_norm(x[:, nn_idx : nn_idx + m]) in order to compare how similar/dissimilar they are?

You can also take the subsequence located at idx and see if there are any other subsequences that are closer than nn_idx by performing a fast similarity search.

What should I do here to improve the outcome?

Remember that scrump is an approximation that gets better and better as you call approx.update() until you reach 100% so that would be your best option. After calling approx.update(), it's possible that nothing changes in terms of the discovered motif (i.e., you got lucky and found the global minima early) but, fundamentally, you can't "find" something that isn't there. scrump updated enough times to reach percentage=1.0 is exact (i.e., is no longer an approximation) and produces the same output as stump.

[seanlaw]

For example, visually, it seems that the segments in red rectangles are more similar than what stumpy has found:

Thank you for providing the viz. Doing what I described above, can you verify that the z-normalized distances that you calculate by hand from the red subsequences are smaller than the global minima approximated by scrump? Then, how many updates are needed to get to the right answer? It might be necessary to use a GPU if you expect less approximated results.