Keywords: . . . . . . . . . 43 linearly separable clusters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . index . . . . . . . . . . . . . . . . . . 73 consumption . . . . . . . . . . . . . . 74 Temporal data, clustering, time warp, temporal alignment kernel, k-means, kernel k-means

Last, but not least, I would like to thank my wife, Zohreh, for standing beside me throughout my career and writing this thesis, and also for her understanding and love during the past few years. She has been my inspiration and motivation

First and foremost, I would like to express my gratitude to my dear adviser, Ahlame DOUZAL, for her support, patience, and encouragement throughout my graduate studies. It has been an honor to be her Ph.D. student. It is not often that one finds an adviser and colleague that always finds the time for listening to the little problems and roadblocks that unavoidably crop up in the course of performing research. Her technical and editorial advice was essential to the completion of this dissertation and has taught me innumerable lessons and insights on the workings of academic research in general. She has been great teacher who has prepared me to get to this place in my academic life.

I have been extremely lucky to have a co-supervisor, Eric GAUSSIER, who cared about my work, and who responded to my questions and queries so promptly. I appreciate his wide knowledge and his assistance in writing reports (i.e., conference papers, journals and this thesis). He provided me with direction, technical support and became more of a mentor and friend, than a professor. I would like to thank my dissertation committee members Mohamed NADIF, Paul HONEINE, Pierre-François MARTEAU and Jose Antonio VILAR FERNANDEZ for their interest in my work. They have each provided helpful feedback.

Completing this work would have been all the more difficult were it not for the support and friendship provided by the other members of the AMA team in Laboratoire d'Informatique de Grenoble. I am indebted to them for their help. I want to thank my parents, who brought a 80-86 computer into the house when I was still impressionable enough to get sucked into a life of computer science. My parents, Ali and Safa, receive my deepest gratitude and love for their dedication and the many years of support during my undergraduate studies that provided the foundation for this work. I doubt that I will ever be able to convey my appreciation fully, but I owe them my eternal gratitude.

Introduction

Due to rapid increase in data size, the idea of discovering hidden information in datasets has been exploded extensively in the last decade. This discovery has been centralized mainly on data mining, classification and clustering. One major problem that arises during the mining process is handling data with temporal feature. Temporal data naturally arise in various emerging applications as sensor networks, dynamic social networks, human mobility or internet of things. Temporal data refers to data, where changes over time or temporal aspects play a central role or are of interest. Unlike static data, there is high dependency among time series or sequences and the appropriate treatment of data dependency or correlation becomes critical in any temporal data processing.

Temporal data mining has recently attracted great attention in the data mining community [START_REF] Bradley | Refining Initial Points for K-Means Clustering[END_REF]; [START_REF] Halkidi | On Clustering Validation Techniques[END_REF]; [START_REF] Kalpakis | Distance Measures for Effective Clustering of ARIMA Time-Series[END_REF]; [START_REF] Keogh | Clustering of Time Series Subsequences is Meaningless: Implications for Previous and Future Research[END_REF]. Basically temporal data mining is concerned with the analysis of temporal data and for finding temporal patterns and regularities in sets of temporal data. Since temporal data mining brings together techniques from different fields such as statistics, machine learning and databases, the literature is diffused among many different sources. According to the techniques of data mining and theory of statistical time series analysis, the theory of temporal data mining may involve the following areas of investigation since a general theory for this purpose is yet to be developed [START_REF] Lin | An Overview of Temporal Data Mining[END_REF]. Temporal data mining tasks include: characterization and representation, similarity computation and comparison, temporal association rules, clustering, classification, prediction and trend analysis, temporal pattern discovery, etc.

Clustering temporal data is a fundamental and important task, usually applied prior to any temporal data analysis or machine learning tasks, for summarization, extraction groups of time series and highlight the main underlying dynamics; all the more crucial for dimensionality reduction in big data context. Clustering problem is about partitioning a given dataset into groups or clusters such that the data points in a cluster are more similar to each other than data points in different clusters [START_REF] Guha | CURE: An Efficient Clustering Algorithm for Large Databases[END_REF]. It may be found under different names in different contexts, such as "unsupervised learning" in pattern recognition, "numerical taxonomy" in biology, "typology" in social sciences and "partition" in graph theory [START_REF] Theodoridis | Pattern Recognition[END_REF]. There are many applications where a temporal data clustering activity is relevant. For example in web activity logs, clusters can indicate navigation patterns of different user groups. In financial data, it would be of interest to group stocks that exhibit similar trends in price movements. Another example could be clustering of biological sequences like proteins or nucleic acids, so that sequences within a group have similar functional properties [F.88]; [START_REF] Miller | A Comprehensive Approach to Clustering of Expressed Human Gene Sequence: The Sequence Tag Alignment and Consensus Knowledge Base[END_REF]; [START_REF] Osato | A computer-based method of selecting clones for a fulllength cDNA project: Simulataneous collection of negligibly redundant and variant cDNAs[END_REF]. Hense, there are a variety of methods for clustering of temporal data and many of usual clustering methods need to compare the time series. Most of these comparison methods are derived from the Dynamic Time Warping (DTW) and align the observations of pairs of time series to identify the delays that can occur between times. Temporal clustering analysis provides an effective manner to discover the inherent structure and condense information over temporal data by exploring dynamic regularities underlying temporal data in an unsupervised learning way. In particular, recent empirical studies in temporal data mining reveal that most of the existing clustering algorithms do not work well due to their complexity of underlying structure and data dependency [START_REF] Keogh | On the need for time series data mining benchmarks: A survey and empirical[END_REF], which poses a real challenge in clustering temporal data of a high dimensionality, complicated temporal correlation, and a substantial amount of noise.

Earlier, the clustering techniques have been extensively applied in various scientific areas. k-means-based clustering, viz. standard k-means [START_REF] Macqueen | Some methods for classification and analysis of multivariate observations[END_REF], k-means++ and all its variations, is among the most popular clustering algorithms, as it provides a good trade-off between the quality of obtained solution and its computational complexity [START_REF] Arthur | K-means++: The Advantages of Careful Seeding[END_REF]. However, time series k-means clustering, under the commonly used dynamic time warping (DTW) [START_REF] Sakoe | Dynamic programming algorithm optimization for spoken word recognition[END_REF]; [START_REF] Kruskall | The symmetric time warping algorithm: From continuous to discrete[END_REF] or several well-established temporal kernels (e.g. κ DTAK , κ GA ) [START_REF] Bahlmann | Online handwriting recognition with support vector machines-a kernel approach[END_REF]; [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF]; [START_REF] Cuturi | A kernel for time series based on global alignments[END_REF]; [START_REF] Cuturi | Fast global alignment kernels[END_REF], is challenging as estimating cluster centroids requires aligning multiple time series simultaneously. Under temporal proximity measures, one standard way to estimate the centroid of two time series is to embed the series into a new Euclidean space according to the optimal mapping between them. To average more than two time series, the problem becomes more complex as one needs to determine a multiple alignment that link simultaneously all the time series. That multiple alignment defines, similarly, a new Euclidean embedding space where time series can be projected and the global centroid estimated. Each dimension identifies, this time, a hyper link that connects more than two elements. Such alignments, referred to as multiple sequence alignment, become computationally prohibitive and impractical when the number of time series and their length increase [START_REF] Thompson | CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice[END_REF]; [START_REF] Wang | On the complexity of multiple sequence alignment[END_REF]; [START_REF] Higgins | T-Coffee: A novel method for fast and accurate multiple sequence alignment[END_REF].

To bypass the centroid estimation problem, costly k-medoids [START_REF] Kaufman | Clustering by means of Medoids, in Statistical Data Analysis Based on the L1-Norm and Related Methods[END_REF] and kernel k-means [START_REF] Girolami | Mercer kernel-based clustering in feature space[END_REF]; [START_REF] Inderjit | Kernel k-means: spectral clustering and normalized cuts[END_REF] are generally used for time series clustering [START_REF] Liao | Clustering of time series data -A survey[END_REF]. For k-medoids, a medoid (i.e. the element that minimizes the distance to the other elements of the cluster), is a good representative of the cluster mainly when time series have similar global dynamics within the class; it is however a bad representative for time series that share only local features [START_REF] Frambourg | Learning Multiple Temporal Matching for Time Series Classification[END_REF]. For kernel k-means [START_REF] Girolami | Mercer kernel-based clustering in feature space[END_REF]; [START_REF] Inderjit | Kernel k-means: spectral clustering and normalized cuts[END_REF], as centroids can not be estimated in the Hilbert space, pairwise comparisons are used to assign a given time series to a cluster. While k-means, of linear complexity, remains a fast algorithm, k-medoids and kernel k-means have a quadratic complexity due to the pairwise comparisons involved. This work proposes a fast, efficient and accurate approach that generalizes the k-means based clustering algorithm for temporal data based on i) an extension of the standard time warp measures to consider both global and local temporal differences and ii) a tractable and fast estimation of the cluster representatives based on the extended time warp measures. Temporal data centroid estimation is formalized as a non-convex quadratic constrained optimization problem, while all the popular existing approaches are of heuristic nature, with no guarantee of optimality. The developed solutions allow for estimating not only the temporal data centroid but also its weighting vector, which indicates the representativeness of the centroid elements. The solutions are particularly studied under the standard Dynamic Time Warping (DTW) [START_REF] Kruskall | The symmetric time warping algorithm: From continuous to discrete[END_REF], Dynamic Temporal Alignment Kernel (κ DTAK ) [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF] and Global Alignment kernels (κ GA and κ TGA ) [START_REF] Cuturi | A kernel for time series based on global alignments[END_REF]; [START_REF] Cuturi | Fast global alignment kernels[END_REF]. A wide range of public and challenging datasets, which are non-isotropic (i.e., non-spherical), not well-isolated (and thus non linearly separable), are used to compare the proposed generalized k-means with k-medoids and kernel k-means. The results of this comparison illustrate the benefits of the proposed method, which outperforms the alternative ones on all datasets. The impact of isotropy and isolation of clusters on the effectiveness of the clustering methods is also discussed.

The main contributions of the thesis are:

• Extend commonly used time warp measures: dynamic time warping, which is a dissimilarity measure, dynamic time warping kernel, temporal alignment kernel, and global alignment kernel, which are three similarity measures, to capture both global and local temporal differences.

• Formalize the temporal centroid estimation issue as an optimization problem and propose a fast and tractable solution under extended time warp measures.

• Propose a generalization of the k-means based clustering for temporal data under the extended time warp measures.

• Show through a deep analysis on a wide range of non-isotropic, linearly non-separable public and challenging datasets that the proposed solutions are faster and outperforms alternative methods, through (a) their efficiency on clustering and (b) the relevance of the centroids estimated.

In the remainder of the thesis, we use bold, lower-case letters for vectors, time series and alignments, the context being clear to differentiate between these elements.

Organisation du manuscrit

This manuscript is organized in four chapters as follows. The first chapter presents the major time series metrics, notations and definitions. As the most conventional metrics between time series are based on the concept of alignment, we then present the definition of temporal alignment, prior to introduce kernels for times series. In the second chapter, we will discuss about time series averaging and state-of-the-art of centroid estimation approaches. To do so, we present the consensus sequence and multiple alignment problem for time series averaging under time warp and review the related progressive and iterative approaches for centroid estimation. We show the importance of proposing a tractable and fast centroid estimation method. In the next chapter, we first study the k-means based clustering algorithms and its variations, their efficiency and complexities. While k-means based clustering is among the most popular clustering algorithms, it is a challenging task because of centroid estimation that needs to deal with the tricky multiple alignment problem under time warp. We then proposes the generalized k-means based clustering for time series under time warp measures. For this, we introduce motivation of this generalization, prior to the problem formalization. We give then an extension of the standard time warp measures and discuss about their properties. In the following, we describe the fully formalized centroid estimation procedure based on the extended time warp measures. We present the generalized centroid-based clustering for temporal data in the context of k-means under both dissimilarity and similarity measures, but as shown later, the proposed solutions being directly applicable to any other centroid-based clustering algorithms. Lastly, the conducted experimentation and results obtained are discussed in the chapter four and we show through a deep analysis on a wide range of non-isotropic, linearly non-separable public data that the proposed solutions outperforms the alternative methods. The quantitative evaluation is finally completed by the qualitative comparison of the obtained centroids to the one of the time series they represent (i.e., ground truth). In the conclusion, we summarize the contributions of this thesis and discuss possible perspectives obtained from the thesis. The temporal sequences and time series are two frequently discussed topics in the literature. The time series can be seen as a special case of sequences, where the order criteria is time. In this chapter, we initially review the state-of-the-art of comparison measurements between pairs of time series. As the most conventional metrics between time series are based on the concept of alignment, we present the definition of temporal alignment, prior to introducing the kernels for times series.

Notations

Introduction

What is a time series? A time series is a kind of sequence with an ordered set of observations where the order criteria is time. A large variety of real world applications, such as meteorology, marketing, geophysics and astrophysics, collect observations that can be represented as time series. The most obvious example for a time series is probably the stock prices over successive trading days or hourly household electric power consumption over successive hours. But not only the industry and financial sector produce a large amount of such data. Social media platforms and messaging services record up to a billion daily interactions [START_REF] Piro | Chat reaches 1 billion messages sent per day[END_REF], which can be treated as time series. Besides the high dimensions of these data, the medical and biological sector provide a great variety of time series, as gene expression data, electrocardiograms, growth development charts and many more. In general, an increasingly large part of worlds data is in the form of time series [START_REF] Maimon | Data Mining and Knowledge Discovery Handbook[END_REF]. Although statisticians have worked with time series for a century, the increasing use of temporal data and its special nature have attracted the interest of many researchers in the field of data mining and analysis [START_REF] Maimon | Data Mining and Knowledge Discovery Handbook[END_REF]; [START_REF] Fu | A Review on Time Series Data Mining[END_REF].

Time series analysis includes methods for analyzing temporal data in order to discover the inherent structure and condense information over temporal data, as well as the meaningful statistics and other characteristics of the data. In the context of statistics, finance, geophysics and meteorology the primary goal of time series analysis is forecasting. In the context of signal processing, control engineering and communication engineering it is used for signal detection and estimation, while in the context of data mining, pattern recognition and machine learning time series analysis can be used for clustering, classification, anomaly detection as well as forecasting. One major problem that arises during time series analysis is comparison of time series or sequences. For this comparison, the most studied approaches rely on proximity measures and metrics. Thus, a crucial question is establishing what we mean by "(dis)similar" data objects (i.e. determining a suitable (dis)similarity measure between two objects). In the specific context of time series data, the concept of (dis)similarity measure is particularly complex due to the dynamic character of the series. Hence, different approaches to define the (dis)similarity between time series have been proposed in the literature and a short overview of well-established ones is presented below.

We start this section by defining the general notion of proximity measures. Next, we discuss major basic proximity measures of time series. Later, we introduce definition of different kernels proximity measures for time series, as well as the explanation of temporal alignment, where the most conventional proximity measures between time series and sequences are based on the concept of alignment.

Comparison of time series

Mining and comparing time series address a large range of challenges, among them: the meaningfulness of the distance, similarity and dissimilarity measures. They are widely used in many research areas and applications. Distances or (dis)similarity measures are essential to solve many signal processing, image processing and pattern recognition problems such as classification, clustering, and retrieval problems. Various distances, similarity and dissimilarity measures that are applicable to compare time series, are presented over the past decade. We start here by defining the general notion of proximity measures.

Similarity measure

Let X be a set of data. A function s: X × X → R is called a similarity on X, it obeys to the following properties, ∀ x, y ∈ X:

• non-negativity: s(x, y) ≥ 0

• symmetry: s(x, y) = s(y, x) • if x ̸ = y, s(x, x) = s(y, y) > s(x, y)
A similarity measure takes on large values for similar objects and zero for very dissimilar objects (e.g. in the context of cluster analysis). In general, the similarity interval is [0, 1], where 1 indicates the maximum of similarity measure.

Dissimilarity measure

A function d: X × X → R is called a dissimilarity on X if, for all x, y ∈ X, it holds the three fundamentals following properties:

• non-negativity: d(x, y) ≥ 0 • symmetry: d(x, y) = d(y, x) • reflexivity: d(x, x) = 0
In some cases, the proximity measurement takes minimum values when the time series are similar together. For this purpose, we introduce a new measure of proximity opposite of the similarity index, called dissimilarity measure. The main transforms used to obtain a dissimilarity d from a similarity s are:

• d(x, y) = (s(x, x)-s(x, y)) • d(x, y) = (s(x, x)-s(x, y)) / s(x, y) • d(x, y) = (s(x, x)-s(x, y)) 1/2

Distance metric

Distance measures play an important role for (dis)similarity problem, in data mining tasks. Concerning a distance measure, it is important to understand if it can be considered metric. A function D: X × X → R is called a distance metric on X if, for all x, y and z ∈ X, it holds the properties:

• non-negativity: D(x, y) ≥ 0

• symmetry: D(x, y) = D(x, y)

• reflexivity: D(x, y) = 0, if and only if x = y

• triangle inequality: D(x, y) ≤ D(x, z)+ D(z, y) It's easy to prove that the non-negativity property is included in the last three properties. If any of these is not obeyed, then the distance is non-metric. Although having the properties of a metric is desirable, a (dis)similarity measure can be quite effective without being a metric. In the following, to simplify, we use broadly the term metric to reference both.

Times series basic metrics

Time series analysis is an active research area with applications in a wide range of fields. One key component in temporal analysis is determining a proper (dis)similarity metric between two time series. This section is dedicated to briefly describe the major well-known measures in a nutshell, which have been grouped into two categories: value-based and behavior-based metrics.

Value-based metrics

A first way to compare time series data involves concept of values and distance metrics, where time series are compared according to their values. This subsection relies on two standard well-known division: (a) without delays (e.g. Minkowski distance) and (b) with delays (e.g. Dynamic Time Warping).

Without delays

Euclidean-(or Euclidian) The most used distance function in many applications, which is commonly accepted as the simplest distances between sequences. The Euclidean distance d E (L 2 norm) between two time series x i = (x i1 , ..., x iT ) and x j = (x j1 , ..., x jT ) of length T , is defined as:

d E (x i , x j ) = T t=1 (x it -x jt ) 2
Minkowski Distance-The generalization of Euclidean Distance is Minkowski Distance, called L p norm. It is defined as:

d Lp (x i , x j ) = p T t=1 (x it -x jt ) p
where p is called the Minkowski order. In fact, for Manhattan distance p = 1, for the Euclidean distance p = 2, while for the Maximum distance p = ∞. All L p -norm distances do not consider the time warp. Unfortunately, they don't correspond to the common understanding of what a time series sequences is, and can not capture flexible similarities.

With delays

Dynamic time warping (DTW)-Searching the best alignment that matches two time series is an important task for many researcher. One of the eminent techniques to execute this task is Dynamic Time Warping (DTW) was introduced in [RD62]; [START_REF] Itakura | Minimum prediction residual principle applied to speech recognition[END_REF]; [START_REF] Kruskall | The symmetric time warping algorithm: From continuous to discrete[END_REF] with application in speech recognition. It finds the optimal alignment between two time series, and captures flexible similarities by aligning the elements inside both sequences. Intuitively, the time series are warped non-linearly in the time dimension to match each other. Simply, it is a generalization of Euclidean distance which allows a non-linear mapping of one time series to another one by minimizing the distance between both. DTW does not require that the two time series data have the same length, and it can handle local time shifting by duplicating (or re-sampling) the previous element of the time sequence.

Let X = {x i } N i=1 be a set of time series x i = (x i1 , ..., x iT ) assumed of length T1 . An alignment π of length |π| = m between two time series x i and x j is defined as the set of m (T ≤ m ≤ 2T -1) couples of aligned elements of x i to elements of x j :

π = ((π 1 (1), π 2 (1)), (π 1 (2), π 2 (2)), ..., (π 1 (m), π 2 (m)))
where π defines a warping function that realizes a mapping from time axis of x i onto time axis of x j , and the applications π 1 and π 2 defined from {1, ..., m} to {1, .., T } obey to the following boundary and monotonicity conditions:

1 = π 1 (1) ≤ π 1 (2) ≤ ... ≤ π 1 (m) = T 1 = π 2 (1) ≤ π 2 (2) ≤ ... ≤ π 2 (m) = T and ∀ l ∈ {1, ..., m}, π 1 (l + 1) ≤ π 1 (l) + 1 and π 2 (l + 1) ≤ π 2 (l) + 1, (π 1 (l + 1) -π 1 (l)) + (π 2 (l + 1) -π 2 (l)) ≥ 1
Intuitively, an alignment π defines a way to associate all elements of two time series. The alignments can be described by paths in the T × T grid as displayed in Figure 1.1, that crosses the elements of x i and x j . For instance, the green path aligns the two time series as: ((x i1 ,x j1 ),(x i2 ,x j2 ),(x i2 ,x j3 ),(x i3 ,x j4 ),(x i4 ,x j4 ),(x i5 ,x j5 ),(x i6 ,x j6 ),(x i7 ,x j7 )). We will denote A as the set of all possible alignments between two time series. Figure 1.1: Three possible alignments (paths) between x i and x j Dynamic time warping is currently a well-known dissimilarity measure on time series and sequences, since it makes them possible to capture temporal distortions. The Dynamic Time Warping (DTW) between time series x i and time series x j , with the aim of minimization of the mapping cost, is defined by:

dtw(x i , x j ) = min π∈A 1 |π| (t ′ ,t)∈π ϕ(x it ′ , x jt ) (1.1)
where ϕ : R × R → R + is a positive, real-valued, divergence function (generally Euclidean norm). Figure 1.2 shows the alignment between two sample time series with delay (dynamic time warping) in comparison with the alignment without delay (Euclidean). While Dynamic Time Warping alignments deal with delays or time shifting (see Figure 1.2-left), the Euclidean alignment π between x i and x j aligns elements observed at the same time:

