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

Issue

I have the following map boundaries in this .gdb folder

and here I have a csv which contains the variables that I want to plot and the coordinates of the points that need to be displayed on the map. My final goal is to create a map with polygons and inside every polygon there should be points according to the coordinates. Every polygon should be colored according to the count of studentid (students) for the year 2019. Any alternative is accepted

I believe that the 1st code chunk below is correct:

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")

Here Im trying to include the student data as well but Im concluding in adataframe with zero data

#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)

If everything above is correct i should do something like:

# 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)

Solution

Your file still doesn't download.

I can give you a generic guide to use ggplot2 to make a map. This will draw the polygons and the points.

You need to modify the Spatial_DataFrames with fortify() to get them into a format ggplot2 can use.

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"   

You can customize the plot with 'color' or 'fill' in the aes().

Answered By - shea