Hierarchical Clustering using Auto-encoded Compact Representation for Time-series Analysis

In our method we first learn a compact representation of the time-series by using an Seq2Seq LSTM Hochreiter and Schmidhuber (1997) auto-encoder. Let $X=\{X_{1},X_{2},....,X_{M}\}$ be a set of $M$ time-series, where $X_{i}=\{x_{i,1},x_{i,2},..,x_{i,k},..x_{i,n}\}$, $x_{i,k}$ gives the value of $i^{th}$ time-series’s $k^{th}$ timestep, and n is the total number of timesteps or length of the time-series. In case of univariate, $x_{i}$ is scalar and is vector in case of multivariate time-series. We learn compact representation of time-series using multi layer Seq2Seq LSTM auto-encoder. The encoder block of the auto-encoder consists of an input layer $h_{0}$ having same number of units as number of timesteps and two hidden layers $h_{1}$ and $h_{2}$ having lengths $h_{l1}$ and $h_{l2}$ respectively. We propose Auto-Encoded Compact Sequence $AECS\in R^{M\times h_{l2}}$, which is the latent representation of layer $h_{2}$ of encoder. The length of AECS is $(h_{l2})$ where $h_{l2}<h_{l1}<n$. A schematic illustration of the auto-encoder model is presented in Figure 2.

We aim to learn a compact representation having a length much less than the original time-series to reduce the computing time as well as to capture important characteristics of time-series. AECS is an undercomplete representation hence bound to capture the important features of the input time-series.
Encoder:
$H_{i}=f(X_{i})\leftarrow h_{i,k},s_{i,k}\leftarrow f_{LSTM}(h_{i,k-1},s_{i,k-1}x_{i,k})$; Decoder:
$X^{\prime}_{i}=g(H_{i})\leftarrow h^{\prime}_{i,k},s^{\prime}_{i,k}\leftarrow f_{LSTM}(h^{\prime}_{i,k-1},s^{\prime}_{i,k-1},h_{i,k})$ ….. (1), where $X\in R^{M\times n},H\in R^{M\times l_{ev}}$ and $X^{\prime}\in R^{M\times n}$; $l_{ev}$: length of the encoder vector.
Here $l_{ev}=length(AECS)=h_{l2}$. The reconstruction loss is Mean Square Error(mse): $l_{recon}=\dfrac{1}{Mn}\sum\limits_{i=1}^{M}\sum\limits_{j=1}^{n}(x_{i,j}-x^{\prime}_{i,j})^{2}$ , where $X_{i}^{\prime}=\{x^{\prime}_{i,1},x^{\prime}_{i,2},.....,x^{\prime}_{i,n}\}$ is the reconstructed output for $i^{th}$ instance.

Hierarchical clustering (HC) produces a hierarchy of clusters formed either by merging or dividing existing clusters recursively. The iterative process of agglomeration/division of the clusters can be visualised in form of a tree like structure known as dendrogram. Each node of the tree/dendrogram contains one or more time-series ($X_{i}$) representing a cluster. Here we use agglomerative hierarchical clustering which takes each input time-series as an individual cluster and then start successively merging pair of clusters having highest similarity between them until a single cluster is formed. Hence it measures similarities among the pair of time-series to put them into same cluster and the dissimilarities between the pair of groups of time-series to take the decision of further fusing the clusters. We have used average linkage Friedman et al. (2001) to obtain pairwise inter-cluster dissimilarities. In average linkage, distance $a_{d}$ between clusters $C_{0}$ and $C_{1}$ is : $a_{d}(C_{0},C_{1})=\dfrac{1}{|C_{0}|\cdot|C_{1}|}\sum\limits_{X_{i}\in C_{0}}\sum\limits_{X_{j}\in C_{1}}dist(X_{i},X_{j})$. ….. (2)
$|C_{i}|$ denotes cardinality or number of members in cluster i. Average linkage aims to form compact and relatively far apart clusters than other linkage methods like single or complete by keeping average pairwise dissimilarity. Usage of linkage methods like single or complete are scope of future study.

