clustMixType : User-Friendly Clustering of Mixed-Type Data in R by

Clustering algorithms are designed to identify groups in data where the traditional emphasis has been on numeric data. In consequence, many existing algorithms are devoted to this kind of data even though a combination of numeric and categorical data is more common in most business applications. Recently, new algorithms for clustering mixed-type data have been proposed based on Huang’s k-prototypes algorithm. This paper describes the R package clustMixType which provides an implementation of k-prototypes in R.


Introduction
Clustering algorithms are designed to identify groups in data where the traditional emphasis has been on numeric data.In consequence, many existing algorithms are devoted to this kind of data even though a combination of numeric and categorical data is more common in most business applications.For an example in the context of credit scoring, see, e.g.Szepannek (2017).The standard way to tackle mixed-type data clustering problems in R is to use either (1) Gower distance (Gower, 1971) via the gower package (van der Loo, 2017) or the daisy(method = "gower") in the cluster package (Maechler et al., 2018); or (2) Hierarchical clustering through hclust() or the agnes() function in cluster.Recent innovations include the package CluMix (Hummel et al., 2017), which combines both Gower distance and hierarchical clustering with some functions for visualization.As this approach requires computation of distances between any two observations, it is not feasible for large data sets.The package flexclust (Leisch, 2006) offers a flexible framework for k-centroids clustering through the function kcca() which allows for arbitrary distance functions.Among the currently pre-implemented kccaFamilies, there is no distance measure for mixed-type data.Alternative approaches based on expectation-maximization are given by the function flexmixedruns() in the fpc package (Hennig, 2018) and the package clustMD (McParland, 2017).Both require the variables in the data set to be ordered according to their data type, and that categorical variables to be preprocessed into integers.The clustMD algorithm (McParland and Gormley, 2016) also allows ordinal variables but is quite computationally intensive.The kamila package (Foss and Markatou, 2018) implements the KAMILA clustering algorithm which uses a kernel density estimation technique in the continuous domain, a multinomial model in the categorical domain, and the Modha-Spangler weighting of variables in which categorical variables have to be transformed into indicator variables in advance (Modha and Spangler, 2003).
Recently, more algorithms for clustering mixed-type data have been proposed in the literature (Amir and Dey, 2007;Dutta et al., 2012;Foss et al., 2016;He et al., 2005;HajKacem et al., 2016;Ji et al., 2012Ji et al., , 2013Ji et al., , 2015;;Lim et al., 2012;Liu et al., 2017;Pham et al., 2011).Many of these are based on the idea of Huang's k-prototypes algorithm (Huang, 1998).The rest of this paper describes the R package clustMixType (Szepannek, 2018), which provides up to the author's knowledge the first implementation of this algorithm in R. The k-modes algorithm (Huang, 1997a) has been implemented in the package klaR (Weihs et al., 2005;Roever et al., 2018) for purely categorical data, but not for the mixed-data case.The rest of the paper is organized as follows: A brief description of the algorithm is followed by the functions in the clustMixType package.Some extensions to the original algorithm are discussed and as well as a worked example application.

k-prototypes clustering
The k-prototypes algorithm belongs to the family of partitional cluster algorithms.Its objective function is given by: where x i , i = 1, . . ., n are the observations in the sample, µ j , j = 1, . . ., k are the cluster prototype observations, and u ij are the elements of the binary partition matrix The R Journal Vol.10/2, December 2018 ISSN 2073-4859 The distance function is given by: where m is an index over all variables in the data set where the first q variables are numeric and the remaining p − q variables are categorical.Note that δ(a, b) = 0 for a = b and δ(a, b) = 1 for a = b, and d() corresponds to weighted sum of Euclidean distance between two points in the metric space and simple matching distance for categorical variables (i.e. the count of mismatches).The trade off between both terms can be controlled by the parameter λ which has to be specified in advance as well as the number of clusters k.For larger values of λ, the impact of the categorical variables increases.
For λ = 0, the impact of the categorical variables vanishes and only numeric variables are taken into account, just as in traditional k-means clustering.
The algorithm iterates in a manner similar to the k-means algorithm (MacQueen, 1967) where for the numeric variables the mean and the categorical variables the mode minimizes the total within cluster distance.The steps of the algorithm are: 1. Initialization with random cluster prototypes.

For each observation do:
(a) Assign observations to its closest prototype according to d().
(b) Update cluster prototypes by cluster-specific means/modes for all variables.
3. As long as any observations have swapped their cluster assignment in 2 or the maximum number of iterations has not been reached: repeat from 2.

