Data and Packages

You’ll need the following packages:

  • rgdal
  • rgeos
  • dismo

You’ll need the California counties and California ozone data for this exercise. We’ll work with the sp paradigm for this activity, although it’s possible to do all of these procedures int he sf world.

Tesselations

A tessellation is a set of non-overlapping polygons that completely tiles a surface. A useful tessellation is the Voronoi tessellation, also known as Thiessen polygons.

Tessellations can be built from a SpatialPoints** object using voronoi() in the dismo package. You should install dismo if you haven’t already (check first).

My tessellation of the California ozone data looked like:

That looks cool, but it would be more helpful if it were cropped to the CA border.

Crop to CA Border

Use crop() (from raster) to crop your tesselation to the california border. My final version looked like:

Now try plotting the ozone points over your tessellation:

Centroids

It is often useful to have a measure of the center point for polygons. One such concept is the centroid, or center of gravity, of a polygon.

You can calculate polygon centroids with gCentroid() from rgeos. Try calculating the centroids of the California counties.

Check out the help entry for gCentroid() to see which arguments you need to use.

Mine looked like this:

I can plot them over the county boundaries for more context:

Notice that some of the counties have compact shapes, while others are very elongated and irregular. Let’s see how different it would look if we tesselated the centroids.

County Tesselation

Now let’s create a tessellation form the centroids.

My tessellated county centroids looked like this:

Areal Neighborhoods

In anticipation of spatial regression, let’s practice building polygon neighborhoods!

Use poly2nb() to create a neighbor list object from your county centroids SPDF. Save it to a variable called nb_cnty.

Neighborhood Summary

You can get a summary by the neighbor using summary()

## Neighbour list object:
## Number of regions: 58 
## Number of nonzero links: 288 
## Percentage nonzero weights: 8.561237 
## Average number of links: 4.965517 
## Link number distribution:
## 
##  2  3  4  5  6  7  8 
##  2  5 15 15 16  2  3 
## 2 least connected regions:
## 34 48 with 2 links
## 3 most connected regions:
## 1 32 36 with 8 links

Neighbor count histogram

You can examine the structure of the neighborhood object and see that it contains a list of the number of neighbors for each polygon. You can plot a histogram of the distribution of neighbor counts

hist(
  unlist(lapply(nb_cnty, length)) - 0.1,
  xlab = "Count of Neighbors",
  main = "Histogram of Neighbor Counts")

Plot the neighborhood network

You can make a plot of the neighborhood network.

First you need to extract the coordinates of the centroids:

ca_cen_coords = coordinates(cen_cnty)

Next I can plot the neighborhood graph and overplot the counties:

par(mar = c(0, 0, 0, 0))
plot.nb(nb_cnty, ca_cen_coords)
plot(ca_cnty, add = T)

That’s a little hard to read, try tinkering with the plot arguments to make it more readable. Here’s my attempt:

Centroid neighorhood network

Try building a polygon network from the Voronoi tessellation of the the counties. Does it have the same structure?

## Neighbour list object:
## Number of regions: 58 
## Number of nonzero links: 286 
## Percentage nonzero weights: 8.501784 
## Average number of links: 4.931034 
## Link number distribution:
## 
##  2  3  4  5  6  7  8 
##  1  8 16  9 19  3  2 
## 1 least connected region:
## 49 with 2 links
## 2 most connected regions:
## 2 33 with 8 links