Identifying k_means clusters number programmaticaly

I am stuck on finding a way to determine the elbow point (the optimal number of clusters to be used) programmaticaly. I need to run k-means on a set of 2D points obtained from an image and group points in regions of interest. There are a different number of region of interest in different images and I have to process a lot of data so manually going and assigning k for every image is not a good solution.

For now I have this code:

from sklearn.cluster import KMeans
from keras import backend as keras

results = model.predict(img, verbose=1, steps=1)

for x in range(0,results.shape):
for y in range(0,results.shape):
if(results[x][y] >= 0.48):
k_means_coords.append([x,y])
else:
results[x][y] = 0

k_means_coords = np.asarray(k_means_coords)
print(k_means_coords)
Sum_of_squared_distances = []

K = range(1, 15)
for k in K:
km = KMeans(n_clusters=k)
km = km.fit(k_means_coords)
Sum_of_squared_distances.append(km.inertia_)

k = get_min_k(Sum_of_squared_distances, K)
plt.plot(K, Sum_of_squared_distances, 'bx-')
plt.xlabel('k')
plt.ylabel('Sum_of_squared_distances')
plt.title('Elbow Method For Optimal k')
plt.show()


It generates this kind of plots:

here a good answer would be 4 and here it would be 3 and this is how results images(matrix) looks • As currently written, the code you've given doesn't suffice to be a minimum working example (MWE). I can guess that you're probably using scikit-learn or maybe keras for some of your functions, but you need to provide some explicit context and detail about the type of data you're working with and the processing you're doing beyond the current code snippet. – Don Walpola Jan 12 at 0:57
• @DonWalpola it is kind of clear that I am using scikit and keras with tf, nevertheless I will add imports and the result image. – Dan Butmalai Jan 12 at 9:30

Just to give idea, we actually try to catch the elbow by looking at slope changes. If the slope of the curve stops abruptly changing, we set it as elbow point by visual inspection. You can simulate such visual inspection by looking at the second-derivative, which captures the changes in first-derivatives, i.e. slopes. Of course, discrete approximations to derivatives (e.g. in this case first derivative is $$x_i-x_{i-1}$$) can be sometimes quite misleading, but in this case, I believe it'll provide some insight. And, with this (discrete approximation) method it might not be possible to select $$K=2$$.

• It is possibly to have clearly 2 regions of clustering, but having something that selects the elbow point automatically is better than not having anything, can you point me to a resource that explains more in depth the solution that uses second derivatives (how i would calculate them from what I have). P.S I added the input image for kmeans. – Dan Butmalai Jan 12 at 9:37
• 2nd order derivative can be approximated as difference of first order derivatives when x-axis step is $1$: $x_i-x_{i-1}-(x_{i-1}-x_{i-2})=x_i+x_{i-2}-2x_{i-1}$. In order to have $x_{i-2}$, i needs to start from $2$, assuming zero-indexing. – gunes Jan 12 at 11:54

This worked out for me: In a for split array of kmeans errors for k = 1->max_number_of_clusters into left and right half like this: (left half)[0:i], (right half)[i-1:len(errors)], find two lines that fit left and right part with np.polyfit, calculate error for both halves. The ellbow point should match the i with minimum error.

Here is the code that performs above described technique:

k_means_data_x = []
k_means_data_y = []

for i in  range(0,len(seq)):
k_means_data_x.append(float(i)+1)
k_means_data_y.append(seq[i])

print( k_means_data_x, k_means_data_y)
min_error = 100500
min_error_index = -1
for i in range(1,len(seq)- 1):
m1, b1 = np.polyfit(k_means_data_x[0:i], k_means_data_y[0:i],1)
m2, b2 = np.polyfit(k_means_data_x[i-1:len(k_means_data_x)], k_means_data_y[i-1:len(k_means_data_x)],1)
print(m1,b1,m2,b2)
plt_y1 = []
plt_y2 = []
error = 0;
for j,poz in zip(k_means_data_x[0:i], range(0,i)):
plt_y1.append(j*m1 + b1)
error += abs(j*m1 + b1 - k_means_data_y[poz])

for j,poz in zip(k_means_data_x[i-1:len(k_means_data_x)],range(i-1,len(k_means_data_x))):
plt_y2.append(j * m2 + b2)
error += abs(j*m2 + b2 - k_means_data_y[poz])
if(min_error > error):
min_error = error
min_error_index = i
print(error, i)
plt.plot(k_means_data_x[0:i], plt_y1, linestyle='solid')
plt.plot(k_means_data_x[i-1:len(k_means_data_x)], plt_y2, linestyle='solid')
print(min_error, min_error_index)


for this image it outputs the following lines and errors(the answer "elbow point" is in the last lnie):

1.3672511428441319 1.3672511428441323 -0.2613457080345315 2.102245992257679
3.145611252164474 1
-1.1445939857833 3.8790962714715653 -0.1656026364489833 1.4639588483540245
1.5498933924053357 2
-0.9193897561561183 3.578823965301989 -0.10029481867056794 1.0068041239051169
0.9236155539718396 3
-0.7355901524672226 3.272491292487162 -0.062334197574044625 0.7284262358639461
0.903667759176304 4
-0.5715963253985074 2.9445036383497336 -0.051881945352224346 0.6482923021633239
1.4737313403253316 5
-0.45611774931386856 2.6750536274855756 -0.03906840030601211 0.545783941793626
1.882568731573826 6
-0.3730698412790054 2.4535925393926084 -0.02716622171240527 0.4465991201802356
2.358273198845083 7 And here is a visualization of lines produced by algorithm: 