In this work we have used following three different distance measures to find the dissimilarity between two time-series (say $X_{i}$ and $X_{j}$, where $X_{i}=\{x_{i,1},x_{i,2},..,x_{i,k},..x_{i,n}\}$ and $X_{j}=\{x_{j,1},x_{j,2},..,x_{j,k},..x_{j,n}\}$, $x_{i,k}$ and $x_{j,k}$ indicates $k^{th}$ timestep of $i^{th}$ and $j^{th}$ time-series respectively) as well as to obtain average linkage among the clusters as follows.
Chebyshev Distance: $max_{k}(|x_{i,k}-x_{j,k}|)$ …..(3): The maximum distance between two data points in any single dimension.
Manhattan Distance: $\sum\limits_{k=1}^{T}{|x_{i,k}-x_{j,k}|}$…..(4): This belongs to Minkowski family which computes distance travelled to get from one data point to the other if a grid-like path is followed.
Mahalanobis Distance: $\sqrt{(X_{i}-X_{j})^{T}\cdot C^{-1}\cdot{(X_{i}-X_{j})}}$…(5): Finds the distance between two data points in multidimensional space. Here $X_{i}$ and $X_{j}$ are two time-series and $C$ is the co-variance matrix between $X_{i}$ and $X_{j}$.

We use Modified Hubert statistic ($\mathcal{T}$), an internal clustering measure which evaluates sum of distance between each pair of instances weighted by distance between the centers of clusters which they belong to. We rank clustering formed by different distance measures based on $\mathcal{T}$. Best clustering has highest $\mathcal{T}$. In case two or more distance measures produces the highest value of $\mathcal{T}$, we report all of them as the best measures. We have measured disagreements between pairs of time-series and separation between clusters using ML (Mahalanobis) to evaluate $\mathcal{T}$, which is one of our prime contributions.
$\mathcal{T}=\frac{2}{n(n-1)}\sum\limits_{X_{i}\in X}\sum\limits_{X_{j}\in X}d(X_{i},X_{j})d(c_{i},c_{j})$ …..(6)
$d(X_{i},X_{j})=d_{ML}(X_{i},X_{j})$ ; $d(c_{i},c_{j})=d_{ML}(c_{i},c_{j})$…..(7) where $C_{i}$ represents $i^{th}$ cluster and $c_{i}$ is the center of cluster $C_{i}$, $d_{ML}(X_{i},X_{j})$ is Mahalanobis distance between time-series $X_{i}$ and $X_{j}$ and $d_{ML}(c_{i},c_{j})$ is Mahalanobis distance between centers of clusters to which $X_{i}$ and $X_{j}$ belongs.

We compute Rand Index (RI) Rand (1971), one of the external clustering validation measures to compare the performance of our proposed mechanism with the benchmark results. RI is a cluster validation measure which considers all pairs of time-series and counts number of pairs assigned in same or different clusters formed by our algorithm to the true clusters based on the given labels. $RI:\frac{TP+TN}{TP+FP+FN+TN}$ …..(8) , where the symbols denote cardinalities of sets of pairs: TP (true positive) denotes elements which are assigned to same cluster that belongs to the same class; TN (true negative) denotes elements which are assigned to different clusters that belongs to different classes; FP (false positive) denotes elements which are assigned to different clusters that belongs to the same class; FN (false negative) denotes elements which are assigned to same cluster that belongs to different classes.

We use univariate and multivariate time-series from UCR Time Series Classification Archive Bagnall et al. (2018b) to evaluate the performance of our proposed algorithm. Additionally we perform our analysis on 2 multivariate datasets from UCI Machine Learning Repository Dua and Graff (2017). Description of univariate and multivariate time-series are presented in Tables 1 and 2 respectively. Each dataset has default train-test split. In our approach we merge the train and test sets and perform Z-normalisation on the merged data.

We first compute Auto-Encoded Compact Representation for time-series following algorithm HC-AECS.AECS() stated in proposed method. We have used Mean-Squared Error (MSE) as loss function and Stochastic Gradient Descent (SGD) Robbins and Monro (1951) as optimizer with learning rate ($lr$) as 0.004 and momentum as 0. We use batch size 32.