π = ((π 1 (1), π 2 (1)), (π 1 (2), π 2 (2)), ..., (π 1 (T ), π 2 (T )))
where ∀t = 1, ..., T, π 1 (t) = π 2 (t) = t, |π| = T (see Figure 1.2-right). According to the alignment definition, the euclidean distance (d E ) between the time series x i and x j is given by:

d E (x i , x j ) def = 1 |π| |π| k=1 ϕ(x i π 1 (k) , x j π 2 (k) ) = 1 T T t=1 ϕ(x it , x jt )
ϕ taken as the Euclidean norm.

Finally, Figure 1.3 shows the optimal alignment path between two sample time series with and without considering time warp. To recall, boundary, monotonicity and continuity are three important properties which are applied to the construction of the path. The boundary condition imposes that the first elements of the two time series are aligned to each other, as well as the last sequences elements. In other word, the alignment refers to the entire both time series sequences. Monotonicity preserve the time-ordering of elements. It means, the alignment path doesn't go back in "time" index. Continuity limits the warping path from long jumps and it guarantees that alignment does not omit important features. A dynamic programming approach is used to find the minimum distance for alignment path. Let x i and x j two time series. A two-dimensional |x i | by |x j | cost matrix D is built up, where the value at D(t, t ′ ) is the minimun distance warp path that can be constructed from the two time series x i = (x i1 , ..., x it ) and x j = (x j1 , ..., x jt ′ ). Finally, the value at D(|x i |, |x j |) will contain the minimum distance between the times series x i and time series x j under time warp. Note that DTW is not a metric as not satisfy the triangle inequality.

Property 1.1 DTW does not follow triangle inequality.

As we can see, for instance, given three time series data,

x i = [0], x j = [1, 2] and x k = [2,3,3], then: dtw(x i , x j ) = 3 , dtw(x j , x k ) = 3 , dtw(x i , x k ) = 8. Hense, dtw(x i , x j ) + dtw(x j , x k ) dtw(x i , x k ).
Time and space complexity of the dynamic time warping is straightforward to define. Each cell in the cost matrix is filled once in constant time. The cost of the optimal alignment can be recursively computed by:

D(t, t ′ ) = d(t, t ′ ) + min ⎧ ⎨ ⎩ D(t -1, t ′ ) D(t -1, t ′ -1) D(t, t ′ -1) ⎫ ⎬ ⎭
where d is a distance function between the elements of time series (e.g. Euclidean). This yields complexity of O(N 2 ). Dynamic time warping is able to find the optimal global alignment path between time series and is probably the most commonly used measure to assess the dissimilarity between them. Dynamic time warping (DTW) has enjoyed success in many areas where its time complexity is not an issue. However, conventional DTW is much too slow for searching an alignment path for large datasets. For this problem, the idea of the DTW technique for aligning time series is changing with the aim of improvement time and space complexity and accuracy. In general, the methods which make DTW faster divide in three main different categories: 1) constraints (e.g., Sakoe-Chiba band), 2) indexing (e.g., piecewise), and 3) data abstraction (e.g., multiscale). Here we briefly mentioned some major ideas that try to enhance this technique. The cost matrix are filled by the DTW algorithm in the shaded areas around the diagonal. In this case, the algorithm finds the optimal alignment warp math through the constraints window. However, the global optimal alignment path will not be found, if it is not fully inside the window. The Sakoe-Chiba dynamic time warping (DTW sc ) between two time series, which is the most well-used constrained DTW, is defined by:

dtw sc (x i , x j ) def = min π∈A C(π) C(π) def = 1 |π| |π| k=1 w π 1 (k),π 2 (k) ϕ(x i π 1 (k) , x j π 2 (k) ) = 1 |π| (t,t ′ )∈π w t,t ′ ϕ(x it , x jt ′ )
where if |tt ′ | < c then w t,t ′ = 1, and else ∞. Function ϕ taken as the euclidean norm and parameter c is the Sakoe-Chiba band width. Figure 1.5 shows examples of the cost matrix and the two constrained optimal alignment path between two sample time series. Piecewise DTW

Piecewise Dynamic Time Warping (PDTW) [START_REF] Keogh | Scaling Up Dynamic Time Warping for Data Mining Applications[END_REF], takes advantage of the fact that one can efficiently approximate most of the time series by a Piecewise Aggregate Approximation (PAA). Let x i = (x i1 , x i2 , . . . , x iT ) be a time series sequence that we wish to reduce its dimensionality to T ′ , where (1 ≤ T ′ ≤ T ). A time series x * i of length T ′ is represented by x * i = (x i1 , ..., x iT ′ ), which t th element of x * i is calculated by the following equation:

x * it = T ′ T T T ′ t t ′ = T T ′ (t-1)+1
x it ′ Simply, to reduce the data from T dimensions to T ′ dimensions, the data is divided into T ′ frames with equal size. The average value of the data falling within each frame is computed and a vector of these values becomes the data reduced representation. The compression rate c (c = T T ′ ) is equal to the ratio of the length of the original time series to the length of its PAA representation. Choosing a value for c, is a trade-off between memory savings and accuracy [START_REF] Keogh | Scaling Up Dynamic Time Warping for Data Mining Applications[END_REF]. Figure 1.6 illustrates a time series and its PAA approximation. Therefore, piecewise method by finding a time series which is most similar to a given one, significantly speed up many DTW applications by reducing the number of times which the DTW should be run.

MultiScale DTW (FastDTW)

To obtain an efficient as well as robust algorithm to compute DTW-based alignments, one can combine the global constraints and dimensionality reduction strategies in some iterative approaches to generate data-dependent constraint regions. The general strategy of multiscale DTW (or FastDTW) [START_REF] Salvador | FastDTW: Toward Accurate Dynamic Time Warping in Linear Time and Space[END_REF] is to recursively project a warping path computed at a coarse resolution level to the next level and then refine the projected warping path. Multiscale DTW or fastDTW speeds up the dynamic time warping algorithm by running DTW on a reduced representation of the data.

Let x i = (x i1 , ..., x iT ) and x j = (x j1 , ..., x jT ) be the sequences to be aligned, having lengths T . The aim is to compute an optimal warping path between two time series x i and x j . The highest resolution level will be called as Level 1. By reducing the feature sampling rate by a factor of f , one achieves time series of length T f . Next, one computes an optimal alignment warping path between two dimension-reduced time series on the resulting resolution level (Level 2). The obtained path is projected onto highest level (Level 1) and it defines a constraint region R. Lastly, an optimal alignment warping path relative to the restricted area R is computed. In general, the overall number of cells to be computed in this process is much smaller than the total number of cells on Level 1 (T 2 ). The constrained warping path may not go along with the optimal alignment warping path. To mitigate this problem, one can increase the constraint region R by adding δ cells to the left, right, top and bottom of every cell in R. The resulting region R δ will be called as δ-neighborhood of R [START_REF] Zinke | Applied Dynamic Programming[END_REF]. Finally, in spite of the great success of DTW and its variants in a diversity of domains, there are still several persistent myths about it, as finding ways to speed up DTW with no (or relaxed) constraints.

Behavior-based metrics

The second category of time series metrics concerns behavior-based metrics, where time series are compared according to their behaviors regardless of their values. That is the case when time series of a same class exhibit similar shape or behavior, and time series of different classes have different shapes. Hence, comparing the time series on the base of their value may not be valid assumption. In this context, we should define which time series are more similar together and which ones are different. The definition of similar, opposite and different behaviors are given as following.

Similar Behavior. Time series x i , x j are of similar behavior, if for each period [t i , t i+1 ], i ∈ {1, ..., T }, they increase or decrease simultaneously with the same growth rate.

Opposite Behavior. Time series x i , x j are of opposite behavior, if for each period [t i , t i+1 ], i ∈ {1, ..., T }, when one time series increases, the other decreases with the same growth rate in absolute value and vice-versa.

Different Behavior. Time series x i , x j are of different behavior, if they are not similar nor opposite (linearly and stochastically independent). Main techniques to recover time series behaviors are: slopes and derivatives comparison, ranks comparison, Pearson and temporal correlation coefficient, and difference between auto-correlation operators. In the following, we briefly describe some well-used behavior-based metrics.

Pearson CORrelation coefficient (COR)

Let X = {x 1 , ..., x N } be a set of time series x i = (x i1 , ..., x iT ) , i ∈ {1, ..., N}. The correlation coefficient between sequences x i and x j is defined by:

cor(x i , x j ) = T t=1 (x it -x i )(x jt -x j ) T t=1 (x it -x i ) 2 T t=1 (x jt -x j ) 2
(1.2) and

x i = 1 T T t=1 x it
Correlation coefficient was first introduced by Bravais and later shown by Pearson [START_REF] Pearson | Contributions to the mathematical theory of evolution[END_REF] to be the best possible correlation between two time series. Till now, many applications in different domains such as speech recognition, system design control, functional MRI and gene expression analysis have used the Pearson correlation coefficient as a behavior proximity measure between time series sequences [START_REF] Macarthur | GATE: Software for the analysis and visualization of high-dimensional time series expression data[END_REF]; [START_REF] Ernst | Clustering short time series gene expression data[END_REF]; [START_REF] Abraham | An integrated framework for simultaneous classification and regression of time-series Data[END_REF]; [START_REF] Cabestaing | Classification of evoked potentials by Pearsonís correlation in a Brain-Computer Interface[END_REF]; [START_REF] Rydell | Robust Correlation Analysis with an Application to Functional MRI[END_REF]. The Pearson correlation coefficient changes between -1 and +1. The case COR = +1, called perfect positive correlation, occurs when two time series sequences perfectly coincide, and the case COR= -1, called the perfect negative correlation, occurs when two time series behave completely opposite. COR= 0 shows that the time series sequences have different behavior. Note that the higher correlation doesn't conclude the similar dynamics.

Temporal CORrelation coefficient (CORT)

To cope with temporal data, a variant of Pearson correlation coefficient, considering temporal dependency within r, is proposed in [START_REF] Douzal-Chouakria | Classification trees for time series[END_REF], called as Temporal CORrelation coefficient (CORT). The authors considered an equivalent formula for the correlation coefficient relying on pairwise values differences as:

cor(x i , x j ) = t,t ′ (x it -x i t ′ )(x jt -x j t ′ ) t,t ′ (x it -x i t ′ ) 2 t,t ′ (x jt -x j t ′ ) 2 (1.3)
One can see that the Pearson correlation coefficient assumes the independence of data as based on the differences between all pairs of observations at [t, t ′ ]; unlike, the behavior proximity needs only to capture how they behave at [t, t + r]. Therefore, the correlation coefficient is biased by all of the remaining pairs of values observed at interval

[t, t ′ ] with |t -t ′ | > r, r ∈ [1, ..., T -1]. Temporal correlation coefficient defined as: cort(x i , x j ) = t,t ′ m tt ′ (x it -x i t ′ )(x jt -x j t ′ ) t,t ′ m tt ′ (x it -x i t ′ ) 2 t,t ′ m tt ′ (x jt -x j t ′ ) 2 (1.4)
where m tt ′ = 1 if |tt ′ | ≤ r, otherwise 0. CORT belongs to the interval [-1, +1]. The value CORT(x i , x j ) = 1 means that in any observed period [t, t ′ ], x i and x j have similar behaviors. The value CORT(x i , x j ) = -1 means that in any observed period [t, t ′ ], x i and x j have opposite behaviors. Lastly, CORT(x i , x j ) = 0 indicates that the time series are stochastically linearly independent. Temporal correlation is sensitive to noise. So, the parameter r can be learned or fixed a priori, and for noisy time series higher value of r is advised. Figure 1.9 illustrates the comparison of the Pearson correlation coefficient and the temporal correlation coefficient between three samples time series. Figure 1.9: Example of comparison of COR vs CORT (r = 1). COR(x 1 , x 2 ) = -0.90, COR(x 1 , x 3 ) = -0.65, COR(x 2 , x 3 ) = 0.87, CORT(x 1 , x 2 ) = -0.86, CORT(x 1 , x 3 ) = -0.71, CORT(x 2 , x 3 ) = 0.93

Difference between Auto-Correlation Operators (DACO)

Auto-correlation is a representation of the degree of similarity which measures the dependency between a given time series and a shifted version of itself over successive time intervals. Let X = {x 1 , ..., x N } be a set of time series x i = (x i1 , ..., x iT ) i ∈ {1, ..., N}. The DACO between the two time series x i and x j defined as [START_REF] Gaidon | A time series kernel for action recognition[END_REF]:

daco(x i , x j ) = ∥x i -xj ∥ 2 (1.5)
where,

xi = (ρ 1 (x i ), ..., ρ k (x i ))
and

ρ τ (x i ) = T -τ t=1 (x it -x i )(x i (t+τ ) -x i ) T t=1 (x it -x i ) 2
and τ is the time lag. Therefore, DACO compares time series by computing the distance between their dynamics, modeled by the auto-correlation operators. Note that the lower DACO doesn't represent the similar behavior.

Combined (Behavior-value-based) metrics

Considering the definition of time series basic metrics, the value-based one is based on the differences between the values of the time series and does not consider the dynamics and behaviors within the time series. Figure 1.10 shows two configurations: the left side, two time series x 1 and x 2 are definitely opposed, a decrease of one corresponding to a growth of the other and vice-versa. While the right side, the two series x 1 and x 3 are in the same direction, but have variations in intensity and scale. However, the overall behavior is similar. But, according to the definition of a value-based metrics, the time series x 1 and x 2 in left, have the same distance (in Euclidean) with the two time series x 1 and x 3 in right. Most of the time, two time series sharing same configuration are considered close, even though they have very different values, and vice-versa. Hence, choosing a proper metric can be crucial and very important for the time series comparison. In some cases, several behavior and value-based metrics may be implied. Some propositions show the benefit of involving both behavior and value-based metrics through a combination function. Therefore, to compare the time series and define a proximity measure covering both the behaviors and values components, a weighted linear (or geometric) function combines behavior and value-based metrics. In this context, they could be shown as:

BV Lin (x i , x j ) = α.V(x i , x j ) + (1 -α).B(x i , x j ) BV Geom (x i , x j ) = (V(x i , x j )) α . (B(x i , x j )) (1-α)
where V and B is a value-based and behavior-based metrics, respectively. The parameter α ∈ [0, 1] is a trade-off between the behavior and value-based components.

On the other hand, a time series may be considered in the spectral representations, which means the time series may be similar because they share the same frequency characteristics. Hence, in some application, the frequential-based metrics (e.g. wavelet transforms, fourier transforms) will be used. More specific works to combine two different metrics through a combination function proposed in [START_REF] Douzal | Adaptive dissimilarity index for measuring time series proximity[END_REF]; [START_REF] Douzal-Chouakria | Classification trees for time series[END_REF].

Kernels for time series

Over the last ten years estimation and learning methods using kernels have become rather popular to cope with non-linearities, particularly in machine learning. Kernel methods [START_REF] Hofmann | Kernel methods in machine learning[END_REF] have been proved useful to handle and analyze structured data such as images, graphs and texts. They map the data from the original space (i.e. input space) via a nonlinear mapping function Φ(.) to a higher dimensional feature space (i.e. Hilbert space), to discover nonlinear patterns (see Figure 1.11). These methods formulate learning problems, in a reproducing kernel Hilbert space, of functions defined on the data domain expanded in terms of a kernel.

< Φ(x i ), Φ(x j ) > = κ(x i , x j ) < x i , x j > Linearly separable Φ(x)

Feature Space Original Space

Linearly non-separable Usually, a kernel function κ(x i , x j ) is used to directly provide the inner products in a feature space without explicitly defining transformation. The kernel corresponds to the dot product in a (usually high-dimensional) feature space. In this space, the estimation methods are linear, but as long as we can formulate everything in terms of kernel evaluations, we explicitly never have to compute in the high-dimensional feature space. Such inner products can be viewed as measuring the similarity between samples. Let X = {x 1 , ..., x N } a set of samples.

Gram matrix (or kernel matrix) of function κ with respect to the {x 1 , ..., x N } is a N × N matrix defined by:

K := [ κ(x i , x j ) ] ij
where κ is a kernel function. If we use algorithms that only depend on the Gram matrix, K, then we never have to know (or compute) the actual features Φ. This is the crucial point of kernel methods.

Positive definite matrix is a real N × N symmetric matrix, that K ij satisfying:

i,j c i c j K ij > 0 for all c i ∈ R.
Positive semidefinite matrix is a real N × N symmetric matrix, that the Gram matrix K ij satisfying: i,j c i c j K ij ≥ 0 for all c i ∈ R. In the following, to simplify, will call them definite. if we really need > 0, we will say strictly positive definite.

A symmetric function κ : X × X → R which for all x i ∈ X, ∀i ∈ {1, ..., N} gives rise to a positive definite Gram matrix is called a Positive definite kernel.

For time series, several kernels are proposed in the last years, properties of which will be discussed later. Ideally, a useful kernel for time series should be both positive definite and able to handle time series of structured data. Here, we refer positive definite kernels as kernels. Note that, for simplicity, we have restricted ourselves to the case of real valued kernels. The linear kernel is the simplest kernel function. It is given by the inner product < x i , x j > plus an optional constant γ. A polynomial kernel allows us to model feature conjunctions up to the order of the polynomial, whereas the Gaussian kernel is an example of Radial Basis Function (RBF) kernel. In Table 1.1, some kernel function examples are given.

Linear Kernel κ(x i , x j ) = (x T i x j + γ) Polynomial Kernel κ(x i , x j ) = (α.x T i x j + γ) δ Gaussian Kernel κ(x i , x j ) = exp( ||xi-xj || 2 -2σ 2 )
Table 1.

1: Kernel function examples

The adjustable parameters σ, α, γ play a major role in the performance of the kernel, and should be carefully tuned to the problem at hand. In the following, different methods to compute the kernel similarity Gram matrix presented.

Temporal correlation coefficient kernel (κ CORT )

cort(x i , x j ), which is a behavior-based metric between two time series x i and x j and described in Section 1.3.2.2, is a linear kernel and κ CORT defined as a positive definite linear kernel as following:

κ CORT (x i , x j ) = cort(x i , x j ) = < △x i , △x j >
It is easy to prove that κ CORT (x i , x j ) can be constructed from inner products in appropriate Hilbert spaces.

Auto-correlation kernel (κ DACO )

To compare the dynamics of two actions Gaidon [START_REF] Gaidon | A time series kernel for action recognition[END_REF] proposed to compute the distance between their respective auto-correlations. This distance is defined as the Hilbert-Schmidt norm of the difference between auto-correlations detailed in Section 1.3.2.3, and called as associated Gaussian RBF kernel "Difference between Auto-Correlation Operators Kernel". κ DACO is defined by:

κ DACO (x i , x j ) = e -1 σ 2 daco(x i , x j )
κ DACO is a positive definite Gaussian kernel, designed specifically for action recognition. Notice that, both κ CORT and κ DACO are temporal kernels under Euclidean alignments. Under time warp, time series should be first synchronized by dtw, then kernels applied. In the following, we introduce the major well-known temporal kernels under dynamic time warping alignment.

Gaussian dynamic time warping kernel (κ GDTW )

Following [START_REF] Bahlmann | Online handwriting recognition with support vector machines-a kernel approach[END_REF], the dynamic time warping distance can be used as a pseudo negative definite kernel to define a pseudo positive definite kernel. κ GDTW is determined by:

κ GDTW (x i , x j ) = e -1 t dtw(x i , x j )
and using Sakoe-Chiba constrained dynamic time warping distance, κ GDTWsc is defined by:

κ GDTWsc (x i , x j ) = e -1 t dtw sc (x i , x j )
where t is a normalization parameter, and sc is a Sakoe-Chiba band. As the DTW is not a metric (invalid triangle inequality), one could fear that the resulting kernel lacks some necessary properties. Hence, general positive definiteness can not be proven for κ GDTW , as simple counterexamples can be found. However, this kernel can produce good results in some cases [START_REF] Haasdonk | Tangent distance kernels for support vector machines[END_REF]; [START_REF] Deoste | Trainging invariant support vector machines[END_REF], but a problem of this approach is the quadratic computational complexity [START_REF] Bahlmann | Online handwriting recognition with support vector machines-a kernel approach[END_REF].

Dynamic time-alignment kernel (κ DTAK )

The Dynamic Time-Alignment Kernel (κ DTAK ) proposed in [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF] adjusts another similarity measure between two time series by considering the arithmetic mean of the kernel values along the alignment path. κ DTAK , a pseudo positive definite kernel , is then defined by:

κ DTAK (x i , x j ) = max π∈A C(π)
where,

C(π) = 1 |π| |π| l=1 φ(x iπ 1 (l) , x jπ 2 (l) ) = 1 |π| |π| (t,t ′ )∈π φ(x it , x jt ′ ) and φ(x it , x jt ′ ) = κ σ (x it , x jt ′ ) = e -1 σ 2 ∥x it -x jt ′ ∥ 2
Dynamic time-alignment kernel function is complicated and difficult to analyze because the input data is a vector sequence with variable length and non-linear time normalization is embedded in the function. It performs DTW in the transformed feature space and finds the optimal path that maximizes the accumulated similarity. Indeed, the local similarity used for the κ GDTW is the Euclidian distance, whereas the one used in κ DTAK is the Gaussian kernel value. κ DTAK is a symmetric kernel function, however, in general, it may still be a non positive definite kernel. Note that, in practice, several ad-hoc methods (e.g. perturb the whole diagonal by the absolute of the smallest eigenvalue) are used to ensure the positive definiteness, when the Gram matrix is not positive definite.

Global alignment kernel (κ GA )

In kernel methods, both large and small similarities matter, they all contribute to the Gram matrix. Global Alignment kernel (κ GA ) [START_REF] Cuturi | A kernel for time series based on global alignments[END_REF] is not based on an optimal path chosen given a criterion, but takes advantage of all score values spanned by all possible alignments. The κ GA which is positive definite under mild conditions2 , seems to do a better job of quantifying all similarities coherently, because it considers all possible alignments. κ GA is defined by:

