title: Build an approximately uniform grid from random sample (python)
body: <p>I want to build a grid from sampled data. I could use a machine learning - clustering algorithm, like k-means, but I want to restrict the centres to be roughly uniformly distributed.</p>
  
  <p>I have come up with an approach using the scikit-learn nearest neighbours search: pick a point at random, delete all points within radius r then repeat. This works well, but wondering if anyone has a better (faster) way of doing this.</p>
  
  <p>In response to comments I have tried two alternate methods, one turns out much slower the other is about the same...</p>
  
  <p>Method 0 (my first attempt):</p>
  
  <pre><code>def get_centers0(X, r): 
  
      N = X.shape[0]
      D = X.shape[1]
      grid = np.zeros([0,D])
      nearest = near.NearestNeighbors(radius = r, algorithm = 'auto')
  
      while N &gt; 0:
          nearest.fit(X)
          x = X[int(random()*N), :]
          _, del_x = nearest.radius_neighbors(x)
          X = np.delete(X, del_x[0], axis = 0)
          grid = np.vstack([grid, x])
          N = X.shape[0]
  
      return grid
  </code></pre>
  
  <p>Method 1 (using the precomputed graph):</p>
  
  <pre><code>def get_centers1(X, r): 
  
      N = X.shape[0]
      D = X.shape[1]
      grid = np.zeros([0,D])
      nearest = near.NearestNeighbors(radius = r, algorithm = 'auto')
      nearest.fit(X)
      graph = nearest.radius_neighbors_graph(X)
  
      #This method is very slow even before doing any 'pruning'
  </code></pre>
  
  <p>Method 2:</p>
  
  <pre><code>def get_centers2(X, r, k): 
  
      N = X.shape[0]
      D = X.shape[1]
      k = k
      grid = np.zeros([0,D])
      nearest = near.NearestNeighbors(radius = r, algorithm = 'auto')
  
      while N &gt; 0:
          nearest.fit(X)
          x = X[np.random.randint(0,N,k), :]
  
          #min_dist = near.NearestNeighbors().fit(x).kneighbors(x, n_neighbors = 1, return_distance = True)
          min_dist = dist(x, k, 2, np.ones(k)) # where dist is a cython compiled function
          x = x[min_dist &lt; 0.1,:]
  
          _, del_x = nearest.radius_neighbors(x)
          X = np.delete(X, del_x[0], axis = 0)
          grid = np.vstack([grid, x])
          N = X.shape[0]
  
      return grid
  </code></pre>
  
  <p>Running these as follows:</p>
  
  <pre><code>N = 50000
  r = 0.1
  x1 = np.random.rand(N)
  x2 = np.random.rand(N)
  X = np.vstack([x1, x2]).T
  
  tic = time.time()
  grid0 = get_centers0(X, r)
  toc = time.time()
  print 'Method 0: ' + str(toc - tic)
  
  tic = time.time()
  get_centers1(X, r)
  toc = time.time()
  print 'Method 1: ' + str(toc - tic)
  
  tic = time.time()
  grid2 = get_centers2(X, r)
  toc = time.time()
  print 'Method 1: ' + str(toc - tic)
  </code></pre>
  
  <p>Method 0 and 2 are about the same...</p>
  
  <pre><code>Method 0: 0.840130090714
  Method 1: 2.23365592957
  Method 2: 0.774812936783
  </code></pre>