The Perl Toolchain Summit needs more sponsors. If your company depends on Perl, please support this very important event.

NAME

PDL::Cluster - PDL interface to the C Clustering Library

SYNOPSIS

 use PDL::Cluster;

 ##-----------------------------------------------------
 ## Data Format
 $d =   42;                     ##-- number of features
 $n = 1024;                     ##-- number of data elements

 $data = random($d,$n);         ##-- data matrix
 $elt  = $data->slice(",($i)"); ##-- element data vector
 $ftr  = $data->slice("($j),"); ##-- feature vector over all elements

 $wts  = ones($d)/$d;           ##-- feature weights
 $msk  = ones($d,$n);           ##-- missing-datum mask (1=ok)

 ##-----------------------------------------------------
 ## Library Utilties

 $mean = $ftr->cmean();
 $median = $ftr->cmedian();

 calculate_weights($data,$msk,$wts, $cutoff,$expnt,
                   $weights);

 ##-----------------------------------------------------
 ## Distance Functions

 clusterdistance($data,$msk,$wts, $n1,$n2,$idx1,$idx2,
                 $dist,
                 $distFlag, $methodFlag2);

 distancematrix($data,$msk,$wts, $distmat, $distFlag);

 ##-----------------------------------------------------
 ## Partitioning Algorithms

 getclustermean($data,$msk,$clusterids,
                $ctrdata, $ctrmask);

 getclustermedian($data,$msk,$clusterids,
                  $ctrdata, $ctrmask);

 getclustermedoid($distmat,$clusterids,$centroids,
                  $errorsums);

 kcluster($k, $data,$msk,$wts, $npass,
          $clusterids, $error, $nfound,
          $distFlag, $methodFlag);

 kmedoids($k, $distmat,$npass,
          $clusterids, $error, $nfound);

 ##-----------------------------------------------------
 ## Hierarchical Algorithms

 treecluster($data,$msk,$wts,
             $tree, $lnkdist,
             $distFlag, $methodFlag);

 treeclusterd($data,$msk,$wts, $distmat,
              $tree, $lnkdist,
              $distFlag, $methodFlag);

 cuttree($tree, $nclusters,
         $clusterids);

 ##-----------------------------------------------------
 ## Self-Organizing Maps

 somcluster($data,$msk,$wts, $nx,$ny,$tau,$niter,
            $clusterids,
            $distFlag);

 ##-----------------------------------------------------
 ## Principal Component Analysis

 pca($U, $S, $V);

 ##-----------------------------------------------------
 ## Extensions

 rowdistances($data,$msk,$wts, $rowids1,$rowids2, $distvec, $distFlag);
 clusterdistances($data,$msk,$wts, $rowids, $index2,
                  $dist,
                  $distFlag, $methodFlag);

 clustersizes($clusterids, $clustersizes);
 clusterelements($clustierids, $clustersizes, $eltids);
 clusterelementmask($clusterids, $eltmask);

 clusterdistancematrix($data,$msk,$wts,
                       $rowids, $clustersizes, $eltids,
                       $dist,
                       $distFlag, $methodFlag);

 clusterenc($clusterids, $clens,$cvals,$crowids, $k);
 clusterdec($clens,$cvals,$crowids, $clusterids, $k);
 clusteroffsets($clusterids, $coffsets,$cvals,$crowids, $k);
 clusterdistancematrixenc($data,$msk,$wts,
                          $clens1,$crowids1, $clens2,$crowids2,
                          $dist,
                          $distFlag, $methodFlag);

FUNCTIONS

cmean

  Signature: (double a(n); double [o]b())

Computes arithmetic mean of the vector $a(). See also PDL::Primitive::avg().

cmean does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

cmedian

  Signature: (double a(n); double [o]b())

Computes median of the vector $a(). See also PDL::Primitive::median().

cmedian does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

calculate_weights

  Signature: (
   double    data(d,n);
   int       mask(d,n);
   double    weight(d);
   double    cutoff();
   double    exponent();
   double [o]oweights(d);
   ; char *distFlag;
)

This function calculates weights for the features using the weighting scheme proposed by Michael Eisen:

 w[i] = 1.0 / sum_{j where dist(i,j)<cutoff} (1 - dist(i,j)/cutoff)^exponent

where the cutoff and the exponent are specified by the user.

calculate_weights does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

clusterdistance

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   int    n1();
   int    n2();
   int    index1(n1);
   int    index2(n2);
   double [o]dist();
   ; 
   char *distFlag;
   char *methodFlag;
   )

Computes distance between two clusters $index1() and $index2(). Each of the $index() vectors represents a single cluster whose values are the row-indices in the $data() matrix of the elements assigned to the respective cluster. $n1() and $n2() are the number of elements in $index1() and $index2(), respectively. Each $index$i() must have at least $n$i() elements allocated.

CAVEAT: the $methodFlag argument is interpreted differently than by the treecluster() method, namely:

a

Distance between the arithmetic means of the two clusters, as for treecluster() "f".

m

Distance between the medians of the two clusters, as for treecluster() "c".

s

Minimum pairwise distance between members of the two clusters, as for treecluster() "s".

x

Maximum pairwise distance between members of the two clusters as for treecluster() "m".

v

Average of the pairwise distances between members of the two clusters, as for treecluster() "a".

clusterdistance does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

distancematrix

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   double [o]dists(n,n);
   ; char *distFlag;
)

Compute triangular distance matrix over all data points.

distancematrix does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

getclustercentroids

  Signature: (
   double data(d,n);
   int    mask(d,n);
   int    clusterids(n);
   double [o]cdata(d,k);
   int    [o]cmask(d,k);
   ; char *ctrMethodFlag;
)

Find cluster centroids by arithmetic mean (ctrMethodFlag="a") or median over each dimension (ctrMethodFlag="m").

getclustercentroids does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

getclustermean

  Signature: (
   double data(d,n);
   int    mask(d,n);
   int    clusterids(n);
   double [o]cdata(d,k);
   int    [o]cmask(d,k);
   )

Really just a wrapper for getclustercentroids(...,"a").

getclustermedian

  Signature: (
   double data(d,n);
   int    mask(d,n);
   int    clusterids(n);
   double [o]cdata(d,k);
   int    [o]cmask(d,k);
   )

Really just a wrapper for getclustercentroids(...,"m").

getclustermedoids

  Signature: (
   double distance(n,n);
   int    clusterids(n);
   int    [o]centroids(k);
   double [o]errors(k);
   )

The getclustermedoid routine calculates the cluster centroids, given to which cluster each element belongs. The centroid is defined as the element with the smallest sum of distances to the other elements.

getclustermedoids does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

kcluster

  Signature: (
   int    nclusters();
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   int    npass();
   int    [o]clusterids(n);
   double [o]error();
   int    [o]nfound();
   ; 
   char *distFlag;
   char *ctrMethodFlag;
   )

K-Means clustering algorithm. The "ctrMethodFlag" determines how clusters centroids are to be computed; see getclustercentroids() for details.

See also: kmedoids().

kcluster does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

kmedoids

  Signature: (
   int    nclusters();
   double distance(n,n);
   int    npass();
   int    [o]clusterids(n);
   double [o]error();
   int    [o]nfound();
   )

K-Medoids clustering algorithm (uses distance matrix).

See also: kcluster().

kmedoids does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

treecluster

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   int    [o]tree(2,n);
   double [o]lnkdist(n);
   ; 
   char *distFlag;
   char *methodFlag;
   )

Hierachical agglomerative clustering.

$tree(2,n) represents the clustering solution. Each row in the matrix describes one linking event, with the two columns containing the name of the nodes that were joined. The original genes are numbered 0..(n-1), nodes are numbered -1..-(n-1). $tree(2,n) thus actually uses only (2,n-1) cells.

$lnkdist(n) represents the distance between the two subnodes that were joined. As for $tree(), $lnkdist() uses only (n-1) cells.

treecluster does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

treeclusterd

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   double distances(n,n);
   int    [o]tree(2,n);
   double [o]lnkdist(n);
   ; 
   char *distFlag;
   char *methodFlag;
   )

Hierachical agglomerative clustering using given distance matrix.

See distancematrix() and treecluster(), above.

treeclusterd does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

cuttree

  Signature: (
   int    tree(2,n);
   int    nclusters();
   int [o]clusterids(n);
   )

Cluster selection for hierarchical clustering trees.

cuttree does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

somcluster

  Signature: (
   double  data(d,n);
   int     mask(d,n);
   double  weight(d);
   int     nxnodes();
   int     nynodes();
   double  inittau();
   int     niter();
   int     [o]clusterids(2,n);
   ; char *distFlag;
)

Self-Organizing Map clustering, does not return centroid data.

somcluster does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