κ GA (x i , x j ) = π∈A |π| l=1 k(x iπ 1 (l) , x jπ 2 (l) )
where k(x, y) = e -λΦσ(x,y) , factor λ > 0, and

Φ σ (x, y) = 1 2σ 2 ∥x -y∥ 2 + log(2 -e -1 2σ 2 ∥x-y∥ 2 )
In the sense of kernel κ GA , two sequences are similar not only if they have one single alignment with high score, which results in a small DTW distance, but also share numerous efficient and suitable alignments. Hence, the function Φ σ is a negative definite kernel, it can be scaled by a factor λ to define a local positive definite kernel e -λΦσ . Global alignment kernels have been relatively successful in different application fields [START_REF] Joder | Temporal integration for audio classification with application to musical instrument classification[END_REF]; [START_REF] Ricci | Learning Pedestrian Trajectories with Kernels[END_REF]; [START_REF] De | Clustering Vessel Trajectories with Alignment Kernels under Trajectory Compression[END_REF] and shown to be competitive when compared with other kernels. However, similar to the κ GDTW and κ DTAK kernels, it has quadratic complexity, O(pT 2 ), where T is the length of time series and p denotes the complexity of the measure between aligned instants.

Triangular global alignment kernel (κ TGA )

The idea of Sakoe-Chiba [START_REF] Sakoe | Dynamic programming algorithm optimization for spoken word recognition[END_REF], to speed up the computation of DTW, applied to the κ GA to introduce fast global alignment kernels, named as Triangular Global Alignment kernels. κ TGA kernel [START_REF] Cuturi | Fast global alignment kernels[END_REF] considers a smaller subset of such alignments rather than κ GA . It is faster to compute and positive definite, and can be seen as trade-off between the full κ GA (accurate but slow) and a Sakoe-Chiba Gaussian kernel (fast but limited). This kernel is defined by:

κ TGA (x i , x j ) = π∈A |π| l=1 k(x iπ 1 (l) , x jπ 2 (l) ) k(x iπ 1 (l) , x jπ 2 (l) ) = w π 1 (l),π 2 (l) k σ (x iπ 1 (l) , x jπ 2 (l) ) 2 -w π 1 (l),π 2 (l) k σ (x iπ 1 (l) , x jπ 2 (l) )
where, w is a radial basis kernel in R (a triangular kernel for integers):

w(t, t ′ ) = 1 - |t -t ′ | c
which guarantees that only the alignments π that are close to the diagonal are considered and the parameters c is the Sakoe-Chiba bandwidth. Note that as c increases the κ TGA converges to the κ GA .

In fact, many generic kernels (e.g. Gaussian kernels), as well as very specific ones, describe different notions of similarity between time series, which do not correspond to any intuitive or easily interpretable high-dimensional representation. They can be used to compute similarities (or distances) in feature spaces using the mapping functions. Precisely, the class of kernels that can be used is larger than those commonly used in the kernel methods which known as positive definite kernels.

Conclusion

A large variety of real world applications, such as meteorology, geophysics and astrophysics, collect observations that can be represented as time series. Time series data mining and time series analysis can be exploited from research areas dealing with sequences and signals, such as pattern recognition, image and signal processing. The main problem that arises during time series analysis is comparison of time series.

There has been active research on quantification of the "similarity", the "dissimilarity" or the "distance" between time series. Even, looking for the patterns and dependencies in the visualizations of time series can be a very exhausting task and the aim of this work is to find the more similar behavior time series, while the value-based proximity measures are the most studied approaches to compare the time series. Most of the time, two time series sharing same configuration are considered close, even they have very different values, their appearances are similar in terms of form. In both comparison metrics, a crucial question is establishing what we mean by "similar" or "dissimilar" time series due to the dynamic of the series (with or without considering delay).

To handle and analysis non-linearly structured data, a simple way is to treat the time series as vectors and simply employ a linear kernel or Radial basis function kernel. For this, different kernels have been proposed. To discover non-linear patterns, they map the data from the original input space to a higher dimensional feature space, called Hilbert space. They can be used to compute similarity between the time series. Ideally, a useful kernel for time series should be positive definite and able to handle the temporal data. Notice that, the computational complexity of kernel construction and evaluation can play a critical role in applying kernel methods to time series data.

In summary, finding a suitable proximity measure is a crucial aspect when dealing with the time series (e.g. a kernel similarity, a similarity or dissimilarity measure, a distance), that captures the essence of the time series according to the domain of application. For example, Euclidean distance is commonly used due to its computational efficiency; however, it is very brittle for time series and small shifts of one time series can result in huge distance changes. Therefore, more sophisticated distances have been devised to be more robust to small fluctuations of the input time series. Notably, Dynamic Time Warping (DTW) has enjoyed success in many areas where its time complexity is not an issue. Using the discussed distances and (dis)similarity proximity measures, we will be able to compare time series and analyze them. This can provide useful insights on the time series, before the averaging and centroid estimation step. In the next chapter, we will discuss about averaging time series problems, difficulties and complexities. Averaging a set of time series, under the frequently used dynamic time warping, needs to address the problem of multiple temporal alignments, a challenging issue in various domains. Under temporal metrics, one standard way to average two time series is to synchronize them and average each pairwise temporal warping alignment. To average more than two time series, the problem becomes more complex as one needs to determine a multiple alignment that link simultaneously all the time series on their commonly shared similar elements. In this chapter, we present initially the definition of a consensus sequence and multiple temporal alignment problems, then we review major progressive and iterative centroid estimation approaches under time warp, and discuss their properties and complexities.

Introduction

Estimating the centroid of a set of time series is an essential task of many data analysis and mining processes, as summarizing a set of time series, extracting temporal prototypes, or clustering time series. Averaging a set of time series, under the frequently used dynamic time warping [START_REF] Kruskall | The symmetric time warping algorithm: From continuous to discrete[END_REF]; [START_REF] Sankoff | Time warps, string edits, and macromolecules: the theory and practice of sequence comparison[END_REF] metric or its variants [START_REF] Itakura | Minimum prediction residual principle applied to speech recognition[END_REF]; [START_REF] Rabiner | A Tutorial on Hidden Markov Models and selected applications in speech recognition[END_REF]; [START_REF] Sakoe | Dynamic programming algorithm optimization for spoken word recognition[END_REF], needs to address the tricky multiple temporal alignments problem, a challenging issue in various domains.

Under time warp, one standard way to estimate the centroid of two time series is first to synchronize them (to find the pairwise links), and then average each pairwise temporal warping alignment. In this way, the centroid can be estimated as the average of the linked elements as illustrated in Figure 2.1. Averaging more than two time series under time warp is a very complex problem, because we need to determine a multiple alignment that links simultaneously all the time series on their commonly shared similar elements. That multiple alignment defines, similarly, a new Euclidean embedding space where time series can be projected and the global centroid estimated. Note that each dimension identifies, this time, a hyper link that connects more than two elements. Finding the multiple temporal alignment of a set of time series or its average sequence is a chicken-or-egg problem: knowing the average sequence under time warp provides a multiple temporal alignment and vice-versa.

A first direct approach to determine a multiple alignment is to search, by dynamic programming, the optimal path within an N-dimensional grid that crosses the N time series. The complexity of this approach prevents its use, as it constitutes an NP-complete problem with a complexity of O(T N ) that increases exponentially with the number of time series N and the time series length T [CL88]; [START_REF] Just | Computational complexity of multiple sequence alignment with SP-score[END_REF]. A second manner, that characterizes progressive approaches, is based on combining progressively pairs of time series centroids to estimate the global centroid. Progressive approaches may however suffer from the error propagation problem: Early errors propagate in all subsequent centroids through the pairwise combinations. The third approach is iterative. It works similarly to the progressive approach but reduces the error propagation by repeatedly refining the centroid and realigning it to the initial time series. The progressive and iterative approaches for centroid estimation are mainly derived from the multiple sequence alignment methods to address the challenging problem of aligning more than two time series [START_REF] Higgins | T-Coffee: A novel method for fast and accurate multiple sequence alignment[END_REF]; [START_REF] Sze | A polynomial time solvable formulation of multiple sequence alignment[END_REF]; [START_REF] Thompson | CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice[END_REF]. In general, the main progressive and iterative approaches are of heuristic nature and are limited to the dynamic time warping metric [START_REF] Soheily-Khah | A Comparison of Progressive and Iterative Centroid Estimation Approaches Under Time Warp[END_REF]. Here, we first introduce the definition of pairwise time series averaging, prior to study the consensus sequence. Next, we review some major progressive and iterative approaches for time series averaging under the dynamic time warping to handle the tricky multiple temporal alignment problem.

Definition 2.1

Under the dynamic time warping, given x i and x j and their dynamic time warping alignment π * , the centroid c(x i , x j ) is defined by averaging their aligned elements as follows:

c(x i , x j ) = (Avg(x iπ * 1 (1) , x jπ * 2 (1) ), ..., Avg(x iπ * 1 (m) , x jπ * 2 (m) ))
where Avg is a standard numerical averaging function.

Consensus sequence

The time series comparison is a key point for many problems, as finding a consensus sequence in molecular biology, bioinformatics or data mining. In the context of sequences, this term is used with two meanings: i) the medoid sequence, ii) the average sequence of the set of time series.

Medoid sequence

The concept refers to a formal definition, corresponding to a sequence in the center of a set of time series. A medoid is a representative time series in a set that minimizes the sum of distances (i.e. inertia) to all other time series within the same set. Figure 2.3 shows a medoid as a consensus between some time series sequences.

In the context of (un)supervised learning, many algorithms require a method to represent information from a set of time series in one and only one sequence. Algorithms like k-medoids are using the medoid of a set of objects, which is better in handling noises and outliers than k-means. But, mostly in other algorithms we need to average a set of time series as an optimal consensus.

Average sequence

Since our purpose is to define an average sequence minimizing the sum of distances to all time series of a set (i.e. inertia), we give the definition of an average sequence when the corresponding distance is DTW. Definition 2.2 Let X = {x 1 , ..., x N } a set of time series x i = (x i1 , ..., x iT ) i ∈ {1, ..., N}, and S is the space of all time series sequences. ∀s ∈ S, the average sequence c should satisfy the following:

N i=1 DTW(x i , c) ≤ N i=1 DTW(x i , s)
Since no information on the length of the average sequence c is available, the search cannot be limited to sequences of a given length, so all possible length values for averages have to be considered. In the following, we will review some techniques used to determine the multiple temporal alignment (and average) of a set of time series.

Multiple temporal alignments

Multiple temporal alignment of time series is an important problem with applications in many scientific domains such as speech recognition [START_REF] Rabiner | Fundamentals of speech recognition[END_REF], computer graphics [START_REF] Bruderlin | Motion signal processing[END_REF], astronomy [START_REF] Listgarten | Multiple alignment of continuous time series[END_REF], computer vision [START_REF] Caspi | Aligning non-overlapping sequences[END_REF], and bio-informatics [START_REF] Aach | Aligning gene expression time series with time warping algorithms[END_REF].

The multiple temporal alignment is computable by extending DTW for aligning N time series. This idea can be generalized by computing DTW in a N -dimensional hypercube. Given this global alignment, the consensus can be found by averaging column by column of the multiple alignments. Thus, the first direct approach to determine a multiple temporal alignment is to search, by dynamic programming, the optimal path within an N -dimensional grid that crosses the N time series.

Dynamic programming

Dynamic programming is a direct method to determine a multiple temporal alignment by search the optimal path within an N-dimensional grid that crosses the N time series. The complexity of this approach prevents its use, as it constitutes an NP-complete problem with a complexity of O(T N ) that increases exponentially with N the number of time series, and T the time series length [START_REF] Carrillo | The multiple sequence alignment problem in biology[END_REF]; [START_REF] Just | Computational complexity of multiple sequence alignment with SP-score[END_REF].

The progressive and iterative approaches for averaging a set of time series are mostly derived from the multiple sequence alignment methods to address the challenging problem of aligning more than two time series [START_REF] Higgins | T-Coffee: A novel method for fast and accurate multiple sequence alignment[END_REF]; [START_REF] Sze | A polynomial time solvable formulation of multiple sequence alignment[END_REF]; [START_REF] Thompson | CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice[END_REF]. To estimate the centroid of more than two time series, several heuristic approaches have been proposed. In the following, we review the main studies related to time series multiple temporal alignments as well as the major progressive and iterative approaches for time series averaging under time warp.

Progressive approaches

The progressive approaches estimate the global centroid c by combining pairwise time series centroids through different strategies. Here, we present the most well-known ones.

NonLinear Alignment and Averaging Filters (NLAAF)

In the past decades, many averaging appraoches have been introduced, but only a few of them have been adapted to time series averaging, mining and clustering. For instance, Gupta et al. in [START_REF] Gupta | Nonlinear alignment and averaging for estimating the evoked potential[END_REF], proposed a time series averaging method based on a tournament scheme, called "NonLinear Alignment and Averaging Filters (NLAAF)". First, pairs of time series are selected randomly, and then aligned according to the DTW. That way, (N/2) 1 averaged sequences are created. The same process is iterated on the estimated centroids , until only one sequence is obtained as a global centroid. In this approach, the averaging method between two time series is applied (N -1) times, as illustrated in Figure 2.4, where c(x i , x j ) refers to the estimated centroid of time series x i and time series x j . 1 N denotes the number of total times series to be averaged.

c(c(c(x 1 ,x 5 ),c(x 3 ,x 2 )),c(c(x 8 ,x 4 ),c(x 6 ,x 7 ))) c(c(x 8 ,x 4 ),c(x 6 ,x 7 )) c(x 6 ,x 7 ) x 7 x 6 c(x 8 ,x 4 ) x 4 x 8 c(c(x 1 ,x 5 ),c(x 3 ,x 2 )) c(x 3 ,x 2 ) x 2 x 3 c(x 1 ,x 5 ) x 5 x 1

Global centroid

Intermediate centroids

Randomize selection pairs

In this method, each element of a centroid is computed as the mean of each linked elements in the DTW alignment. The major drawback of NLAAF lies in the growth of its resulting length, because each use of the averaging method can almost double the length of the average sequence. As classical datasets comprise thousands of time series, with each one including hundreds of data points, simply storing the resulting average may be impossible. This length problem is moreover worsened by the quadratic computational complexity of DTW, that grows bi-linearly with the lengths of the sequences. That is why NLAAF method is generally used in conjunction with a process reducing length of the average, unfortunately leading to information loss and unsatisfactory approximation. Additionally, the average strongly depends on the random selection of time series and different choices lead to different results.

Prioritized Shape Averaging (PSA)

To avoid the bias induced by random selection, Niennattrakul et al. among others [START_REF] Thompson | CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice[END_REF]; [START_REF] Niennattrakul | On clustering multimedia time series data using k-means and dynamic time warping[END_REF]; [START_REF] Niennattrakul | Shape averaging under time warping[END_REF] proposed a framework of shape averaging called "Prioritized Shape Averaging (PSA)" based on the hierarchical clustering. The pairwaise time series centering is guided by the dendrogram obtained through the hierarchical clustering.

The PSA uses hierarchical clustering as a manner to identify priorities between time series. In particular, to estimate the global centroid, the set is first clustered using the agglomerative clustering to get a hierarchical relationship among the whole time series. The simple linkage or complete linkage is considered in general to fasten the dendrogram build, where almost the average linkage or centroids are the best-performed methods in result. Subsequently, the pairwise time series centroids are combined respectively to their clustering order in the dendrogram (see Figure 2.5). Each parent node is averaged in a bottom-up manner using a weighted DTW averaging, called SDTW 2 . So, the most similar time series are averaged first. Note that the weight of an averaged sequence is calculated from the number of the time series upon which the averaged sequence is formed. Initially, all the time series have the same weight of one.

x 1 x 5 x 3 x 2 x 4 x 1 x 5 x 3 x 2 x 4 c(x 1 ,x 5 ) c(x 2 ,x 4 ) c(c(x 1 ,x 5 ),x 3 ) c(c(c(x 1 ,x 5 ),x 3 ),c(x 2 ,x 4 ))
An example of averaging six sample time series using PSA is described in Figure 2.6. According to the dendrogram obtained thorough hierarchical clustering, first the time series x 2 and x 3 are averaged. The average sequence denoted by c(x 2 , x 3 ), has the weight of two. Then, the intermediate centroid c(x 1 , c(x 2 , x 3 )) is computed by averaging the time series x 1 and the average sequence c(x 2 , x 3 ). The intermediate centroid c(x 1 , c(x 2 , x 3 )) will have the weight of three, since the time series x 1 and c(x 2 , x 3 ) have weight of one and two, respectively. The process continues till one obtains a global centroid. Although this hierarchical averaging approach aims to remove the bias induced by random selection, the growth length of the average sequence remains a problem. Additionally, local averaging strategies like NLAAF or PSA may let an initial approximation error propagate throughout the averaging process. If the averaging process has to be repeated (e.g. during k-means clustering iterations), the effects may dramatically alter the quality of the result. This is why a global averaging approach is desirable, where time series would be averaged all together, with no sensitivity to their order of consideration.

x 6 x 5 x 4 x 3 x 2 x 1 c(x 2 ,x 3 ) c(x 1 ,c(x 2 ,x 3 )) c(c(x 1 ,c(x 2 ,x 3 )),x 4 ) c(x 5 ,x

Cross-Words Reference Template (CWRT)

A direct way to estimate the centroid is proposed by Abdulla [START_REF] Abdulla | Cross-words reference template for DTW-based speech recognition sys-tems[END_REF], where a dynamic time warping between each time series and a reference one, generally the time series medoid, is first performed. Each time series is then described in a representation space defined by the reference medoid by resampling, stretching and shortening operations (see Figure 2.7). Finally the global centroid c is computed by averaging the time-aligned time series across each point. The approach is called "Cross-Words Reference Template (CWRT)".

The global estimated centroid c has the same length as the reference time series (e.g. medoid), and the result does not depend on the order in which time series are processed. But the method is a heuristic approach, with no guarantee of optimality.

x 21 x 22

x 23 x 24 x 25 x 26 x 27

x 11 x 12 x 13 x 14 x 15 x 16 x 17

x 31 x 32 x 33 x 34 x 35

x 36 x 37

Iterative approaches

The progressive approaches described above suffer from the early error propagation through the set of pairwise centering combinations. To avoid that, iterative approaches proceed similarly to the progressive ones, but reduce the error propagation by repeatedly refining the centroid and realigning it to the initial time series, until its stabilization [START_REF] Hautamaki | Time-series clustering by approximate prototypes[END_REF]; [START_REF] Petitjean | A global averaging method for dynamic time warping, with applications to clustering[END_REF].

Dtw Barycenter Averaging (DBA)

Petitjean et al. in [START_REF] Petitjean | A global averaging method for dynamic time warping, with applications to clustering[END_REF] proposed a global averaging method, called "Dtw Barycenter Averaging (DBA)". The method consists in iteratively refining an initially average sequence (potentially arbitrary), in order to minimize its distance to the set of time series. The aim is to minimize inertia as the sum of squared DTW distances from the average sequence to the set of time series. Technically, for each refinement iteration DBA works in two steps:

• Computing the DTW between each time series and the temporary average sequence to find the temporal alignments.

• Updating each element of the average sequence with the barycenter of the elements aligned to it during the first step.

Figure 2.8 shows three iterations of DBA on an example with two time series x i and x j , while c is the global estimated centroid.

(iteration 1) (iteration 2) (iteration n)

x i x j x i x j x i x j c
Figure 2.8: DBA iteratively adjusting the average of two time series x i and x j

In a nutshell, the DBA under time warp is a global approach that can average a set of time series all together. The global estimated centroid has the same length as the initial average sequence, and like CWRT, the result is not depending on the order of time series. However the time complexity of DBA approach is smaller than NLAAF and PSA [START_REF] Petitjean | A global averaging method for dynamic time warping, with applications to clustering[END_REF], but the time complexity problem remains.

All the progressive and iterative methods discussed in this chapter are heuristic, with no guarantee of optimality [START_REF] Soheily-Khah | Progressive and Iterative Approaches for Time Series Averaging[END_REF]; [START_REF] Soheily-Khah | A Comparison of Progressive and Iterative Centroid Estimation Approaches Under Time Warp[END_REF]. Lastly, if the provided approximations are accurate for globally similar time series, they are in general poor for time series which share local characteristics with distinctive global behaviors. To circumvent these problems, we proposed in [START_REF] Soheily-Khah | Generalized k-means-based clustering for temporal data under weighted and kernel time warp[END_REF] a tractable fast centroid estimation that captures both global and local temporal features under time warp measures. It formalizes the multiple time series averaging problem as an optimization problem and proposes a solution yielding a local optimum. We will discuss the details of that in next sections.

Conclusion

We have seen that most of the works on averaging a set of time series can be analyzed along two dimensions: first, the way they consider the individual time series for averaging, and second, the way they compute the elements of the resulting average sequences. These two characteristics have proved useful to classify the existing time series averaging techniques. They are also useful to expand new solutions. The main shortcoming of most existing methods is their use of pairwise averaging. To compute the mean of N time series by pairwise averaging, the order of time series influences the quality of the result. This is why a global averaging approach is desirable, with no sensitivity to their order of consideration.

Let us summarize the main characteristics of the described above approaches. In both NLAAF and PSA, the length of the global centroid increases with the number of time series to be averaged. The length of the centroids estimated by CWRT or DBA is however the same as the reference time series length. Furthermore, all the progressive and iterative approaches are heuristic, with no guarantee of optimality. Lastly, even if the provided approximations are accurate for globally similar time series, they are in general poor for time series that share local characteristics with distinctive global behaviors.

In next chapter, we present the k-means clustering approach, the most popular clustering algorithm, and its variations, prior to introduce the proposed fast and accurate centroid estimation approach under time warp to capture both global and local temporal features and its application in clustering. In this chapter, the above questions are addressed through the following points: i) an extension of the standard time warp measures to consider both global and local temporal differences, ii) a tractable and fast estimation of the cluster representatives based on the extended time warp measures, and iii) a generalization of k-means clustering for temporal data under the extended time warp measures.

