random_cluster
Signature: (short [o]cluster(o,c); int obs=>o; int clu=>c)
Creates masks for random mutually exclusive clusters. Accepts two parameters, num_obs and num_cluster.
Extra parameter turns into extra dim in mask. May loop a long time if num_cluster approaches num_obs
because empty cluster is not allowed.
my $cluster = random_cluster( $num_obs, $num_cluster );
which_cluster
Signature: (short a(o,c); indx [o]b(o))
Given cluster mask dim [obs x clu], returns the cluster index to which an obs belong.
Does not support overlapping clusters. If an obs has TRUE value for multiple clusters, the returned index
is the first cluster the obs belongs to. If an obs has no TRUE value for any cluster, the return val is
set to -1 or BAD if the input mask has badflag set.
Usage:
# create a cluster mask dim [obs x clu]
perldl> p $c_mask = iv_cluster [qw(a a b b c c)]
[
[1 1 0 0 0 0]
[0 0 1 1 0 0]
[0 0 0 0 1 1]
]
# get cluster membership list dim [obs]
perldl> p $ic = $c_mask->which_cluster
[0 0 1 1 2 2]
which_cluster processes bad values. It will set the bad-value flag of all output ndarrays if the flag is
set for any of the input ndarrays.
assign
Signature: (data(o,v); centroid(c,v); short [o]cluster(o,c))
Takes data pdl dim [obs x var] and centroid pdl dim [cluster x var] and returns mask dim [obs x cluster]
to cluster membership. An obs is assigned to the first cluster with the smallest distance (ie sum squared
error) to cluster centroid. With bad value, obs is assigned by smallest mean squared error across
variables.
perldl> $centroid = ones 2, 3
perldl> $centroid(0,) .= 0
perldl> p $centroid
[
[0 1]
[0 1]
[0 1]
]
perldl> $b = qsort( random 4, 3 )
perldl> p $b
[
[0.022774068 0.032513883 0.13890034 0.30942479]
[ 0.16943853 0.50262636 0.56251531 0.7152271]
[ 0.23964483 0.59932745 0.60967495 0.78452117]
]
# notice that 1st 3 obs in $b are on average closer to 0
# and last obs closer to 1
perldl> p $b->assign( $centroid )
[
[1 1 1 0] # cluster 0 membership
[0 0 0 1] # cluster 1 membership
]
assign processes bad values. It will set the bad-value flag of all output ndarrays if the flag is set
for any of the input ndarrays.
centroid
Signature: (data(o,v); cluster(o,c); float+ [o]m(c,v); float+ [o]ss(c,v))
Takes data dim [obs x var] and mask dim [obs x cluster], returns mean and ss (ms when data contains bad
values) dim [cluster x var], using data where mask == 1. Multiple cluster membership for an obs is okay.
If a cluster is empty all means and ss are set to zero for that cluster.
# data is 10 obs x 3 var
perldl> p $d = sequence 10, 3
[
[ 0 1 2 3 4 5 6 7 8 9]
[10 11 12 13 14 15 16 17 18 19]
[20 21 22 23 24 25 26 27 28 29]
]
# create two clusters by value on 1st var
perldl> p $a = $d( ,(0)) <= 5
[1 1 1 1 1 1 0 0 0 0]
perldl> p $b = $d( ,(0)) > 5
[0 0 0 0 0 0 1 1 1 1]
perldl> p $c = cat $a, $b
[
[1 1 1 1 1 1 0 0 0 0]
[0 0 0 0 0 0 1 1 1 1]
]
perldl> p $d->centroid($c)
# mean for 2 cluster x 3 var
[
[ 2.5 7.5]
[12.5 17.5]
[22.5 27.5]
]
# ss for 2 cluster x 3 var
[
[17.5 5]
[17.5 5]
[17.5 5]
]
centroid processes bad values. It will set the bad-value flag of all output ndarrays if the flag is set
for any of the input ndarrays.
kmeans
Implements classic k-means cluster analysis. Given a number of observations with values on a set of
variables, kmeans puts the observations into clusters that maximizes within-cluster similarity with
respect to the variables. Tries several different random seeding and clustering in parallel. Stops when
cluster assignment of the observations no longer changes. Returns the best result in terms of R2 from the
random-seeding trials.
Instead of random seeding, kmeans also accepts manual seeding. This is done by providing a centroid to
the function, in which case clustering will proceed from the centroid and there is no multiple tries.
There are two distinct advantages from seeding with a centroid compared to seeding with predefined
cluster membership of a subset of the observations ie "seeds",
(1) a centroid could come from a previous study with a different set of observations;
(2) a centroid could even be "fictional", or in more proper parlance, an idealized prototype with respect
to the actual data. For example, if there are 10 person's ratings of 1 to 5 on 4 movies, ie a ratings pdl
of dim [10 obs x 4 var], providing a centroid like
[
[5 0 0 0]
[0 5 0 0]
[0 0 5 0]
[0 0 0 5]
]
will produce 4 clusters of people with each cluster favoring a different one of the 4 movies. Clusters
from an idealized centroid may not give the best result in terms of R2, but they sure are a lot more
interpretable.
If clustering has to be done from predefined clusters of seeds, simply calculate the centroid using the
centroid function and feed it to kmeans,
my ($centroid, $ss) = $rating($iseeds, )->centroid( $seeds_cluster );
my %k = $rating->kmeans( { CNTRD=>$centroid } );
kmeans supports bad value*.
Default options (case insensitive):
V => 1, # prints simple status
FULL => 0, # returns results for all seeding trials
CNTRD => PDL->null, # optional. pdl [clu x var]. disables next 3 opts
NTRY => 5, # num of random seeding trials
NSEED => 1000, # num of initial seeds, use NSEED up to max obs
NCLUS => 3, # num of clusters
Usage:
# suppose we have 4 person's ratings on 5 movies
perldl> p $rating = ceil( random(4, 5) * 5 )
[
[3 2 2 3]
[2 4 5 4]
[5 3 2 3]
[3 3 1 5]
[4 3 3 2]
]
# we want to put the 4 persons into 2 groups
perldl> %k = $rating->kmeans( {NCLUS=>2} )
# by default prints back options used
# as well as info for all tries and iterations
CNTRD => Null
FULL => 0
NCLUS => 2
NSEED => 4
NTRY => 5
V => 1
ss total: 20.5
iter 0 R2 [0.024390244 0.024390244 0.26829268 0.4796748 0.4796748]
iter 1 R2 [0.46341463 0.46341463 0.4796748 0.4796748 0.4796748]
perldl> p "$_\t$k{$_}\n" for (sort keys %k)
R2 0.479674796747968
centroid # mean ratings for 2 group x 5 movies
[
[ 3 2.3333333]
[ 2 4.3333333]
[ 5 2.6666667]
[ 3 3]
[ 4 2.6666667]
]
cluster # 4 persons' membership in two groups
[
[1 0 0 0]
[0 1 1 1]
]
n [1 3] # cluster size
ss
[
[ 0 0.66666667]
[ 0 0.66666667]
[ 0 0.66666667]
[ 0 8]
[ 0 0.66666667]
]
Now, for the valiant, kmeans is threadable. Say you gathered 10 persons' ratings on 5 movies from 2
countries, so the data is dim [10,5,2], and you want to put the 10 persons from each country into 3
clusters, just specify NCLUS => [3,1], and there you have it. The key is for NCLUS to include
$data->ndims - 1 numbers. The 1 in [3,1] turns into a dummy dim, so the 3-cluster operation is repeated
on both countries. Similarly, when seeding, CNTRD needs to have ndims that at least match the data ndims.
Extra dims in CNTRD will lead to threading (convenient if you want to try out different centroid
locations, for example, but you will have to hand pick the best result). See stats_kmeans.t for examples
w 3D and 4D data.
*With bad value, R2 is based on average of variances instead of sum squared error.