pca

  Signature: (
   double  [o]U(d,n);
   double  [o]S(d);
   double  [o]V(d,d);
   )

Principal Component Analysis (SVD), operates in-place on $U() and requires ($SIZE(n) >= $SIZE(d)).

pca does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

rowdistances

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   int    rowids1(ncmps);
   int    rowids2(ncmps);
   double [o]dist(ncmps);
   ; char *distFlag;
)

Computes pairwise distances between rows of $data(). $rowids1() contains the row-indices of the left (first) comparison operand, and $rowids2() the row-indices of the right (second) comparison operand. Since each of these are assumed to be indices into the first dimension $data(), it should be the case that:

 0 <= $rowids1(i),rowids2(i) < $SIZE(n)    for 0 <= i < $SIZE(ncmps)

See also clusterdistance(), clusterdistances(), clusterdistancematrixenc(), distancematrix().

rowdistances does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

clusterdistances

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   int    rowids(nr);
   int    index2(n2);
   double [o]dist(nr);
   ; 
   char *distFlag;
   char *methodFlag;
   )

Computes pairwise distance(s) from each of $rowids() as a singleton cluster with the cluster represented by $index2(), which should be an index vector as for clusterdistance(). See also clusterdistancematrixenc().

clusterdistances does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

clustersizes

  Signature: (int clusterids(n); int [o]clustersizes(k))

Computes the size (number of elements) of each cluster in $clusterids(). Useful for allocating less than maximmal space for $clusterelements().

The output piddle should never be marked BAD.

clusterelements

  Signature: (int clusterids(n); int [o]clustersizes(k); int [o]eltids(mcsize,k))

Converts the vector $clusterids() to a matrix $eltids() of element (row) indices indexed by cluster-id. $mcsize() is the maximum number of elements per cluster, at most $n. The output PDLs $clustersizes() and $eltids() can be passed to clusterdistancematrix().

clusterelements does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

clusterelementmask

  Signature: (int clusterids(n); byte [o]eltmask(k,n))

Get boolean membership mask $eltmask() based on cluster assignment in $clusterids(). No value in $clusterids() may be greater than or equal to $k. On completion, $eltmask(k,n) is a true value iff $clusterids(n)=$k.

clusterelementmask does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

clusterdistancematrix

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   int    rowids(nr);
   int    clustersizes(k);
   int    eltids(mcsize,k);
   double [o]dist(k,nr);
   ; 
   char *distFlag;
   char *methodFlag;
   )

DEPRECATED in favor of clusterdistancematrixenc(). In the future, this method is expected to become a wrapper for clusterdistancematrixenc().

Computes distance between each row index in $rowids() considered as a singleton cluster and each of the $k clusters whose elements are given by a single row of $eltids(). $clustersizes() and $eltids() are as output by the clusterelements() method.

See also clusterdistance(), clusterdistances(), clustersizes(), clusterelements(), clusterdistancematrixenc().

clusterdistancematrix does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

clusterenc

  Signature: (
   int    clusterids(n);
   int [o]clusterlens(k1);
   int [o]clustervals(k1);
   int [o]clusterrows(n);
   ;
   int k1;
   )

Encodes datum-to-cluster vector $clusterids() for efficiently mapping clusters-to-data. Returned PDL $clusterlens() holds the lengths of each cluster containing at least one element. $clustervals() holds the IDs of such clusters as they appear as values in $clusterids(). $clusterrows() is such that:

 all( rld($clusterlens, $clustervals) == $clusterids )

... if all available cluster-ids are in use.

If specified, $k1 is a perl scalar holding the number of clusters (maximum cluster index + 1); an appropriate value will guessed from $clusterids() otherwise.

Really just a wrapper for some lower-level PDL and PDL::Cluster calls.

clusterdec

  Signature: (
   int    clusterlens(k1);
   int    clustervals(k1);
   int    clusterrows(n);
   int [o]clusterids(n);
   )

Decodes cluster-to-datum vectors ($clusterlens,$clustervals,$clusterrows) into a single datum-to-cluster vector $clusterids(). $(clusterlens,$clustervals,$clusterrows) are as returned by the clusterenc() method.

Un-addressed row-index values in $clusterrows() will be assigned the pseudo-cluster (-1) in $clusterids().

Really just a wrapper for some lower-level PDL calls.

clusteroffsets

  Signature: (
   int    clusterids(n);
   int [o]clusteroffsets(k1+1);
   int [o]clustervals(k1);
   int [o]clusterrows(n);
   ;
   int k1;
   )

