Perform point-in-polygon operations

Given a set of point locations (e.g., field sightings of some organism) and polygons circumscribing areas of interest (e.g., U. S. National Parks), you might find yourself asking whether each point does or does not lie within any polygons. And for each point that does, you might then want to know which polygon it lies within. These kinds of containment questions fall in the domain of what are often referred to as point-in-polygon operations. In some cases you might even have points that are contained by multiple polygons, although here we'll restrict our analysis to a polygon dataset with non-overlapping areas.

As with many tasks, R provides several ways to determine which points (if any) are contained within a given polygon, and conversely which polygons (if any) contain a given point. The approach demonstrated here takes advantage of the powerful over method provided by the sp package.

Input Data

  • CSV table of (fictionalized) brown bear sightings in Alaska, each containing an arbitrary ID and spatial location specified as a lat-lon coordinate pair. These are derived from occurrence data obtained via a GBIF search, but with query results subsequently anonymized, associated spatial coordinates randomly fuzzed, and records arbitrarily shuffled and subset for demonstration purposes.
  • Polygon shapefile containing the boundaries of US National Parks greater than 100,000 acres in size. This is version 1.3.0 of the Parks and Protected Lands dataset obtained from naturalearthdata.com.

Workflow

  1. Read in the bear sightings data and tell R to treat it as a set of spatial points.
  2. Read in the National Park polygons.
  3. Identify which points lie within one of the National Parks.
  4. For each point, get the name of the containing park (if any), and add it to the bear sighting data table.
  5. Write results to file, in both CSV and ESRI Shapefile formats, and draw a map of the bear sightings and parks.

Output

  • CSV table similar to the input dataset, but with an additional column specifying the park (if any) in which the bear was sighted.
  • Point shapefile identical to the CSV, but in a format more amenable to direct manipulation in a GIS.
  • Map visualization:

Map of bear sightings

Notes

  • The 'over' method used here has replaced the now-deprecated 'overlay' method used in previous incarnations of this demo script. This new method is quite similar in spirit, but differs somewhat in the combinations of arguments it accepts, and in the form of the answer returned in each case.

R script


require(sp)
require(rgdal)
require(maps)

# read in bear data, and turn it into a SpatialPointsDataFrame
bears 

<p class="rteindent2"><em><a href="/files/scicomp/Dloads/demo-point-in-polygon.zip">Download R script and data</a> for this example</em></p>&#13;
&#13;
<h3>Learning More:</h3>&#13;
<ul><li><a href="/files/The_GIS_Primer_Buckley.pdf"><strong>The GIS Primer: The Nature of Spatial Information</strong></a></li>&#13;
</ul>

scicomp: