Week 7: Combining Data From Different Sources

PPOL 6805: GIS for Spatial Data Science

Workshop Sessions
Author

Christy Hsu

Published

October 9, 2026

R Coding Workshop: 5th Meeting

Outline

  • resources
  • Combining Data Frames
  • more on df to sf, sf to df
  • debugging

Working with Multiple Data Frames

library(tidyverse)
── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.2.1     ✔ readr     2.2.0
✔ forcats   1.0.1     ✔ stringr   1.6.0
✔ ggplot2   4.0.2     ✔ tibble    3.3.1
✔ lubridate 1.9.5     ✔ tidyr     1.3.2
✔ purrr     1.2.1     
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(sf)
Linking to GEOS 3.13.0, GDAL 3.8.5, PROJ 9.5.1; sf_use_s2() is TRUE
library(terra)
terra 1.9.34

Attaching package: 'terra'

The following object is masked from 'package:tidyr':

    extract
library(mapview)
library(rnaturalearth)

Combining data frames of the same shape

  • bind_cols() and bind_rows() can combine two or more data frames
  • cbind() and rbind()
division_df <- tibble(
    division_str = c(
        'New England', 'Middle Atlantic', 'East North Central',
        'West North Central', 'South Atlantic', 'East South Central',
        'West South Central', 'Mountain', 'Pacific'
    )
)
division_code_df <- tibble(
    division_code = 1:9
)
paste(
    'Dimensions of division_df:', dim(division_df),
    'Dimensions of division_code_df', dim(division_code_df)) |> print()
[1] "Dimensions of division_df: 9 Dimensions of division_code_df 9"
[2] "Dimensions of division_df: 1 Dimensions of division_code_df 1"
division_lookup <- bind_cols(
    division_code_df, division_df
)
division_lookup
division_code division_str
1 New England
2 Middle Atlantic
3 East North Central
4 West North Central
5 South Atlantic
6 East South Central
7 West South Central
8 Mountain
9 Pacific

bind_rows():

  • check if colnames align
book_df <- read_csv('data/book-challenge11.csv')
Rows: 931 Columns: 17
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr   (3): title, author, state
dbl  (13): book_id, year, removed, explicit, antifamily, occult, language, l...
date  (1): date

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
removed_tb <- book_df |> count(state, removed) |> filter(removed == 1)
challenged_tb <- book_df |> count(state)
removed_tb |> dim() |> print()
[1] 39  3
challenged_tb |> dim() |> print()
[1] 47  2
not_removed_tb <- book_df |>
  group_by(state) |>
  summarize(n = sum(removed == 0)) |>
  mutate(var = 'not removed')
not_removed_tb |> tail()
state n var
VA 16 not removed
VT 13 not removed
WA 6 not removed
WI 7 not removed
WV 3 not removed
WY 3 not removed
removed_tb <- book_df |>
  group_by(state) |>
  summarize(n = sum(removed == 1)) |>
  mutate(var = 'removed')
removed_tb |> tail()
state n var
VA 18 removed
VT 0 removed
WA 3 removed
WI 3 removed
WV 0 removed
WY 1 removed
all_tb <- book_df |>
  count(state) |>
  mutate(var = 'all')
all_tb |> tail()
state n var
VA 34 all
VT 13 all
WA 9 all
WI 10 all
WV 3 all
WY 4 all
count_df <- bind_rows(all_tb, removed_tb, not_removed_tb)
count_df |> dim()
[1] 141   3
count_df |> arrange(state) |> tail()
state n var
WV 3 all
WV 0 removed
WV 3 not removed
WY 4 all
WY 1 removed
WY 3 not removed
removed_df <- book_df |>
  mutate(removed = as.factor(removed)) |>
  count(state, removed, .drop = FALSE)
removed_df |> dim()
[1] 94  3

Joins

  • left_join(), inner_join(), right_join(), full_join(), semi_join(), and anti_join().
  • Binary operation on dataframes
    • takes pair of data frames (x, y)
region_lookup <- tibble(
    region_code = 1:4,
    region_str = c('Northeast', 'Midwest', 'South','West')
)
regions <- state.region |> 
  recode("North Central" = "Midwest")
regions
 [1] South     West      West      South     West      West      Northeast
 [8] South     South     South     West      West      Midwest   Midwest  
[15] Midwest   Midwest   South     South     Northeast South     Northeast
[22] Midwest   Midwest   South     Midwest   West      Midwest   West     
[29] Northeast Northeast West      Northeast South     Midwest   Midwest  
[36] South     West      Northeast Northeast South     Midwest   South    
[43] South     West      Northeast South     West      South     Midwest  
[50] West     
Levels: Northeast South Midwest West
state_region_df <- tibble(
  state_abb = state.abb,
  region_str = regions
)

