I have a locus L of points (lat, long). And I would like to find N=10 points (let's call them warehouses) such that:

$$loss = \sum_{l \in L} maximum_{w \in W}(distance(l, w))^2 $$

is minimized.

Is there a documented algorithm or approach that solves this problem? Right now I am thinking Excel may be able to handle this task. However I have too much data for Excel and will need to implement this in Python / Pandas.

  • $\begingroup$ You haven't really defined what you mean by distance. If the points are close together then you may be alright using Euclidean distance but technically distance on a sphere (as with lat,lon) should use Haversine Distance. $\endgroup$
    – Eumenedies
    Aug 14, 2017 at 8:18

2 Answers 2


I can tell you how I would do it, but there is almost certainly a faster implementation.

Assuming you start with, for each point in $L$, the distances to each warehouse, $w \in W$. These distances should be calculated by the haversine formula.

You can find the distance to the $N$th closest point in $w$ by using the quickselect algorithm. This is very similar to the quicksort algorithm but only sorting the parts that you care about. The average case for quickselect is $O(N)$ but you'll need to repeat for each $l \in L$.

Note that, since the square is monotonic for positive distances, you only need to minimise

$$\sum_{l \in L} maximum_{w \in W}(distance(l, w)) $$

I found a handy implementation of the quickselect algorithm on KoderDojo

  • $\begingroup$ "you only need to minimise..." - but this is the main problem here - how to minimize that expression? $\endgroup$
    – VividD
    Aug 14, 2017 at 13:53
  • $\begingroup$ Just pointing out that it is not necessary to square the haversine distances which might save a bit of time. If you were using Euclidean distance then it would be less complex to leave them as squared distance but I believe haversine gives you the distance rather than squared distance. $\endgroup$
    – Eumenedies
    Aug 14, 2017 at 14:23

Scipy has tools most most of this already.

locations = train[['latitude', 'longitude']].values
center = locations.mean(axis=0)
warehouses = np.repeat(np.expand_dims(center, 0), 20, 0)
warehouses = warehouses.flatten()

def Distances(warehouses):
    warehouses = np.reshape(warehouses, [20, 2])
    distances = scipy.spatial.distance.cdist(locations, warehouses)
    closests = distances.min(axis=1)
    other_way = distances.min(axis=0)
    return np.append(closests, other_way)

x = scipy.optimize.least_squares(Distances, warehouses, verbose=2)
  • $\begingroup$ If quickselect can be used here it will be a lot faster. $\endgroup$ Aug 14, 2017 at 12:52
  • $\begingroup$ Nice find. However, I think minimax and least square optimization do not yield the same result. $\endgroup$
    – VividD
    Aug 14, 2017 at 13:51

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.