A small R script to see whether you should go beyond 50th Street in West Philly.

west-philly-incidents.R

library(data.table)
library(lubridate)
library(dplyr)
library(ggmap)
library(stargazer)

# Data: https://www.opendataphilly.org/dataset/crime-incidents

philly <- fread("incidents.csv",
                select = c("location_block","text_general_code","dispatch_time","dispatch_date", "lat","lng"))

# Subset incidents in West Philly / University City
uc <- philly[philly$lat < 39.966 & philly$lat > 39.942 & philly$lng > -75.25 & philly$lng < -75.19,]


# Create variable with year
uc$year <- year(uc$dispatch_date)

# Subset year 2018
uc <- uc[uc$year == 2018,]


# Table 1
t1 <- as.data.frame(head(sort(table(uc$text_general_code),decreasing = T),10))
stargazer(t1, out = "tab_1.tex",summary = F)

# Table 2
t2 <- as.data.frame(head(sort(table(uc$location_block),decreasing = T),10))
stargazer(t2, out = "tab_2.tex",summary = F)

# Plot 1
qmplot(lng, lat, data = uc, geom = "blank", 
       zoom = 15, maptype = "toner", darken = 0.1) %>% 
  ggsave(filename = "plot_1.png")

# Plot 2
uc %>%
  # Subsetting
  filter(text_general_code %in% c("Other Assaults",
                                       "Aggravated Assault No Firearm",
                                       "Robbery Firearm")) %>% 
  # Plotting
qmplot(lng, lat, data = ., geom = "blank", 
       zoom = 15, maptype = "toner", darken = 0.1, legend = "topright") +
  stat_density_2d(aes(fill = ..level..), geom = "polygon", alpha = .3, color = NA) +
  scale_fill_gradient2("Assault and Robbery\nIncidents in 2018:", low = "green", mid = "yellow", high = "red", midpoint = 700)

ggsave(filename = "plot_2.png")