Speeding-up dynamic time warping

4 6 8

Figure 2.8. 3D representation of cumulative distance matrixD, where height is given^c by the value of the considered matrix value. The optimal warping path is marked by black circles.

most comfortable way to travel fromD(1,^c 1m)toD(n, m).^c

2.4.3 Speeding-up dynamic time warping

Although the calculation time of DTW could be reduced by dynamic programming toO(nm), it is still computationally intensive especially in the case of long time series. In order to speed up its calculation, different approaches are used, often simultaneously.

Global constraints

The most popular method to speed up DTW is the application of a global constraint that limits the warping window — the walkable area — in matrixD. The required^c time to compute DTW decreases toO(nw), wherewis the average width of the warping window. Two of the most frequently used global constraints, the Sakoe-Chiba band [57] and the Itakura-parallelogram [58] are shown in Figure 2.9.

Ratanamahatana and Keogh [10] presented the R-K band using a heuristic search algorithm that automatically learns the constraint from the data and locally shrinks the Sakoe-Chiba band making it even narrower. Yuet al.[59] replaced this heuristic search algorithm and utilized a large margin criterion to have a better generalization ability on unseen test data. Although the R-K band usually provides better results than the Sakoe-Chiba band, it still suffers from overfitting issues [60] and in many cases the effectiveness of constrained and unconstrained DTW is the same [33].

1 18 20

1

(a) Constraint of Sakoe-Chiba band

1 18

20

1

(b) Constraint of Itakura-parallelogram

Figure 2.9. The effects of the Sakeo-Chiba band and the Itakura-parallelogram on the warping path

Moreover, Kurbalijaet al.[61] empirically proved that DTW is very sensitive to the introduction of any global constraint and even a small change in an already narrow band changes the relation between the time series. Last but not least, application of a global constraint can make DTW rigid for real-time applications where there is no chance to do proper preprocessing and/or compensate the initial/ending shifts due to time or hardware limit. On the other hand, it has to be underlined that constraining is a must when fast DTW computation is prioritized, as the most effecting lower bounding functions and indexing methods are based on constraining [62, 63].

An additional reason alongside the provided speed up for using global constraints is that they lower the chance of unnatural warpings. Such unwanted warpings can be produced where a small part of one time series maps onto a large part of the other time series. While global constraining may seem adequate for this purpose, its shortcomings mentioned above have to be considered. Thus, several other methods were proposed to have better control on the unwanted warpings. These methods are reviewed in Section 5.1 and a novel a technique is introduced that considers multiple features of the time series to address the problem of unwanted warpings.

Dimensionality reduction

Another popular method to speed up DTW calculation is to reduce the length of the time series. Such length shrinkage speeds up the computation time of unconstrained DTW quadratically and by the factor of abandoned elements in case of constrained DTW. The dimensionality reduction is usually done by resempling or by selecting one of the popular time series representations, from which four were presented in

Section 2.2.

Some representations provide another advantage, such that the focus of the dissimilarity measure can be put to other feature(s) that are not directly available. For example, using DFT, underlying processes can be compared based on their spectral properties.

Lower bounding and indexing

Several problems — such as nearest neighbor classification and similarity search — require comparison of the query time series to a — usually large — set of time series.

In such cases, even if all of the aforementioned speed up techniques are utilized, the computation time is usually still too high to compare the query directly to all items in the dataset. Such sequential scanning can work on a smaller dataset, but another solution has to be found for large databases. Obvious choice would be the usage of any off-the-shelf spatial access method like the R-trees [64]; however, the time series usually contains samples more than spatial access methods can handle and they would degrade to linear scanning eventually.

Faloutsos et al. [65] introduced the generic multimedia indexing (GEMINI) approach that provided a framework for the fast access of basically any kind of data stored in that large datasets. It is based on the assumption that if the original, n-dimensional data can be represented by a lower,N < n-dimensional feature set

— called indexes — , where isN small enough to be searched efficiently with any spatial access method, degradation to linear scan and false dismissal can be avoided when the following lower bounding lemma is satisfied:

d_N(X_N, Y_N)≤d_n(X_n, Y_n), (2.32) whereX_n,Y_n andX_N,Y_N aren andN-dimensional data,d_n andd_N denotes distance measures defined innandN-dimensional space, respectively.

Steps of the GEMINI framework were defined as follows:

1. SelectN < n-dimensional features (indexes) of the originaln-dimensional data.

2. Map the selected features into a point in theN-dimensional feature space.

3. Organize the points with any existing spatial access method.

4. Retrieve candidates from theN-dimensional space.

5. Discard the false alarms from the candidates.

For better understanding, an example is presented on how fast time series retrieval can be done with the help of GEMINI framework. Let us consider a database that contains time series with length of n and the task is to find the k nearest (most similar) item of them for query item C_n. Utilizing the GEMINI framework and discrete Fourier transform (DFT), the task can be solved trough the following steps:

