Posted to Tutorials  /  Tags: , ,

# Voronoi Diagram and Delaunay Triangulation in R

The `deldir` package by Rolf Turner makes the calculations and plotting straightforward, with a few lines of code.

A Voronoi diagram splits up a plane based on a set of original points. Each polygon, or Voronoi cell, contains an original point and all areas that are closer to that point than any other. The Delaunay triangulation of the same set of original points is the dual graph of the Voronoi diagram.

We see them around these parts occasionally.

The main thing you need to know for this R tutorial though is that the `deldir` package by Rolf Turner does most of the work for you. Calculations and plotting only takes a couple lines of code.

### Setup

R is an open source statistical computing language useful for analysis. You can also use it as a visualization tool.This tutorial makes use of R and a handful of packages, so install these first if you haven’t already.

To install R, go to the homepage, click on the link to download R, and select the file that corresponds to your system. It’s typically a straightforward installation with only a few clickety clicks.

With R installed, open the console (i.e. the app) and install the necessary packages using the `install.packages()` function. You definitely need the `deldir` package:

```install.packages("deldir", dependencies=TRUE)
```

The second part of this tutorial does some stuff with spatial data and you’ll need some packages for that.

```install.packages(c("sp", "maps", "maptools", "mapproj"), dependencies=TRUE)
```

If you don’t care about applying the Voronoi and Delaunay to spatial data though, you can skip these packages if you like.

Unzip the file and set your working directory in R to the resulting directory.Finally, download the source for this tutorial, and use it to follow along. The data also includes data that you will plot later on.

### Calculations and Plotting

```# Generate points
x <- rnorm(500, 0, 1.5)
y <- rnorm(500, 0, 1)
```

This gives you 500 x-coordinates and 500 y-coordinates sampled from normal distributions with mean 0 and standard deviation 1.5, and mean 0 and standard deviation 1, respectively.

To calculate the Voronoi tessellation and the Delaunay triangulation, just call the `deldir()` function on `x` and `y`.

```# Calculate tessellation and triangulation
vtess <- deldir(x, y)
```

Then plot on `x` and `y` to setup a plotting window, without any points.

```plot(x, y, type="n", asp=1)
```

You should see a blank plot space. Add the points:

```points(x, y, pch=20, col="red", cex=0.5)
```

Then plot the results of `deldir()`, specifying wlines as “tess” to to plot a Voronoi diagram and add to `TRUE` too add the diagram to the current plotting window.

```plot(vtess, wlines="tess", wpoints="none", number=FALSE, add=TRUE, lty=1)
```

The result:

Nice.

Similarly, you can plot the triangulation without any new calculations. Just set wlines to “triang” instead of “tess” like in the previous function call.

```# Delaunay triangulation
plot(x, y, type="n", asp=1)
plot(vtess, wlines="triang", wpoints="none", number=FALSE, add=TRUE, lty=1)
points(x, y, pch=20, col="black", cex=0.5)
```

The result:

The process isn’t much different with real data. Let’s try it out with airport locations in the conterminous United States like in Mike Bostock’s d3.js example.

Be sure to set your working directory to where you saved the tutorial download.First, load the data with `read.csv()`.

```# Load data
```

Sorry non-conterminous United States.This dataset includes all US locations, which means it includes Hawaii, Alaska, and Puerto Rico. For this example, we’re only interested in the conterminous United States, so filter accordingly.

```source("lib/latlong2state.R")
airports\$state <- latlong2state(airports[,c("longitude", "latitude")])
airports_contig <- na.omit(airports)
```

I covered this in more detail in the gridded maps tutorial.The `latlong2state()` function, which came directly from a StackOverflow answer. it takes a data frame with longitude and latitude, in that order, and returns the the state that a longitude-latitude point is in. If the location is not in the conterminous United States the function returns `NA`. So run `latlong2state()` on the airports data frame, and then omit any rows that have `NA` state.

This gives you a data frame that looks like this:

```   latitude  longitude       state
1  31.95376  -89.23450 mississippi
2  30.68586  -95.01793       texas
4  42.74135  -78.05208    new york
5  30.68801  -81.90594     florida
6  34.49167  -88.20111 mississippi
7  32.85049  -86.61145     alabama
8  43.08751  -88.17787   wisconsin
...
```

We’ve looked at map projections quite a few times in previous tutorials. You can find more here in the geographic paths tutorial. A projection isn’t necessary for this example, but we’re trying to duplicate the Bostock example. Plus, it looks nice.There are 2,988 locations total. Project the points to Albers:

```# Projection
library(mapproj)
airports_projected <- mapproject(airports_contig\$longitude, airports_contig\$latitude, "albers", param=c(39,45))
```

Like before, you can plot the (projected) points.

```# Plot
par(mar=c(0,0,0,0))
plot(airports_projected, asp=1, type="n", bty="n", xlab="", ylab="", axes=FALSE)
points(airports_projected, pch=20, cex=0.1, col="red")
```

Think of the data as a standard x-y data frame. It just happens to be geographic, where x is longitude and y is latitude.

With this in mind, you calculate the Voronoi tessellation and Delaunay triangulation in the same way with `deldir()`.

```# Voronoi
vtess <- deldir(airports_projected\$x, airports_projected\$y)
plot(vtess, wlines="tess", wpoints="none", number=FALSE, add=TRUE, lty=1)
```

Note that add in `deldir()` is set to `TRUE` to add to the current plot (a map in this case). You get what you expect.

Same deal with the Delaunay triangulation.

```# Delaunay triangulation
vtess <- deldir(airports_projected\$x, airports_projected\$y)
plot(vtess, wlines="triang", wpoints="none", number=FALSE, add=TRUE, lty=1)
```

Yay.

### Wrapping up

If you need a Voronoi tessellation or a Delaunay triangulation based on a set of points in a two-dimensional space, the `deldir` package has you covered in a couple of lines of code. The package provides other functions too like `plot.tile.list()` which lets you color cells on a per-cell basis. Get all the details on the package page.

You might also be interested in Paul Murrell’s implementation of Voronoi treemaps.

Nathan Yau is a statistician who works primarily with visualization. He earned his PhD in statistics from UCLA, is the author of two best-selling books — Data Points and Visualize This — and runs FlowingData. Introvert. Likes food. Likes beer. Follow him @flowingdata.

Become a member. Learn to visualize your data. Support FlowingData.

### Membership

This is for people who want to learn to make and design data graphics. Your support goes directly to FlowingData, an independently run site.

### Benefits of Membership

• Instant access to tutorials on how to make and design data graphics
• Source code and files to use with your own data
• Four-week course on visualization in R
• Hand-picked links and resources from around the web