Skip to Content
Course content

111: Coordinate Reference Systems in R

Click on the "Edit" button in the top corner of the screen to edit your slide content.

It's a classic mistake. I've seen it happen in production code more times than I'd like to admit: someone loads a CSV of latitude and longitude coordinates, treats them as simple X and Y coordinates in a Cartesian plane, and runs a distance calculation. On the surface, the code runs without an error. The numbers look reasonable. But the underlying math is fundamentally broken because the Earth isn't a flat piece of graph paper.

The danger of Euclidean math on a sphere

Imagine you're trying to calculate the distance between the Empire State Building and Central Park. If you just take the longitude and latitude as numbers and plug them into the Pythagorean theorem—which is essentially what happens if you use the base dist() function in R on a data frame of coordinates—you are performing "Euclidean" distance. You're treating degrees as if they were meters.

# The naive, wrong way
nyc_points <- data.frame(
  lon = c(-73.9857, -73.9654), 
  lat = c(40.7484, 40.7829)
)
# This calculates distance in "decimal degrees", which is meaningless for measurement
naive_dist <- sqrt(sum((nyc_points[1,] - nyc_points[2,])^2))

The problem here is that a degree of latitude is always roughly the same distance, but a degree of longitude shrinks as you move toward the poles. If you do this in New York, you get one kind of error; if you do it in Oslo, your results will be wildly different. You aren't measuring distance; you're measuring an arc on a sphere using a ruler meant for a tabletop. It's a silent failure, which is the most dangerous kind of bug in data science.

Moving from degrees to meters

To do this right, we need to use the sf (simple features) package to tell R exactly what these numbers represent. This is where Coordinate Reference Systems (CRS) come in. Most of the world uses WGS84 (EPSG: 4326) for GPS coordinates. It's a global standard, but it's a geographic coordinate system, meaning it's defined by angles, not linear distance.

The "better way" involves two steps: defining the starting CRS and then transforming those coordinates into a projected coordinate system—one that flattens a specific part of the earth into a 2D plane using meters. For New York, we'd use a UTM (Universal Transverse Mercator) zone.

library(sf)

# 1. Define the points and assign the WGS84 CRS
nyc_sf <- st_as_sf(nyc_points, coords = c("lon", "lat"), crs = 4326)

# 2. Transform to a projected system (UTM Zone 18N for NYC)
nyc_projected <- st_transform(nyc_sf, crs = 32618)

# Now the distance is calculated in meters
actual_dist <- st_distance(nyc_projected[1,], nyc_projected[2,])

The trade-off of precision

You might ask why we don't just stay in WGS84 and use st_distance() directly. You actually can; sf is smart enough to realize it's dealing with a sphere and will use "great circle" calculations (Haversine formula). However, as soon as you start doing more complex spatial operations—like calculating the area of a park, creating a 500-meter buffer around a subway station, or intersecting two polygons—spherical math becomes computationally expensive and conceptually messy.

By transforming to a projected CRS (like UTM), you're essentially choosing a specific "map projection" that minimizes distortion for your specific area of interest. The cost is that you have to look up the correct EPSG code for your region. If you use a projection meant for Maine to analyze data in California, you'll introduce new distortions. It's a trade-off between the convenience of a global system and the mathematical accuracy of a local one.




📋 Practical Task

Calculating Accurate Distances Between NYC Landmarks

You have been given a dataset of three landmarks in New York City. Your goal is to prove why the naive approach fails and implement the correct spatial workflow.

# Setup data
landmarks <- data.frame(
  name = c("Empire State Building", "Central Park", "Statue of Liberty"),
  lon = c(-73.9857, -73.9654, -74.0445),
  lat = c(40.7484, 40.7829, 40.6892)
)

Your Task:

  • Calculate the "naive" distance between the Empire State Building and the Statue of Liberty using the basic distance formula on the raw longitude and latitude.
  • Convert the landmarks data frame into an sf object using the WGS84 coordinate system (EPSG: 4326).
  • Transform that object into the UTM Zone 18N projection (EPSG: 32618).
  • Calculate the actual distance in meters between the Empire State Building and the Statue of Liberty.
  • Print both the naive result and the actual result to see the magnitude of the error.
Rating
0 0

There are no comments for now.

to be the first to leave a comment.