{'fraction_non_alphanumeric': 0.039832285115303984, 'fraction_numerical': 0.01195745011258638, 'mean_word_length': 4.553830674141153, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 7, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Model-based fuzzy time series clustering of conditional higher moments

This paper develops a new time series clustering procedure allowing for heteroskedasticity, non-normality and model's non-linearity. At this aim, we follow a fuzzy approach. Specifically, considering a Dynamic Conditional Score (DCS) model, we propose to cluster time series according to their estimated conditional moments via the Autocorrelation-based fuzzy C-means (A-FCM) algorithm. The DCS parametric modelling is appealing because of its generality and computational feasibility. The usefulness of the proposed procedure is illustrated using an experiment with simulated data and several empirical applications with financial time series assuming both linear and nonlinear models' specification and under several assumptions about time series density function.


Introduction
Time-series data are of interest due to their ubiquity in various areas ranging from science and engineering to finance and economics. While each time series consists of many data points, it can also be seen as a single object. Clustering -i.e., grouping objects with maximum similarity with other objects within the same group -is a useful tool for exploratory analysis as it identifies common structures within data. Moreover, clustering procedures for time series have interesting applications in finance. For example, trough clustering methods, it is possible to build portfolios of similar stocks (e.g. [35]). Clustering complex objects like time series is a difficult task because this kind of data has particular features like serial correlation and usually are both noisy and heteroskedastic with shifts [2]. From a methodological point of view, time series clustering methods can be divided into three main classes [40]: observation-based, feature-based and model-based clustering. The first group classifies time series on the ground of the raw data (see, e.g. [1,17,23,20,26,37]). In the series with different lengths, the observationbased clustering methods are usually built upon the so-called Dynamic Time Warping (DTW). This technique allows finding an optimal alignment between two given sequences of different length. (e.g. [12], [28]). The methods in the second group cluster objects according to some of their features, like the autocorrelation function [3,11,22] or the parameters of the wavelet decomposition (e.g. [29,43]). The models of the third group assume that the time series are generated from the same statistical model, and group objects according to the estimated parameters. In this context, we mention [48] that proposed a distance measure between ARMA processes. This study has been extended by several authors (e.g. [18,46]). Considering heteroskedasticity, [46] proposed a distance measure between GARCH processes. Then, he clustered time series according to their estimated conditional and unconditional variances. Paper [47] extended [46] for multivariate financial time series. Moreover, [27] proposed a new robust approach for clustering heteroskedastic time series by introducing the fuzzy logic in a GARCH-based approach. In this method, each time series could be in more than one cluster with a certain level of probability following a partitioning around medoids (PAM) procedure. This paper presents a clustering procedure of model-based type. Recent papers that follow a model-based approach (e.g. [46,47,27]) propose to cluster time series according to conditional variance estimates. Time-varying variance is a crucial aspect of heteroskedastic time series, but it is not the only one we should look at to cluster dynamic objects. Indeed, it is reasonable to assume that the entire conditional distribution is time-varying in a time series framework. Some authors (e.g. [23]) showed the usefulness of considering higher moments for clustering. In particular, [23] developed a clustering procedure to group time series with similar tails. There are several statistical models able to capture the dynamic behaviour of conditional higher moments. Examples are [33] and [39]. Recently, [19,31] developed the Dynamic Conditional Score (DCS also called Generalized Autoregressive Score), a very general statistical model that considers the score function of the predictive model density as the driving mechanism for timevarying parameters. A wide class of GARCH-type processes are special cases of the DCS model (see, e.g. [19,5]). Since the DCS model is based on the score, it exploits the data's complete density structure rather than just a few moments.
Similarly to the proposal of [46,47,27], which classifies time series according to the estimated conditional variances, it is possible to cluster time series according to the estimated conditional moments obtained via a Dynamic Conditional Score. One of the crucial points is to determine a similarity measure for the estimated conditional higher moments. Note that the estimated conditional moments are essentially time series themselves. Various distances have been proposed for clustering evolutive objects [40,2,25]. We take inspiration from the works carried out by [22,27] since we adopt a fuzzy logic within a model-based clustering procedure. In particular, we use the Autocorrelation-based fuzzy C-means algorithm (the A-FCM of [22]) to cluster time series with similar conditional moments. Such conditional moments are estimated through a DCS procedure. In general, the fuzzy approach for clustering is very appealing for time series data, as shown by [27] that demonstrated its validity in model-based procedures. The time series's dynamic behaviour plays an essential role in the clustering process, as argued by several authors (e.g. [43,42,28,44]). The reasons for that are different. The fuzzy approach allows one time series to be allocated to two or more clusters with a membership degree that represents the uncertainty related to each time series assignment to each cluster. The membership degree implicitly indicates the presence of a second-best cluster almost as good as the first best. Traditional clustering methods are not able to highlight such conclusions. Moreover, in the real world, identifying a clear boundary between clusters is not easy, so a fuzzy approach appears more attractive than a deterministic one. The literature also highlights that fuzzy clustering procedures are computationally more efficient [27]. In the end, the fuzzy logic has been proved very useful for clustering time series with particular properties [16], like the financial one ( [41,21,34]). Furthermore, thanks to the generality of the underlying statistical model, we are able to consider both linear and nonlinear model specifications and to cluster also non-Gaussian time series. Importantly, our approach considers time-varying higher moments and does not restrict only to the conditional variance. The present paper is structured as follows. In the next section, the statistical model underlying our model-based procedure is presented. Then, in Section 3, we discuss the fuzzy clustering algorithm in detail. To this aim, we first specify the employed distance measure and then the procedure used for estimating and clustering the conditional moments. The proposed method's usefulness is illustrated using several empirical applications with financial time series assuming both linear and nonlinear models' specification and under several assumptions about time series density function. In particular, Section 4 presents the used dataset and outlines the methodologies, while Section 5 collects and discusses the obtained empirical results. The figures related to the estimated conditional moments and their autocorrelation functions are postponed in a devoted Appendix for easy reading.
2 Dynamic Conditional Score model

Basic statistical model
Let us consider a time series (y t : t = 1, . . . , T ) with f t the time-varying parameter vector observed at time t. In our context, the available information set at time t is F t and -as we will see below -it is given by a collection of some previous realizations of the f 's and the y's. Then, define θ a vector of static parameters. In the Dynamic Conditional Score model we assume that y t is generated by the observation density: where the dependence on the parameter θ -which is conveniently let be explicit in (1) -is due to the dependence of f t on θ.
Indeed, given two integers 0 ≤ n, m ≤ T −1, we formalize the DCS(n, m) model for the t-th realization f t of the time-varying parameter vector as: where ω is a real vector and the A's and the B's are real matrices with an appropriate dimension -all the scalar parameters in ω, A 1 , . . . , A n , B 1 , . . . , B m are collected in the vector θ introduced in (1); s t is a type of score of the conditional distribution defined in (1), and it is a function of the data and the parameters, so that s t = s t (y t , f t , F t ; θ). More specifically, this class of models suppose that the evolution of the time-varying parameter vector f t is driven by a vector proportional to the score of density (1) together with an autoregressive component. Indeed, s t is defined as: is a positive definite scaling matrix known at time t and ∇ t (y t , f t , F t ; θ) is the score of y t evaluated with respect to f t , i.e.: Since the score depends on the complete density and not only on some moments of y t , the DCS(n, m) model uses the full density structure for updating f t . We have to highlight that we could get different DCS(n, m) specifications depending by the choice about scaling S t we make. A very appealing feature of DCS(n, m) model is that the vector of parameters θ could be estimated by maximum likelihood as the solution of the following problem (see e.g. [19]):θ = arg max θ T t=1 t (5) where t = log p(y t |f t , F t ; θ) for a realization of y t .

Some special cases
As we said, the Dynamic Conditional Score model is very general since we could get different model specifications considering different assumptions about the scaling (3), the observation density (1) and f t . A common approach, proposed by [19], is to scale the score in order to account for its variance. Particularly, the authors proposed to scale using a matrix S t equal to the inverse of the information matrix of f t to a power γ ≥ 0: with I: where E denotes the expectation at time t − 1 and the conditional score ∇ t is defined as (4). The parameter γ is fixed by the user and usually takes value in the set {0, 1 2 , 1}. When γ = 0, then S t is the identity matrix -i.e., S t = I -and there is no scaling. Differently, if γ = 1, then the conditional score ∇ t is premultiplied by the inverse of (7) to obtain (3). Finally, γ = 1 2 leads to ∇ t scaled to its square-root to obtain (3). We notice that setting γ = 0 one has that the DCS(n, m) model contains the Autoregressive Conditional Multinomial (ACM) model of [49] and the GARMA models of [7] as particular cases.
If we suppose f t = log(σ 2 t ) from a DCS(1, 1) we get Exponential GARCH (E-GARCH) of [45]. Differently, if we suppose t-student density with v degrees of freedom for y t we get the Beta-t-EGARCH of [32] and [30]. Interestingly, the Beta-t-EGARCH is different from the t-GARCH(1,1) model of [10]. Indeed, differently [10], the t − DCS(1, 1) does not results in a substantial increase in the variance after a large absolute change in y t . Another interesting example to mention is when we have y t ∼ G.E.D.(µ t , φ t , v t ). Indeed, assuming in this case f t = σ 2 t we obtain the GED-GARCH(m,n) pro-cess [53], while letting all parameters f t = (µ t , φ t , v t ) to be time vary-varying we get the GED − DCS(n, m) model [5].
A very important feature of this family of distributions is that they include also other common distributions, for different values of shape parameter v. This assumption has been revealed as very useful for modeling financial time series (e.g. [14], [15], [45], [51]). As we can see, the DCS(n, m) could be in principle nonlinear and able to model conditional heteroskedastic time series as well as time series with conditional time-varying higher moments.

Fuzzy clustering of DCS-based conditional moments
We adopt a fuzzy approach in order to account for the uncertainty in time series clustering. The time clustering criterion is grounded on the conditional higher moments. The outline of the procedure is as follows. We first estimate R conditional moments according to the Dynamic Conditional Score model presented in the previous section, hence adopting a model-based procedure. Then, we cluster time series according to the features of the estimated conditional moments. In details, time series with similar r-th conditional moment -where similarity is measured through a predefined concept of distance, see the next subsection -are placed in the same cluster, for each r = 1, . . . , R. Hence, by estimating R moments, we can define R levels of clustering ( [46]). Within each level r = 1, . . . , R, we cluster time series according to their conditional moment over time path similarity. Statistical reasoning suggests that once two series belong to the same cluster for all the R conditional moments, then they share similar conditional distribution. This remark is also in the light of authoritative studies (see, e.g. [46,47]). The adopted distance measure is of Euclidean type. Specifically, we introduce a similarity criterion based on the Euclidean distance between a sufficiently large set of estimated autocorrelations.
Let defineρ (r) l,k the estimated autocorrelation at lag l of the conditional rth moment for a given k-th time series. Estimated autocorrelationρ (r) l,k could be obtained with the usual estimator: whereȳ (r) is the process' mean. A simple distance between two time series X and Y of r-th conditional moment that is based on estimated conditional moments' autocorrelation can be defined as: that is equal to the squared Euclidean distance between the estimated autocorrelation functions. In using the distance defined in (9), we account for the autocorrelation structure of the processes, hence being in line with the fuzzy clustering approach proposed in [22]. The procedure we propose is based on the fuzzy C-means (FCM) algorithm [38,8]. In particular, we adopt a fuzzy clustering algorithm for time series with similar autocorrelation structure of conditional moments using the autocorrelation-based fuzzy C-means (A-FCM) algorithm of [22]. The proposed clustering model could be formalized as follows: under the constraints: where u k,c denotes the membership degree of the r-th conditional moments of the k-th time series to the c-th cluster, the parameter m > 1 controls for the fuzziness of the partition,ρ (r) l,k is the k-th unconditional moments' autocorrelation at l-th lag andρ (r) l,c represents the c-th centroid autocorrelation. Following [22], the optimal solutions of the model (10) are equal to: and:ρ To solve the constrained optimization problem in (10), let us consider the following Lagrangian function: Then it follows that: By substitutions of the two equations above we get the solution for the membership degree u k,c . Similarly it could be shown that by taking the derivative of the Lagrangian function above with respect toρ (r) l,c we get the solution showed before. Moreover, we have to highlight that the solution ofρ (r) l,c satisfies for the internality property the autocorrelation constraints −1 ≤ρ (r) l,c ≤ 1. Summarising, the procedure we propose involves the following steps: 1. For each k-th time series (y k,t : t ≥ 0), we estimate the DCS model with its parameters.
At the end of the procedure, we could define r-levels of clusters depending on how many conditional moments we estimate. Then, we could be able to understand also which time series are located in the same cluster for more than one conditional moment. In doing so, we obtain a clustering of the considered time series.

Empirical experiments: data and methods
This section outlines the main ingredients of the proposed empirical application.

Data
We consider the dataset of [4], consisting in daily log-returns of the Dow Jones 30 constituents from 3/03/1987 to 3/02/2009 for a total of T = 5, 521 observations for each k-th series 2 . The constituents are labelled according to the standard nomenclature. Time series plots are shown in Fig. 1.

Figure 1: Stock returns' time series
All the considered stock returns are stationary. Indeed, the absence of a deterministic trend is clear from Fig. 1. Moreover, an unreported analysis shows also the absence of a unit root, hence getting additional evidence of the stationarity. We notice that the Dow Jones 30 constituents do not share the same variability since some of them have more a volatile pattern. Then, most of them are far from the origin on the y-axis. However, we could not say anything more about data distribution. A common approach is to assume normality. However, real data are far to be normally distributed, especially financial returns. Indeed, as we can see from stock returns' empirical densities in Fig. 2, data are not normally distributed. This consideration is valid for almost all of them. Specifically, they show heavy tails and different degrees of skewness. As further evidence, in the following Tab. 1 are reported the results of Jarque-Bera (JB) normality test [36] for each k-th time series. Considering the empirical kurtosis and skewness, under the null hypothesis of the test, we have that observed data are generated by a normal distribution. According to [50], both [36] and its adjusted version developed by [52] are on the whole, the most powerful normality tests.  For all stock returns, we reject the null. If we suppose that higher moments of the data distribution could be time-varying, these findings justify different modelling approaches that are discussed in the following subsection.

Methodology; employed classes of DCS models
First, we consider a linear and Gaussian model specification. This kind of assumption is in line with other classical papers of model-based financial time series clustering (e.g. [46,27]). Then, we relax such a condition by considering the t-Dynamic Conditional Score model -also called t-Beta-EGARCH -that is both non-Gaussian and nonlinear, and the skew-t-Dynamic Conditional Score (skew-t-DCS), that allows also for time-varying skewness in asset returns. In the first case, we have two time-varying moments; in the second one, we have three conditional moments (location, scale and shape). For comparison purposes, in the third and last experiment, we estimate four conditional moments (location, scale, skewness and shape).

Gaussian Dynamic Conditional Score
In the first application, we specify a Gaussian-DCS(1,1) model for all time series. In particular, we suppose y t ∼ N (µ t , σ 2 t ) where: The updating mechanism for time varying parameters µ t and σ 2 t could be specified as follows: where θ is the vector containing all the scalar parameters in ω, A, B and s t is scaled by (6) setting γ = 1. In particular, the conditional score vectors are in this case given by: Summarizing, the model's variables and parameters are: In order to estimate conditional moments we consider the following equations for each k-th time series: Once parameters are estimated thorough the usual maximum likelihood approach, conditional moments estimatesμ t andσ t are obtained as:

t-Dynamic Conditional Score
As we have seen, stock returns are not normally distributed. In order to account for non-normality we now assume that y t follows a t-student distribution with location µ t , scale φ t and degrees of freedom v t > 2 with its density given by: Now we have that f t = (µ t , φ t , v t ). Assuming a DCS(1,1), the uptading mechanism for f t is again given by: In the case of t-student density, the conditional score vectors are equal to: where ψ(·) is the Digamma function. By scaling the conditional score with I −1 in (6), the updating mechanism for conditional moments could be written as: In compact form: This model has also been defined as Beta-t-EGARCH model by [31,30].
To estimate conditional moments, we estimate parameters with MLE and replace them within the conditional moments' equations.

Skew-t Dynamic Conditional Score
In the last experiment we account for dynamic conditional skewness in asset returns. The density associated to the skew-t distribution is [6]: where the vector of the four time varying parameters is f t = (µ t , φ t , λ t , v t ), the term T is given by: and: In such a case, the conditional scores are equal to (see [24]): Written in compact form, we estimate the following skew-t-DCS(1,1) model: where now we have: The conditional moments equation could be written as: Where the third one allows for time varying skewness. Replacing MLE estimatesω,Â,B we get estimated conditional moments.

Empirical experiments: results
In what follows, we provide a detailed explanation of clustering results.

Fuzzy clustering under normality assumption
The first step is to estimate the model presented in the Section 4.2.1. Conditional moments time series are shown, for the sake of brevity, within the Appendix section. Before the identification of the clusters through the application of the A-FCM clustering algorithm of [22], we have found the optimal number of clusters. Following a fuzzy approach, in the Tab. 2 are reported the values of the fuzzy silhouette criterion proposed by [13]. The optimal number of clusters is the one with the highest value.  Table 2: Fuzzy silhouette of [13] Therefore, we have C = 2 number of clusters for all R = 2 conditional moments. In the first, we have two groups with the same conditional mean while in the second we have two groups with the same conditional variance. Particularly the DCS-AFCM procedure results in the following assignment: Further, in the Tab. 3 we reported the membership degrees of the assignment or, alternatively, the level of probability that each time series has to be in each C = 2 clusters. In the last columns, instead, are reported the final assignment.  Table 3: Clustering based on conditional mean (left) and variance (right) estimates. Membership degrees for each cluster are reported, while the group's assignment (assigned ) after the DCS-based AFCM procedure is in the last column.
The 86.7% of the assets are clustered in the same group about conditional mean estimated trough a Gaussian Dynamic Conditional Score model while, as the other side, 76.7% of them are within the second group according to the conditional variances' pattern. Moreover, cluster assignment is much more uncertain for some assets when the conditional variance is considered, rather than the conditional mean. For example, some assets like AA or AIG have almost the same membership degree for both clusters 1 and 2 when we take the conditional variance. This never happens for conditional mean clustering. Therefore, here we have a clear example of the fuzzy clustering importance; we could be able to assign asset like AA and AIG to both clusters with almost the same probability. Moreover, starting from these R = 2 levels of clusters, we are also able to group time series with exactly the same conditional distribution over time. This could be done by checking for assets that are within the same group for both conditional mean and variance. Particularly, in this case, the assets BAC, GM and AIG are always in the second group, so we say that they share the same conditional distribution assuming only two time-varying moments. On the other side, the assets AXP, CAT, CVX, KO, MMM and MRK are always in the first group according to both moments. Hence, we could generate an additional cluster with those assets that share the same conditional distribution. Nevertheless, those clusters could be formed according to certain probabilities. Suppose we assign, about conditional variance, with the 46% of probability both AA and AIG to the cluster 1 instead of 2. In this case AIG does not have anymore the same conditional distribution of BAC and GM, while AA will have the same conditional distribution of AXP, CAT, CVX, KO, MMM and MRK with the 46% of probability.

Fuzzy clustering under t-student density
Considering stock retuns' non normality, we estimated the t − DCS(1, 1) model presented in the Section 4.2.2. The optimal number of clusters has been decided according to the fuzzy silhouette criterion whose results are showed in the Tab  Table 4: Fuzzy silhouette of [13] The clusters obtained according to the t-DCS-AFCM procedure are the following: Moreover, in the Tab. 5 we reported the uncertainty about clustering assignment in such scenario. Differently from Gaussian scenario, here we have R = 3 levels of clusters. In the first, we have two groups with the same conditional location, in the second we have two groups with the same conditional scale, and in the last one we have two groups with the same conditional shape.
In the case, we have, as the previous one, almost no uncertainty in clustering time series according to their first conditional moment (see Tab. 5). Instead, a bit of uncertainty rise in clustering two assets (BAC and T ) about their second conditional moment (Tab. 5). Instead, we have also low uncertainty in clustering time series according to conditional shape (Tab. 5). About clusters composition, the first level has the 86.7% of the assets within the same group for conditional mean, even if now most of them are located in the group 2 instead of 1. The cluster level r = 2 is more balanced since the 43.3% of the assets are in the first group and the rest in the second one. At this point, we have that most of the assets are located in the same group. Indeed, we have that 11 of the 16 assets within the same cluster for conditional scale are also in the same group for the conditional mean. Therefore, we could conclude that 11 assets (AXP, BAC, HD, IBM, JNJ, KO, MCD, MMM, MSFT, PG, VTX, WMT ) are similar in terms of conditional location and scale. Nevertheless, not all of them share the same conditional distribution. Indeed, only AXP and BAC are always in the same group, and we say that they share a similar conditional distribution. However, we have to highlight a probabilistic argument that is very important is such a scenario. Indeed, we should account for the uncertainty about the conditional second-moment group assignment for the stock BAC. With a probability of 45% its conditional scale could be assigned to group 1 instead of 2. In this case, with such level of probability, no stock could be considered with a similar conditional distribution.  Table 5: Clustering based on conditional location (left), scale (center) and shape (right) estimates. Membership degrees for each cluster are reported, while the group's assignment (assigned ) after the DCS-based AFCM procedure is in the last columns.

Fuzzy clustering under skew-t density
The last experiment allows for time-varying skewness with the additional equation (17) that has to be estimated. Then we cluster time series according to their moments' estimate via the skew-t-DCS model. Optimal number of clusters is shown in Tab. 6.  Table 6: Fuzzy silhouette of [13] With this last experiment we do not get the same optimal number of clusters for all considered higher conditional moments. Indeed, while we consider C = 3 clusters for conditional location, we have C = 2 for both conditional scale and skewness and C = 5 for conditional shape. In this case, it is evidently more difficult than in the previous models to have time series within the same clusters for all conditional moments. The results of the proposed skew-t-DCS model-based A-FCM algorithm are summarized by the following scheme: Once again, clustering made according to conditional mean is the one with the lowest level of uncertainty since the membership degree of any asset to the alternative groups is tiny (Tab. 7). Differently, as we have also seen for the other two experiments, the clusters based on the conditional scale are the ones with the highest level of uncertainty (Tab. 8). For example, some assets -like PFE and BAC -have 35% probability to be clustered in a different group; for PG such a probability is 40%, while others like C and BA also a higher probability of 47% to be in another cluster. In the end, while also the conditional skewness-based clusters show some relevant cases of uncertainty (e.g. HPQ, JPM, MRK, T or UTX ); differently, shape-based clusters present a low level of uncertainty (see Tab. 9 and Tab. 10). Overall, while we can build R clustering levels, in this context, accounting for cluster uncertainty evident for some conditional moments like scale and skewness, we are not able to find stocks always within the same clusters. In other words, all stocks have different conditional distributions regarding higher moments. If we look at just some of them, we could easily get in error.

Experiments with simulated data
As robustness study, in this section we illustrate the results of the proposed clustering procedure with simulated data. At this aim, we generate N = 10 independent time series of different lengths T = 50, 200, 1000 belonging to different Data Generating Processes (D.G.P.). The consideration of different lengths is useful to evaluate the performances of the proposed algorithm for either short or long time series. Moreover, to evaluate how the selection of the number of lags L affects the results, we consider as in [25] different choice L = 10, 25, 50. The performance of the proposed clustering procedure is evaluated, following [25], as the correct classification rate over M = 300 replications. We suppose the presence of C = 2 clusters with N = 5 time series in each group, where the N = 10 time series are generated by the following D.G.P.: that is the Beta-Skew-t-EGARCH of [30], hence a non-Gaussian and highly non-linear process, with ω be the log-scale intercept, φ is the persistency parameter, α the ARCH parameter, u k,t the conditional score with respect to η k,t and β the leverage parameter. We suppose that the time series in c = 1 are generated by different D.G.P. (17) with respect to those in c = 2. In particular, we first suppose that the time series in c = 1 are generated by a D.G.P. with leverage (so β = 0) with the same parameters of those in c = 2 that, instead, is supposed to be generated by a process without leverage. The Fig. 3 shows the difference of four time series simulated by the process (17) with exactly the same parameters. However, while for two we suppose the leverage (β = 0), for the other two we assume no leverage (β = 0). As evident, the time series generated with the process augmented for the leverage are characterized by similar fluctuations of those without leverage but with shocks (e.g. stock returns immediately after a financial crisis) that dramatically changes the scale of the figures (see Fig. 3). In a second scenario, instead, we suppose that the time series in the two groups are generated by the D.G.P. in (17) both without leverage but with different parameters. The proposed DCS-A-FCM algorithm is compared with respect the Autocorrelation based Fuzzy C-Means of [22] that is close in the spirit to our proposal. For a fair comparison, besides we know the true D.G.P. (17), we perform the classification with a Dynamic Conditional Score model based on a wrong density specification. Indeed, we suppose a Gaussian distribution of the data while they are Skew-t distributed, so we estimate the conditional mean and the conditional variance with a Gaussian-DCS(1,1) as presented in the Section 4.  ACF is computed at different lags L. Both the groups are generated with the D.G.P. (17) replacing ω = 1e −6 , φ = 0.6 and α = 0.1. However, while in the first group β = 0 in the second group β = 1. The classification rate refers to the percentage of the average correct classification over the M = 300 trials.
First of all it is evident the superior classification ability of the DCS-based clustering procedures with respect the A-FCM for all the considered time series lengths, even if we specified a wrong density function. In particular, the classification based on the conditional variance provides better results (up to 72% in the case of long time series) than the classification done according to the conditional mean. The time series length is an important element to evaluate the robustness of the result for either short and long time series. On one hand, the A-FCM algorithm ensures good performances with long time series (T = 500) with an accuracy close to 65% in the best case but it is not able to perform better then our proposal also in this case (both with 67% and 72% of correct classification rate). On the other hand, for short time series, the DCS-based procedures show robust results with the classification accuracy about the 70%, while the A-FCM has an accuracy level a bit higher than 50%. Moreover, it is not clear what is the effect of the lag's length L. Indeed, in the case of short time series length, for the A-FCM algorithm it seems that increasing L reduces the classification accuracy, while for the DCS-A-FCM with the conditional mean similar are robust to the selection of the lag length. In the end, the DCS-A-FCM with the conditional variance ensures higher classification accuracy with increasing L. In the case of long time series, instead, for all the classification algorithms it seems that with increasing lag length L the miss-classification rate decreases, with the highest accuracy obtained for all of them specifying L = 50. The Tab. 12 shows, instead, the classification accuracy for the second scenario where the time series are generated with the D.G.P. (17) without leverage and different parameters. As in the previous case, the time series are simulated for different lengths T and the classification is done according to different lags' selection L.  ACF is computed at different lags L. Both the groups have β = 0 but with different parameters. The first group is generated with the D.G.P. (17) replacing ω = 2e −6 , φ = 0.4 and α = 0.2. The second group is generated by the D.G.P. (17) assuming ω = 1e −6 , φ = 0.8 and α = 0.05. The classification rate refers to the percentage of the average correct classification over the M = 300 trials.
In this second scenario, less influenced by the presence of outliers, the classification accuracy of the DCS-based clustering algorithms is more evident. The benchmark clustering model, the A-FCM, shows in most of the considered cases an accuracy close to 50%, while the proposed algorithms reach also the 82% of correct classification rate. As in the previous example, clustering according to conditional variance provide a much higher classification accuracy with respect the clustering based on the conditional mean. More in details, for the case of short time series, the over performances of the proposed DCS-A-FCM procedures are higher up to the 10% of greater classification accuracy, with percentage always much above the 60% threshold, that is the highest rate for the benchmark A-FCM. In the case of long time series, the conditional moment-based clustering approaches provide over performances up to 26% of greater classification accuracy rate. In the end, in this second scenario the results seems more robust to the choice of the lag length L for all the studied clustering algorithms. However, according to the simulations' results, with long time series a lag length L = 50 it appears the best choice.

Conclusions
In this paper, we develop a new model-based fuzzy clustering algorithm able to classify time series according to their conditional higher moments, estimated via a Dynamic Conditional Score. Since the implemented statistical model is based on the score, it exploits the complete density structure of the data rather than just a few moments. The similarity between conditional moments' estimates is measured according to a Euclidian distance of the estimated sample autocorrelations. In other words, we propose to cluster time series according to some features of conditional moments. Thanks to the generality of the underlying statistical model, the framework presented here generalizes most of the model-based procedures existing in the literature. In this respect, it is worth to mention that a wide class of GARCH-type processes are just special cases of the Dynamic Conditional Score. By specifying different distribution densities for the considered time series, we can cluster challenging instances like non-Gaussian heteroskedastic time series that also show time-varying higher moments. Moreover, the Dynamic Conditional Score allows us to consider also nonlinear models underlying the clustering procedure. The Beta-t-EGARCH is a clear example of a possible nonlinear and non-Gaussian model specification. Interestingly, we adopt a fuzzy clustering perspective. In doing so, we admit that each time series can be in more than one cluster with a certain level of probability. Indeed, the fuzzy approach implicitly indicates the presence of a second-best cluster -sometimes, almost as good as the first best; this is a property that is missing in the traditional clustering methods. Moreover, in the real world, the identification of a clear boundary between clusters is not an easy task, so a fuzzy approach appears more attractive than a deterministic one. We present applications to real financial data based on stock returns. In this respect, we show that it is possible to obtain a number of reference conditional moments for clustering, which depends on the number of timevarying parameters in the specified Dynamic Conditional Score. In the case of a Gaussian-DCS model, we get R = 2 time-varying moments (mean and variance) with two clusters. In the case of a t-DCS model, instead, we have R = 3 time-varying parameters and the same number of clusters. Trough the DCS-based fuzzy clustering we could obtain also thinner clusters with time series that share a similar conditional distribution if those series are placed in the same cluster R times. Future researches could be devoted to the possible applications of the proposed clustering procedure in other real-world frameworks. An interesting proposal might be the application of the presented clustering procedure to the construction of portfolios of stocks when there is a huge amount of available assets. Nevertheless, the clustering procedure we propose is general, and it could also be applied in other disciplines such as engineering, like the relevant case of signal processes clustering.