Let us take a representative dataset Adiac to explain thoroughly. It is an univariate time-series with 781 instances (M = 781) where each time-series is of length 176 (n=176) i.e $X\in R^{781\times 176}$. The input layer $h_{0}$ produces a representation of the same length as that of the number of time steps i.e $h_{0}\in R^{781\times 176}$. This representation is passed to first hidden layer $h_{1}$ which scales down the length of each instance to $h_{l1}=16$ ($h_{1}\in R^{781\times 16}$). This compressed sequence is fed to layer $h_{2}$ of length $h_{l2}=12$ where we get a further compressed representation for each $X_{i}$ of Adiac as $AECS_{Adiac}\in R^{781\times 12}$.

Also let us consider a multivariate dataset ERing having 60 instances (M=60) of time-series length 65 (n=65) with 4 dimensions i.e $X\in R^{60\times 65\times 4}$. Each timestep consists of a vector of size 4 corresponding to the dimensions. Here after passing the dataset through the input layer $h_{0}$ we get a representation for each time-series of length equal to the number of timesteps. Correspondingly layers $h_{1}$ and $h_{2}$ produces representations of length $h_{l1}$ and $h_{l2}$ for each time-series respectively.

We started our analysis using 7 different distance measures - Chebyshev, Cosine, Euclidean, Canberra, Manhattan, Mahalanobis and Cross-correlation. After extensive analysis, we have concluded that Chebyshev, Manhattan and Mahalanobis performs better than the other measures on raw time-series as well as on proposed compact representation of time-series - AECS. Figure 3 depicts RI measures of HC on raw time-series using the above mentioned seven distance measures on four representative datasets. We can see either Chebyshev (CH), Manhattan (MA), Mahalanobis (ML) performs best in each one of these datasets. Hence, we use Chebyshev, Manhattan and Mahalanobis as the distance measures in our algorithm.

Next we calculate $\mathcal{T}$ based on Eq.6 and Eq.7 and HC-AECS.BestCluster() on the clusterings formed by above concluded three distance measures, $\mathcal{T}_{CH}$ for CH, $\mathcal{T}_{MA}$ for MA and $\mathcal{T}_{ML}$ for ML. Extending the example presented in section 4.2, for Adiac we obtained $\mathcal{T}_{CH}=1.38$, $\mathcal{T}_{MA}=1.56$, $\mathcal{T}_{ML}=1.90$. Hence we choose clusters formed by ML as it has highest $\mathcal{T}$. While validation we see RI of clusters formed by ML is indeed the highest among the three RI measures ( RI (CH) = 0.698, RI (MA) = 0.780, RI (ML) = 0.926).

We illustrate further in Figure 4 the plot of $\mathcal{T}$ vs RI for the said 3 distance measures on 3 univariate (DistalPhalanxOAG (DPOAG), CricketX, and SwedishLeaf) and 3 multivariate (Libras, Handwriting and NATOPS) time-series. In this figure we can observe that for each time-series, $\mathcal{T}$ and RI are linearly dependent signifying the clustering having the best $\mathcal{T}$ also has the highest RI using either CH, MA, ML. As depicted in the figure for the univariate time-series CH produces the best clustering for DPOAG and ML produces best clustering for CricketX and SwedishLeaf. Whereas for multivariate time-series, MA performs best for NATOPS and ML performs best for Libras and Handwriting.

We compare HC-AECS with important State-of-the-Art clustering mechanisms on time-series like HC with parametric Dynamic Time Warping (DTW)-derivative(DDTW) (HC-DDTW)Łuczak (2016) and K-shapePaparrizos and Gravano (2015). Authors of Łuczak (2016) have already established that they perform better than K-means hence we don’t include K-means in our analysis.

We have used RI as external clustering measure to perform unbiased comparison with the prior work as it has presented performance measure using RI. Further we have extended our analysis by using Normalized Mutual Information (NMI) Strehl and Ghosh (2002) as external clustering measure across a subset of considered 52 time-series. We have observed that NMI reflects also the same trend as obtained using RI as depicted in Table 3. As an example DistalPhalanxOAG, where NMI(HC-DDTW) = 0.414, NMI(K-Shape) = 0.252 and NMI(HC-AECS) = {0.458, 0.321, 0.017} for the distance measures CH, MA and ML respectively. Here HC-AECS using CH gives the best performance as obtained using RI.