Bergabung dengan poligon di R

Saya bertanya-tanya bagaimana cara menggabungkan poligon spasial menggunakan kode R?

Saya bekerja dengan data sensus di mana area tertentu berubah dari waktu ke waktu dan saya ingin bergabung dengan poligon dan data yang sesuai dan hanya melaporkan area yang bergabung. Saya mempertahankan daftar poligon yang memiliki perubahan sensus ke sensus dan bahwa saya berencana untuk bergabung. Saya ingin menggunakan daftar nama area ini sebagai daftar pencarian untuk diterapkan pada data sensus dari tahun yang berbeda.

Saya bertanya-tanya apa fungsi R yang digunakan untuk menggabungkan poligon yang dipilih dan data masing-masing. Saya telah meng-google-nya tetapi hanya menjadi bingung oleh hasil.

Solusi berikut didasarkan pada posting oleh Roger Bivand di R-sig-Geo . Saya mengambil contohnya mengganti shapefile Jerman dengan beberapa data sensus dari Oregon yang dapat Anda unduh dari sini (ambil semua komponen shapefile dari 'kabupaten Oregon dan data sensus').

Mari kita mulai dengan memuat paket-paket yang diperlukan dan mengimpor shapefile ke R.

# Required packages
libs <- c("rgdal", "maptools", "gridExtra")
lapply(libs, require, character.only = TRUE)

# Import Oregon census data
oregon <- readOGR(dsn = "path/to/data", layer = "orcounty")
oregon.coords <- coordinates(oregon)

Selanjutnya, Anda memerlukan beberapa variabel pengelompokan untuk mengumpulkan data. Dalam contoh kami, pengelompokan hanya didasarkan pada koordinat county tunggal. Lihat gambar di bawah, batas hitam menunjukkan poligon asli, sedangkan batas merah mewakili poligon yang dikumpulkan oregon.id.

# Generate IDs for grouping
oregon.id <- cut(oregon.coords[,1], quantile(oregon.coords[,1]), include.lowest=TRUE)

# Merge polygons by ID
oregon.union <- unionSpatialPolygons(oregon, oregon.id)

# Plotting
plot(oregon)
plot(oregon.union, add = TRUE, border = "red", lwd = 2)

Sejauh ini baik. Namun, atribut data yang terkait dengan subregional shapefile asli (mis. Kepadatan populasi, luas, dll.) Hilang saat beraksi unionSpatialPolygons. Saya kira Anda juga ingin menggabungkan data sensus Anda yang terkait dengan shapefile, sehingga Anda memerlukan langkah menengah.

Pertama-tama Anda harus mengonversi poligon Anda ke bingkai data untuk melakukan agregasi. Sekarang mari kita ambil kolom atribut data enam hingga delapan ("AREA", "POP1990", "POP1997") dan agregat sesuai dengan fungsi penerapan ID di atas sum.

# Convert SpatialPolygons to data frame
oregon.df <- as(oregon, "data.frame")

# Aggregate and sum desired data attributes by ID list
oregon.df.agg <- aggregate(oregon.df[, 6:8], list(oregon.id), sum)
row.names(oregon.df.agg) <- as.character(oregon.df.agg$Group.1)

Akhirnya, konversikan kembali kerangka data Anda ke SpatialPolygonsDataFramepenyediaan shapefile yang sebelumnya tidak disatukan oregon.uniondan Anda memperoleh poligon umum dan data sensus Anda yang diperoleh dari langkah agregasi perangkuman di atas.

# Reconvert data frame to SpatialPolygons
oregon.shp.agg <- SpatialPolygonsDataFrame(oregon.union, oregon.df.agg)

# Plotting
grid.arrange(spplot(oregon, "AREA", main = "Oregon: original county area"), 
             spplot(oregon.shp.agg, "AREA", main = "Oregon: aggregated county area"), ncol = 1)

Berikut adalah solusi menggunakan paket sf:

library(tidycensus)
library(dplyr)
library(sf)
library(ggplot2)

# get data from tindycensus for demonstration (note you need an API key, folow instructions here: https://walkerke.github.io/tidycensus/articles/basic-usage.html)
census <- tidycensus::get_acs(geography = "tract", variables = "B19013_001",
                           state = "TX", county = "Tarrant", geometry = TRUE) %>% 
  arrange(NAME)

# reduce dataset size
census <- census[1:8,]

# create grouping variable
group_1 <- census$GEOID[1:2]
group_2 <- census$GEOID[6:8]

census <- census %>% mutate(group = case_when(GEOID %in% group_1 ~ 'newgroup1',
                                              GEOID %in% group_2 ~ 'newgroup2',
                                              TRUE ~ GEOID))

# summarise by grouping variable (performs a union on grouped polygons and sums 'estimate')
census2 <- group_by(census, group) %>% 
  summarise(estimate = sum(estimate), do_union = TRUE)

# visualise using ggplot2 development version and facet by merged/unmerged datasets
plot_data <- rbind(census %>% select(group, estimate) %>%
                     mutate(facet = "unmerged"), 
                   census2 %>% mutate(facet = "merged"))

gp <- ggplot() + 
      geom_sf(data = plot_data, aes(fill = estimate), color = 'white') + 
      scale_fill_viridis_c() + 
      facet_wrap(~facet, ncol = 1)