left_join() retains the number of rows of the data frame on the left (x)

state_region_df <- state_region_df |>
  left_join(region_lookup, join_by(region_str))
state_region_df
state_abb region_str region_code
AL South 3
AK West 4
AZ West 4
AR South 3
CA West 4
CO West 4
CT Northeast 1
DE South 3
FL South 3
GA South 3
HI West 4
ID West 4
IL Midwest 2
IN Midwest 2
IA Midwest 2
KS Midwest 2
KY South 3
LA South 3
ME Northeast 1
MD South 3
MA Northeast 1
MI Midwest 2
MN Midwest 2
MS South 3
MO Midwest 2
MT West 4
NE Midwest 2
NV West 4
NH Northeast 1
NJ Northeast 1
NM West 4
NY Northeast 1
NC South 3
ND Midwest 2
OH Midwest 2
OK South 3
OR West 4
PA Northeast 1
RI Northeast 1
SC South 3
SD Midwest 2
TN South 3
TX South 3
UT West 4
VT Northeast 1
VA South 3
WA West 4
WV South 3
WI Midwest 2
WY West 4

Joins sf and df

state_vars <- c('REGION10', 'DIVISION10', 'GEOID10', 'STUSPS10', 'NAME10', 'ALAND10', 
'AWATER10', 'INTPTLAT10', 'INTPTLON10', 'geometry')
st_layers('data/tl_2010_us_state10')
name geomtype driver features fields crs
tl_2010_us_state10 Polygon ESRI Shapefile 52 14 NAD83 , GEOGCRS[“NAD83”,
DATUM["North American Datum 1983",
    ELLIPSOID["GRS 1980",6378137,298.257222101,
        LENGTHUNIT["metre",1]]],
PRIMEM["Greenwich",0,
    ANGLEUNIT["degree",0.0174532925199433]],
CS[ellipsoidal,2],
    AXIS["latitude",north,
        ORDER[1],
        ANGLEUNIT["degree",0.0174532925199433]],
    AXIS["longitude",east,
        ORDER[2],
        ANGLEUNIT["degree",0.0174532925199433]],
ID["EPSG",4269]] |
us_state_sf <- st_read('data/tl_2010_us_state10')|>
  select(any_of(state_vars))
Reading layer `tl_2010_us_state10' from data source 
  `/Users/jpj/gtown-local/ppol6805-workshop/w07/data/tl_2010_us_state10' 
  using driver `ESRI Shapefile'
Simple feature collection with 52 features and 14 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -179.2311 ymin: 17.83151 xmax: 179.8597 ymax: 71.44106
Geodetic CRS:  NAD83
us_state_sf |> head(4)
REGION10 DIVISION10 GEOID10 STUSPS10 NAME10 ALAND10 AWATER10 INTPTLAT10 INTPTLON10 geometry
4 8 56 WY Wyoming 251470069067 1864445306 +42.9918024 -107.5419255 MULTIPOLYGON (((-108.6213 4…
1 2 42 PA Pennsylvania 115883064314 3397122731 +40.9042486 -077.8280624 MULTIPOLYGON (((-80.51909 3…
2 3 39 OH Ohio 105828706692 10269012119 +40.4149297 -082.7119975 MULTIPOLYGON (((-84.05271 3…
4 8 35 NM New Mexico 314160748240 756659673 +34.4391265 -106.1261511 MULTIPOLYGON (((-109.0462 3…
us_state_sf |> st_crs()
Coordinate Reference System:
  User input: NAD83 
  wkt:
GEOGCRS["NAD83",
    DATUM["North American Datum 1983",
        ELLIPSOID["GRS 1980",6378137,298.257222101,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["latitude",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["longitude",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    ID["EPSG",4269]]
lib_df <- read_csv('data/public-libraries-2010.csv')
Rows: 55 Columns: 10
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (1): STABR
dbl (9): POPU_LSA, POPU_ST, CENTLIB, BRANLIB, BKMOB, BKVOL, LIBRARIA, VISITS...

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
lib_df |> summary()
    STABR              POPU_LSA           POPU_ST            CENTLIB     
 Length:55          Min.   :   53883   Min.   :   53883   Min.   :  0.0  
 Class :character   1st Qu.: 1409268   1st Qu.: 1322416   1st Qu.: 54.5  
 Mode  :character   Median : 4141264   Median : 3751351   Median :100.0  
                    Mean   : 5595942   Mean   : 5663819   Mean   :165.7  
                    3rd Qu.: 6540910   3rd Qu.: 6497622   3rd Qu.:234.0  
                    Max.   :38647288   Max.   :38648090   Max.   :755.0  
    BRANLIB          BKMOB           BKVOL             LIBRARIA     
 Min.   :  1.0   Min.   : 0.00   Min.   :      -1   Min.   :  -1.0  
 1st Qu.: 28.5   1st Qu.: 2.00   1st Qu.: 4423164   1st Qu.: 206.1  
 Median : 83.0   Median : 7.00   Median : 9575400   Median : 587.0  
 Mean   :140.1   Mean   :13.44   Mean   :14772714   Mean   : 855.7  
 3rd Qu.:192.0   3rd Qu.:17.50   3rd Qu.:17074348   3rd Qu.:1063.4  
 Max.   :942.0   Max.   :81.00   Max.   :74753745   Max.   :4062.5  
     VISITS              REGBOR        
 Min.   :       -1   Min.   :      -1  
 1st Qu.:  6216945   1st Qu.:  763359  
 Median : 18282842   Median : 1938999  
 Mean   : 28628594   Mean   : 3118459  
 3rd Qu.: 39948818   3rd Qu.: 4129365  
 Max.   :178979300   Max.   :22275618  
lib_sf <- lib_df |> left_join(us_state_sf, join_by(STABR == STUSPS10))
lib_sf |> mapview()
Error:
! 
oops! Arguments xcol and/or ycol are missing!
You probably expected lib_sf to be a spatial object. 
However it is of class spec_tbl_df. 
Either convert lib_sf to a spatial object or provide xcol and ycol.
oops! Arguments xcol and/or ycol are missing!
You probably expected lib_sf to be a spatial object. 
However it is of class tbl_df. 
Either convert lib_sf to a spatial object or provide xcol and ycol.
oops! Arguments xcol and/or ycol are missing!
You probably expected lib_sf to be a spatial object. 
However it is of class tbl. 
Either convert lib_sf to a spatial object or provide xcol and ycol.
oops! Arguments xcol and/or ycol are missing!
You probably expected lib_sf to be a spatial object. 
However it is of class data.frame. 
Either convert lib_sf to a spatial object or provide xcol and ycol.

Joins df to sf

  • check missing values in geometry column
lib_sf <- lib_df |> filter(!(STABR %in% c('GU', 'MP', 'VI'))) |> left_join(us_state_sf, join_by(STABR == STUSPS10))
lib_sf <- lib_sf |>
  st_as_sf(
    geometry = lib_sf$geometry,
    crs = 4269
)
lib_sf |> mapview(zcol = 'CENTLIB')
lib_sf |> class()
[1] "sf"          "spec_tbl_df" "tbl_df"      "tbl"         "data.frame" 
# result of subsetting
lib_sf2 <- us_state_sf |> left_join(lib_df, join_by(STUSPS10 == STABR))
lib_sf2 |> mapview(zcol = 'CENTLIB')
lib_sf2 |> head()
REGION10 DIVISION10 GEOID10 STUSPS10 NAME10 ALAND10 AWATER10 INTPTLAT10 INTPTLON10 POPU_LSA POPU_ST CENTLIB BRANLIB BKMOB BKVOL LIBRARIA VISITS REGBOR geometry
4 8 56 WY Wyoming 251470069067 1864445306 +42.9918024 -107.5419255 544270 544270 23 53 2 2497545 189.16 3872783 383126 MULTIPOLYGON (((-108.6213 4…
1 2 42 PA Pennsylvania 115883064314 3397122731 +40.9042486 -077.8280624 12054201 12284183 451 183 29 27790282 1445.81 47188171 5549200 MULTIPOLYGON (((-80.51909 3…
2 3 39 OH Ohio 105828706692 10269012119 +40.4149297 -082.7119975 11551941 11551941 240 480 58 45224425 2725.65 88255852 8767349 MULTIPOLYGON (((-84.05271 3…
4 8 35 NM New Mexico 314160748240 756659673 +34.4391265 -106.1261511 1588981 2009671 91 27 1 4518130 290.37 8324986 1140637 MULTIPOLYGON (((-109.0462 3…
3 5 24 MD Maryland 25141638381 6989579585 +38.9466584 -076.6744939 5633514 5618344 15 169 18 13954140 1297.29 33662473 3271760 MULTIPOLYGON (((-75.74776 3…
1 1 44 RI Rhode Island 2677566454 1323668539 +41.5978358 -071.5252895 1430385 1052587 47 26 2 4636067 232.11 6279148 574776 MULTIPOLYGON (((-71.65321 4…

Debugging