Clustering

import math

# Shop customers: (visits per month, average basket in hundred taka). No labels.
PTS = [(2, 3), (3, 2), (2, 2), (3, 3), (4, 2),          # rare, small baskets
       (10, 3), (11, 4), (12, 3), (10, 4), (11, 2),     # frequent, small baskets
       (3, 12), (2, 13), (4, 14), (3, 13), (4, 12),     # rare, big baskets
       (20, 20)]                                         # one unusual customer

def dist(a, b):
    return math.hypot(a[0] - b[0], a[1] - b[1])

def kmeans(k, iters=20):
    # start from spread-out points: the first point, then repeatedly the farthest from those chosen
    centers = [PTS[0]]
    while len(centers) < k:
        centers.append(max(PTS, key=lambda p: min(dist(p, c) for c in centers)))
    for _ in range(iters):
        groups = [[] for _ in range(k)]
        for p in PTS:                                   # step 1: assign each point to its nearest center
            groups[min(range(k), key=lambda i: dist(p, centers[i]))].append(p)
        centers = [(sum(x for x, _ in g) / len(g), sum(y for _, y in g) / len(g)) for g in groups]  # step 2: move centers
    inertia = sum(min(dist(p, c) ** 2 for c in centers) for p in PTS)
    return centers, groups, inertia

print("k  inertia   (sum of squared distances to the nearest center)")
for k in range(1, 6):
    print(f"{k}  {kmeans(k)[2]:8.1f}")

for k in (3, 4):
    centers, groups, _ = kmeans(k)
    print(f"\nk = {k}:")
    for c, g in sorted(zip(centers, groups)):
        print(f"  center ({c[0]:5.2f}, {c[1]:5.2f})  {len(g)} customers")

# DBSCAN: a point with at least MIN neighbours within EPS starts or grows a cluster; others are noise.
EPS, MIN = 2.0, 3
labels = [None] * len(PTS)
cluster = 0
for i, p in enumerate(PTS):
    if labels[i] is not None:
        continue
    near = [j for j, q in enumerate(PTS) if dist(p, q) <= EPS]
    if len(near) < MIN:
        labels[i] = -1                                  # noise for now
        continue
    cluster += 1
    queue = near
    while queue:
        j = queue.pop()
        if labels[j] in (None, -1):
            labels[j] = cluster
            more = [m for m, q in enumerate(PTS) if dist(PTS[j], q) <= EPS]
            if len(more) >= MIN:
                queue.extend(m for m in more if labels[m] in (None, -1))
print(f"\nDBSCAN eps={EPS} min={MIN}: {cluster} clusters")
for c in sorted(set(labels)):
    members = [PTS[i] for i in range(len(PTS)) if labels[i] == c]
    print(f"  {'noise' if c == -1 else 'cluster ' + str(c)}: {len(members)} customers, first {members[0]}")

package main

import (
	"fmt"
	"math"
	"sort"
)

type pt struct{ x, y float64 }

// Shop customers: (visits per month, average basket in hundred taka). No labels.
var pts = []pt{{2, 3}, {3, 2}, {2, 2}, {3, 3}, {4, 2}, // rare, small baskets
	{10, 3}, {11, 4}, {12, 3}, {10, 4}, {11, 2}, // frequent, small baskets
	{3, 12}, {2, 13}, {4, 14}, {3, 13}, {4, 12}, // rare, big baskets
	{20, 20}} // one unusual customer

func dist(a, b pt) float64 { return math.Hypot(a.x-b.x, a.y-b.y) }

func nearest(p pt, cs []pt) int {
	best := 0
	for i := range cs {
		if dist(p, cs[i]) < dist(p, cs[best]) {
			best = i
		}
	}
	return best
}

func kmeans(k, iters int) ([]pt, [][]pt, float64) {
	// start from spread-out points: the first point, then repeatedly the farthest from those chosen
	centers := []pt{pts[0]}
	for len(centers) < k {
		far, farD := pts[0], -1.0
		for _, p := range pts {
			if d := dist(p, centers[nearest(p, centers)]); d > farD {
				far, farD = p, d
			}
		}
		centers = append(centers, far)
	}
	var groups [][]pt
	for it := 0; it < iters; it++ {
		groups = make([][]pt, k)
		for _, p := range pts { // step 1: assign each point to its nearest center
			i := nearest(p, centers)
			groups[i] = append(groups[i], p)
		}
		for i, g := range groups { // step 2: move centers
			sx, sy := 0.0, 0.0
			for _, p := range g {
				sx, sy = sx+p.x, sy+p.y
			}
			centers[i] = pt{sx / float64(len(g)), sy / float64(len(g))}
		}
	}
	inertia := 0.0
	for _, p := range pts {
		d := dist(p, centers[nearest(p, centers)])
		inertia += d * d
	}
	return centers, groups, inertia
}

func main() {
	fmt.Println("k  inertia   (sum of squared distances to the nearest center)")
	for k := 1; k <= 5; k++ {
		_, _, in := kmeans(k, 20)
		fmt.Printf("%d  %8.1f\n", k, in)
	}
	for _, k := range []int{3, 4} {
		centers, groups, _ := kmeans(k, 20)
		idx := make([]int, k)
		for i := range idx {
			idx[i] = i
		}
		sort.Slice(idx, func(a, b int) bool {
			ca, cb := centers[idx[a]], centers[idx[b]]
			return ca.x < cb.x || (ca.x == cb.x && ca.y < cb.y)
		})
		fmt.Printf("\nk = %d:\n", k)
		for _, i := range idx {
			fmt.Printf("  center (%5.2f, %5.2f)  %d customers\n", centers[i].x, centers[i].y, len(groups[i]))
		}
	}

	// DBSCAN: a point with at least MIN neighbours within EPS starts or grows a cluster; others are noise.
	const eps, minPts = 2.0, 3
	labels := make([]int, len(pts)) // 0 = unvisited, -1 = noise, 1.. = cluster
	neighbours := func(p pt) []int {
		var out []int
		for j, q := range pts {
			if dist(p, q) <= eps {
				out = append(out, j)
			}
		}
		return out
	}
	cluster := 0
	for i, p := range pts {
		if labels[i] != 0 {
			continue
		}
		near := neighbours(p)
		if len(near) < minPts {
			labels[i] = -1 // noise for now
			continue
		}
		cluster++
		queue := near
		for len(queue) > 0 {
			j := queue[len(queue)-1]
			queue = queue[:len(queue)-1]
			if labels[j] == 0 || labels[j] == -1 {
				labels[j] = cluster
				if more := neighbours(pts[j]); len(more) >= minPts {
					for _, m := range more {
						if labels[m] == 0 || labels[m] == -1 {
							queue = append(queue, m)
						}
					}
				}
			}
		}
	}
	fmt.Printf("\nDBSCAN eps=%.1f min=%d: %d clusters\n", eps, minPts, cluster)
	for c := -1; c <= cluster; c++ {
		var members []pt
		for i, l := range labels {
			if l == c {
				members = append(members, pts[i])
			}
		}
		if len(members) == 0 {
			continue
		}
		name := fmt.Sprintf("cluster %d", c)
		if c == -1 {
			name = "noise"
		}
		fmt.Printf("  %s: %d customers, first (%v, %v)\n", name, len(members), members[0].x, members[0].y)
	}
}