k-prototypes in R
An implementation of the k-prototypes algorithm is given by the function kproto(x, k, lambda = NULL, iter.max= 100, nstart = 1, na.rm = TRUE) where • x is a data frame with both numeric and factor variables.As opposed to other existing R packages, the factor variables do not need to be preprocessed in advance and the order of the variables does not matter.
• k is the number of clusters which has to be pre-specified.Alternatively, it can also be a vector of observation indices or a data frame of prototypes with the same columns as x.If ever at the initialization or during the iteration process identical prototypes do occur, the number of clusters will be reduced accordingly.
• lambda > 0 is a real valued parameter that controls the trade off between Euclidean distance for numeric variables and simple matching distance for factor variables for cluster assignment.
If no λ is specified the parameter is set automatically based on the data and a heuristic using the function lambdaest().Alternatively, a vector of length ncol(x) can be passed to lambda (cf.Section on Extensions to the original algorithm).
• iter.maxsets the maximum number of iterations, just as in kmeans().The algorithm may stop prior to iter.max if no observations swap clusters.
• nstart may be set to a value > 1 to run k-prototypes multiple times.Similar to k-means, the result of k-prototypes depends on its initialization.If nstart > 1, the best solution (i.e. the one that minimizes E) is returned.
• Generally, the algorithm can deal with missing data but as a default NAs are removed by na.rm = TRUE.
Two additional arguments, verbose and keep.data,can control whether information on missing values should be printed and whether the original data should be stored in the output object.The keep.data=TRUE option is required for the default call to the summary() function, but in case of large x, it can be set to FALSE to save memory.
The output is an object of class "kproto".For convenience, the elements are designed to be compatible with those of class "kmeans": • cluster is an integer vector of cluster assignments • centers stores the prototypes in a data frame The R Journal Vol.10/2, December 2018 ISSN 2073-4859 • size is a vector of corresponding cluster sizes • withinss returns the sum over all within cluster distances to the prototype for each cluster • tot.withinss is their sum over all clusters which corresponds to the objective function E • dists returns a matrix of all observations' distances to all prototypes in order to investigate the crispness of the clustering • lambda and iter store the specified arguments of the function call • trace lists the objective function E as well as the number of swapped observations during the iteration process Unlike "kmeans", the "kproto" class is accompanied corresponding predict() and summary() methods.The predict.kproto() method can be used to assign clusters to new data.Like many of its cousins, it is called by predict(object, newdata) The output again consists of two elements: a vector cluster of cluster assignments and a matrix dists of all observations' distances to all prototypes.
The investigation resulting from a cluster analysis typically consists of identifying the differences between the clusters, or in this specific case, those of the k prototypes.For practical applications besides the cluster sizes, it is of further interest to take into account the homogeneity of the clusters.For numeric variables, this can be done by calling the R function summary().For categorical variables the representativity of the prototypes is given their frequency distribution obtained by prop.table().The summary.kproto() method applies these methods to the variables conditional to the resulting clusters and returns a comparative results table of the clusters for each variable.
The summary is not restricted to the training data but it can further be applied to new data by calling summary(object,data) where data are new data that will be internally passed to the predict() method on the object of class "kproto".If no new data is specified (default: data = NULL), the function requires object to contain the original data (argument keep.data= TRUE).In addition, a function clprofiles(object, x, vars = NULL) supports the analysis of the clusters by visualization of cluster profiles based on an object of class "kproto" and data x.Note that the latter may also have different variables compared to the original data, such as for profiling variables that were not used for clustering.As opposed to summary.kproto(),no new prediction is done but the cluster assignments of the "kproto" object given by object$cluster are used.For this reason, the observations in x must be the same as in the original data.Via the vars argument, a subset of variables can be specified either by a vector of indices or variable names.Note that clprofiles() is not restricted to objects of class "kproto" but can also be applied to other cluster objects as long as they are of a "kmeans"-like structure with elements cluster and size.