Clustering: An introduction

Clustering is a fundamental task in data mining and can be considered as one of the most important unsupervised learning problems. There is an increasing interest in the use of clustering methods in machine learning, pattern recognition, image processing and information retrieval. It has also a rich history in other disciplines such as biology, psychology, archaeology, geology and marketing. The goal of clustering is to partition data points into homogeneous groups (or clusters) in such a manner that similar data points are grouped together, while different data points belong to different groups. Formally, the clustering structure is represented as a set of clusters C = {C 1 , ..., C k } of X, the set of data, such that:

k i=1 C i = X
and for all i ̸ = j;

C i C j = ∅
The clustering methods can be classified according to: i) the type of input data to the algorithm, ii) the clustering criteria defining the similarity, dissimilarity or distance between data points, and iii) the theory and fundamental concepts. Consequently many clustering algorithms have been proposed in the literature, each one uses a different scientific discipline. Fraley and Raftery [START_REF] Fraley | How many clusters? Which clustering methods? Answers via model-based cluster analysis[END_REF] suggest dividing the clustering algorithms into two main groups: hierarchical and partitioning. Han and Kamber [START_REF] Han | Data Mining: Concepts and Techniques[END_REF] propose categorizing the algorithms into additional three main categories: density-based, model-based clustering and grid-based methods. Hierarchical methods make the clusters by recursively partitioning the data points in either a top-down or bottom-up manner. For example, in agglomerative hierarchical clustering, each data point initially represents a cluster of its own. Then clusters are merged, according to some similarity measure, until the desired cluster structure is obtained. The result of this clustering method is a dendrogram. Density-based methods assume the data points that belong to each cluster are drawn from a specific probability distribution [START_REF] Banfield | Model-based Gaussian and non-Gaussian clustering[END_REF]. The idea is to continue growing the given cluster as long as the density (or the number of data points) in the neighborhood exceeds some pre-defined threshold. The density-based methods are designed for discovering clusters of arbitrary shape which are not necessarily convex. Model-based clustering methods attempt to optimize the fit between the given data and some mathematical models. These methods find characteristic descriptions for each group. The most frequently used model-based clustering methods are decision trees and neural networks. Grid-based methods partition the space into a finite number of cells that form a grid structure. All of the operations for clustering are performed on the grid structure. The grid-based have the fastest processing time that typically depends on the number of the grid instead of the data points [START_REF] Han | Data Mining: Concepts and Techniques[END_REF].

In this thesis, we mainly focus on partitioning methods, which use an iterative way to create the clusters by moving data points from one cluster to another, based on a distance measure, starting from an initial partitioning. Such clustering methods typically require that the number of clusters k will be pre-set by the user (e.g. elbow method1 , see Figure 3.1).

k-means

k-means algorithm is one of the most popular and simplest unsupervised clustering techniques and it is commonly used in medical imaging, bio-metrics, computer vision, pattern recognition and related fields. Even though k-means was first proposed over 50 years ago [START_REF] Ball | ISODATA: A novel method of data analysis and pattern classification[END_REF]; [START_REF] Macqueen | Some methods for classification and analysis of multivariate observations[END_REF], it is still one of the most widely used algorithms for clustering.

Generally, k-means is a clustering method that aims to find k centroids, one for each cluster, that minimize the sum of distance of each data point from its respective cluster centroid. It solves, for x i ∈ X:

argmin {C 1 ,...,C k } k j=1 x i ∈C j d(x i , c j ) (3.1)
where (C 1 , ..., C k ) are k non-overlapping clusters, c j is the representative of cluster C j , and d is a distance function (e.g., Euclidean, DTW).

Algorithm 1 k-means clustering (X, k)

Input: X = (x 1 , x 2 , ..., x N ) Input: k the number of clusters Output: {c 1 , c 2 , ..., c k } (set of cluster centroids) p = 0 Randomly choose k objects and make them as initial centroids ({c

(0) 1 , c (0) 2 , ..., c (0)
k }) repeat

Assign each data point to the cluster with the nearest centroid p ← p + 1 // Centroid update for j:=1 to k do Update the centroid c (p) j of each cluster using Eq. 3.2 end for until c

(p) j ≈ c (p-1) j ∀j = 1, 2, ..., k Return c 1 , c 2 , ... , c k
The k-means clustering is presented in Algorithm 1. To initialize the k-means, one needs to specify the number of clusters k. The algorithm starts with an initial set of cluster centers (i.e. centroids) {c 1 , c 2 , ..., c k }, chosen at random or according to some heuristic procedure. The clustering algorithm uses an iterative refinement technique. In each iteration, each data point is assigned to its nearest cluster centroid. Then the cluster centroids are re-calculated. The centroid c, the center of each cluster, is calculated as the mean of all the data points belonging to that cluster:

c j = 1 N j x i ∈cluster j x i (3.2)
where N j = |cluster j | is the number of data points belonging to cluster j.

The refinement steps are repeated until the centroids no longer move. The complexity of each iteration of the k-means algorithm performed on N data points is O(k × N ). This linear complexity is one of the reasons for the popularity of the k-means clustering algorithms. Even if the number of data points is substantially large, this algorithm is computationally attractive. Other reasons for the k-means clustering algorithm's popularity are simplicity of implementation and speed of convergence. A proof of the finite convergence of the k-means algorithms is given in [START_REF] Shokri | K-Means-Type Algorithms: A Generalized Convergence Theorem and Characterization of Local Optimality[END_REF]. However, it is shown that under certain conditions the algorithm may fail to converge to a local minimum.

The standard k-means algorithm has some limitations. Firstly, the empty clusters can be obtained if no data points are allocated to a cluster during the assignment step. There exists several strategies to handle the empty clusters; such as choosing a data point from the cluster with the highest within inertia, instead of the centroid of the empty cluster. Secondly, the solution depends heavily on the initial cluster centroids, since the k-means cost function is non-convex. Sometimes different initializations lead to very different final clustering results (see Figure 3.2). A simple but very popular solution for this problem is the use of multiple runs (restarts), where the clusters centroids are randomly placed at different initial positions, hence better local minima can be found. Still we have to decide on the number of restarts and also we are never sure if the initializations are sufficient or not. In addition, the algorithm is sensitive to the noises and outliers. Data points with extremely large values as outliers may substantially distort the distribution of the data. k-medoids clustering which is more robust than k-means in the presence of outliers will be discuss in the following.

k-medoids

k-medoids clustering [KR87] is very similar to the k-means clustering algorithm. The major difference between them is that: while a cluster is represented with its center in the k-means algorithm, it is represented with the most centrally located data point in a cluster (i.e. medoid) in the k-medoids clustering. k-medoids uses representative data points as reference points instead of taking the mean value of the data points in each cluster [START_REF] Zhong | Improved K-Means Clustering Algorithm for Exploring Local Protein Sequence Motifs Representing Common Structural Property[END_REF]. The k-medoid minimizes the sum of distance of each data point from its respective cluster medoid, which solves:

argmin {C 1 ,...,C k } k j=1 x i ∈C j d(x i , m j ) (3.3)
where (C 1 , ..., C k ) are k non-overlapping clusters, m j is the medoid of cluster C j , and d is a distance function. The most common realization of k-medoid is Partitioning Around Medoids (PAM) [START_REF] Kaufman | Finding Groups in Data. An Introduction to Cluster Analysis[END_REF], which is known to be the most powerful. The algorithm has two phases: i) in the first phase, "Build ", a collection of k data points are selected for an initial set of cluster medoids, ii) in the second phase, "Swap", one tries to improve the quality of clustering by exchanging the selected data points with non-selected data points. PAM clustering works effectively for small datasets, but does not scale well for large datasets due to its complexity [START_REF] Ng | Effecient and Effictive Clustering Methods for Spatial Data Mining[END_REF]; [START_REF] Kamber | Spatial clustering methods in data mining: A survey[END_REF]. To deal with very large datasets, a sampling-based method, called CLARA (Clustering LARge Applications) can be used. The CLARA chooses multiple samples of data, applies the PAM clustering on each sample, and returns its best clustering as output [START_REF] Kaufman | Finding Groups in Data. An Introduction to Cluster Analysis[END_REF]. Finally, the choice of clustering algorithm depends both on the type of data and on the particular purpose and application.

Kernel k-means

A major drawback to k-means clustering algorithm is that it does not work for non-linearly separable clusters (see Figure 3.4). Kernel-based clustering techniques address this limitation by introducing a kernel function to capture the non-linear structure in data. Kernel k-means is an extension of standard k-means algorithm that maps data points from the input space to a higher dimensional feature space through a non-linear transformation Φ and minimizes the clustering error in the feature space [START_REF] Lee | Evaluation of the performance of clustering algorithm in kernel-induced feature space[END_REF]. When k-means is applied in this feature space, the linear separators in the feature space correspond to nonlinear separators in the input space. So, separation and clustering may be easier in higher dimensional feature space (see Figure 3.5). The objective function that kernel k-means tries to minimize is the clustering error in feature space, and can be written as:

argmin {C 1 ,...,C k } k j=1 x i ∈C j ∥Φ(x i ) -m j ∥ 2 (3.4)
while (C 1 , ..., C k ) are k clusters in feature space, Φ is a nonlinear transformation, and the centroid or pre-image of cluster C j is denoted by m j , defined as:

m j = 1 |C j | x i ∈C j Φ(x i ) (3.5)
Note that, cluster centroids in the feature space cannot be computed. We don't have explicit knowledge of representations of transformation function Φ(.) and the inverse of Φ to map from feature space back to input space typically does not exist (the pre-image problem). Therefore, on can expand the minimization problem to solve it. If we expand the distance computation ∥Φ(x i )m j ∥2 in the objective function (Eq. 3.4) by using the scalar product and the kernel trick 2 , we obtain the following:

∥Φ(x i ) -m j ∥ 2 = ⎛ ⎜ ⎜ ⎜ ⎝ Φ(x i ).Φ(x i ) - 2 x i ′ ∈C j Φ(x i ).Φ(x i ′ ) |C j | + x i ′ ,x i ′′ ∈C j Φ(x i ′ ).Φ(x i ′′ ) |C j | 2 ⎞ ⎟ ⎟ ⎟ ⎠ (3.6)
Hence only inner products are used in the computation of the Euclidean distance between a data point and a centroid of the cluster in the high-dimensional feature space. As a result, if we are given a kernel matrix K, where K ij = κ(x i , x i ′ ) = Φ(x i ).Φ(x i ′ ), the clustering can be performed without knowing explicit representations of Φ(x i ) and Φ(x i ′ ).

While, kernel-based clustering algorithm can capture the nonlinear structure in data, they do not scale well in term of speed and memory requirements when the number of data points to be clustered exceeds of thousands. It requires computation and storage of N × N kernel matrix in memory, rendering it non-scalable to large dataset. But, what is the largest number of data points, N , that can be handled? Using an intel xeon e7-8837 processor machine, octa-core, 2.8GHz, 4TB max memory, with N greater than one million data points, may take several days to compute the kernel matrix alone. Thus, distributed and approximate versions of kernel k-means clustering algorithms recommended to handle large datasets.

Since a large fraction of attention from the data mining community has focuses on temporal data, of all the techniques applied clustering is perhaps the most frequently used, we are interested in clustering temporal data in this thesis.

Motivation

Clustering time series is an important task in providing useful information in various domain, and usually applied prior to any temporal data analysis or machine learning tasks, to extract groups of time series and highlight the main underlying dynamics. Recently, it has attracted huge concentration in the data mining community [START_REF] Bradley | Refining Initial Points for K-Means Clustering[END_REF]; [WW00]; [START_REF] Kalpakis | Distance Measures for Effective Clustering of ARIMA Time-Series[END_REF]; [START_REF] Keogh | Clustering of Time Series Subsequences is Meaningless: Implications for Previous and Future Research[END_REF]. k-means-based clustering, viz. standard k-means, k-means++, and all its variations, is among the most popular clustering algorithms, because it provides a good trade-off between the quality of the solution obtained and its computational complexity [START_REF] Arthur | K-means++: The Advantages of Careful Seeding[END_REF]. However, the k-means clustering of temporal data under the widely used dynamic time warping (DTW) or the several well-established temporal kernels [START_REF] Bahlmann | Online handwriting recognition with support vector machines-a kernel approach[END_REF]; [START_REF] Cuturi | Fast global alignment kernels[END_REF]; [START_REF] Cuturi | A kernel for time series based on global alignments[END_REF]; [START_REF] Kruskall | The symmetric time warping algorithm: From continuous to discrete[END_REF]; [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF]; [START_REF] Marteau | On Recursive Edit Distance Kernels with Application to Time Series Classification[END_REF] is challenging, because estimating the cluster centroids requires aligning multiple temporal data simultaneously. Such alignments, referred to as multiple sequence alignments [START_REF] Higgins | T-Coffee: A novel method for fast and accurate multiple sequence alignment[END_REF]; [START_REF] Thompson | CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice[END_REF], become computationally prohibitive, costly and impractical when the data size increases. To do so, some progressive and iterative heuristics have been proposed but are limited to the standard DTW and to the temporal data of the same global behavior [ACS03]; [CNH00]; [START_REF] Petitjean | A global averaging method for dynamic time warping, with applications to clustering[END_REF]; [START_REF] Thompson | CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, position-specific gap penalties and weight matrix choice[END_REF]. For temporal data clustering, to bypass the centroid estimation problem, the k-medoids and the kernel k-means [START_REF] Inderjit | Kernel k-means: spectral clustering and normalized cuts[END_REF]; [START_REF] Girolami | Mercer kernel-based clustering in feature space[END_REF] are generally used [START_REF] Liao | Clustering of time series data -A survey[END_REF]. For the k-medoids, a medoid is a good representative of data that has the same global dynamic (or behavior) but inappropriate for capturing local temporal features [START_REF] Frambourg | Learning Multiple Temporal Matching for Time Series Classification[END_REF]. For the kernel k-means, although efficient for non-linearly separable clusters since centroids cannot be estimated in the Hilbert space, it refers to pairwise comparisons for the cluster assignment step. While k-means-based clustering, of linear complexity, remains a fast algorithm, the k-medoids and the kernel k-means clustering have a quadratic complexity due to the pairwise comparisons involved.

In this chapter, we propose a fast and accurate approach that generalizes the k-means based clustering algorithm for temporal data based on i) an extension of the standard time warp measures to consider both global and local temporal differences and ii) a tractable and fast estimation of the cluster representatives based on the extended time warp measures. Temporal data centroid estimation is formalized as a non-convex quadratic constrained optimization problem. The developed solutions allow to estimate not only the temporal data centroid but also its weighting vector, which indicates the representativeness of the centroid elements. The solutions are particularly studied on four time warp measures that are commonly used in practice: the standard dynamic time warping [START_REF] Kruskall | The symmetric time warping algorithm: From continuous to discrete[END_REF], the dynamic temporal alignment kernel [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF], the Gaussian dynamic time warping kernel [START_REF] Bahlmann | Online handwriting recognition with support vector machines-a kernel approach[END_REF], and the global alignment kernel [START_REF] Cuturi | A kernel for time series based on global alignments[END_REF]; [START_REF] Cuturi | Fast global alignment kernels[END_REF].

Problem statement

The k-means algorithm aims at providing a partition of a set of data points in distinct clusters such that the inertia3 within each cluster is minimized. The k-means clustering algorithm was originally developed with the Euclidean distance, the representative of each cluster being defined as the center of gravity of the cluster. This algorithm can, however, be generalized to arbitrary dissimilarities (resp. similarities) by replacing the representative update step with an explicit optimization problem that yields, for each cluster, the representative that minimizes (resp. maximizes) the total dissimilarity (resp. similarity) to all the data points of that cluster [START_REF] Shokri | K-Means-Type Algorithms: A Generalized Convergence Theorem and Characterization of Local Optimality[END_REF]; [START_REF] Hastie | The elements of statistical learning: data mining, inference and prediction[END_REF]; [START_REF] Bock Hans-Hermann | Origins and extensions of the k-means algorithm in cluster analysis[END_REF]. Here we focus on (dis)similarity measures between time series that are based on time alignments because such measures (which encompass the standard dynamic time warping and its variants/extensions) are the ones most commonly used in the literature.

In the following, X denotes a set of univariate discrete time series x = (x 1 , ..., x T ) of assumed length T 4 . An alignment π of length |π| between two time series is defined as the sequence of couples of aligned elements of the two time series. Intuitively, an alignment π between two time series describes a way to associate each element of one time series to one or more elements of another time series and vice-versa. Such an alignment can be simply represented by a path in the T × T grid, where the monotonicity conditions ensure that the alignment path is neither going back nor jumping. We will denote A as the set of all the alignments between two time series. For measures based on the time warp alignments, the integration of a weighting vector allows one to differently weigh the different time stamps of the series under consideration. This notion is used for the centroid (or representative) of each cluster c = (c 1 , ..., c T ), which is no longer a simple data point (i.e. a simple time series) but rather a time series with an associated weighting vector w = (w 1 , ..., w T ). The role of the weight vector is to indicate, for each cluster centroid (or representative), the importance of each time stamp, while this importance varying from cluster to cluster. Let d be a dissimilarity measure defined on R T ×(R T ) 2 , such that d provides a measure of dissimilarity between a given time series x and a weighted centroid (c, w). The generalized k-means clustering algorithm aims to find a partition of the data points in k clusters (C 1 , ..., C k ) such that the intra-cluster dissimilarity is minimized. The associated minimization problem for x ∈ X can be written as:

argmin {C 1 ,...,C k } k i=1 x∈C i d(x, (c i , w i )) (3.7)
where (c i , w i ) ∈ R T × R T , the weighted representative of cluster C i , is defined by:

⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ (c i , w i ) = argmin c , w x∈C i d(x, (c, w))
subject to :

T t=1 w t = 1, w t > 0, ∀t (3.8)
where the given constraints guarantee that the problem is not degenerate. The generalized k-means clustering algorithm, given in Algorithm 2, is a direct generalization of the standard k-means algorithm in which the representative (or centroid) update step has been replaced with a more general representative update step that does not rely on the Euclidean distance.

Algorithm 2 Generalized k-means clustering (X, k)

1: Input: X, k 2: Output: {C 1 , ..., C k } 3: p = 0 4: Select k representatives {c (p) 1 , ..., c (p)
k }, either randomly or through strategies à la k-means++ [AV07] 5: repeat 6:

Cluster assignment C

(p) i ← {x ∈ X | c (p) i = argmin cj ,1≤j≤k d(x, (c (p) j , w (p) j ))}, 1 ≤ i ≤ k 7:
Representative update (c

(p+1) i , w (p+1) i ) ← h(C (p) i )
where the function h : P(X) → R T × R T satisfies:

x∈C (p) i d(x, h(C (p) i )) ≤ x∈C (p) i d(x, (c (p) i , w (p) i ))
(3.9)

8: p ← p + 1 9: until all clusters {C 1 , ..., C k } are stable 10: Return {C 1 , ..., C k }
The function h (line 7 of Algorithm2), referred to as the representative function, provides, from a set of time series points, a point in R T with an associated weighting vector (h : P(X) → R T × R T ). As soon as the function h satisfies the condition expressed in Inequality (3.9), the intra-cluster dissimilarity decreases in the representative update step. This decrease is maximum when the function h is defined using Eq. (3.8)5 . The intra-cluster dissimilarity also decreases in the cluster assignment step, which is identical to the one used in the standard k-means clustering algorithm. Thus the generalized k-means clustering algorithm provided above converges [START_REF] Shokri | K-Means-Type Algorithms: A Generalized Convergence Theorem and Characterization of Local Optimality[END_REF]; [START_REF] Hastie | The elements of statistical learning: data mining, inference and prediction[END_REF]; [START_REF] Bock Hans-Hermann | Origins and extensions of the k-means algorithm in cluster analysis[END_REF]. The above formulation can of course be adapted to the similarity measures, by substituting the dissimilarity d with a similarity s and arg min with arg max in Eqs. (3.7) and (3.8), and by considering representative estimation functions h satisfying:

x∈C (p) i s(x, h(C (p) i )) ≥ x∈C (p) i s(x, (c (p) i , w (p) i ))
This condition indicates that the solution obtained at iteration p, h(C (p) i ), should be better than the solution which obtained at the previous iteration, (c

(p) i , w (p) i ).
In the following, we first give the definition of the extended time warp measures, and then under these metrics we formalize the problem of centroid estimation in the representative update step.

Extended time warp measures

We focus in this thesis on four time warp measures that are commonly used in practice: the Dynamic Time Warping (DTW) [START_REF] Kruskall | The symmetric time warping algorithm: From continuous to discrete[END_REF], which is a dissimilarity measure, the Dynamic Temporal Alignment Kernel (κ DTAK ) [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF], the Gaussian Dynamic Time Warping (κ GDTW ) [START_REF] Bahlmann | Online handwriting recognition with support vector machines-a kernel approach[END_REF], and the Global Alignment kernel (κ GA ) [Cut+07]; [START_REF] Cuturi | Fast global alignment kernels[END_REF], which are three similarity measures and constitute a reference in kernel machines in several domains, such as computer vision, speech recognition or machine learning [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF]; [START_REF] Bailer | Sequence kernels for clustering and visualizing near duplicate video segments[END_REF]; [START_REF] Zhou | Hierarchical aligned cluster analysis for temporal clustering of human motion[END_REF]. To take into account both global and local temporal differences, we propose an extension of the four measures with a weighted centroid, as discussed above. The extensions mainly introduce a weighted warping function that guides the learned alignments according to the representative (or centroid) elements importance to capture both global and local temporal features [START_REF] Soheily-Khah | Generalized k-means-based clustering for temporal data under weighted and kernel time warp[END_REF].

Weighted dynamic time warping (WDTW)

Definition 3.1 The Weighted Dynamic Time Warping (WDTW) between the time series x and the weighted time series (c, w) is defined by:

