# Cluster tendency using Hopkins statistic implementation in Python

The Hopkins statistic, is a statistic which gives a value which indicates the cluster tendency, in other words: how well the data can be clustered.

If the value is between {0.01, ...,0.3}, the data is regularly spaced.

If the value is around 0.5, it is random.

If the value is between {0.7, ..., 0.99}, it has a high tendency to cluster.

I have a question about my implementation of the Hopkins statistic.

Is it correct? If so, other people can use it :)

X is the data with shape (n,m).

d = len(vars) # columns
n = len(X) # rows
m = int(0.1 * n) # heuristic from article [1]

from sklearn.neighbors import NearestNeighbors
nbrs = NearestNeighbors(n_neighbors=1, algorithm='brute').fit(X)

from random import sample
rand_X = sample(range(0, n, 1), m)

ujd = []
wjd = []
for j in range(0, m):
u_dist, _ = nbrs.kneighbors(np.random.normal(size=(1, d)).reshape(1, -1), 2, return_distance=True)
ujd.append(u_dist[0][1])
w_dist, _ = nbrs.kneighbors(X[rand_X[j]].reshape(1, -1), 2, return_distance=True)
wjd.append(w_dist[0][1])

H = sum(ujd) / (sum(ujd) + sum(wjd))
print H


Any recommendations are much appreciated.

[1] Validating Clusters using the Hopkins Statistic from IEEE 2004.

• what does variable m means here? Jul 7, 2018 at 13:58
• The amount of samples. You take a subset of your data to decrease computing time. Jul 9, 2018 at 5:48
• Why do you take the second nearest neighbor for u_dist?
– sds
Mar 30, 2020 at 17:42
• Good question. Perhaps because the first nearest neighbor is the considered data point itself. The second nearest neighbor is then the closest point to the considered point. To be honest, might be worth checking. Apr 1, 2020 at 7:11
• When I look at the source document [1], the distances are raised to the dth power. I don't know what the python function to calculate distance is giving us exactly, but it looks like the distances are raised to the 1st power instead of the dth power. Any comments? Jan 28, 2022 at 22:39

Almost correct, the only problem is that you are using np.random to randomly generate points in feature space. This means that you assume that the points will be in the range of [0,1]. I've tweaked your code to check the min/max range of input data and adjust accordingly. The you can find the updated version here (with link to this post of course)

• Actually it is good practice in machine learning to normalize your features. Normalization is done to not have a group of features influence the clustering just because these features have higher values than the other. . Please vote for the question, I need the points :) Jun 24, 2017 at 8:49

There are a few errors in this implementation:

1. In your u_dist you are sampling from a normal distribution. You should be sampling from a uniform distribution .
2. You are sampling also from a standard normal which is centered around 0, while you should be sampling in the entire space of the data (so to use min and max on your data to find the edges).
3. In the u_dist you are taking the 2nd closest neighbor, instead of the 1st.
4. You should exponentiate the distances to the dimension of your data.

In addition, the implementation could be improved:

• You are using for loops instead of vector operations.
• You are mixing numpy arrays with python lists, which is IMO not very elegant.
• The results are not reproducible, as they depend on the internal seed of numpy/python.

Here is how your code would look like after the changes:

import numpy as np
from sklearn.neighbors import NearestNeighbors

def hopkins(X, portion=0.1, seed=247):
# X: numpy array of shape (n_samples, n_features)
n = X.shape[0]
d = X.shape[1]
m = int(portion * n)

np.random.seed(seed)
nbrs = NearestNeighbors(n_neighbors=1).fit(X)
# u_dist
rand_X = np.random.uniform(X.min(axis=0), X.max(axis=0), size=(m,d))
u_dist = nbrs.kneighbors(rand_X, return_distance=True)[0]
# w_dist
idx = np.random.choice(n, size=m, replace=False)
w_dist = nbrs.kneighbors(X[idx,:], 2, return_distance=True)[0][:,1]

U = (u_dist**d).sum()
W = (w_dist**d).sum()
H = U / (U + W)
return H


# Clustered Data

For clustered data (and here clustered means - not uniformly spread across the space; so even a sample from a normal distribution is considered extremely clustered) you would get something like this:

c = 2 # how much clusters are far from each other
means = c*torch.tensor([[0.,0], [1.,1], [1,-1], [-1,1], [-1,-1]])  # draw from 5 clusters
cov = torch.tensor([[1.,0.],[0,3]]) # 2d covariance matrix
# draw the samples
n = 500  # actual samples are x5 this
samples = D.MultivariateNormal(means, cov).sample((n,))
# plot the clusters
X = samples.reshape(-1,2).numpy()
y = np.tile(range(1,6), n)
plt.scatter(X[:,0], X[:,1], c=y, s=5)


# calculate the hopkins statistic
hopkins(X)
> 0.8983462887476471


# Uniform

For uniform data you should get something around 0.5:

X = D.Uniform(torch.tensor([-3.,-3]), torch.tensor([3,3.])).sample((5*n,))
plt.scatter(X[:,0], X[:,1], s=5)


hopkins(X.numpy())
>0.47382462294058925


# Grid

For a grid you should get a very low statistic, even around 0:

x = np.linspace(-3, 3, 50)
y = np.linspace(-3, 3, 50)
xv, yv = np.meshgrid(x, y)
X = np.c_[xv.reshape(-1,1), yv.reshape(-1,1)]
plt.scatter(X[:,0], X[:,1], s=5)


hopkins(X)
> 0.13711754341114704

• Thank you for the improvement and clear explanation. Much appreciated! Apr 20, 2023 at 11:09