SSB - Sum of squares between clusters

Question

I got a little confused with the squares and the sums. As far as I know, the variance or total sum of squares (TSS) is smth like

$\sum_{i}^{n} (x_i - \bar x)^2$

and the sum of squares within (SSW) is

$\sum_{j}^{K} \sum_{i}^{n} (x_i - c_j)^2 \qquad i \in C_j$ where k ist the number of clusters

and that

$TSS = SSW + SSB$

Correct so far?

I therefore can do $TSS - SSW = SSB$

But what is the direct way to get $SSB$ from a given codebook? Preferablly I'd like to know, how to get there in numpy/scipy cause reading the equations is kind of hard for me..

import numpy as np
from scipy.cluster.vq import vq

# X.shape[0] -> observations
# X.shape[1] -> features
partition, euc_distance_to_centroids = vq(X, codebook)

TSS = np.sum((X-X.mean(0))**2)      
SSW = np.sum(euc_distance_to_centroids**2)
SSB = TSS - SSW

I think I'm missing the number of observations per cluster, when doing the SSB

>>> X
array([[ 2.,  4.,  2.],
       [ 1.,  3.,  1.],
       [ 3.,  4.,  2.],
       [ 2.,  3.,  2.],
       [ 1.,  5.,  5.]])
>>> codebook
array([[ 1.  ,  3.  ,  1.  ],
       [ 2.33,  3.67,  2.  ],
       [ 1.  ,  5.  ,  5.  ]])
>>> TSS
14.800000000000001
>>> SSW
1.3333333333333333
>>> SSB
13.466666666666667
>>> ((X.mean(0)-codebook)**2).sum() # How do I put the "num_clust_obs" in here?
12.542222222222223                  # Obviously not correct..

Answer 1

Could not find a shorter way, than this one

import numpy as np
from scipy.cluster.vq import vq

X = np.array([[ 2.,  4.,  2.],
              [ 1.,  3.,  1.],
              [ 3.,  4.,  2.],
              [ 2.,  3.,  2.],
              [ 1.,  5.,  5.]])

codebook = np.array([[ 1.  ,  3.  ,  1.  ],
                     [ 2.33,  3.67,  2.  ],
                     [ 1.  ,  5.  ,  5.  ]])

partition, euc_distance_to_centroids = vq(X, codebook)

TSS = np.sum((X-X.mean(0))**2)
SSW = np.sum(euc_distance_to_centroids**2)
SSB = TSS - SSW

# The 'direct' way
B = []
c = X.mean(0)
for i in range(partition.max()+1):
    ci = X[partition == i].mean(0)
    B.append(np.bincount(partition)[i]*np.sum((ci - c)**2))
SSB_ = np.sum(B)

print(TSS, SSW, SSB, SSB_)

prints

14.8 1.3334 13.4666 13.4666666667

Answer 2

The direct way of computing SSB:

$$ \begin{array}{} SSB = \sum_{k}^Kn_k\sum_{j}^{p}(\bar x_{jk} - \bar x_{j.})^2& \text{where }\begin{array}{}n_k = \text{number of elements in cluster k}\\\bar x_{jk}=\text{mean for variable j in cluster k}\\ \bar x_{j.} = \text{Mean for variable x}\end{array} \end{array} $$

In computation we will show using R first then python:

Using faithful dataset:

R

result <- kmeans(faithful, 2)
centers <- result$center
partition <- result$cluster
SSB <- colSums((t(centers) - colMeans(faithful))**2) %*% table(partition)

# compare the computed SSB and the one returned by kmeans
c(SSB, result$betweenss)
[1] 41538.39 41538.39

python

from statsmodels.datasets import get_rdataset
from sklearn.cluster import k_means
X = datasets.get_rdataset('faithful').data.to_numpy()
centers, partition, SSW = k_means(X, 2, n_init = 'auto')
((X.mean(0) - centers)**2).sum(1).dot(np.unique(partition, return_counts = True)[1])
41538.38830431382

Note that using the data you provided, the results were abit different. But the code should definitely work. Also you could consider using matrix multiplication to do the same thing.