WDTW(x, (c, w)) = min π∈A 1 |π| (t ′ ,t)∈π f (w t ) ϕ(x t ′ , c t ) C(π) (3.10) = C(π * )
where A (as before) is the set of all alignments possible between two time series, f : (0, 1] → R + is a non-increasing function and ϕ : R × R → R + is a positive, real-valued, dissimilarity function.

The cost function C computes the sum of the weighted dissimilarities ϕ between x and (c, w) through the alignment π. When the weights are uniform (or f is a constant function) and when ϕ is the Euclidean distance, Eq. 3.10 corresponds to the well known Dynamic Time Warping (DTW) [START_REF] Kruskall | The symmetric time warping algorithm: From continuous to discrete[END_REF]; [START_REF] Sankoff | Time warps, string edits, and macromolecules: the theory and practice of sequence comparison[END_REF]. The above definition thus generalizes the DTW to the case where the different instants are weighted. The fact that f is non-increasing guarantees that the most important instants (i.e. the instants with the higher weights) of the centroid c should be privileged in the optimal alignment that minimizes the cost function C. Lastly, the optimal alignment π * is obtained through the same dynamic programming procedure as the one used for the standard dynamic time warping. and w = (w 1 , ..., w 7 ) in the 7 × 7 grid. The value of each cell is the weighted divergence f (w t ) ϕ t ′ t = f (w t ) ϕ(x it ′ , c t ) between the aligned elements x t ′ and c t . The optimal path π * (the green one) that minimizes the average weighted divergence is given by π 1 = (1, 2, 2, 3, 4, 5, 6, 7) and π 2 = (1, 2, 3, 4, 4, 5, 6, 7).

The complexity of computing the extended time warp measure WDTW is O(T 2 ), because the optimal alignment is obtained by dynamic programming6 . However, one can speed it up by considering instead of A, a subset of alignments usually around the diagonal [SC71]; [START_REF] Itakura | Minimum prediction residual principle applied to speech recognition[END_REF]; [START_REF] Sakoe | Dynamic programming algorithm optimization for spoken word recognition[END_REF].

Weighted dynamic temporal alignment kernel (WK DTAK )

More recently, several DTW-based temporal kernels that allow one to process time series with kernel machines have been introduced. The most common is the Dynamic Temporal Alignment Kernel (κ DTAK ) [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF]. Here, we propose an extension of the pseudo-positive definite kernel κ DTAK .

Definition 3.2

The Weighted Dynamic Temporal Alignment Kernel (WK DTAK ) between the time series x and the weighted time series (c, w) is defined by:

WK DTAK (x, (c, w)) = max π∈A 1 |π| (t ′ ,t)∈π f (w t ) κ(x t ′ , c t ) Cκ(π) (3.11) = C κ (π * )
where f : (0, 1] → R + is a non-decreasing function and κ is a Gaussian kernel (e -1 2σ 2 ∥x-y∥ 2 ) with an associated free parameter σ corresponding to the standard deviation which used to measure the similarity between aligned elements.

When the weights are uniform (or the warping function f is a constant function) and when κ is the standard Gaussian kernel, the Eq. 3.11 corresponds to the Dynamic Time Alignment Kernel (κ DTAK ) introduced in [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF]. The non-decreasing property of f ensures that the most important instants are privilaged in the optimal alignment that maximizes the cost function C κ . The complexity of WK DTAK , similar to the WDTW, is quadratic in T , since the optimal alignment is again obtained by dynamic programming.

Weighted Gaussian dynamic time warping (WK GDTW ) Definition 3.3

The Weighted Gaussian Dynamic Time Warping Kernel (WK GDTW ) between the time series x and the weighted time series (c, w) is defined by:

WK GDTW (x, (c, w)) = max π∈A exp ⎡ ⎣ -1 λ 1 |π| (t ′ ,t)∈π f (w t ) ϕ(x t ′ , c t ) ⎤ ⎦ Cκ(π) (3.12) = C κ (π * )
where f : (0, 1] → R + is a non-increasing function and ϕ : R × R → R + is an Euclidean norm.

The cost function C k computes the exponential of the weighted dynamic time warping distance between x and (c, w) through the alignment π. When the weights are uniform (or f is a constant function), Eq. 3.11 corresponds to the Gaussian Dynamic Time Warping kernel (κ GDTW ) introduced in [START_REF] Bahlmann | Online handwriting recognition with support vector machines-a kernel approach[END_REF].

As before, the optimal alignment π * is obtained through the same dynamic programming procedure as the one used for the standard dynamic time warping.

Weighted kernel global alignment (WK GA )

Similarly, the extension of the Global Alignment kernel [START_REF] Cuturi | A kernel for time series based on global alignments[END_REF], which defines a true positive definite kernel, under fair conditions, on the basis of all of the alignments π ∈ A, can be defined as follows.

Definition 3.4

The Weighted Global Alignment Kernel (WK GA ) between the time series x and the weighted time series (c, w) is defined by:

WK GA (x, (c, w)) = π∈A (t ′ ,t)∈π f (w t ) κ(x t ′ , c t ) Cκ(π) (3.13) = C κ (π * )
where f : (0, 1] → R + is a non-decreasing function and κ is a local similarity kernel induced from the divergence ϕ as κ = e -λϕ .

The centroid c, with its weight vector w, of a set of time series X corresponds to the weighted time series that minimizes (or maximizes) a dissimilarity (or similarity) measure with all the elements of X. Considering respectively WDTW, WK DTAK , WK GDTW and WK GA metrics defined above as dissimilarity and similarity measures, we obtain following optimization problems, solutions of which correspond to the weighted centroids: where π * i corresponds to the optimal alignment (i.e. the alignment of minimal cost for WDTW, or maximal cost for WK DTAK , WK GDTW or WK GA ) between the time series x i and the weighted time series (c, w), obtained through dynamic programming, subject to: T t=1 w t = 1 and w t > 0, ∀t Depending on the choice of f and ϕ, the restriction of g to w with c and Π * fixed, denoted here g(w/c, Π * ) with Π * = (π * 1 , ..., π * N ), as well as the restriction of g to c with w and Π * fixed (denoted g(c/w, Π * )), are not necessarily pseudo-convex7 . This makes the WDTW problem difficult to solve and even ill-posed. The standard alternating procedure, that would amount here for the WDTW problem to alternately minimizing g(w/c, Π * ) and g(c/w, Π * ), is not guaranteed to converge and to find a local minimum. The same is true for the WK DT AK , WK GDT W or WK GA problems, for which the pseudo-concavity of the functions g κ (w/c, Π * ) and g κ (c/w, Π * ) depends on the choice of f and κ. We establish here a theorem, the proof of which is given in Appendix A, that provides sufficient conditions on f and ϕ for the functions g(w/c, Π * ) and g(c/w, Π * ) to be pseudo-convex.

Theorem 3.1 Let f : [0, 1] → R + be the non-increasing function used in g, and let ϕ : R × R → R + be the positive, real-valued function used in g. Let us furthermore define:

Condition 1: ∀a ∈ [-1, 1], ∀x ∈ [0, 1], ∞ k=2 a k k! f (k) (x) ≥ 0 Condition 2: ∀(c, c ′ ) ∈ R 2 , ∀x ∈ R, ∞ k=2 (c ′ -c) k k! ∂ k ϕ(x,c) ∂c k ≥ 0 Then: If Condition 1 holds, then g(w/c, Π * ) is pseudo-convex; If Condition 2 holds, then g(c/w, Π * ) is pseudo-convex.
A similar theorem, also proven in Appendix A, can be established for f and κ, with conditions on their form so that the functions g κ (w/c, Π * ) and g κ (c/w, Π * ) are pseudoconcave.

Theorem 3.2 Let f : [0, 1] → R + be the non-decreasing function used in g κ , let κ be the positive definite symmetric kernel used in g κ and let Conditions 3 and 4 be defined as:

Condition 3: ∀a ∈ [-1, 1], ∀x ∈ [0, 1], ∞ k=2 a k k! f (k) (x) ≤ 0 Condition 4: ∀(c, c ′ ) ∈ R 2 , ∀x ∈ R, ∞ k=2 (c ′ -c) k k! ∂ k κ(x,c) ∂c k ≤ 0 Then: If Condition 3 holds, then g κ (w/c, Π * ) is pseudo-concave; If Condition 4 holds, then g κ (c/w, Π * ) is pseudo-concave.
These theorems allow one to choose functions that guarantee that the alternating approach to solve the WDTW, WK DTAK , WK GDTW or WK GA problems is well-defined and converges towards local minima (for WDTW) and maxima (for WK DTAK , WK GDTW or WK GA ). We now proceed to the presentation of this alternating approach as well as of the functions we have retained.

Centroid estimation formalization and solution

We describe here the general strategy followed to estimate the representatives (or centroids) of a cluster of data points (X), prior to study the solution. This strategy leads to the four extended time warp measures introduced above.

Representative update through alternative optimization

The problem given in Eq. 3.8 (and its maximization counterpart for similarity measures) can be solved by alternatively minimizing (or maximizing) for c and w, with the other one being fixed, and by recomputing the optimal alignment warping path when necessary. This strategy is based on the following steps:

1. Set initial values for c and w and compute for these values Π * = {π * x / x ∈ X}, where Π * denotes the optimal alignments warping path between c and the time series in X. Note that for WK GA , computing Π * is not needed because this kernel is based on the sum of all alignments. Steps 2, 3 and 4 lead to a decrease (resp. increase) in the overall dissimilarity (resp. similarity) between the weighted centroid (c, w) considered at each iteration and the data points in X. Indeed, each step respectively aims at finding values of c, w and π that minimize (resp. maximize) the objective function, while the other quantities being fixed. Thus, the above strategy defines an algorithm that converges, and the associated function, which from X produces (c, w), is a valid representative function. The minimum (resp. maximum) in steps 2 and 3 can be obtained by computing the partial derivatives w.r.t c and w and either setting those partial derivatives to 0 and solving for c and w or using gradient descent approaches. Depending on the convexity (resp. concavity) properties of the measures, the minimum (resp. maximum) in step 2 and 3 may be local. This is not a problem per se because the procedure still converges while decreasing (resp. increasing) the inertia between the current pair (c, w) and the data points in X. Now, we study the application of this strategy to each of the time warp measures retained by proposing the corresponding solutions.

Compute: argmin

Solution for WDTW

The extended form of the dynamic time warping measure is defined as:

WDTW(x, (c, w)) = min π∈A ⎛ ⎝ 1 |π| (t ′ ,t)∈π f (w t ) ϕ(x t ′ , c t ) ⎞ ⎠ (3.14)
in which f (w t ) is a non increasing warping function. The constraints, T t=1 w t = 1 and w t > 0, in the WDTW problem are based on quasi-concave functions (respectively an hyperplan and identity functions). Thus provided, that the functions f and ϕ in g satisfy Conditions 1 and 2 of Theorem 3.1, the Karush-Kuhn-Tucker conditions hold and the optimization problems in steps 2 and 3 can be solved with the method of Lagrange multipliers. At each step (2, 3 and 4), the function g is thus non-increasing and the overall procedure converges to a local minimum of g. For the function ϕ, we use here the standard Euclidean distance. In this case, Condition 2 of Theorem 3.1 amounts to

∞ k=2 (c ′ -c) k k! ∂ k ϕ(x, c) ∂c k = (c -c ′ ) 2 > 0
as all the partial derivatives of order greater than 2 are null. ϕ thus satisfies condition 2. For the function f , we consider here two different families:

f (x) = x -α e -αx
where α ∈ R + controls the influence of the weighting scheme. The negative exponent guarantees that the most important instants (i.e., those highly weighted) are privileged in the optimal alignment that minimizes the WDTW. For α = 0, Eq. 3.14 leads to the standard DTW.

The dominating term, for sufficiently smooth functions, in

∞ k=2 a k k! f (k) (x)
used in Condition 1 of Theorem 3.1 is obtained with the second order derivative (as a ∈ [-1, 1]). Condition 1 thus amounts, for sufficiently smooth functions, to a convexity condition that is satisfied by the two smooth families of functions above. The solution of step 2 can thus be obtained by equating the partial derivative of g(c/w, Π * ) with respect to c to 0 and solving for c. The following proposition, the proof of which is given in Appendix B, provides the form of the solution obtained by doing so.

Given w = (w 1 , ..., w T ) and Π * = {π *

x / x ∈ X}, the function defined in Eq. 3.14 is convex in c and the centroid (or representative) c that minimizes the sum of within cluster dissimilarities is obtained by solving the partial derivative equation, leading to, ∀t, 1 ≤ t ≤ T :

c t = x∈X 1 |π * x | (t ′ ,t)∈π * x x t ′ x∈X |N (t, x)| |π * x | (3.15)
where π * x denotes the optimal alignment path for x ∈ X and |N (t, x)| = {t ′ / (t ′ , t) ∈ π * x } denotes the number of time instants of x aligned to time t of centroid c.

For step 3, the solution for weight estimation of WDTW is obtained by equating the partial derivative of the Lagrangian of g(w/c, Π * ), subject to T t=1 w t = 1 and w t > 0, ∀t, with respect to w to 0 and solving for w. The following propositions, also proven in Appendix B, provide the form of solutions obtained. Note that closed form solutions are obtained here for both step 2 and 3. Given c = (c 1 , ..., c T ) and Π * = {π *

x / x ∈ X}, the weight vector w that minimizes g(w/c, Π * ), ∀t, 1 ≤ t ≤ T is given by: for f (x) = x -α ;

w t = A 1 1+α t T t=1 A 1 1+α t for f (x) = e -αx ; w t = 1 α log ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ A t ( T t=1 A t ) 1/T ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ + 1 T with A t = x∈X 1 |π * x | (t ′ ,t)∈π * x (x t ′ -c t ) 2 (3.16)
The solution to Eq. 3.8, corresponding to the representative update step of the generalized k-means clustering algorithm (Algorithm 2), is thus finally obtained through Algorithm 3.

Algorithm 3 Centroid wdtw 1: Input: X = {x 1 , ..., x N } 2: Output: (c,w) 3: Initialization:

• c (0) randomly selected from X

• w (0) = ( 1 T , ..., 1 T ) • Π * (0) : optimal alignments between X and (c (0) , w (0) )

• p = 0 4: repeat 5: p ← p + 1 6: Update c (p) t using Eq. 3.15, 1 ≤ t ≤ T 7: Update w (p) t using Eq. 3.16, 1 ≤ t ≤ T 8:
Update Π * (p) : optimal alignments between X and (c (p) , w (p) ) 9: until (c (p) , w (p) ) ≈ (c (p-1) , w (p-1) ) 10: Return (c,w)

Solution for WK DTAK

The extension of the κ DTAK is defined as:

WK DTAK (x, (c,w)) = max π∈A ⎛ ⎝ 1 |π| (t ′ ,t)∈π f (w t ) κ(x t ′ , c t ) ⎞ ⎠ (3.17)
where a Gaussian kernel (e -1 2σ 2 ∥x-y∥ 2

) with an associated free parameter σ corresponding to the standard deviation is used to measure the similarity between aligned elements. Similarly, the constraints, T t=1 w t = 1 and w t > 0, ∀t, provided that the functions f and κ in g κ satisfy Conditions 3 and 4 of Theorem 3.2, the Karush-Kuhn-Tucker conditions hold and the optimization problems in steps 2 and 3 can be solved with the method of Lagrange multipliers. At each step (2, 3 and 4), the function g κ is thus non-decreasing and the overall procedure converges to a local maximum of g κ . We consider here a simple gaussian kernel for κ, that is at the basis of several DTW kernels, as in [START_REF] Shimodaira | Dynamic time-alignment kernel in support vector machine[END_REF]; [START_REF] Cuturi | A kernel for time series based on global alignments[END_REF]; [START_REF] Zhou | Aligned cluster analysis for temporal segmentation of human motion[END_REF]; [START_REF] Cuturi | Fast global alignment kernels[END_REF]; [START_REF] Zhou | Generalized time warping for multimodal alignment of human motion[END_REF].

For the function f , we consider the two families of functions, which can be seen as the counterpart of the functions used for WDTW:

f (x) = x α α < 1 log(αx) α > 0
In this case, f is a non-decreasing warping function that ensures that the most important instants are privileged in the optimal alignment. The standard κ DTAK formulation is obtained with α = 0. Note also that Eq. 3.17 can be viewed as a pre-image problem [START_REF] Honeine | Preimage problem in kernel-based machine learning[END_REF] that aims at estimating the pre-image, in the input space, of the barycenter of the data points in the feature space induced by the WK DTAK kernel.

As before, the dominating term, for sufficiently smooth functions, in ∞ k=2 a k k! f (k) (x) used in Condition 3 of Theorem 3.2 is obtained with the second order derivative (as a ∈ [-1, 1]).

Condition 3 thus amounts, for sufficiently smooth functions, to a concavity condition that is satisfied by the two smooth families of functions above.

Given w = (w 1 , ..., w T ) and Π * = {π *

x / x ∈ X}, the function defined in Eq. 3.17 is concave in c. The solution of step 2 can be thus be obtained by equating the partial derivative of g κ (c/w, Π * ) with respect to c to 0 and solving to c. We however know of no closed-form solution in this case, and resort to a stochastics gradient ascent method based on the following update rules at iteration p:

c (p+1) t = c (p) t + η (p) ∂L ∂c (p) t and
η (p+1) = η (p) p (η (0) = 1)
with the partial derivative equation ∀t, 1 ≤ t ≤ T :

∂L ∂c t = x∈X 1 |π * x | (t ′ ,t)∈π * x f (w t ) (x t ′ -c t ) σ 2 e (-(x t ′ -c t ) 2 2σ 2 ) (3.18)
where L is given by Eq. 3.17.

For the weights estimation, for step 3, the function defined in Eq. 3.17 is concave in w as α ∈ [0; 1]. The solution of the partial derivative equation, obtained by equating the partial derivative of the Lagrangian, integrating the constraints on w to 0, can easily be obtained (proof is given in Appendix C). Given c = (c 1 , ..., c T ) and Π * = {π *

x / x ∈ X}, the weight vector w that maximizes the sum of intra-cluster similarities subject to T t=1 w t = 1 and w t > 0, ∀t, is defined by: for f (x) = x α ;

w t = A 1 1-α t T t=1 A 1 1-α t for f (x) = log (αx); w t = A t T t=1 A t with A t = x∈X 1 |π * x | (t ′ ,t)∈π * x e (- (x t ′ -c t ) 2 2σ 2 ) (3.19)
Algorithm 4 Centroid wk dtak 1: Input: X = {x 1 , ..., x N } 2: Output: (c,w) 3: Initialization:

• c (0) randomly selected from X

• w (0) = ( 1 T , ..., 1 T ) • Π * (0) : optimal alignments between X and (c (0) , w (0) ) Update Π * (p) : optimal alignments between X and (c (p) , w (p) ) 14: until (c (p) , w (p) ) ≈ (c (p-1) , w (p-1) ) 15: Return (c,w) Algorithm 4 summarizes the steps required to solve the centroid estimation problem (Eq. 3.8) for WK DTAK , corresponding to the representative update step of the generalized k-means algorithm. The gradient ascent steps for estimating centroid c correspond to lines 4 to 11.

• p = 0 4: repeat 5: p ← p + 1 6: // Update c (p) t , 1 ≤ t ≤ T : 7: q = 0, η = 1 8: repeat 9: q ← (q + 1), η ← η q 10: c (q+1) ← c (q) + η ∂L ∂c ,

Solution for WK GDTW

The extended form of the Gaussian dynamic time warping kernel is defined as:

WK DTAK (x, (c, w)) = max π∈A ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ exp ⎡ ⎣ -1 λ 1 |π * | (t ′ ,t)∈π * f (w t ) (x t ′ -c t ) 2 ⎤ ⎦ A(w,c) ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ (3.20)
Given w = (w 1 , ..., w T ) and Π * = {π *

x / x ∈ X}, considering f (w t ) a non-increasing function, the solution can thus be obtained with the method of Lagrange multipliers. The solutions for both centroid and weight estimation lead however to no closed form solutions, and resort to a stochastic gradient ascent method based on the following update rules at iteration p, for the centroid c:

c (p+1) t = c (p) t + η (p) ∂L ∂c (p) t and η (p+1) = η (p) p (η (0) = 1)
with

∂L ∂c t = -2 λ ⎡ ⎣ x∈X 1 |π * x | ⎛ ⎝ t ′′ ∈T,(t ′′ ,t)∈π * x f (w t )(c t -x t ′′ ) ⎞ ⎠ . A(w, c) ⎤ ⎦ (3.21)
and for the weights vector w:

w (p+1) t = w (p) t + η (p) ∂L w ∂c (p) t and η (p+1) = η (p) p (η (0) = 1)
with

∂L w ∂w t = α λ w -(α+1) t ⎡ ⎣ x∈X 1 |π * x | ⎛ ⎝ t ′′ ∈T,(t ′′ ,t)∈π * x (x t ′′ -c t ) 2 ⎞ ⎠ . A(w, c) ⎤ ⎦ for f (x) = x -α , and
∂L w ∂w t = α λ e (-α wt) ⎡ ⎣ x∈X 1 |π * x | ⎛ ⎝ t ′′ ∈T,(t ′′ ,t)∈π * x (x t ′′ -c t ) 2 ⎞ ⎠ . A(w, c) ⎤ ⎦ (3.22) for f (x) = e -αx .
The solutions, for step 2 and 3, are obtained by equating the partial derivative of the Lagrangian (integrating the constraints) of g κ (c/w, Π * ) with respect to c to 0 and solving for c and g κ (w/c, Π * ) with respect to w to 0 and solving for w. The proofs are given in Appendix D. Algorithm 5 summarizes the steps to solve the centroid estimation problem (Eq. 3.8) for WK GDTW . Both centroid and weight estimated by a gradiant ascent method.

Algorithm 5 Centroid wk gdtw 1: Input: X = {x 1 , ..., x N } 2: Output: (c,w) 3: Initialization:

• c (0) randomly selected from X

• w (0) = ( 1 T , ..., 1 T ) • Π * (0) : optimal alignments between X and (c (0) , w (0) )