1. Create truncated DFT representation of all time series stored in the database by storing their firstN Fourier coefficients only. Also keep the original n-dimensional sequences.

2. Store the truncated DFT representations in any spatial access method for fast access.

3. As soon as query seriesC_nis available, create its truncated DFT representation and retrieve thek nearest items from the indexes.

4. Retrieve the original sequences of thesek items and compute the actual n-dimensional distances (d_n) between them and C_n. Record the maximum, max(d_n).

5. Retrieve all items from the indexes that are closer than max(d_n) to C_n’s truncated DFT representation.

6. Retrieve the original sequences of the lastly retrieved items and compute the actualn-dimensional distances between them andCn. Select the time series with theksmallest distances.

Note how the last three steps take advantage of the lower bounding lemma. As the distances are smaller in the indexing space, it is ensured for every item being closer to queryC_nthanmax(d_n)in the indexing space that its real distance will be bigger. Thus, in the last step, only the false alarms have to be discarded. It is also worth to mention that the tighter lower bounding results less unnecessary calculations in this last step, i.e. the closer the distances in indexing space to the distances in the original space is, the better the performance will be.

Several indexing approach have been introduced for different time series repre-sentations based on the GEMINI framework. The first such indexing schemes were based on the fact that the Euclidean distance is preserved under transformation such

as Fourier and Haar (wavelet) transform [65, 26] and lower bounding can be ensured by truncating the coefficients.

However, many of the time series representations are not based on such trans-formations and new lower bounding distances were introduced for PAA [34], PLA [66], etc. The same logic was fallowed by Kimet al.[67] when indexing of DTW was allowed. They, however, created only a four-dimensional index space, thus it could nottake advantage of multi-dimensional index structures that scale well to higher dimensions. In addition, although four features are extracted, only one of them (determined at query time) is actually used in the lower bounding function, thus the lower bound is very loose, and many false alarms are generated, each of which will require evaluation with the quadratic-time DTW algorithm[62]. To correct these drawbacks, Keogh and Ratanamahatana [62] definedLBKeoghlower bounding distance using global constrains as envelops around the time series, and indexed the original time series data using a special PAA representation subject to these bounds. Based on their work, several other lower bounding measures were introduced [56, 68, 69]. Although these measures often provided a tighter lower than LB_Keogh, they were altogether slower due to the increased time in index calculation or in similarity search in the indexed space [33]. Lemire [63] offered a real speed up factor of 1-1.3 due to the applied tighter lower bound and the fact that his similarity search algorithm was not required to run significantly more times than forLB_Keogh [33].

Recent approaches

As multicore systems are not limited to mainframes anymore, CPUs with more than one core have became standard in almost all areas from notebooks to small embedded systems. Following this trend, several methods have been introduced to utilize all available cores for calculating DTW. Srikanthanet al.[70] reported that search of time series under DTW in their tests took3 hours and 2 minutes on a single core. The [8-core version]was able to complete the computation in 23 minutes[70].

Graphics cards can also be utilized to speed up DTW calculation, here a speed up of one [71] to two [72, 73] orders of magnitude was reached.

Another intuitive approach was presented by Rakthanmanon et al. [20] who optimized every step of sequential scan under DTW. By introducing several early abandoning steps and cascading the lower bounding measures, better speed was reached than with any existing index based approach.

As of the time of writing this thesis, the fastest way to execute search under DTW is to utilize a field-programmable gate array (FPGA). Best to the author’s knowledge, Sartet al.[74] were the first who utilized FPGAs for the calculation of DTW and a speed up magnitude order of four was reached. Wanget al.[18] raised the bar even higher making possible to search trough a whole year of electrocardiograph (ECG) data with DTW less than one minute only.

Considering the application of speed up techniques

Although it may seem desirable to apply all the aforementioned speed up techniques alongside DTW, their application is by far not as obvious as one would expect.

Global constraints have clear advantages but they will also disable DTW’s unique property, which enables it to compensate the initial/ending shifts of time series in real-time applications. Similarly, dimensionality reduction in most cases lead to loss of fine details. Application of any lower bounding measure with indexing requires storage of the indexes, on-line computation of the lower bounding representations and searching in the index space. The application will decide whether these efforts pay off or not. For example, it makes no sense to use any of these techniques when the query time series has to be classified using a small number of classes, when time series have a short length or when severe initial/ending shifts can occur.

Another important aspect to be considered for embedded systems is constant run-time. Several applications — such as automotive, aviation, aerospace, etc. — requires a precise run-time measurement of the given algorithm, and when indexing is applied, it can be complicated to always provide worst case scenarios. Instead, constrained or un-constrained DTW can provide predictable run-time.