Reference Map: Philadelphia Shootings

# Read GeoJSON file and convert to simple feature object
shootings <- st_read("https://raw.githubusercontent.com/kyleknox3/STA553/main/week9/shootings.geojson", quiet = TRUE)
phillyNeighbor  <- st_read("https://pengdsci.github.io/STA553VIZ/w08/Neighborhoods_Philadelphia.geojson")

## Reading layer `Neighborhoods_Philadelphia' from data source 
##   `https://pengdsci.github.io/STA553VIZ/w08/Neighborhoods_Philadelphia.geojson' 
##   using driver `GeoJSON'
## Simple feature collection with 158 features and 8 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -75.28027 ymin: 39.867 xmax: -74.95576 ymax: 40.13799
## Geodetic CRS:  WGS 84

philly  <- st_read("https://pengdsci.github.io/STA553VIZ/w08/PhillyNeighborhood-blocks.geojson") # block level data

## Reading layer `PhillyNeighborhood-blocks' from data source 
##   `https://pengdsci.github.io/STA553VIZ/w08/PhillyNeighborhood-blocks.geojson' 
##   using driver `GeoJSON'
## Simple feature collection with 17555 features and 7 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -75.28027 ymin: 39.867 xmax: -74.95576 ymax: 40.13799
## Geodetic CRS:  WGS 84

# Extract latitude and longitude from the geometry column
shootings_coords <- st_coordinates(shootings)

# Rename the columns to 'longitude' and 'latitude'
colnames(shootings_coords) <- c("longitude", "latitude")

# Bind the coordinates data frame to the original data frame
shootings <- cbind(shootings, shootings_coords)

# Convert age to numeric and handle missing values
shootings$age <- as.numeric(as.character(shootings$age))
shootings$age <- ifelse(is.na(shootings$age), 0, shootings$age)

# Convert fatal to "Yes" or "No"
shootings$fatal <- ifelse(shootings$fatal == 1, "Yes", "No")

# Create age group variable
shootings$age_group <- cut(shootings$age, breaks = c(0, 18, 30, 40, 50, Inf),
                           labels = c("0-18", "19-30", "31-40", "41-50", "50+"))

# Create text for hover or click
shootings$text <- paste("Fatal: ", shootings$fatal, "<br>",
                        "Age Group: ", shootings$age_group, "<br>",
                        "Year: ", shootings$year, "<br>",
                        "Race: ", shootings$race, "<br>",
                        "Sex: ", shootings$sex, "<br>",
                        "Victims Wound: ", shootings$wound, "<br>",
                        "Offender Deceased: ", shootings$offender_deceased, "<br>",
                        "Location: ", shootings$location, "<br>")

# Define color palette for age groups
color_palette <- c("0-18" = "red", "19-30" = "blue", "31-40" = "green", "41-50" = "orange", "50+" = "purple")

# Create the interactive map
interactive_map <- plot_ly(data = shootings, lat = ~latitude, lon = ~longitude, type = "scattermapbox",
                           mode = "markers", 
                           marker = list(size = 9, color = ~age_group, colors = color_palette),
                           text = ~text, # Display text on hover or click
                           hoverinfo = "text") %>%
  layout(mapbox = list(style = "open-street-map", 
                       center = list(lon = mean(shootings$longitude, na.rm = TRUE), lat = mean(shootings$latitude, na.rm = TRUE)), 
                       zoom = 10),
         title = "Philadelphia Shootings in 2023",
         legend = list(itemsizing = "constant"))

# Display the map
interactive_map