• p = 0 4: repeat 5: p ← p + 1 6: // Centroid update 7: c (p) t ← Gradient-Ascent(c (p-1) t , w (p-1) t , Π * (p-1) ), 1 ≤ t ≤ T 8: // Weight update 9: w (p) t ← Gradient-Ascent(c (p) t , w (p-1) t , Π * (p-1) ), 1 ≤ t ≤ T 10: // Alignment update 11:
Π * (p) : optimal alignments between X and (c (p) , w (p) ) 12: until g κ (c (p) , w (p) ) ≈ g κ (c (p-1) , w (p-1) ) 13: Return (c,w)

The case of WK GA

The extension of the κ GA , for f (w t ) = w α t , is defined as:

WK GA (x, (c,w)) = π∈A (t ′ ,t)∈π w α t e -λΦ(x t ′ ,ct) (3.23) with Φ(x t ′ , c t ) = 1 2σ 2 ∥x t ′ -c t ∥ 2 + log(2 -e -1 2σ 2 ∥x t ′ -ct∥ 2 )
where, similar to the WK DTAK , a Gaussian kernel is used as the similarity measure between aligned elements. α plays the same role as before and lies again in interval [0; 1]. The parameter λ is used to attenuate the diagonal dominance problem which arises mainly for time series of significantly different lengths. Although diagonal dominance leads to the positive definite Gram matrices, it reduces the performance of the measure in practice. Lastly, the complexity of WK GA is similar to that of WDTW, WK DTAK and WK GDTW as one can rely on a recurrence formula to estimate the cost of a given pair (t ′ , t) on the basis of the costs of the three previous pairs (t ′ -1, t), (t ′ -1, t -1) and (t ′ , t -1) [START_REF] Cuturi | A kernel for time series based on global alignments[END_REF]; [START_REF] Cuturi | Fast global alignment kernels[END_REF]. As for WK DTAK and WK GDTW the time computation can be reduced by considering only a subset of alignments; for WK GA , a fastest version exists based on a pairing of Gaussian and triangular kernels [START_REF] Cuturi | Fast global alignment kernels[END_REF].

Given w = (w 1 , ..., w T ), the partial derivative of WK GA for updating the centroid c is, ∀t, 1 ≤ t ≤ T :

∂L ∂c t = x∈X π∈A ⎛ ⎝ (t ′ ,t)∈π - (x t ′ -c t ) σ 2 (2 -e -1 2σ 2 (x t ′ -ct) 2 ) ⎞ ⎠ × (t ′ ,t)∈π w α t e -λΦ(x t ′ , ct) (3.24)
where Φ(x t ′ , c t ) and L correspond to Eq. 3.23. We know that there is no closed-form solution to the equation ∂L ∂ct = 0. More importantly, the computation of the above derivative given in Eq. 3.24 cannot benefit from the recurrence formulas of WK GA , which ensure an efficient computation of the measure (see proofs in Appendix E). Thus, if we want to use gradient ascend methods, we need to compute the scores associated with each alignment and sum them, which is impractical in situations where T is large (more than a few tens). For this reason, we will not use WK GA in generalized k-means.

In the remainder of this thesis, for the experiments, we thus focus on the following extended (dis)similarity measures for time series: 1) WDTW as defined by Eq. 3.14, 2) WK DTAK as defined by Eq. 3.17, and 3) WK GDTW as defined by Eq. 3.20.

Conclusion

The k-means-based clustering is among the most popular clustering algorithms, because it provides a good trade-off between the quality of solution obtained and its computational complexity. However, the k-means clustering of temporal data under the widely used dynamic time warping or temporal kernels is challenging, because estimating cluster centroids requires aligning multiple temporal data simultaneously. Such alignments, become computationally prohibitive, costly and impractical when the data size increases. For temporal data clustering, to bypass the centroid estimation problem, the k-medoids and the kernel k-means are generally used. But, while k-means-based clustering of linear complexity, remains a fast algorithm, the k-medoids and kernel k-means clustering have a quadratic complexity due to the pairwise comparisons involved.

Here, we introduced a generalized centroid-based clustering algorithm for temporal data under time warp measures. For this, we proposed an extension of the common time warp measures to capture both global and local differences, and a tractable, fast and efficient estimation of the cluster centroids, under the extended time warp measures. This cluster centroid estimation is formalized as a quadratic constrained optimization problem. Finally, the developed solutions allow for estimating not only the temporal data centroid but also its weight vector, which indicates the importance of centroid's elements. The solutions generalize the averaging problem to time series of distinct behaviors that share local characteristics, through several warping functions.

Chapter 4

Experimental study In this chapter, we first describe the datasets retained to conduct our experiments to evaluate the effectiveness of the proposed approach, then specify the validation process, prior to present and discuss about the results obtained. Later, we mainly compare the proposed generalized k-means algorithm based on the extended time warp measures, to the alternative k-means clustering, standard k-medoids and kernel k-means clustering approaches. A wide range of datasets is used to evaluate the efficiency of the proposed approach, in the context of classification. The quantitative evaluation is finally completed by a visualization and qualitative comparison of the centroids obtained.

Data description

Our experiments are conducted on 24 public datasets1 . The classes composing the datasets are known in advance and define the ground truth partitions. The experimentations are first carried out on standard well-known datasets, which define a favorable case for the averaging and clustering task as time series behave similarly within classes (e.g. cbf, cc and gunpoint), as illustrated in Figure 4.1. We then consider more complex datasets (e.g. bme, umd and consseason). The second type of datasets is more challenging: a) it is related to time series that have different global behaviors within classes, while still sharing local characteristics, b) the commonly shared characteristics may appear at varying time stamps within the classes, and c) be hidden by noise (see Figure 4.2). For example, bme includes two challenging classes begin and end characterized by a small bell arising at the initial and final period respectively. The overall behavior may be different depending on whether the large bell is pointing upward or downward. umd introduces more complexity with the classes up and down characterized by a small bell that may occur at different time stamps. To characterize these datasets, and because clustering methods such as k-means algorithms are known to be more efficient when the clusters are isotropic (i.e. spherical) and well isolated, for each dataset we computed the Bartlett's test of sphericity2 [CL71]; [START_REF] Donald | Stopping rules in principal components analysis: a comparison of heuristical and statistical approaches[END_REF] as well as the isolation ratio. The isotropy is measured as the average of Bartlett's p-values for each cluster, as the p-value corresponds to the probability of being spherical. The cluster isolation ratio corresponds to the ratio of the sum of the dynamic time warping dissimilarities within clusters to the total DTW dissimilarities. The lower isolation ratio indicates the more isolated clusters. Table 4.1 describes the datasets considered, with their main characteristics: number of clusters (k), number of time series as dataset size (N ), time series length (T ), p-values for the isotropy and an isolation ratio. In particular, Bartlett's p-values ≤ 0.15 indicative of non-isotropic datasets, and isolation ratio ≥ 0.2, indicative of non-well-isolated datasets, are in bold. For instance, twopatterns is composed of clusters that are globally non-spherical but reasonably well-isolated, whereas yoga clusters are relatively spherical but not well-isolated. As one can note from the Bartlett's p-values and isolation ratios displayed in Table 4.1, most of the datasets considered are composed of clusters that are either non-isotropic or not well-isolated (or both) and are thus challenging for the k-means and kernel k-means clustering algorithms [START_REF] Girolami | Mercer kernel-based clustering in feature space[END_REF]. In addition, to visualize the underlying structure of the datasets, we performed a multidimensional scaling 3 [CC01] on the DTW pairwise dissimilarity matrix. Figure 4.3 displays the cluster structures on the first plan; with a stress4 lower than 20%, while the representations obtained can be considered as accurate images of the underlying structures. Finally, as one can note, the datasets retained have very diverse structures (and shapes) and many of them are not linearly separable.

Validation process

We compare here the proposed generalized k-means-based clustering algorithm (denoted as Gk-means) with the k-means based on the four well-known averaging methods (i.e., NLAAF, PSA, CWRT and DBA), k-medoids (k-med) under Euclidean and the standard DTW, and kernel k-means5 (Kk-means) based on Euclidean, the standard κ GDTW , κ DTAK , κ GA and κ TGA . We focus here on the above approaches because they constitute the most frequently used variants of the k-means for time series clustering. As mentioned earlier, the k-medoids and the kernel k-means clustering algorithms are mainly used to bypass the centroid estimation process required by k-means with time warp measures. For our comparison, we rely on the Rand index 6 [Ran71], which measures the agreement between the obtained clusters and the ground truth ones, to evaluate each method. The Rand index lies in [0, 1] and the higher index, the better the agreement is. In particular, the maximal value of Rand index, 1, is reached when the clusters (or partitions) are identical. For all the methods, the parameters taken within a validation set through a standard line/grid search process, as described in Table 4.2, and the parameters retained are the ones that maximize the Rand index on the validation set. Finally, the results reported hereafter are averaged after 10 repetitions of the corresponding algorithm. Because κ GDTW and κ DTAK are not strictly definite positive kernels [START_REF] Cuturi | A kernel for time series based on global alignments[END_REF]; [START_REF] Cuturi | Fast global alignment kernels[END_REF], we systematically added a scaled identity matrix δ I to the Gram matrices with negative eigenvalues, while δ is the absolute value of the smallest negative eigenvalue and

Method Metric Line / Grid values Kernel k gdtw t ∈ {0.2, 0.5, 1, 2, 5} ⊙ med(dtw(x, y)) k-means k dtak σ ∈ {0.2, 0.5, 1, 2, 5} ⊙ med(∥x -y∥) k ga σ ∈ {0.2, 0.5, 1, 2, 5} ⊙ med(∥x -y∥) , λ ∈ {10 -4 , 10 -3 , 10 -2 , 10 -1 , 1, 10} k tga σ ∈ {0.2, 0.5, 1, 2, 5} ⊙ med(∥x -y∥) , c ∈ {0.25, 0.5} ⊙ med(∥x∥) Generelized wdtw α ∈ {10 -2 , 10 -1 , 10 0 } k-means wk gdtw α ∈ {10 -2 , 10 -1 , 10 0 } , t ∈ {0.2, 0.5, 1, 2, 5} ⊙ med(dtw(x, y)) wk dtak α ∈ {10 -2 , 10 -1 , 10 0 } , σ ∈ {0.2, 0.5, 1, 2, 5} ⊙ med(∥x -y∥)
Table 4.2: Parameter Line / Grid: The ⊙ multiplication operator is element-wise (e.g. {1, 2, 3} ⊙ med = {med, 2med, 3med}). med(x) stands for the empirical median of x evaluated on the validation set and x and y are vectors sampled randomly within time series in the validation set.

I is the identity matrix. For κ GA method, the estimated value of the regularization parameter λ ≈ 1 demonstrates that the problem of diagonal dominance has not been encountered for the retained datasets.

Experimental results

Results for time series clustering

In the context of clustering, the Rand indices and time consumptions for each method, on each dataset, are reported in Tables 4.3 and 4.4, respectively. Results in bold correspond to the best values, while the statistically significantly lower or higher ones (two-sided t-test at 5% risk),respectively, are shown with the down or up arrows.

Based on the Rand indices displayed in Table 4.3, one can note that Gk-means clustering method with WDTW leads to the best clustering results overall (23 datasets out of 24), followed by Gk-means method with WK DT AK . Furthermore, the clustering approaches based on Gk-means introduced here outperform k-medoids and kernel k-means clustering almost on all datasets. Table 4.4 shows that Gk-means with WDTW leads almost to the fastest clustering method (18 datasets out of 24).

Impact of isotropy and isolation on clustering quality

Two most challenging datasets are apparently umd and bme, with almost the lowest Rand indices varying in [0.52, 0.62] for k-means, k-medoids and kernel k-means and in [0.60, 0.69] for Gk-means clustering methods. In addition, umd and bme are composed of weakly isolated and weakly spherical classes (see Table 4.1), include time series of distinct global behaviors within clusters [START_REF] Frambourg | Learning Multiple Temporal Matching for Time Series Classification[END_REF]. Figure 4.4 shows the time series behaviors within the three classes begin, middle and end of dataset bme, and Figure 4.5 illustrates the time series behaviors within the three classes up, middle and down of dataset umd. Thus they require a similarity or dissimilarity measure able to capture the commonly shared features within clusters. Gk-means clustering method, by relying on weighted centroids, can capture the characteristics shared locally within the clusters and improves the performances over the standard clustering methods, which fail to capture such local features. The second category of challenging datasets consists of coffee, ecg200 ,gunpoint, lighting2 and yoga and is composed of weakly isolated but highly (or mildly) spherical clusters (see Table 4 From the above results, we can obtain that all the k-means clustering, k-medoids and kernel k-means methods provide relatively good results even when the clusters are non-spherical (or non-isotropic); however, they are limited when the clusters are non-isolated. By contrast, Gk-means clustering with WDTW obtains very good clustering performances on all datasets, particularly on challenging datasets composed of clusters that are both non-spherical and not well-isolated.

Relationships between methods and datasets

To compare globally the different clustering methods, we rely on a multiple correspondence analysis (MCA7 ), to analyze the 17 methods (considered as individuals) and the 24 datasets (considered as categorical variables).

Multiple correspondence analysis is a data analysis technique for nominal categorical data and can be viewed as an extension of correspondence analysis (CA) which allows one to analyze the pattern of relationships of several categorical dependent variables. It can also incorporate quantitative variables. MCA is concerned with relationships within a set of variables, which usually are homogeneous, and allows the direct representation of individuals as points in geometric space. To do so, each method is described by a vector ("-", "+","++",...), with as many dimensions as there are datasets, in which the modalities "++", "+" and "-" indicate whether the Rand index of a method on a dataset is respectively highly greater, greater or lower than the mean Rand index obtained for that dataset over all the methods. Distinct groups of methods, corresponding to distinct ways to perform on the different datasets, can be distinguished. The first group is defined by the k-means using NLAAF, k-medoids and kernel k-means (with Euclidean) and is opposed to the other methods as it yields the lowest performances (corresponding to modality "-") particularly on cc, fish and symboles. The second group is defined by the k-means using PSA and the three kernel k-means approaches (κ GA , κ TGA and κ GDTWsc ) that yield the lowest performances on beef, bme, medicalimages, oliveoil, trace, swedishleaf and twopatterns.

To understand this limitation, we have studied the cluster structure induced by the dissimilarities based on the κ DTAK and κ GA Gram matrices. Lastly, the last group includes the Gk-means approaches which yield the best performances (corresponding to modality "++") almost for all datasets. We can in particular observe that Gk-means with WDTW outperforms the other methods particularly on lighting2, lighting7, medicalimages, osuleaf, and yoga.

The blue arrow in Figure 4.7 is known as the Guttman effect [START_REF] Le | Geometric data analysis: from correspondence analysis to structured data analysis[END_REF]; it starts from the upper left corner, that contains mainly datasets of modalities "-", then goes through the lower left corner, containing mixed modalities ("-" and "+"), through the lower right corner, again containing mixed modalities ("+" and "++") , and finally ends in the upper right corner containing datasets of modalities "++". The blue arrow suggests a global ranking of the methods based on their performance on all the datasets: Gk-means with WDTW is the best method, followed by Gk-means with WK DTAK , followed by Gk-means with WK GDTW , k-means using DBA, κ DTAK and k-medoids (with DTW), which yield similar results, and finally followed by the rest (see Figure 4.7). 4.4 shows the time consumptions of the different clustering methods. As a reminder, results in bold correspond to the fastest methods with lowest time consumption, while the significantly higher ones are shown with the red up arrows. As one can note, Gk-means clustering is the fastest approach, with a approximately constant time requirement for all of the methods on small datasets (e.g. beef, coffee, oliveoil). For large datasets (e.g. cbf, cc, medicalimages, swedishleaf, symboles, yoga), both k-medoids and kernel k-means clustering are significantly slower mainly due to the involved pairwise comparisons, whereas the Gk-means clustering, particularly with WDTW, remains remarkably fast. The standard k-means based on NLAAF and PSA averaging methods are highly time-consuming, largely because of the progressive increase of the centroid length during the pairwise combination process, which makes these approaches unusable for large time series datasets.

Figure 4.8 shows the multiple correspondence analysis of clustering time consumption to compare globally the clustering time consumption of different methods. The blue arrow starts from the lower left corner, that contains mainly the fastest methods, then goes through the upper left corner and finally ends in the lower right corner containing the lowest clustering methods. The first group is defined by the k-means using CWRT, the standard Gaussian kernel k-means and k-medoids (with Euclidean), while all of these methods use Euclidean and therefore are not comparable with the other methods under time warp. Note that in clustering time consumption, Table 4.4, we ignored these three methods. The Gk-means approaches yield the fastest clustering results overall (18 datasets out of 24), followed by k-means with DBA. The third group is defined by kernel k-means (κ DTAK , κ GDTW , κ GA and κ TGA ) and k-medoids (with DTW). Lastly, as mentioned, k-means with NLAAF and PSA with huge difference yield the lowest performances. The obtained empirical results are in agreement with the theoretical complexities of each method. Let N be the number of time series to be clustered in k clusters, with the time series length T , and p the complexity of the measure between aligned instants (p = 5 for WDTW, p = 8 for WK DTAK , p = 9 for WK GDTW , and p = 18 for WK GA -see Eqs. 3.14, 3.17, 3.20 and 3.23). Therefore, the complexity of Gk-means clustering is O(pN kT 2 ), the one of k-medoids is O(pN (Nk)T 2 ) and the complexity of kernel k-means clustering is O(pN 2 T 2 ), showing that kernel k-means clustering is slower than k-medoids, which is in turn slower than the Gk-means clustering. Lastly, the difference between Gk-means clustering with WDTW and WK DTAK (or WK GDTW ) is mainly due to the ascent gradient part, while the difference between kernel k-means clustering with κ DTAK and κ GA is explained by the different values of p for each method. (c i , w i ) = Centroid Estimation(X train |Y train == i) using Algo. 3, 4 or 5 8: end for 9: for ∀x j ∈ X test do 10: Classify x j by nearest weighted centroid using Eq. 3.10, 3.11 or 3.12 11: end for 12: Return class labels for ∀ x ∈ X test For the experiment, we here rely on the classification time consumption, presented in Table 4.7, while the main contribution is to speed up temporal classification. To verify the accuracy of the method, we have confidence in the classification error rate (see Table 4.5 and 4.6). The lower the error rate, the better the agreement is. Lastly, we consider a goodness criteria which is the ratio of "fastness" to "accuracy". For this, we define "fastness" as the normalized time consumption (so that all values are in (0,1)), and (1 -error rate) is the "accuracy". The lower value is the better classifier. Similar to the previous tables in clustering section, results in bold correspond to the best values, while the statistically significantly lower or higher ones (t-test at 5% risk) are represented by the down or up arrows, respectively.

Classification error rate

Based on the classification error rate displayed in Table 4.5, for k = 1, nearest neighbor classification methods lead to the best results in overall (23 datasets out of 24), followed by the proposed weighted centroid estimation methods. Using a small number of weighted centroids per class obtained by generalized k-means based clustering within each class, we can beat the nearest neighbor classifier. For example, for dataset beef, if we assume k=2 and we consider two weighted centroids per class, by generalized k-means based clustering, we outperform all the nearest neighbors classifiers by error rate = 0.433. A second example, for dataset osuleaf, if we assume k=2, generalized k-means based clustering in each class, obtains two weighted centroids per class. Thus, the classification error rate using WDTW will be 0.295, which outperforms all the other classifiers.

The experiments show that the proposed nearest weighted centroid classification methods significantly outperform the alternative nearest centroid classification approaches. To validate the claim, Table 4.6 presents all the nearest centroid classification error rate results. As one can note, classification using WDTW yields the best classification error rate overall (19 datasets out of 24), followed by WK DTAK .

To compare globally the different nearest centroid classification methods, as before, we rely on MCA9 . To do so, each method is described by a vector ("-", "+","++"), in which the modalities "++", "+" and "-" indicate whether the classification error rate of a method on a dataset is respectively highly greater, greater or lower than the mean error rate obtained for that dataset over all the methods. Lower error rate indicated better classification. According to the Figure 4.9, nearest centroid classification using PSA and NLAAF (to estimate the centroid) have the higher error rates almost for all the datasets, as they yield the lowest classification performances (corresponding to modality "++"). The group includes the the proposed nearest weighted centroid classification approaches which yield the best performances (corresponding to modality "-") almost for all datasets (23 datasets out of 24), followed by nearest neighbor classification using DBA. We can in particular observe that nearest weighted centroid classifier based on WDTW outperforms the other methods particularly on bme, consseason, fish, symboles, trace, and yoga.

The blue arrow suggests a global ranking of the methods based on their performance on all the datasets: it starts from upper left corner, that contains mainly datasets of modalities "-", then goes through the lower left corner, containing mixed modalities ("-" and "+"), through the lower right corner, containing mixed modalities ("+" and "++") , and finally ends in the upper right corner, again containing mixed modalities ("+" and "++"). The projection of the methods on the blue arrow suggest a global ranking: classification with nearest weighted centroid using WDTW is the best method, followed by WK DTAK , DBA and WK GDTW and finally followed by the rest.

Classification goodness criteria

Lastly, for a global comparison considering both classification error rate and run time index, we define here a goodness criteria as a ratio of "fastness" to "accuracy". The "fastness" (normalized time consumption), should be lower and the "accuracy" (1 -error rate), should be higher.

classification goodness = fastness accuracy = normalized run time 1error rate Table 4.8 presents the classification average ranking on different methods, according to the defined goodness criteria. Note that, we ignore the k-NN classification based on Euclidean and kernel RBF, whereas they are not comparable with the classification methods which use time warp measures. The best performance belongs to the proposed nearest weighted centroid classifier with WDTW with the average rank of 1.81, followed by the nearest centroid classifier using using DBA, with the average rank of 2.48, followed by the nearest weighted centroid classifier using WK DTAK . The classifiers which use nearest centroid with NLAAF or PSA, and nearest neighbor using kernels κ GDTW or κ GA yield the lowest performances according to the global average ranking.

A closer look at the centroids

We finally visualize here some of the centroids obtained by the different methods and compare their shape to the one of the time series they represent. Figures 4.14 to 4.19 display the centroids obtained by the different methods. The first line shows some randomly selected time series of the averaged sets, the second line the centroids provided by the alternative methods and the last line the centroids obtained by WDTW, WK DTAK and WK GDTW . Finally, the estimated weights are shown below the corresponding centroids. To assess the weights obtained, we rely on the following observation: a uniform weight distribution is expected for time series that behave similarly within classes, as the centroid's elements are equivalently shared within the class, whereas a non-uniform distribution is awaited for time series that behave differently within classes, as higher weights are assigned to shared centroid's elements. Thus, it is possible to assess the relevance of the estimated weights by examining in which cases their distribution is uniform, and in which cases it is not. To do so, we perform for each times series collection a X 2 -test (at 5% risk) measuring the distance between the estimated weights and a uniform model. A uniform weights distribution is accepted for a p-value greater than 0.05 (e.g. cc-"Cyclic"), otherwise rejected (e.g. umd-"Down").

As one can note, for standard datasets, almost all the approaches succeed in obtaining centroids more or less similar to the initial time series (given at the first line). For instance, in Figures 4.15 for the classes "Cyclic" of the cc dataset. However, we observe generally less representative centroids for NLAAF and PSA with a drastically large centroid's length of about 10 3 elements vs. 10 2 for the other methods. For umd data, Figure 4.16 shows that the WDTW, WK DTAK and WK GDTW provide the most representative centroid with the highest weights successfully assigned to the begin and end regions that characterize the umd dataset.

To have a closer look at the ability of the proposed method in centroid estimation, here we justify our claim on two complex and completely noisy datasets. The spiral data, proposed in [START_REF] Zhou | Canonical time warping for alignment of human behavior[END_REF], consists of 3-D spatio-temporal time series (2-D spatial and 1-D temporal) generated from latent time series:

x i = U T i (z + b i 1 T l )M i e T i ∈ R 3 * n i
where the canonical time series z ∈ R 2 * l is a curve in two dimensions (x, y). U i ∈ R 2 * 2 and b i ∈ R 2 are randomly generated projection matrix and translation vector respectively. The binary matrix M i ∈ {0, 1} l * n i is generated by randomly choosing n i ≤ l columns from I l for temporal distortion. The spatial dimension e i ∈ R n i is generated with zero-mean Gaussian noise. The latent time series z and three sample generated time series are visualized in Figure 4.11. According to Figure 4.17, on spiral data, only WDTW, WK DTAK and WK GDTW succeed in estimating centroids with an exact spiral shape, as well as providing an effective uniform weighting, as spiral relies on time series of a similar global behavior. For the more complex spiral2 collection, Figure 4.18 shows the ability of the proposed methods (WDTW, WK DTAK and WK GDTW ), to circumvent the noise problem and to reveal the locally shared signature. For example, in Figure 4.19, we can see clearly for the dimension x of spiral2 that (a) the initial hidden signature is well revealed by the estimated centroids, (b) the noise is removed, and (c) the corresponding region is highly weighted as expected.

Discussion

From Table 4.3, we can see that the generalized k-means clustering leads to the highest Rand indexes, where the best scores (indicated in bold) are reached by WDTW for almost all datasets. To have a global view of clustering results, Figure 4.22 shows the average Rand index values of all datasets, for each clustering method. The black color illustrates the best and white color presents the worst one. From Table 4.4, the results reveal WDTW the fastest methods under time warp and k-means using PSA the slowest one. The Generalized k-means using WK DTAK and WK GDTW appears slightly slower that WDTW; this is due to the gradient ascent research part used in the former. To evaluate the relevance of the estimated weights, we performed a X 2 -test10 (at 5% risk, the significance level of the hypothesis test) measuring the distance between the estimated weights and a uniform model. The p-value is the probability of observing a test statistic more extreme than the observed value under the null hypothesis. Table 4.9 allows one to corroborate the characteristics of some sample datasets for WDTW, WK DTAK and WK GDTW : the uniform hypothesis is accepted for all standard datasets as well as for spiral dataset (i.e. on all the datasets in which time series behave similarly within the classes) as p-values are greater than 0.05, and strongly rejected for the complex datasets like bme, umd and spiral2, in which time series behave differently with in the classes, as p-values are lower than 0.05.

Regarding the classification results (

Conclusion and perspectives

The DTW and kernel-DTW are among the most frequently used metrics for time series in several domains as signal processing, data mining or machine learning. However, for clustering time series, approaches are generally limited to the k-medoids or kernel k-means to circumvent centroid estimation under DTW and the tricky multiple temporal alignment problem. This work introduces a generalized k-means-based clustering algorithm for temporal data under extended time warp measures. For this, we propose i) an extension of the common time warp measures and ii) an accurate, fast and efficient solution to the problem of time series averaging, under the extended time warp measures, that captures local temporal features. Additionally, we introduce an extension of time series averaging to kernel-DTW metric as it constitutes a crucial metric for many kernel approaches; Yields a fast method to compute the centroid of time series. Accordingly, we generalize the averaging problem to time series of distinct behaviors that share local characteristics by estimating a weight vector, through several warping function. The integration of a weight vector permits us to weigh time stamps differently and indicate the importance of each time stamp. Lastly, we use the proposed centroid estimation and generalized k-means clustering approach in the classification context. In this way, to speed up the nearest neighbor classification, we suggest a fast accurate nearest weighted centroid classifier.