Encodes datum-to-cluster vector $clusterids() for efficiently mapping clusters-to-data. Like clusterenc(), but returns cumulative offsets instead of lengths.

Really just a wrapper for clusterenc(), cumusumover(), and append().

clusterdistancematrixenc

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   int    clens1(k1);
   int    crowids1(nc1);
   int    clens2(k2);
   int    crowids2(nc2);
   double [o]dist(k1,k2);
   ; 
   char *distFlag;
   char *methodFlag;
   )

Computes cluster-distance between each pair of clusters in (sequence($k1) x sequence($k2)), where 'x' is the Cartesian product. Cluster contents are passed as pairs ($clens(),$crowids()) as returned by the clusterenc() function (assuming that the $cvals() vector returned by clusterenc() is a flat sequence).

The deprecated method clusterdistancematrix() can be simulated by this function in the following manner: if a clusterdistancematrix() call was:

 clustersizes   ($cids, $csizes=zeroes(long,$k));
 clusterelements($cids, $celts =zeroes(long,$csizes->max)-1);
 clusterdistancematrix($data,$msk,$wt, $rowids, $csizes,$celts,
                       $cdmat=zeroes(double,$k,$rowids->dim(0)),
                       $distFlag, $methodFlag
                      );

Then the corresponding use of clusterdistancematrixenc() would be:

 ($clens,$cvals,$crows) = clusterenc($cids);
 clusterdistancematrixenc($data,$msk,$wt,
                          $clens,        $crows,   ##-- "real" clusters in output dim 0
                          $rowids->ones, $rowids,  ##-- $rowids as singleton clusters in output dim 1
                          $cdmat=zeroes(double,$clens->dim(0),$rowids->dim(0)),
                          $distFlag, $methodFlag);

If your $cvals() are not a flat sequence, you will probably need to do some index-twiddling to get things into the proper shape:

 if ( !all($cvals==$cvals->sequence) || $cvals->dim(0) != $k )
 {
   my $cdmat0 = $cdmat;
   my $nr     = $rowids->dim(0);
   $cdmat     = pdl(double,"inf")->slice("*$k,*$nr")->make_physical(); ##-- "missing" distances are infinite
   $cdmat->dice_axis(0,$cvals) .= $cdmat0;
 }

$distFlag and $methodFlag are interpreted as for clusterdistance().

See also clusterenc(), clusterdistancematrix().

clusterdistancematrixenc does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

clusterdistancesenc

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   int    coffsets1(k1);
   int    crowids1(nc1);
   int    cwhich1(ncmps);
   int    coffsets2(k2);
   int    crowids2(nc2);
   int    cwhich2(ncmps);
   double [o]dists(ncmps);
   ; 
   char *distFlag;
   char *methodFlag;
   )

Computes cluster-distance between selected pairs of co-indexed clusters in ($cwhich1,$cwhich2). Cluster contents are passed as pairs ($coffsetsX(),$crowidsX()) as returned by the clusteroffsets() function.

$distFlag and $methodFlag are interpreted as for clusterdistance().

See also clusterenc(), clusterdistancematrixenc().

clusterdistancesenc does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

getclusterwsum

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double clusterwts(k,n);
   double [o]cdata(d,k);
   int    [o]cmask(d,k);
   )

Find cluster centroids by weighted sum. This can be considered an expensive generalization of the getclustermean() and getclustermedian() functions. Here, the input PDLs $data() and $mask(), as well as the output PDL $cdata() are as for getclustermean(). The matrix $clusterwts() determines the relative weight of each data row in determining the centroid of each cluster, potentially useful for "fuzzy" clustering. The equation used to compute cluster means is:

 $cdata(d,k) = sum_{n} $clusterwts(k,n) * $data(d,n) * $mask(d,n)

For centroids in the same range as data elements, $clusterwts() should sum to 1 over each column (k):

 all($clusterwts->xchg(0,1)->sumover == 1)

getclustermean() can be simulated by instantiating $clusterwts() with a uniform distribution over cluster elements:

 $clusterwts = zeroes($k,$n);
 $clusterwts->indexND(cat($clusterids, xvals($clusterids))->xchg(0,1)) .= 1;
 $clusterwts /= $clusterwts->xchg(0,1)->sumover;
 getclusterwsum($data,$mask, $clusterwts, $cdata=zeroes($d,$k));

