2

I am researching Lloyd iteration, an iterative algorithm that distributes points in a space. In each iteration, Lloyd's algorithm builds a Voronoi map with each input point in its own Voronoi cell, then centers each point in its Voronoi cell.

I'm seeing some strange behavior with Scipy's Voronoi implementation though: it seems that certain points are materializing out of nowhere in some iterations. The image below captures this behavior. If you look closely, you'll see two new points and cells appear in the center of the map a few iterations in:

enter image description here

Here's the code used to generate the Voronoi distributions:

from scipy.spatial import Voronoi, voronoi_plot_2d
import matplotlib.pyplot as plt
import numpy as np
import umap, os

def find_centroid(verts):
  '''Return the centroid of a polygon described by `verts`'''
  area = 0
  x = 0
  y = 0
  for i in range(len(verts)-1):
    step = (verts[i, 0] * verts[i+1, 1]) - (verts[i+1, 0] * verts[i, 1])
    area += step
    x += (verts[i, 0] + verts[i+1, 0]) * step
    y += (verts[i, 1] + verts[i+1, 1]) * step
  if area == 0: area += 0.01
  return np.array([  (1/(3*area))*x,  (1/(3*area))*y  ])

def lloyd_iterate(X):
  voronoi = Voronoi(X, qhull_options='Qbb Qc Qx')
  centroids = []
  for i in voronoi.regions:
    region = voronoi.vertices[i + [i[0]]]
    centroids.append( find_centroid( region ) )
  return np.array(centroids)

def plot(X, name):
  '''Plot the Voronoi map of 2D numpy array X'''
  v = Voronoi(X, qhull_options='Qbb Qc Qx')
  plot = voronoi_plot_2d(v, show_vertices=False, line_colors='y', line_alpha=0.5, point_size=5)
  plot.set_figheight(14)
  plot.set_figwidth(20)
  plt.axis([-10, 10, -10, 10])
  if not os.path.exists('plots'): os.makedirs('plots')
  plot.savefig( 'plots/' + str(name) + '.png' )

# get 1000 observations in two dimensions and plot their Voronoi map
X = np.random.rand(1000, 4)
X = umap.UMAP().fit_transform(X)
plot(X, 0)

# run several iterations, plotting each result
for i in range(20):
  X = lloyd_iterate(X)
  plot(X, i)

Am I overlooking something, or is there something funny going on here? Any insight others can offer would be greatly appreciated.

duhaime
  • 25,611
  • 17
  • 169
  • 224
  • Cool animation. Not sure if this is still needed, but a useful debugging technique is to keep information about how each point moves over time, and draw a polyline. If points are jumping around that will be more obvious. – tfinniga Oct 09 '18 at 22:21
  • @tfinniga thanks for your note. Is it possible for points to "jump" in a Voronoi diagram? I'm trying to understand why that would be possible... – duhaime Oct 10 '18 at 00:53
  • 1
    After looking at your code, I'm not seeing any obvious problems. I don't think that voronoi diagrams should have new points jump into existence.. if you think about them as the dual of delauney triangulations, new points shouldn't be created. It might be something weird related to the centroid computation and the point at infinity? Are you getting NaNs or very large numbers? I'd try printing out reams of info (one line per centroid computation) and look for weirdness – tfinniga Oct 10 '18 at 05:27
  • That's the thing, I didn't see any NaN's or infinity's or oddities. I didn't seed the random number generator for that trial, so I can't get the same data again unfortunately... – duhaime Oct 10 '18 at 10:38

0 Answers0