The efficiency of the proposed centroid estimation solution for classification and clustering algorithms is analyzed on a wide range of public standard and challenging datasets, which are non-spherical, non-well-isolated and (/or) linearly non-separable. The experimental validation is based on standard datasets in which time series share similar behaviors within classes, as well as on more complex datasets exhibiting time series that share only local characteristics, that are multidimensional and noisy.

The approaches are compared through several criteria. For quantitative evaluation, we consider a Rand index, a goodness clustering criterion, and the space and time requirements, as well as the classification error rate as a classification goodness. On the other hand, the qualitative one consists in the visualization and comparison of the centroids obtained by different methods, and the relevance of the weights estimated. In addition, we propose here to investigate the performance of the different methods on the different datasets through a multiple correspondence analysis that highlights the relationships between the methods and the datasets, indicates which method performs well or badly on which datasets, and suggests a global ranking of the different methods. Both the quantitative evaluation and the qualitative one demonstrate the effectiveness of the proposed approaches to being faster and more accurate, which outperform the other methods, in the context of averaging, classification and clustering.

In particular, the approaches introduced in this study have the following characteristics:

• The centroids obtained are more representative of the set, than the centroids obtained by other methods;

• The weights associated to the centroids reflect the elements that are shared locally in a set of time series;

• The estimated centroids yield better clusters in the k-means-based clustering than the centroids obtained by other methods;

• The time and space requirements are lower in comparison with the time and space requirements of other methods.

Perspectives

Explore new temporal metrics and warping functions In this thesis, we focus on four temporal metrics: the dynamic time warping, the dynamic temporal alignment kernel, the Gaussian dynamic time warping and the global alignment kernel. It could be interesting to integrate other temporal metrics as well as different warping functions in our proposed framework to improve the results, specially for very challenging datasets. In this work, we consider only two different warping functions for each temporal metric. Even we are globally best method than the rest, one can explore more warping function, specially for challenging and noisy datasets with low values of Rand index.

Furthermore, in the context of time complexity, although we obtain the best run time in comparison with the alternatives, but still one can investigate for new warping functions, which leads to a suitable and closed-form solution, to circumvent the gradient-based optimization process in weight estimation steps and in result decrease the time consumption.

Lastly, the generalized k-means-based clustering on the WK DTAK and WK GDTW suggests a new promising alternative approach to kernel k-means, that encourages to deepen this issue for further temporal kernels and pre-image problem in future works.

Faster clustering/classification process by temporal segmentation For several large and challenging datasets as character trajectories 11 , using the proposed weight and centroid estimation approach, one can simply distinguish the patterns of interest. Figure 4.23 shows the obtained segments of interest for some sample classes of character trajectories data. Pink color shows the segment of interest (instances with higher weights) for each sample character.

In this case, a second development is to decrease time and space complexities using a temporal segmentation of time series. The idea is to localize segments of interest, according Study the potential of this thesis in other machine learning algorithms The next advancement will be evaluating the extended temporal measures in other algorithms. In this work, we studied time series averaging under time warp for clustering and classification. The next step is to evaluate them in an other classifier such as a decision tree or a support vector machine or extension of the proposed weighted temporal measures in the fuzzy c-means clustering (i.e. soft clustering) context, while in some application (e.g. bio-informatics, image analysis or marketing) each data point can belong to more than one cluster.

Other propositions to define the combined metric In this thesis, for each proposition, we considered only one single metric (e.g. WDTW) to compare time series. The idea of combined metrics can also be the other proposition to gain better results for some multidimensional challenging data sets. Investigate the combination of some behavior-based with value-based metrics, to deal with different dimension of data can be an interesting suggestion for future works. where Π * is the set of all known (and fixed) alignments for the time series x in X. As only positive quantities are involved, A(X, c, Π * , t) ≥ 0. We have: (f (w ′ t )f (w t )) A(X, c, Π * , t)

A Taylor expansion yields:

f (w ′ t ) = f (w t ) + (w ′ t -w t )f ′ (w t ) + ∞ k=2 (w ′ t -w t ) k k! f (k) (w t ) Thus, if ∞ k=2 (w ′ t -w t ) k k! f (k) (w t ) ≥ 0 then f (w ′ t ) -f (w t ) ≥ (w ′ t -w t )f ′ (w t )
and: where Π * is the set of all known (and fixed) alignments for the time series x in X. As only positive quantities are involved, B(X, w, Π * , t) ≥ 0. We have: The proof for Theorem 3.2 directly parallels the ones above, the only difference being in the sign of the inequalities manipulated.

(i) g(c
1 |π * x | (t ′ ,t)∈π * x f (w t ) ϕ(x t ′ , c t )
with:

f (x) =
x -α α > 0 e -αx α > 0 and ϕ the Euclidean distance.

Given w = (w 1 , ..., w T ) and Π * = {π * x / x ∈ X}, the function defined in Eq. 3.14 is convex in c and the centroid c that minimizes the sum of cluster dissimilarities is obtained by solving the partial derivative equation, ∀t, 1 ≤ t ≤ T : The solution for weight estimation of WDTW is obtained by equating the partial derivative of the Lagrangian of g(w/c, Π * ), subject to T t=1 w t = 1 and w t > 0, ∀t, with respect to w to 0 and solving for w, ∀t, 1 ≤ t ≤ T : for f (x) = x -α ;

c t =
w t = A 1 1+α t T t=1 A 1 1+α t (B.2)
for f (x) = e -αx ;

w t = 1 α log ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ A t ( T t=1 A t ) 1/T ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ + 1 T (B.3) with A t = x∈X 1 |π * x | (t ′ ,t)∈π * x (x t ′ -c t ) 2
Proof: As f satisfies Condition 1 of Theorem 3.1 and the constraints are pseudo-convex, the Karush-Kuhn-Tucker conditions are satisfied. The Lagrangian of the above problem is defined as:

L = x∈X 1 |π * x | (t,t ′ )∈π * x f (w t ) ϕ(x t ′ , c t ) + λ(1 - T t=1 w t )
and the solution to the above problem is obtained by solving:

∂L ∂w t = 0
For f (x) = x -α with α > 0, one gets:

∂L
w t = α 1/(α+1) (-λ) 1/α+1 (A t ) 1/(α+1)
Summing over t and equating to 1 (the constraint) yields: From Eq. B.4, we have:

w t = 1 -α [log λ -log (-α A t )]
and considering Eq. B.5:

w t = 1 -α log (-α) + 1 T T t=1 log A t - α T -(log (-α) + log A t ) ⇒ w t = 1 -α 1 T T t=1 log A t - α T -log A t = 1 α log A t - 1 T T t=1 log A t + 1 T
and finally

w t = 1 α log ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ A t ( T t=1 A t ) 1/T ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ + 1 T
Note that, the weight w t > 0 of the element t of c that ensures the inertia minimization is defined as the proportion of the divergence induced by the elements t ′ of the series x ∈ X aligned to t.

Appendix C

Proof of Solutions 3.18 and 3.19

We present here the proof for proposed solutions 3.18 and 3.19. Let us first recall that: ) .

Given c = (c 1 , ..., c T ) and Π * = {π * x / x ∈ X}, the weight vector w that maximizes the sum of intra-cluster similarities g κ (c, w, Π * ), subject to T t=1 w t = 1 and w t > 0, ∀t, is defined by: for f (x) = x α ; A

w t = A 1 1-α t T t=1 A 1 1-α t (C.
( 1 1-α ) t (C.5) From Eq. C.4 and replacing the value for λ in the previous Eq. C.5, we have:

w t = α 1 1-α A 1 1-α t α ( 1 1-α ) T t=1 A ( 1 1-α ) t
and finally yields Eq. C.2: Proof: We can not conclude if A 1 ≤ 0 ⇒ A 2 ≤ 0. Therefore, there is no proof for the pseudo-concavity of g k (w, c) and we can not benefit from the recurrence formulas of WK GA .

w t = A 1 1-α t T t=1 A 1 1-α
g k (c/w, Π * ) = x∈X exp ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ -1 λ 1 |π * x | (t ′ ,t)∈π * x f (w t ) (x t ′ -c t ) 2 F (c) ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ = x∈X exp -1 λ F (c) ⇒ ∂g k (c

Introduction1

  Comparison of time series: state-of-the-art 1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Comparison of time series . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Times series basic metrics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4 Kernels for time series . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 Time series averaging and centroid estimation 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.2 Consensus sequence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.3 Multiple temporal alignments . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 Generalized k-means for time series under time warp measures 3.1 Clustering: An introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3 Problem statement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4 Extended time warp measures . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.5 Centroid estimation formalization and solution . . . . . . . . . . . . . . . . . . 3.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 Experimental study iii iv Contents 4.1 Data description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 4.2 Validation process . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 4.3 Experimental results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 4.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89 Conclusion 93 A Proof of Theorems 3.1 and 3.2 97 B Proof of Solutions 3.15 and 3.16 101

1. 3

 3 .3 Combined (Behavior-value-based) metrics . . . . . . . . . . . . . . . . . . 20 1.4 Kernels for time series . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 1.4.1 Temporal correlation coefficient kernel (κ CORT ) . . . . . . . . . . . . . . . 23 1.4.2 Auto-correlation kernel (κ DACO ) . . . . . . . . . . . . . . . . . . . . . . . 23 1.4.3 Gaussian dynamic time warping kernel (κ GDTW ) . . . . . . . . . . . . . . 24 1.4.4 Dynamic time-alignment kernel (κ DTAK ) . . . . . . . . . . . . . . . . . . . 24 1.4.5 Global alignment kernel (κ GA ) . . . . . . . . . . . . . . . . . . . . . . . . 25 1.4.6 Triangular global alignment kernel (κ TGA ) . . . . . . . . . . . . . . . . . . 25 1.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26

Figure 1 . 2 :

 12 Figure 1.2: dynamic time warping alignment (left) vs. euclidean alignment (right)

Figure 1 . 3 :

 13 Figure 1.3: The optimal alignment path between two sample time series with time warp (left), without time warp (right)

  Constrained DTW Constraints are widely used to speed up dynamic time warping algorithm. The most commonly used ones is the Sakoe-Chiba band (symmetric and asymmetric) [SC71]; [SC78] and Itakura parallelogram [Ita75], which are shown in Figure 1.4.

Figure 1 . 4 :

 14 Figure 1.4: Speed up DTW using constraints: (a)Sakoe-Chiba band (b)Asymmetric Sakoe-Chiba band (c)Itakura parallelogram

Figure 1 . 5 :

 15 Figure 1.5: The constrained optimal alignment path between two sample time series (Sakoe-Chiba band = 5 -left), (Itakura parallelogram -right)

Figure 1 . 6 :

 16 Figure 1.6: The time series sequence x i , and its PAA x * i

Figure 1 . 7 :

 17 Figure 1.7: Alignment between two time series by (a)DTW (b)PDTW

  Figure 1.8 demonstrates an example of multiscale DTW for alignment two time series in the direction of reduce time and space complexity.

Figure 1

 1 Figure 1.8: (a) The optimal alignment warp path on level 2. (b) The optimal alignment warp path with respect to the constraint region R by projecting alignment path to level 1. (c) The optimal alignment warp path with δ = 2

Figure 1 .

 1 Figure 1.10: Comparison of time series: similar in value, opposite in behavior (left), similar in behavior (right)

Figure 1 .

 1 Figure 1.11: Kernel trick: embed data into high dimension space (feature space)

y 1 y 2 y 3 y 4 y 5 y 6 y 7 x 1 x 2 x 3 x 4 x 5 x 6 x 7 a 1 a 2 a 3 a 4 a 5 a 6 a 7 a 8 a 9 iFigure 2 . 1 :

 74921 Figure 2.1: Pairwise temporal warping alignment between x and y (left), the estimation of the centroid c of time series x and time series y in the embedding Euclidean space (right). A dimension i in the embedded Euclidean space identifies a pairwise link a i .

Figure 2 .

 2 Figure 2.2 illustrates the pairwise time series centroid estimation under the dynamic time warping. Note that the length T ≤ m ≤ 2T -1 of the centroid c(x i , x j ) is generally higher than the length of the centered time series.

Figure 2 . 2 :

 22 Figure 2.2: Three possible alignments (paths) are displayed between x i and x j , π * (the green one) being the dynamic time warping one. The centroid c(x i , x j ) = (c 1 , ..., c 8 ) is defined as the average of the linked elements through π * , with for instance c 3 = Avg(x i2 , x j3 ).

Figure 2 . 3 :

 23 Figure 2.3: Medoid as a prototype

Figure 2 . 4 :

 24 Figure 2.4: Centroid estimation by random pairwise centroid combination.

Figure 2 . 5 :

 25 Figure 2.5: Centroid estimation by centroid combination through clustering dendrogram.

Figure 2 . 6 :

 26 Figure 2.6: Example of six time series sequence averaging using PSA

Figure 2 . 7 :

 27 Figure 2.7: Centroid estimation based on a reference time series. The DTW is performed between time series x 1 , x 2 and the reference time series x 3 (left). Time series x 1 and x 2 are embedded in the space defined by x 3 (right) where the global centroid is estimated, and 'avg' is the standard mean function.

elbow

  Number of clusters kTotal with-in clusters sum of squares

Figure 3 . 1 :

 31 Figure 3.1: The elbow method suggests k=3 cluster solutions

Figure 3

 3 Figure 3.2: k-means: different initializations lead to different clustering results

Figure 3 . 3 :

 33 Figure 3.3: Outliers effect: k-means clustering (left) vs. k-medoids clustering (right)

Figure 3 . 4 :

 34 Figure 3.4: Non-linearly separable clusters

Figure 3 . 5 :

 35 Figure 3.5: The kernel trick -complex in low dimension (left), simple in higher dimension (right)

Figure 3

 3 Figure 3.6: k-means clustering (left) vs. kernel k-means clustering (right)

Figure 3 . 7 :

 37 Figure 3.7: Three possible alignments are displayed between x = (x 1 , ..., x 7 ) and c = (c 1 , ..., c 7 )

  g κ (c, w) = argmax c , w N i=1wk dtak (x i , (c, w)) ga (x i , (c, w))

  1 and w t > 0, ∀t 4. Compute Π * for the values of c and w obtained in the previous steps 5. Go back to step 2 until c and w are stable.

Sommaire 4 . 1

 41 Data description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 4.2 Validation process . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 4.3 Experimental results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 4.3.1 Results for time series clustering . . . . . . . . . . . . . . . . . . . . . . . 69 4.3.2 Results for time series classification . . . . . . . . . . . . . . . . . . . . . . 76 4.3.3 A closer look at the centroids . . . . . . . . . . . . . . . . . . . . . . . . . 82 4.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89

Figure 4 . 1 :

 41 Figure 4.1: The time series behaviors within the classes "Funnel", "Cyclic" and "Gun" of the datasets cbf, cc and gunpoint, respectively.

Figure 4 . 2 :

 42 Figure 4.2: The time series behaviors within the classes "Begin", "Down" and "Warm" of the datasets bme, umd and consseason, respectively.

Figure 4 . 3 :

 43 Figure 4.3: Structures underlying datasets

Figure 4 . 4 :

 44 Figure 4.4: The time series behaviors within the three classes of dataset bme

Figure 4 . 5 :

 45 Figure 4.5: The time series behaviors within the three classes of dataset umd

  .1). All the k-means clustering, k-medoids and kernel k-means methods lead to weak clustering results with a Rand index lying in [0.48, 0.59], where the Gk-means clustering (with WDTW) leads to noticeably better clustering performances with a Rand index in [0.60, 0.72]. A third category of datasets is defined by car, cbf and medicalimages is composed of mildly isolated but weakly spherical (or isotropic) clusters. Although k-means clustering, k-medoids and kernel k-means methods yield reasonable clustering performances, with a Rand index in [0.50, 0.74], the Gk-means clustering (with WDTW) yields the best clustering results with a Rand index in [0.77, 0.87]. For the remaining datasets, composed of highly isolated and highly (or mildly) spherical clusters, good clustering results (Rand index in [0.54, 0.97]) are obtained using k-means clustering, k-medoids and kernel k-means methods, but are improved almost significantly by using Gk-means clustering (with WDTW) with a Rand index in [0.75, 0.98].

  Figure 4.6 shows the projection obtained through MDS on the first plan for the dataset cc, with the dissimilarities 8 derived from κ DTAK . As one can see, the red and blue clusters, non linearly separable with DTW (Figure 4.3 -second row and second column), are now well separated. The situation for κ GA directly parallels the one for κ DTAK .

Figure 4

 4 Figure 4.6: CC clusters based on κ DTAK dissimilarities

Figure 4

 4 Figure 4.7: Global comparison of the clustering Rand index

Figure 4

 4 Figure 4.8: Global comparison of the clustering Time consumption

Figure 4 .

 4 Figure 4.10 shows the multiple correspondence analysis of classification time consumption to compare globally the classification time consumption of different methods. The first group is defined by the k-NN standard Gaussian kernel (RBF) and k-NN (with Euclidean), followed by the nearest neighbor classifiers using CWRT, DBA, WDTW, WK DTAK and WK GDTW . As you can see in the figure 4.10, k-NN using kernel κ DTAKsc is next in ranking based on the time consumption, followed by κ TGA , k-NN using DTW sc and kernel κ GDTWsc . The third group is contains the k-nearest neighbor classifiers using DTW and kernels. The last group, the classification approaches which use centroid estimation with NLAAF and PSA yield the lowest performances.

Figure 4 .

 4 Figure 4.10: Global comparison of the classification Time consumption (k=1)

Figure 4 . 11 :

 411 Figure 4.11: Latent curve z and three induced instances x 1 , x 2 , x 3 without noise (left), and with noise e i (right)spiral dataset

  Figure 4.13 shows the progression of their x, y and e dimensions over time.

Figure 4

 4 Figure 4.12: Latent curve z and three induced instances x 1 , x 2 , x 3 sharing local characteristics for the spiral2 dataset

Figure 4 .

 4 Figure 4.13: spiral2: Progression of the 2-D spatial coordinates (x, y) and the noise dimension e over time for the spiral2 dataset

Figure 4 .

 4 Figure 4.14: cbf-"Funnel" centroids: (a) Initial time series, (b) NLAAF, (c) PSA, (d) CWRT, (e) DBA, (f) WDTW, (g) WK DTAK , (h) WK GDTW

Figure 4 .

 4 Figure 4.15: cc-"Cyclic" centroids: (a) Initial time series, (b) NLAAF, (c) PSA, (d) CWRT, (e) DBA, (f) WDTW, (g) WK DTAK , (h) WK GDTW

Figure 4 .

 4 Figure 4.16: umd-"Down" centroids: (a) Initial time series, (b) NLAAF, (c) PSA, (d) CWRT, (e) DBA, (f) WDTW, (g) WK DTAK , (h) WK GDTW

Figure

  Figure 4.17: spiral centroids: (a) Initial time series, (b) NLAAF, (c) PSA, (d) CWRT, (e) DBA,

  Figure 4.17: spiral centroids: (a) Initial time series, (b) NLAAF, (c) PSA, (d) CWRT, (e) DBA, (f) WDTW, (g) WK DTAK , (h) WK GDTW

Figure 4 .

 4 Figure 4.18: spiral2 centroids: (a) Initial time series, (b) NLAAF, (c) PSA, (d) CWRT, (e) DBA, (f) WDTW, (g) WK DTAK , (h) WK GDTW

Figure 4 .

 4 Figure 4.19: spiral2-x centroids: (a) Initial time series, (b) NLAAF, (c) PSA, (d) CWRT, (e) DBA, (f) WDTW, (g) WK DTAK , (h) WK GDTW

Figures 4 .

 4 Figures 4.20 and 4.21 complete these results and show the progressive estimation of the centroid and weights through the first iterations of WDTW and WK DTAK , respectively. These figures show the ability of the proposed methods to remove the noise and to capture the shared signature fastly during the first few iterations of the estimation process. The first row shows the initial centroid (left) and the initial weights (right), which are from uniform distribution.Each row shows one iteration, and finally, the last row, represents the estimated centroid and its weight vector.

Figure 4 .

 4 Figure 4.20: spiral2-x centroids: centroid (left) and weight (right) estimation through the first iterations of WDTW

Figure 4 .

 4 Figure 4.21: spiral2-x centroids: centroid (left) and weight (right) estimation through the first iterations of WK DTAK

Figure 4 . 22 :

 422 Figure 4.22: Comparison of average Rand index values

Figure 4 .

 4 Figure 4.23: Three samples classes of character trajectories-"e", "o" and "p": the ground truth (top), the segments of interest (down)

  Let us first recall Theorem 3.1.Theorem 3.1 Let f : [0, 1] → R + be the non-increasing function used in g, and let ϕ : R×R → R + be the positive, real-valued function used in g. Let us furthermore define:Condition 1: ∀a ∈ [-1, 1], ∀x ∈ [0, 1], ∞ k=2 a k k! f (k) (x) ≥ 0 Condition 2: ∀(c, c ′ ) ∈ R 2 , ∀x ∈ R, ∞ k=2 (c ′ -c) k k! ∂ k ϕ(x,c) ∂c k , then g(w/c, Π * ) is pseudo-convex; If Condition 2 holds, then g(c/w, Π * ) is pseudo-convex. Proof: g(w/c, Π * ) is pseudo-convex iff: ∀ w, w ′ ∇g(w/c, Π * ) (w ′w) ≥ 0 =⇒ g(w ′ /c, Π * ) ≥ g(w/c, Π * ) Let A(X, c, Π * , t) = x∈X 1 |π * x | (t ′ ,t)∈π * x ϕ(x t ′ , c t )

  ,t)∈π * x f (w t )ϕ(x t ′ , c t ), (ii) ∂g(w t /c, Π * ) ∂ wt = f ′ (w t ) A(X, c, Π * , t), (iii) ∇g(w/c, Π * )(w ′w) = T t=1 f ′ (w t )(w ′ tw t ) A(X, c, Π * , t), (iv) g(w ′ /c, Π * )g(w/c, Π * ) = T t=1

  w ′ t )f (w t ))A(X, c, Π * , t) ≥ T t=1 f ′ (w t )(w ′ tw t )A(X, c, Π * , t)that is:g(w ′ /c, Π * )g(w/c, Π * ) ≥ ∇g(w/c, Π * )(w ′w) Thus, if ∇g(w/c, Π * )(w ′w) ≥ 0, then g(w ′ /c, Π * ) ≥ g(w/c, Π * ), which proves the implication from Condition 1.Proof: g(c/w, Π * ) is pseudo-convex iff: ∀ c, c ′ ∇g(c/w, Π * ) (c ′c) ≥ 0 =⇒ g(c ′ /w, Π * ) ≥ g(c/w, Π * ) Let B(X, w, Π * , t) = x∈X 1 |π * x | (t ′ ,t)∈π * x f (w t )

  that is:g(c ′ /w, Π * )g(c/w, Π * ) ≥ ∇g(c/w, Π * )(c ′c)Thus, if ∇g(c/w, Π * )(c ′c) ≥ 0, then g(c ′ /w, Π * ) ≥ g(c/w, Π * ), which proves the implication from Condition 2.

  ) 1/α+1 (A t ) 1/(α+1) = 1 so, T t=1 (A t ) 1/(α+1) . α 1/(α+1) (-λ) 1/α+1 = 1 thus λ = -α ( f (x) = e -αx ; ,t)∈π * x -α e (-αwt) ϕ(x t ′ , c t ) -log (-α A t ) = 1 ⇒ T. log λ -T t=1 log (-α A t ) = -α ⇒ T. log λ = T t=1 log (-α A t )α ⇒ T. log λ = T. log (-α) + T t=1 log A tα thus log λ = log (-α) +

  g κ (c, w) = x∈X 1 |π * x | (t ′ ,t)∈π * x f (w t ) κ(x it ′ , c t )with:f (x) = x α α < 1 log(αx) α > 0and κ the Gaussian kernel:κ(x, c) = e (-∥x-c∥ 2 2σ 2 )(x and c are two real numbers). Given w = (w 1 , ..., w T ) and Π * = {π *x / x ∈ X}, for the centroid c:

  2)for f (x) = log (αx); As f satisfies Condition 3 of Theorem 3.2 and the constraints are pseudo-concave, the Karush-Kuhn-Tucker conditions are satisfied. The Lagrangian of the above problem is defined as:t ′ )∈π * x f (w t ) κ(x it ′ , c t ) + λ(1 -T t=1 w t )and the solution to the above problem is obtained by solving ∂L ∂w = 0.For f (x) = x α (0 ≤ α ≤ 1): over t and equating to 1 (the constraint) leads to:

