Find the number of clusters in a list of integers

user_1422096 0
Viewed 478

Let's consider the distance d(a, b) = number of digits which are pairwise different in a and b, e.g.:

d(1003000000, 1000090000) = 2   #   the 4th and 6th digits don't match

(we only work with 10-digit numbers) and this list:

L = [2678888873,
     2678878873,  # distance 1 from L[0]
     1000000000,  
     1000040000,  # distance 1 from L[2]
     1000300000,  # distance 1 from L[2], distance 2 from L[3]
     1000300009,  # distance 1 from L[4], distance 2 from L[2]
    ]

I would like to find the minimal number of points P such that each integer in the list is at a distance <= 1 of a point in P.

Here I think this number is 3: every number in the list is at distance <= 1 of 2678888873, 1000000000, or 1000300009.

I imagine an O(n^2) algorithm is possible by first computing a distance matrix i.e. M[i, j] = d(L[i], L[j]).

Is there a better way to do this, especially using Numpy? (maybe there's a built-in algorithm in Numpy/Scipy?)

PS: If we see these 10-digit integers as strings, we're close to finding a minimal number of clusters in a list of many words with a Levenshtein distance.

PS2: I know realize this distance has a name on strings: Hamming distance.

cluster-analysis nearest-neighbor levenshtein-distance numpy python
1 Answers

Let's see what we know from a the distance metric. Given a number P (not necessarily in L), if two members of L are within distance 1 of P, they each share 9 digits with P, but not necessarily the same ones, so they are only guaranteed to share 8 digits with each other. So any two numbers that have distance 2 are guaranteed to two unique Ps that are distance 1 from each of them (and distance 2 from each other as well). You can use this information to reduce the amount of brute-force effort required to optimize the selection of P.

Let's say you have a distance matrix. You can immediately discard rows (or columns) that don't have entries less than 3: they are their own cluster automatically. For the remaining entries that are equal to 2, construct a list of possible P values. Find the number of elements of L that are within 1 of each element of P (another distance matrix). Sort P by the number of neighbors, and select. You will need to update the matrix at each iteration as you remove members with maximal neighbors to avoid inefficient grouping due to overlap (members of L that are near multiple members of P).

You can compute a distance matrix for L in numpy by first converting it to a 2D array of digits:

L = np.array([2678888873, 2678878873, 1000000000, 1000040000, 1000300000, 1000300009])

z = 10     # Number of digits
n = len(L) # Number of numbers

dec = 10**np.arange(z).reshape(-1, 1).astype(np.int64)
digits = (L // dec) % 10

digits is now a 10xN array:

array([[3, 3, 0, 0, 0, 9],
       [7, 7, 0, 0, 0, 0],
       [8, 8, 0, 0, 0, 0],
       [8, 8, 0, 0, 0, 0],
       [8, 7, 0, 4, 0, 0],
       [8, 8, 0, 0, 3, 3],
       [8, 8, 0, 0, 0, 0],
       [7, 7, 0, 0, 0, 0],
       [6, 6, 0, 0, 0, 0],
       [2, 2, 1, 1, 1, 1]], dtype=int64)

You can compute the distance between digits and itself, or digits and any other 10xM array using != and sum along the right axis:

distance = (digits[:, None, :] != digits[..., None]).sum(axis=0)

The result:

array([[ 0,  1, 10, 10, 10, 10],
       [ 1,  0, 10, 10, 10, 10],
       [10, 10,  0,  1,  1,  2],
       [10, 10,  1,  0,  2,  3],
       [10, 10,  1,  2,  0,  1],
       [10, 10,  2,  3,  1,  0]])

We are only concerned with the upper (or lower) triangle of that matrix, so we can immediately mask out the other triangle:

distance[np.tril_indices(n)] = z + 1

Find all candidate values of P: all elements of L, but also all pairs between elements that have distance 2:

# Find indices of pairs that differ by 2
indices = np.nonzero(distance == 2)
# Extract those numbers as 10xKx2 array
d = digits[:, np.stack(indices, axis=1)]
# Compute where the difference is nonzero (Kx2)
locs = np.diff(d, axis=2).astype(bool).squeeze()
# Find the index of the first digit to replace (K)
s = np.argmax(locs, axis=0)

The extra values of P are constructed from each half of d, with the digits represented by k replaced from the other half:

P0 = digits[:, indices[0]]
P1 = digits[:, indices[1]]
k = np.arange(s.size)
tmp = P0[s, k]
P0[s, k] = P1[s, k]
P1[s, k] = tmp

Pextra = np.unique(np.concatenate((P0, P1), axis=1)

So now you can compute the total set of possibilities for P:

P = np.concatenate((digits, Pextra), axis=1)
distance2 = (P[:, None, :] != digits[..., None]).sum(axis=0)

You can discard any elements of Pextra that match with elements of digits based on the distance:

mask = np.concatenate((np.ones(n, bool), distance2[:, n:].all(axis=0)))
P = P[:, mask]
distance2 = distance2[:, mask]

Now you can iteratively distance P with L, and select the best values of P, removing any values that have been selected from the distance matrix. A greedy selection from P will not necessarily be optimal, since an alternative combination may require fewer elements due to overlaps, but that is a matter for a simple (but somewhat expensive) graph traversal algorithm. The following snippet just shows a simple greedy selection, which will work fine for your toy example:

distMask = distance2 <= 1
quality = distMask.sum(axis=0)

clusters = []
accounted = 0
while accounted < n:
    # Get the cluster location
    best = np.argmax(quality)
    # Get the cluster number
    clusters.append(P[:, best].dot(dec).item())
    # Remove numbers in cluser from consideration
    accounted += quality[best]
    quality -= distMask[distMask[:, best], :].sum(axis=0)

The last couple of steps can be optimized using sets and graphs, but this shows a starting point for a valid approach. This is going to be slow for large data, but probably not prohibitively so. Do some benchmarks to decide how much time you want to spend optimizing vs just running the algorithm.

user_29887300
Related