mstump include with subspace

[JulioV]

Hi! I'm trying to force mstump to include the first column (a) in its results but instead, I get columns e and c when k =1. I managed to reduce this problem to a few columns and the following code, what am I missing?

m = 7
mps, indices = stumpy.mstump(df, m, include=[0])
motifs_idx = np.argsort(mps, axis=1)[:, :2]
k = 1
S = stumpy.subspace(df, m, motifs_idx[k][0], indices[k][motifs_idx[k][0]], k)
print(f"For k = {k}, the {k + 1}-dimensional subspace includes subsequences from {df.columns[S].values}")
for i, s in enumerate(S):
    print(np.min(mps[s]))

I am using stumpy-1.8.0

Test data

a,b,c,d,e
14.04,0.0,456.29999999999995,0.0,1467.1799999999998
15.209999999999988,0.0,519.48,0.0,1435.59
24.57,73.71,368.54999999999995,189.54,1119.6899999999998
3.51,31.589999999999996,277.28999999999996,143.91,912.5999999999999
19.305,101.78999999999999,365.03999999999996,315.9,1804.1399999999999
5.265,45.629999999999995,238.67999999999998,224.64,814.3199999999999
8.189999999999987,28.08,491.4,143.91,1597.05
10.53,14.04,382.59,35.099999999999994,1014.39
3.51,56.16,252.71999999999997,105.3,540.54
10.53,129.87,298.34999999999997,379.08,1288.1699999999998
8.774999999999999,0.0,266.76,0.0,565.11
19.305,45.629999999999995,358.02,115.83,1147.77
3.51,3.51,259.74,38.61,621.27
3.51,17.549999999999997,210.6,45.629999999999995,435.23999999999995
3.51,21.06,214.10999999999999,42.12,410.66999999999996
10.53,0.0,428.21999999999997,0.0,1060.02
3.51,0.0,326.43,0.0,782.7299999999999
12.870000000000012,0.0,442.26,0.0,972.27
7.02,0.0,203.57999999999998,0.0,435.23999999999995
3.51,0.0,203.57999999999998,0.0,431.72999999999996
10.53,24.57,343.97999999999996,140.39999999999998,902.0699999999999
7.02,17.549999999999997,351.0,35.099999999999994,902.0699999999999
10.53,28.08,214.10999999999999,63.17999999999999,505.43999999999994
12.285,28.08,245.7,129.87,610.74
12.285,24.57,207.08999999999997,410.66999999999996,796.77
9.360000000000012,0.0,147.42,0.0,336.96
7.02,0.0,291.33,0.0,673.92
7.02,0.0,143.91,0.0,312.39
5.265,0.0,200.07,0.0,372.06
24.57,0.0,431.72999999999996,0.0,1172.34
14.04,17.549999999999997,403.65,24.57,1081.08
21.06,0.0,403.65,0.0,1031.9399999999998
10.53,21.06,91.25999999999999,200.07,368.54999999999995
10.53,0.0,28.08,0.0,35.099999999999994
5.265,10.53,224.64,28.08,547.56
11.699999999999987,0.0,231.66,0.0,375.57
8.189999999999987,0.0,294.84,0.0,589.68
12.285,21.06,340.46999999999997,59.669999999999995,1175.85
5.265,0.0,140.39999999999998,0.0,228.14999999999998
8.189999999999987,66.69,329.94,196.56,1010.8799999999999
3.51,7.02,301.85999999999996,31.589999999999996,958.2299999999999
4.679999999999988,70.19999999999999,280.79999999999995,256.22999999999996,870.4799999999999
5.850000000000011,0.0,249.20999999999998,0.0,456.29999999999995
10.53,7.02,214.10999999999999,42.12,551.0699999999999
7.02,14.04,161.45999999999998,35.099999999999994,351.0
5.850000000000011,140.39999999999998,196.56,456.29999999999995,1074.06
14.04,101.78999999999999,179.01,477.35999999999996,1158.3
3.51,0.0,238.67999999999998,0.0,568.62
7.02,0.0,263.25,0.0,572.13
15.794999999999998,84.24,221.13,333.45,1140.75
4.679999999999988,0.0,242.19,0.0,473.84999999999997
9.360000000000012,91.25999999999999,256.22999999999996,333.45,1432.08
9.360000000000012,52.65,214.10999999999999,252.71999999999997,1319.76

[seanlaw]

@JulioV Great question and I think there is an opportunity to make this more explicit in the Tutorial. It looks like:

mps, indices = stumpy.mstump(df, m, include=[0])

Is correct. However, if you look at the documentation for subspace, you also need to provide an include parameter as well:

S = stumpy.subspace(df, m, motifs_idx[k][0], indices[k][motifs_idx[k][0]], k, include=[0])

I believe that this is what you are looking for and I've updated the tutorial to reflect this. Will you report back and let us know if this solves your problem?

[JulioV]

Thank you, it works as expected with the include parameter. I have another, perhaps naive, question. Why do all the mstump examples always identify two motifs? Is it possible to extract more than two?

[seanlaw]

Good question!

The first thing to note is that the global motif pair are complementary subsequences to each other but that need not be the case after the global pair. So, technically speaking, mstump doesn't actually identify ANY motifs and, instead, it computes a multi-dimensional matrix profile. Then, you scan this matrix profile along each dimension and look for the index locations of the minimum value(s):

motifs_idx = np.argsort(mps, axis=1)[:, :2]

It is this sorting that reveals the global motif pair (we extracted the top 2 by indexing with [:, :2]). In this case, we are finding subsequences that have the absolute closest match (and it is possible that it is a poor match as well!). Imagine if we had a multi-dimensional time series generated from random sampling, then that means that there is certainly no motifs. However, if you perform the np.argsort(mps, axis=1) then you'll "get an answer" but you should really visualize or compare this top subsequence pair with other potential subsequence pairs. Usually, the top pair is the "most" interesting to look at. Afterward, you may look at np.argsort(mps, axis=1)[:, 3] or np.argsort(mps, axis=1)[:, 4] or np.argsort(mps, axis=1)[:, 5]. Or, more generally, you may "discover" the top-K potential motifs.

[JulioV]

Ah I see, so if I do motifs_idx = np.argsort(mps, axis=1)[:, :3] with k=1. Then I can plot those three motifs between two time series like this:

f,ax = plt.subplots(k+1, 1, figsize=(20,20),sharex=True)

for i, s in enumerate(S):
    tmp = df[df.columns[s]]
    ax[i].plot(tmp)
    ax[i].set_title("Variable:{}".format(df.columns[s]))
    ax[i].plot(tmp[motifs_idx[k, 0] : motifs_idx[k, 0] + m], c='orange',linewidth=4)
    ax[i].plot(tmp[motifs_idx[k, 1] : motifs_idx[k, 1] + m], c='red',linewidth=4)
    ax[i].plot(tmp[motifs_idx[k, 2] : motifs_idx[k, 2] + m], c='blue',linewidth=4)

    ax[i].axvline(x=motifs_idx[k, 0] ,linestyle='dashed',c='k',linewidth=2)
    ax[i].axvline(x=motifs_idx[k, 1] ,linestyle='dashed',c='k',linewidth=2)
    ax[i].axvline(x=motifs_idx[k, 2] ,linestyle='dashed',c='k',linewidth=2)

[seanlaw]

Yes, but note that these three motifs may not have anything to do with each other. In fact, they are likely to be completely independent as you increase the 3 to some higher number.

[JulioV]

Oh ok. Thanks a lot for all your responses!