Similarly, getclustermedian() can be simulated by setting $clusterwts() to 1 for cluster medians and otherwise to 0. More sophisticated centroid discovery methods can be computed by this function by setting $clusterwts(k,n) to some estimate of the conditional probability of the datum at row $n given the cluster with index $k: p(Elt==n|Cluster==k). One way to achieve such an estimate is to use (normalized inverses of) the singleton-row-to-cluster distances as output by clusterdistancematrix().

getclusterwsum does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

attachtonearest

  Signature: (
   double data(d,n);
   int    mask(d,n);
   double weight(d);
   int    rowids(nr);
   double cdata(d,k);
   int    cmask(d,k);
   int    [o]clusterids(nr);
   double [o]cdist(nr);
   ; 
   char *distFlag;
   char *methodFlag;
   )

Assigns each specified data row to the nearest cluster centroid. Data elements are given by $data() and $mask(), feature weights are given by $weight(), as usual. Cluster centroids are defined by by $cdata() and $cmask(), and the indices of rows to be attached are given in the vector $rowids(). The output vector $clusterids() contains for each specified row index the identifier of the nearest cluster centroid. The vector $cdist() contains the distance to the best clusters.

See also: clusterdistancematrix(), attachtonearestd().

attachtonearest does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

attachtonearestd

  Signature: (
   double cdistmat(k,n);
   int rowids(nr);
   int [o]clusterids(nr);
   double [o]dists(nr);
   )

Assigns each specified data row to the nearest cluster centroid, as for attachtonearest(), given the datum-to-cluster distance matrix $cdistmat(). Currently just a wrapper for a few PDL calls. In scalar context returns $clusterids(), in list context returns the list ($clusterids(),$dists()).

checkprototypes

  Signature: (
   protos(k);
   [o]cprotos(k);
   byte [t]otmp(n);
   ; int nsize => n)

(Deterministic)