Extensions to the original algorithm
For unsupervised clustering problems, the user typically has to specify the impact of the specific variables on the desired cluster solution which is controlled by the parameter λ.Generally, small values λ ∼ 0 emphasize numeric variables and will lead to results similar to standard k-means whereas larger values of λ lead to an increased influence of categorical variables.In Huang (1997b), the average standard deviation σ of the numeric variables is the suggested choice for λ and in some practical applications in the paper values of 1 3 σ ≤ λ ≤ 2 3 σ are used.The function lambdaest(x, num.method = 1, fac.method = 1, outtype = "numeric") provides different data based heuristics for the choice of λ: The average variance σ 2 (num.method= 1) or standard deviation σ (num.method= 2) over all numeric variables is related to the average concentration h cat of all categorical variables.We compute h cat by averaging either h m = 1 − ∑ c p 2 mc (fac.method= 1) or h m = 1 − max c p mc (fac.method= 2) over all variables m where c are the categories of the factor variables.We set λ = σ t h cat , t ∈ {1, 2} as a user-friendly default choice to prevent over-emphasizing either numeric or categorical variables.If kproto() is called without specifying lambda, the parameter is automatically set using num.method= fac.method= 1.Note that this should be considered a starting point for further analysis; the explicit choice of λ should be done carefully based on the application context.
Originally, λ is real-valued, but in order to up-or downweight the relevance of single variables in a specific application context, the function kproto() is extended to accept vectors as input where each The R Journal Vol.10/2, December 2018 ISSN 2073-4859 element corresponds to a variable specific weight, λ m .The formula for distance computation changes to: Note that the choice of λ only affects the assignment step for the observations but not the computation of the prototype given a cluster of observations.By changing the outtype argument into a vector, the function lambdaest() returns a vector of λ m s.In order to support a user-specific definition of λ based on the variables' variabilities, outtype = "variation" returns a vector of original values of variability for all variables in terms of the quantities described above.
An issue of great practical relevance is the ability of an algorithm to deal with missing values which can be solved by an intuitive extension of k-prototypes.During the iterations, cluster prototypes can be updated by ignoring NAs: Both means for numeric variables as well as modes for factors can be computed based on the available data.Similarly, distances for cluster assignment can be computed for each observation based on the available variables only.This not only allows cluster assignment for observations with missing values, but already takes these observations into account when the clusters are formed.By using kproto(), this can be obtained by setting the argument na.rm = FALSE.Note that in practice, this option should be handled with care unless the number of missing values is very small.The representativeness of a prototype might become questionable if its means and modes do not represent the major part of the cluster's observations.Finally, a modification of the original algorithm presented in Section k-prototypes clustering allows for vector-wise computation of the iteration steps, reducing computation time.The update of the prototypes is not done after each new cluster assignment, but once each time the whole data set has been reassigned.The modified k-prototypes algorithm consists of the following steps: 1. Initialization with random cluster prototypes.
2. Assign all observations to its closest prototype according to d().

Update cluster prototypes.
4. As long as any observations have swapped their cluster assignment in 2 or the maximum number of iterations has not been reached: repeat from 2.

Example
As an example, data x with two numeric and two categorical variables are simulated according to the documentation in ?kproto:Using set.seed(42), four clusters j = 1, . . ., 4 are designed such that two pairs can be separated only by their numeric variables and the other two pairs only by their categorical variables.The numeric variables are generated as normally distributed random variables with cluster specific means µ 1 = µ 2 = −µ 3 = −µ 4 = Φ −1 (0.9) and the categorical variables have two levels (A and B) each with a cluster specific probability p 1 1 summarizes the clusters.It can be seen that both numeric and categorical variables are needed in order to identify all four clusters.
By construction, taking into account either only numeric (k-means) or only factor variables (kmodes) will not be able to identify the underlying cluster structure without further preprocessing of the data in this example.A performance comparison using the Rand index (Rand, 1971) as computed by the package clusteval (Ramey, 2012) results in rand indices of 0.728 (k-means) and 0.733 (k-modes).As already seen above, the prototypes as identified by clustMixType do represent the true clusters quite well, as the corresponding Rand index improves to 0.870.library(klaR) library(clusteval) kmres <-kmeans(x[,3:4], 4) # kmeans using numerics only kmores <-kmodes(x[,1:2], 4) # kmodes using factors only cluster_similarity(kpres$cluster, clusid, similarity = "rand") cluster_similarity(kmres$cluster, clusid, similarity = "rand") cluster_similarity(kmores$cluster, clusid, similarity = "rand") The runtime of kproto() is linear in the number of observations (Huang, 1997b) and thus it is also applicable to large data sets.half of them categorical).It is possible to run kproto() for more than hundred variables which is far beyond most practical applications where human interpretation of the resulting clusters is desired.Note that clustMixType is written in R and currently no C++ code is used to speed up computations which could be a subject of future work.
In order to determine the appropriate number of clusters for a data set, we can use the standard scree test.In this case, the objective function E is given by the output's tot.withinss element.The kproto() function is run multiple times for varying numbers of clusters (but fixed λ) and the number of clusters is chosen as the minimum k from whereon no strong improvements of E are possible.In Figure 2 (right), an elbow is visible at the correct number of clusters in the sample data; recall that we simulatd from four clusters.Note that from a practitioner's point of view, an appropriate solution requires clusters that are well represented by their prototypes.For this reason, the choice of the number of clusters should further take into account a homogeneity analysis of the clusters as returned by summary(), clprofiles() or the withinss element of the "kproto" output.
Figure 2 (left) shows the behaviour of the runtime for 50 repeated simulations of the example data with an increased number of variables (half of them numeric and The R Journal Vol.10/2, December 2018 ISSN 2073-4859

Table 1 :
Separation of clusters in the example.