slowly getting there
[imago.git] / src / k_means.py
1 """K-means module."""
2
3 import random
4
5 def cluster(k, d, data, i_centers=None):
6     """Find *k* clusters on *d* dimensional *data*."""
7     if i_centers:
8         old_centers = i_centers
9     else:
10         borders = [(min(p[0][i] for p in data), max(p[0][i] for p in data))
11                for i in range(d)]
12         old_centers = [[(h - l) * random.random() + l for (l, h) in borders]
13                for _ in range(k)]
14     clusters, centers = next_step(old_centers, data)
15     while delta(old_centers, centers) > 0:
16         old_centers = centers
17         clusters, centers = next_step(old_centers, data)
18     dst = lambda c, p: sum((a - b) ** 2 for (a, b) in zip(p, c)) ** 0.5
19     score = sum([sum(map(lambda p: dst(c, p[0]), clus)) for clus, c in
20                  zip(clusters, centers)])
21     return clusters, score
22
23 def next_step(centers, data):
24     """Compute new clusters and centers."""
25     clusters = [[] for _ in centers]
26     for point in data:
27         clusters[nearest(centers, point)].append(point)
28     centers = [centroid(c) for c in clusters]
29     return clusters, centers
30
31 def nearest(centers, point):
32     """Find the nearest cluster *center* for *point*."""
33     d, i = min(((sum((p - c) ** 2 for (p, c) in zip(point[0], center)) ** 0.5 ,
34                 index) if center else (float('inf'), len(centers)))
35                for (index, center) in enumerate(centers))
36     return i
37
38 def centroid(cluster):
39     """Find the centroid of the *cluster*."""
40     # TODO is this just a mean of coordinates?
41     # TODO should we try different definitions?
42     l = float(len(cluster))
43     try:
44         d = len(cluster[0][0]) #TODO empty cluster error
45     except IndexError:
46         return None
47     return [sum(c[0][i] for c in cluster) / l for i in range(d)]
48
49 def delta(c1, c2):
50     """Find the absolute distance between two lists of points."""
51     # TODO rewrite this to a sane form
52     return sum((sum(abs(cc1 - cc2) for (cc1, cc2) in zip (ccc1, ccc2)) if ccc2
53                else 0.) for (ccc1, ccc2) in zip(c1, c2))