Ensure that the assignment $protos() from $k objects to integer "prototype" indices in the range [0,$n( contains no repetitions of any of the $n possible prototype values. One use for this function is the restriction of (randomly generated) potential clustering solutions for $k clusters in which each cluster is represented by a "prototypical" element from a data sample of size $n.

Requires: $n >= $k.

checkprototypes does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

checkpartitions

  Signature: (
   part(n);
   [o]cpart(n);
   [t]ptmp(k);
   ; int ksize => k)

(Deterministic)

Ensure that the partitioning $part() of $n objects into $k bins (identified by integer values in the range [0,$k-1]) contains at least one instance of each of the $k possible values. One use for this function is the restriction of (randomly generated) potential clustering solutions for $n elements into $k clusters to those which assign at least one element to each cluster.

Requires: $n >= $k.

checkpartitions does not process bad values. It will set the bad-value flag of all output piddles if the flag is set for any of the input piddles.

randomprototypes

  Signature: (int k; int n; [o]prototypes(k))

Generate a random set of $k prototype indices drawn from $n objects, ensuring that no object is used more than once. Calls checkprototypes().

See also: checkprototypes(), randomassign(), checkpartitions(), randompartition().

randompartition

  Signature: (int k; int n; [o]partition(n))

Generate a partitioning of $n objects into $k clusters, ensuring that every cluster contains at least one object. Calls checkpartitions(). This method is identical in functionality to randomassign(), but may be faster if $k is significantly smaller than $n.

See also: randomassign(), checkpartitions(), checkprototypes(), randomprototypes().

COMMON ARGUMENTS

Many of the functions described above require one or more of the following parameters:

d

The number of features defined for each data element.

n

The number of data elements to be clustered.

k
nclusters

The number of desired clusters.

data(d,n)

A matrix representing the data to be clustered, double-valued.

mask(d,n)

A matrix indicating which data values are missing. If mask(i,j) == 0, then data(i,j) is treated as missing.

weights(d)

The (feature-) weights that are used to calculate the distance.

Warning: Not all distance metrics make use of weights; you must provide some nonetheless.

clusterids(n)

A clustering solution. $clusterids() maps data elements (row indices in $data()) to values in the range [0,$k-1].

Distance Metrics

Distances between data elements (and cluster centroids, where applicable) are computed using one of a number of built-in metrics. Which metric is to be used for a given computation is indicated by a character flag denoted above with $distFlag(). In the following, w[i] represents a weighting factor in the $weights() matrix, and $W represents the total of all weights.

Currently implemented distance metrics and the corresponding flags are:

e

Pseudo-Euclidean distance:

 dist_e(x,y) = 1/W * sum_{i=1..d} w[i] * (x[i] - y[i])^2

Note that this is not the "true" Euclidean distance, which is defined as:

 dist_E(x,y) = sqrt( sum_{i=1..d} (x[i] - y[i])^2 )
b

City-block ("Manhattan") distance:

 dist_b(x,y) = 1/W * sum_{i=1..d} w[i] * |x[i] - y[i]|
c

Pearson correlation distance:

 dist_c(x,y) = 1-r(x,y)

where r is the Pearson correlation coefficient:

 r(x,y) = 1/d * sum_{i=1..d} (x[i]-mean(x))/stddev(x) * (y[i]-mean(y))/stddev(y)
a

Absolute value of the correlation,

 dist_a(x,y) = 1-|r(x,y)|

where r(x,y) is the Pearson correlation coefficient.

u

Uncentered correlation (cosine of the angle):

 dist_u(x,y) = 1-r_u(x,y)

where:

 r_u(x,y) = 1/d * sum_{i=1..d} (x[i]/sigma0(x)) * (y[i]/sigma0(y))

and:

 sigma0(w) = sqrt( 1/d * sum_{i=1..d} w[i]^2 )
x

Absolute uncentered correlation,

 dist_x(x,y) = 1-|r_u(x,y)|
s

Spearman's rank correlation.

 dist_s(x,y) = 1-r_s(x,y) ~= dist_c(ranks(x),ranks(y))

where r_s(x,y) is the Spearman rank correlation. Weights are ignored.

k

Kendall's tau (does not use weights).

 dist_k(x,y) = 1 - tau(x,y)
(other values)

For other values of dist, the default (Euclidean distance) is used.

For hierarchical clustering, the 'link method' must be specified by a character flag, denoted above as $methodFlag. Known link methods are:

s

Pairwise minimum-linkage ("single") clustering.

Defines the distance between two clusters as the least distance between any two of their respective elements.

m

Pairwise maximum-linkage ("complete") clustering.

Defines the distance between two clusters as the greatest distance between any two of their respective elements.

a

Pairwise average-linkage clustering (centroid distance using arithmetic mean).

Defines the distance between two clusters as the distance between their respective centroids, where each cluster centroid is defined as the arithmetic mean of that cluster's elements.

c

Pairwise centroid-linkage clustering (centroid distance using median).

Identifies the distance between two clusters as the distance between their respective centroids, where each cluster centroid is computed as the median of that cluster's elements.

(other values)

Behavior for other values is currently undefined.

For the first three, either the distance matrix or the gene expression data is sufficient to perform the clustering algorithm. For pairwise centroid-linkage clustering, however, the gene expression data are always needed, even if the distance matrix itself is available.

ACKNOWLEDGEMENTS

Perl by Larry Wall.

PDL by Karl Glazebrook, Tuomas J. Lukka, Christian Soeller, and others.

C Clustering Library by Michiel de Hoon, Seiya Imoto, and Satoru Miyano.

Orignal Algorithm::Cluster module by John Nolan and Michiel de Hoon.

KNOWN BUGS

Dimensional requirements are sometimes too strict.

Passing weights to Spearman and Kendall link methods wastes space.

AUTHOR

Bryan Jurish <moocow@cpan.org> wrote and maintains the PDL::Cluster distribution.

Michiel de Hoon wrote the underlying C clustering library for cDNA microarray data.

COPYRIGHT

PDL::Cluster is a set of wrappers around the C Clustering library for cDNA microarray data.

  • The C clustering library for cDNA microarray data. Copyright (C) 2002-2005 Michiel Jan Laurens de Hoon.

    This library was written at the Laboratory of DNA Information Analysis, Human Genome Center, Institute of Medical Science, University of Tokyo, 4-6-1 Shirokanedai, Minato-ku, Tokyo 108-8639, Japan. Contact: michiel.dehoon 'AT' riken.jp

    See the files REAMDE.cluster, cluster.c and cluster.h in the PDL::Cluster distribution for details.

  • PDL::Cluster wrappers copyright (C) Bryan Jurish 2005-2018. All rights reserved. This package is free software, and entirely without warranty. You may redistribute it and/or modify it under the same terms as Perl itself.

SEE ALSO

perl(1), PDL(3perl), Algorithm::Cluster(3perl), cluster(1), http://bonsai.hgc.jp/~mdehoon/software/cluster/software.htm