.t

  Proof of Solutions 3.21 and 3.22Let us first recall that:g κ (c, w) = N i=1 wk gdtw (x i , (c, w)) ,t)∈π * x f (w t ) (x t ′c t ) 2 ⎤ ⎦ A(w,c) with: f (x) = x -α α > 0 e -αx α > 0 Given w = (w 1 , ..., w T ) and Π * = {π * x / x ∈ X}, for the centroid c: ′′ ∈T,(t ′′ ,t)∈π * x f (w t )(c tx t ′′ ) ⎞ ⎠ . A(w, c)

Table 3 .

 3 1.

	k-means	k-medoids
	complexity per iteration O(k × N )	complexity per iteration O(k × (N -k) 2 )
	sensitive to noises and outliers	less sensitive to noises and outliers
	implementation of algorithm is easy implementation of algorithm is complicated
	require to specify k	require to specify k
	comparatively more efficient	comparatively less efficient

Table 3 .

 3

1: Comparison between k-means and k-medoids

  with ∂L ∂c given by Eq. 3.18

	11:	until ∂L ∂c ≈ 0
	12:	Update w	(p) t	using Eq. 3.19, 1 ≤ t ≤ T
	13:

Table 4 .

 4 3: Clustering Rand index

  In this section, we evaluate the relevance of the proposed estimated centroids and weights in a classification context. Among many algorithms to classify temporal data, the k-nearest neighbor with the DTW is performant. However, it is very time consuming as it needs to compare time series query with all training samples. While the idea of data reduction to mitigate this problem has a long history, classification based on nearest centroid could be one solution. To do so, we propose to use k-nearest weighted centroid classifier based on the above extended time warp measures (see Algorithm 6). Using generalized k-means-based clustering within each class, one can obtain a small number of weighted centroids per class.Algorithm 6 Nearest weighted centroid classification (X train ,Y train ,X test ,k) 1: Input: X train = {x 1 , ..., x Ntrain } 2: Input: Y train = labels 3: Input: X test = {x 1 , ..., x Ntest } 4: Input: k = class number 5: Output: class labels for ∀ x ∈ X test 6: for i : 1 to k do

	4.3.2 Results for time series classification
	Kernel k-means				k-medoid		k-means		Generalized k-means
	(Kk-means)				(k-med)	nlaaf	psa	cwrt dba	(Gk-means)
	rbf k gdtw k gdtwsc k dtak k dtaksc	k ga	k tga	e d	dtw	dtwsc	dtw	dtw	dtw	dtw	wdtw wk dtak wk gdtw
	7:

Table 4

 4

	.5: Classification error rate (k=1)

Table 4 .

 4 Complexity and time considerations Table4.7 shows the time consumption for different classification methods. In time comparison, we ignore the classification methods using k-NN with standard Gaussian kernel (RBF) and the k-NN (with Euclidean). While these methods use Euclidean, they are not comparable with the other methods under time warp. As one can note, the nearest centroid classifier using CWRT and DBA are the fastes classification methods, similar to the proposed nearest weighted centroid classification methods. The slight difference is related to the weight computation for each centroid. Whereas we significantly improve the accuracy in the proposed nearest weighted centroid classifiers, one can ignore this non-significant difference in run time. As we mentioned before, the main drawback of NLAAF and PSA averaging methods lie in the growth of their resulting lengths. While nearest centroid classification methods normally should be faster than the nearest neighbor ones, the nearest centroid classifiers using NLAAF or PSA are highly time-consuming (particulary for large datasets), because of the progressive increase of the centroid length during the pairwise combination process.

6: Nearest centroid classification error rate (k=1) Figure 4.9: Global comparison of the nearest centroid classification error rate (k=1)

  Table 4.8), for all datasets, the highest classification goodness scores are obtained by a nearest centroid classifier using WDTW, followed by DBA and WK DTAK . Lastly, the qualitative evaluations (Figures 4.14 to 4.19, consisting in the visualization and comparison of the centroids obtained by different methods) show the effectiveness of the proposed approaches to provide accurate time series averaging for standard and complex datasets.

			p-values
		wdtw	wk DTAK	wk GDTW
	bme	≃ 0.00	≃ 0.00	≃ 0.00
	cbf	1	1	1
	cc	1	1	1
	umd	≃ 0.00	≃ 0.00	≃ 0.00
	spiral	1	1	1
	spiral2	≃ 0.00	≃ 0.00	≃ 0.00
	Table 4.9: p-value (X 2 -test: uniformity distribution test)

  /w, Π * ) = A Taylor expansion yields:ϕ(x t ′ , c ′ t ) = ϕ(x t ′ , c t ) + (c ′ tc t ) ∂ϕ(x t ′ , c t ) ∂c t + ∂ k ϕ(x t ′ , c t ) ∂ k c t ′ , c ′ t )ϕ(x t ′ , c t ) ≥ (c ′ tc t ) ∂ϕ(x t ′ , c t ) ∂c t

	Thus, if then	(ii) (iii)	x∈X ∂ϕ(x t ′ , c t ) 1 |π * x | (t ′ ,t)∈π * x ∂c t B(X, w, Π * , t), f (w t )ϕ(x t ′ , c t ), ∇g(c/w, Π T ∂g(c t /w, Π * ) ∂ ct = t=1 ∂ϕ(x t ′ , c t ) ∂c t (c ′ k=2 (c ′ t -c t ) k k! ∞ k=2 (c ′ t -c t ) k k! ∂ k ϕ(x t ′ , c t ) ∂ k c t ≥ 0 t -c ∞ ϕ(x t

* )(c ′c) = t ) B(X, w, Π * , t), (iv) g(c ′ /w, Π * )g(c/w, Π * ) = T t=1 (ϕ(x t ′ , c ′ t )ϕ(x t ′ , c t )) B(X, w, Π * , t)

and:

T t=1 (ϕ(x t ′ , c ′ t )ϕ(x t ′ , c t ))B(X, w, Π * , t) ≥ T t=1 (c ′ tc t ) ∂ϕ(x t ′ , c t ) ∂c t B(X, w, Π * , t)

  *x denotes the optimal alignment for x ∈ X and |N (t, x)| = {t ′ / (t ′ , t) ∈ π *x } denotes the number of time instants of x aligned to time t of c.Proof: As ϕ satisfies Condition 2 of Theorem 3.1, the Karush-Kuhn-Tucker conditions are satisfied and the solution to the minimization problem above is obtained by solving:

		x∈X	1 |π * x |	(t ′ ,t)∈π * x	c t =	x∈X	1 |π * x |	x (t ′ ,t)∈π *	x t ′
	and to Eq.B.1.
							x∈X	1 |π * x | |N (t, x)| (t ′ ,t)∈π * x	x t ′	(B.1)
								x∈X	|π * x |
	where π ∂g(c/w, Π * ) ∂c	= 0
	Thus	∂g(c/w, Π * ) ∂c	=		x∈X	1 |π * x |	x (t ′ ,t)∈π *	2 f (w t ) (c t -x t ′ ) = 0
	then					1
			x∈X	|π * x |

(t ′ ,t)∈π * x (c tx t ′ ) = 0

which leads to:

  ′′ ∈T,(t ′′ ,t)∈π* x f (w t ) 2 (c tx t ′′ )The solution for the weights vector w, is obtained by equating the partial derivative of the Lagrangian of g κ (c/w, Π * ) with respect to w to 0 and solving for w: (c, w, X, π) . n(π, t)(iii) ∇g k (w).(w ′w) = (w t 2 ) . n(π, t) . k(x t 1 , c t 2 ) and A 2 = x∈X π∈A (t 1 ,t 2 )∈π k(x t 1 , c t 2 ) .

	= ∂g k (c/w, Π * ) x∈X λ |π * -1 x | ∂c t = -2 ⎡ ⎣ λ ⎡ t ⎤ ⎦ .exp ⎣ x∈X 1 |π * x | ⎛ ⎝ t ′′ ∈T,(t ′′ ,t)∈π * x f (w t ) (c t -x t ′′ ) -1 λ and for f (x) = e -αx : So,: ∂g k (w/c, Π * ) ∂w t = α λ e (-αwt) ⎡ ⎣ x∈X 1 |π * x | ⎛ ⎝ t ′′ ∈T,(t ′′ ,t)∈π * x (x t ′′ -c t ) 2 F (c) ⎞ ⎠ . A(w, c) ⎞ ⎠ . A(w, c) ⎤ ⎦ ⎤ ⎦ . = f ′ (w t ) f (w t ) x∈X π∈A t=1 x∈X π∈A f (w t ) (w ′ t -w t ) g k (c, w, X, π) . n(π, t) f ′ (w t ) T t=1 T t=1 g k T f (w ′ t ) n(π,t) -f (w t ) n(π,t)
	.
	(iv) g k (w ′ |c) -g k (w|c)
	=	x∈X	⎡ ⎣	w π∈A (t ′ ,t)∈π (p+1) t f (w ′ = w t ) k(x t ′ , c t ) (p) t + η (p) ∂L w ∂c (p) t ⎤ ⎦ -	and η (p+1) =	η (p) p	(η (0) = 1)
	with
										⎡	⎛	⎤
				∂L w ∂w t	=	α λ	w	-(α+1) t	⎣	x∈X	1 |π * x |	⎝	⎦	(D.2)
	for f (x) = x -α , and
										⎡		⎛	⎞	⎤
				∂L w ∂w t	=	α λ	e (-αwt)	⎣	x∈X	1 |π * x |	⎝ t ′′ ∈T,(t ′′ ,t)∈π * x	t=1 (x t ′′ -c t ) 2 T f (w t ) n(π,t) ⎠ . A(w, c) t ′ ∈π(t) ⎦ k(x t ′ , c t ) (D.3)
	for f (x) = e -αx . Thus,
					⎡						⎤
							T					T
	Proof: ξ =	g x∈X ⎣ f (w ′ t ) n(π,t) x∈X π∈A t=1 ⎛ T x∈X π∈A t=1 t ′ ∈π(t) = ⎝ k(x t ′ , c t ) exp ⎡ ⎣ -1 λ t ′ ∈π(t) ⎞ ⎠ . 1 |π * x |	f (w t ) n(π,t) f (w t ) (x t ′ -c t ) 2 t ′ ∈π(t) ⎤ k(x t ′ , c t ) ⎦ T t=1 A(w,c) t=1 (t ′ ,t)∈π * T x t=1 f (w ′ t ) n(π,t) -f (w t ) n(π,t)	⎦
	for f (x) = x -α : Suppose;
												⎡	⎛	⎞	⎤
		∂g k (w/c, Π * ) ∂w t A 1 = T t=1 f ′ (w t ) = α λ f (w t )	x∈X x∈X |π * -1 λ 1 x | t -w t ) -(α+1) t ⎣ . (w ′ w	∂F (c) ∂c t ⎝ t ′′ ∈T,(t ′′ ,t)∈π * . exp x	λ -1 (x t ′′ -c t ) 2 F (c)	⎠ . A(w, c)	⎦

/w, Π * ) ∂c t = t ′′ ∈T,(t ′′ ,t)∈π * x (x t ′′c t ) 2 ⎞ ⎠ . A(w, c) k (w/c, Π * ) = x∈X ⎡ ⎣ π∈A (t ′ ,t)∈π f (w t ) k(x t ′ , c t ) ⎤ ⎦ = x∈X π∈A ⎡ ⎣ (t ′ ,t)∈π f (w ′ t ) k(x t ′ , c t ) -(t ′ ,t)∈π f (w t ) k(x t ′ , c t ) ⎤ ⎦ ξ

where,

(t ′ ,t)∈π f (w t ) k(x t ′ , c t ) = T t=1 t ′ ∈π(t) f (w t ) k(x t ′ , c t ) = k(x t ′ , c t ) -x∈X π∈A (t 1 ,t 2 )∈π f

One can make this assumption as dynamic time warping can be applied equally on time series of same or different lengths.

κGA is positive definite, if κ/(1 + κ) is positive definite.

The idea is this: run k-means clustering on the dataset for a range of values of k (start with k = 1), and for each value of k compute the inertia as the sum of distances between each member of a cluster and its centroid. Then, plot a line chart of this inertia for each value of k. If the line chart looks like an arm, then the "elbow" on the arm is the value of k that is the best.

∥Φ(xi) -Φ(xj)∥ 2 = κ(xi, xi) -2κ(xi, xj) + κ(xj, xj)

Inertia is defined as the sum of distances between any data point in the cluster and the centroid (or representative) of the cluster

All the results can however be directly extended to multivariate time series, possibly of different lengths, as the temporal alignments, at the core of the (dis)similarity measures we consider, can be defined in those cases as well.

Depending on the dissimilarity measure used, it may not be possible to obtain the point that minimizes Eq.(3.8), and therefore "looser" functions, based on Inequality (3.9), have to be considered.

T is length of each time series.

This is typically the case if ϕ is concave. A pseudo-convex function is a function that behaves like a convex function with respect to finding its local minima, but need not actually be convex.

umd, bme at http://ama.liglab.fr/~douzal/tools.html, the rest of the data at http://www.cs.ucr. edu/~eamonn/time_series_data/

Barspher Matlab function

Stress is defined as the error between the distances induced by the first plan and the original dissimilarities.

kmedoids, kernelkmeans Matlab functions

We have used the MCA function available in the R package FactoMineR.

The dissimilarities are defined as d dtak (i, j) = κDTAK(i, i) + κDTAK(j, j) -2κDTAK(i, j).

MCA function in the R package FactoMineR.

chi2gof function in Matlab 2015b, which returns the p-value of the hypothesis test

Acknowledgments

Proof of Solutions 3.15 and 3.16

We present here the proofs for proposed WDTW solutions 3.15 and 3.16. Let us first recall that:

Proofs for the Case of WK GA

Let us first recall Theorem 3.2. Theorem 3.2 Let f : [0, 1] → R + be the non-decreasing function used in g κ , let κ be the positive definite symmetric kernel used in g κ and let Conditions 3 and 4 be defined as:

We have: