通过将 .shp 文件与另一个具有坐标的数据框结合使用，创建具有彩色多边形和坐标点的地图

Question

中有以下地图边界

这里我有一个 csv，其中包含我要绘制的变量和需要在地图上显示的点的坐标。我的最终目标是创建一个带有多边形的地图，并且在每个多边形内部应该有根据坐标的点。每个多边形都应根据 2019 年的 studentid（学生）数量进行着色。接受任何替代方案

我相信下面的第一个代码块是正确的：

library(sf)
library(tidyverse)
library(data.table)
library(tigris)
library(leaflet)
library(mapview)

options(tigris_use_cache = TRUE)

# To keep enough digits on coords
options(digits = 11)

#coordinate reference system (long-lat system)
cr_sys = 4326
# Shp file for hs boundaries (constitutes overall district bounds)
hs_bounds <- st_read("C:/Users/makis/Documents/school/TPS_schools.shp")


# Read the feature class
#fc <- readOGR(dsn=fgdb )
#fc <- spTransform(fc, CRS("+proj=longlat +datum=WGS84 +no_defs"))

# Convert hs_bounds into longlat coord system
hs_bounds <- hs_bounds %>% 
  st_transform(4326)


tmp <- list.files(pattern = "school_report_data_fake.csv")

raw_master <- lapply(tmp,
                     function(x) read_csv(x,guess_max = 5000)) %>% 
  rbindlist(., fill = TRUE)


# r blocks in tps
tps_blocks <- blocks(state = "OK") %>% 
  st_as_sf() %>% 
  st_transform(crs = 4326) %>% 
  st_intersection(hs_bounds)

tps_bgs <- block_groups(state = "OK") %>% 
  st_as_sf() %>% 
  st_transform(crs = 4326) %>% 
  st_intersection(hs_bounds)

mapview(hs_bounds) 

# Display all tps block groups on interactive map
tps_blocks_map <- mapview(tps_bgs) %>% 
  addFeatures(., hs_bounds)

# convert to df and remove geometry bc its a list col
tps_blocks_df <- tps_blocks %>% 
  as.data.frame() %>% 
  select(-geometry)


# Export blocks in tps. GEOID10 is the unique identifier for the block
write_csv(tps_blocks_df, path = "C:/Users/makis/Documents/school/tps_blocks.csv")

我在这里也尝试包括学生数据，但我在一个数据帧为零的数据框中结束

#r students by geography

student_geos <- raw_master %>% 

  #filter for students active in a given year
  filter(year == 2019) %>% 


  # filter(row_number() %in% sample(length(year), 20000)) %>% 
  # Parse lat/long. I believe that I should do something here with the lat and long 
#and some variable of the csv like the geocode variable that is used here
#a similar should be present in my csv file as well
  #mutate(lat = as.numeric(str_extract(geocode, "[0-9]+.[0-9]+"))) %>% 
  #mutate(lon = as.numeric(str_extract(geocode, "-[0-9]+.[0-9]+"))) %>% 



  # Please don't ask me why this rowwise is necessary
  rowwise() %>% 

  # Create sf point for each set of coords
  mutate(pt = st_sfc(st_point(x = c(lon, lat)), crs = 4326)) %>% 

  # Turn df into sfc then take intersection of pts and blocks
  st_as_sf() %>% 
  st_intersection(tps_blocks)

# convert to df and remove geometry bc its a list col
student_geos_df <- student_geos %>% 
  as.data.frame() %>% 
  select(-pt)

如果上面的一切都是正确的，我应该这样做：

# enrollment by tract
tract_enrol <- student_geos %>% 
  as.data.frame() %>% 
  group_by(year, TRACTCE10) %>% 
  summarize(enrollment = n())

# convert list of tracts into sfc
tracts <- tracts(state = "OK",
                 county = c("Tulsa", "Osage", "Wagoner", "Creek"),
                 year = 2010) %>% 
  st_as_sf() %>% 
  as.data.frame() %>% 
#I guess student id instead of TRACTE10 here
  inner_join(tract_enrol, by = "TRACTCE10") %>% 
  st_as_sf()



mapview(tracts, zcol = "enrollment", legend = TRUE)

Answer 1

您的文件仍未下载。

我可以为您提供使用 ggplot2 制作地图的通用指南。这将绘制多边形和点。

您需要将 Spatial_DataFrames 修改为 fortify() 以使它们成为 ggplot2 可以使用的格式。

library(ggplot2)

hs_b2 <- fortify(hs_bounds) #or instead of "hs_bounds", "tracts" or whatever your polygon 
                            #is called. If that doesn't work you need "<-as.data.frame()". 
                            #Make sure the output has a separate column for x and y.

#repeat for the points (student) object.
student_2 <- fortify(<studentpointsobject>) 

ggplot(data=<student_2>, aes(x=x,y=y)) +  
    geom_polygon(data=hs_b2, aes(x=long, y=lat, group=group) ,        #this will create the polygon
    colour="grey90", alpha=1,               #one of your color options for polygons
    fill="grey40") +                        #one of your color options for polygons
    theme(axis.ticks.y = element_blank(),axis.text.y = element_blank(), # get rid of x ticks/text
      axis.ticks.x = element_blank(),axis.text.x = element_blank()) +  # get rid of y ticks/text
    geom_point(aes(color="grey20"))         #to get the points drawn. mess with "fill" and "color"

您可以在 aes() 中使用 'color' 或 'fill' 自定义绘图。

通过将 .shp 文件与另一个具有坐标的数据框结合使用，创建具有彩色多边形和坐标点的地图

Create a map with colored polygons and coordinate points by using a .shp file in combination with another dataframe with coordinates

r

android-mapview

sf