Rev 8309 | Blame | Compare with Previous | Last modification | View Log | Download | RSS feed
\name{clara}\alias{clara}\title{Clustering Large Applications}\description{Computes a \code{"clara"} object, a \code{\link{list}} representing aclustering of the data into \code{k} clusters.}\usage{clara(x, k, metric = c("euclidean", "manhattan", "jaccard"),stand = FALSE, cluster.only = FALSE, samples = 5,sampsize = min(n, 40 + 2 * k), trace = 0, medoids.x = TRUE,keep.data = medoids.x, rngR = FALSE, pamLike = FALSE, correct.d = TRUE)}\arguments{\item{x}{data matrix or data frame, each row corresponds to an observation,and each column corresponds to a variable. All variables must benumeric (or logical).Missing values (NAs) are allowed.}\item{k}{integer, the number of clusters.It is required that \eqn{0 < k < n} where \eqn{n} is the number ofobservations (i.e., n = \code{nrow(x)}).}\item{metric}{character string specifying the metric to be used for calculatingdissimilarities between observations.The currently available options are "euclidean", "manhattan","jaccard". % , and "gower". For the latter, see \code{\link{daisy}()}.Euclidean distances are root sum-of-squares of differences, andmanhattan distances are the sum of absolute differences.}\item{stand}{logical, indicating if the measurements in \code{x} arestandardized before calculating the dissimilarities. Measurementsare standardized for each variable (column), by subtracting thevariable's mean value and dividing by the variable's mean absolutedeviation.}\item{cluster.only}{logical; if true, only the clustering will becomputed and returned, see details.}\item{samples}{integer, say \eqn{N}, the number of samples to be drawn from thedataset. The default, \code{N = 5}, is rather small for historical (andnow back compatibility) reasons and we \emph{recommend to set\code{samples} an order of magnitude larger}.}\item{sampsize}{integer, say \eqn{j}, the number of observations in eachsample. \code{sampsize} should be higher than the number of clusters(\code{k}) and at most the number of observations (\eqn{n =}\code{nrow(x)}). While computational effort is proportional to \eqn{j^2},see note below, it may still be advisable to set\eqn{j = }\code{sampsize} to a \emph{larger} value than the (historical) default.}\item{trace}{integer indicating a \emph{trace level} for diagnosticoutput during the algorithm.}\item{medoids.x}{logical indicating if the medoids should bereturned, identically to some rows of the input data \code{x}. If\code{FALSE}, \code{keep.data} must be false as well, and the medoidindices, i.e., row numbers of the medoids will still be returned(\code{i.med} component), and the algorithm saves space by needingone copy less of \code{x}.}\item{keep.data}{logical indicating if the (\emph{scaled} if\code{stand} is true) data should be kept in the result.% (\code{keepdata} is equivalent to \code{keep.data} where the former% is deprecated.)Setting this to \code{FALSE} saves memory (and hence time), butdisables \code{\link{clusplot}()}ing of the result. Use\code{medoids.x = FALSE} to save even more memory.}\item{rngR}{logical indicating if \R's random number generator shouldbe used instead of the primitive clara()-builtin one. If true, thisalso means that each call to \code{clara()} returns a different result-- though only slightly different in good situations.}\item{pamLike}{logical indicating if the \dQuote{swap} phase (see\code{\link{pam}}, in C code) should use the same algorithm as\code{\link{pam}()}. Note that from Kaufman and Rousseeuw'sdescription this \emph{should} have been true always, but as theoriginal Fortran code and the subsequent port to C has alwayscontained a small one-letter change (a typo according to Martin Maechler)with respect to PAM, the default, \code{pamLike = FALSE} has been chosen toremain back compatible rather than \dQuote{PAM compatible}.}\item{correct.d}{logical or integer indicating that---only in the caseof \code{NA}s present in \code{x}---the correct distance computationshould be used instead of the wrong formula which has been presentin the original Fortran code and been in use up to early 2016.Because the new correct formula is not back compatible, for the timebeing, a warning is signalled in this case, unless the user explicitlyspecifies \code{correct.d}.}}\value{If \code{cluster.only} is false (as by default),an object of class \code{"clara"} representing the clustering. See\code{\link{clara.object}} for details.If \code{cluster.only} is true, the result is the "clustering", aninteger vector of length \eqn{n} with entries from \code{1:k}.}\details{\code{clara} (for "euclidean" and "manhattan") is fully described inchapter 3 of Kaufman and Rousseeuw (1990).Compared to other partitioning methods such as \code{pam}, it can deal withmuch larger datasets. Internally, this is achieved by consideringsub-datasets of fixed size (\code{sampsize}) such that the time andstorage requirements become linear in \eqn{n} rather than quadratic.Each sub-dataset is partitioned into \code{k} clusters using the samealgorithm as in \code{\link{pam}}.\crOnce \code{k} representative objects have been selected from thesub-dataset, each observation of the entire dataset is assignedto the nearest medoid.The mean (equivalent to the sum) of the dissimilarities of theobservations to their closest medoid is used as a measure of thequality of the clustering. The sub-dataset for which the mean (orsum) is minimal, is retained. A further analysis is carried out onthe final partition.Each sub-dataset is forced to contain the medoids obtained from thebest sub-dataset until then. Randomly drawn observations are added tothis set until \code{sampsize} has been reached.When \code{cluster.only} is true, the result is simply a (possiblynamed) integer vector specifying the clustering, i.e.,\cr\code{clara(x,k, cluster.only=TRUE)} is the same as \cr\code{clara(x,k)$clustering} but computed more efficiently.}\note{%% mostly by Martin Maechler :By default, the random sampling is implemented with a \emph{very}simple scheme (with period \eqn{2^{16} = 65536}) inside the Fortrancode, independently of \R's random number generation, and as a matterof fact, deterministically. Alternatively, we recommend setting\code{rngR = TRUE} which uses \R's random number generators. Then,\code{clara()} results are made reproducible typically by using\code{\link{set.seed}()} before calling \code{clara}.The storage requirement of \code{clara} computation (for small\code{k}) is about\eqn{O(n \times p) + O(j^2)}{O(n * p) + O(j^2)} where\eqn{j = \code{sampsize}}, and \eqn{(n,p) = \code{dim(x)}}.The CPU computing time (again assuming small \code{k}) is about\eqn{O(n \times p \times j^2 \times N)}{O(n * p * j^2 * N)}, where\eqn{N = \code{samples}}.For \dQuote{small} datasets, the function \code{\link{pam}} can be useddirectly. What can be considered \emph{small}, is really a functionof available computing power, both memory (RAM) and speed.Originally (1990), \dQuote{small} meant less than 100 observations;in 1997, the authors said \emph{\dQuote{small (say with fewer than 200observations)}}; as of 2006, you can use \code{\link{pam}} withseveral thousand observations.}\author{Kaufman and Rousseeuw (see \code{\link{agnes}}), originally.Metric \code{"jaccard"}: Kamil Kozlowski (\code{@ownedoutcomes.com})and Kamil Jadeszko.All arguments from \code{trace} on, and most \R documentation and alltests by Martin Maechler.%% Kasper Fischer-Rasmussen provided \R and C code for \code{metric = "gower"}.}\seealso{\code{\link{agnes}} for background and references;\code{\link{clara.object}}, \code{\link{pam}},\code{\link{partition.object}}, \code{\link{plot.partition}}.}\examples{## generate 500 objects, divided into 2 clusters.x <- rbind(cbind(rnorm(200,0,8), rnorm(200,0,8)),cbind(rnorm(300,50,8), rnorm(300,50,8)))clarax <- clara(x, 2, samples=50)claraxclarax$clusinfo## using pamLike=TRUE gives the same (apart from the 'call'):all.equal(clarax[-8],clara(x, 2, samples=50, pamLike = TRUE)[-8])plot(clarax)## cluster.only = TRUE -- save some memory/time :clclus <- clara(x, 2, samples=50, cluster.only = TRUE)stopifnot(identical(clclus, clarax$clustering))## 'xclara' is an artificial data set with 3 clusters of 1000 bivariate## objects each.data(xclara)(clx3 <- clara(xclara, 3))## "better" number of samplescl.3 <- clara(xclara, 3, samples=100)## but that did not change the result here:stopifnot(cl.3$clustering == clx3$clustering)## Plot similar to Figure 5 in Struyf et al (1996)\dontrun{plot(clx3, ask = TRUE)}\dontshow{plot(clx3)}## Try 100 times *different* random samples -- for reliability:nSim <- 100nCl <- 3 # = no.classesset.seed(421)# (reproducibility)cl <- matrix(NA,nrow(xclara), nSim)for(i in 1:nSim)cl[,i] <- clara(xclara, nCl, medoids.x = FALSE, rngR = TRUE)$clusteringtcl <- apply(cl,1, tabulate, nbins = nCl)## those that are not always in same cluster (5 out of 3000 for this seed):(iDoubt <- which(apply(tcl,2, function(n) all(n < nSim))))if(length(iDoubt)) { # (not for all seeds)tabD <- tcl[,iDoubt, drop=FALSE]dimnames(tabD) <- list(cluster = paste(1:nCl), obs = format(iDoubt))t(tabD) # how many times in which clusters}}% end{examples}\keyword{cluster}