Evaluating and Forecasting for Philadelphia’s Bike Share Expansion
Rebekah Adams, Minwook Kang, Aidan Rhianne
2023-04-30
Return to MUSA 801 Projects page
This project was produced as part of the University of Pennsylvania’s Master of Urban Spatial Analytics (MUSA) Spring 2023 Practicum course (MUSA 801) taught by Michael Fichman and Matthew Harris. Along with our instructors, we would like to give special thanks to Ian Smith and Laura Culp of Bicycle Transit Systems, and Waffiyah Murray of OTIS for providing insight and support over the course of the semester.
Ⅰ. Introduction
1-1. Project Overview
Philadelphia’s bike share system, Indego, has provided communities across the city access to public bikes for more than eight years. Today, the entire network consists of more than 200 stations offering more than 2,000 bikes to residents and visitors alike. The service area, however, is heavily concentrated in Center City with pockets of stations in surrounding neighborhoods. To reach more residents, in 2022 the Office of Transportation, Infrastructure, and Sustainability (OTIS) launched Indego’s multi-year expansion plan in partnership with the system’s operator, Bicycle Transit Systems (BTS). Guided by more than 1,300 responses to Indego’s suggest-a-station feature, the expansion plan aims to further densify areas of high ridership and includes proposals to enlarge the service area by fixing bike share docks in outer parts of the city and infilling areas on the fringe.
We created Station Planner for Indego’s planners. As they plan for Indego’s expansion, our tool lets them evaluate proposed station configurations in terms of:
1. Expected ridership: The total estimated annual ridership,
2. Accessibility: The total population with access to at least one station within a 5-minute walk, and;
3. Equity of the system: The socio-economic characteristics of those with access to the network.
This report will present the process and code our team used to create Station Planner’s predictive model. First, we discuss our data exploration and feature engineering process. Through qualitative research including discussions with OTIS and BTS, we engineered to features that correlate with customers’ choice and ability to use the Indego bike share system. Then we detail our model building processes, including the evaluation of our models in terms of accuracy and generalizability. Finally, our model is tested on Indego’s 2023 Q1 ridership data to demonstrate the applicability of our model.
Ⅱ. Exploring the Data & Feature Engineering
2-1. Set up
2-1-1. Analysis Grid
To facilitate our spatial analysis, we began by creating a fishnet grid across Philadelphia, with each cell a thousand feet square, to reflect a 5-minute walk, which is the estimated maximum amount of time someone would walk to get a bike (NACTO, 2015). Dividing Philadelphia’s geography into these grid cell squares allow us to model how spatial and demographic features vary over the city, and will ultimately help us predict future ridership across space.
# ======== Packages and Functions Set Up ================
# Rmarkdown global setting
knitr::opts_chunk$set(echo = TRUE)
knitr::opts_chunk$set(message = FALSE)
knitr::opts_chunk$set(warning = FALSE)
knitr::opts_chunk$set(cache = TRUE)
knitr::opts_chunk$set(fig.align = 'center')
# Modelling packages
library(reticulate)
use_condaenv("keras-tf", required = T)
library(performanceEstimation)
library(keras)
library(tensorflow)
library(viridis)
library(gridExtra)
library(tidyverse)
library(sf)
library(LiblineaR)
library(poissonreg)
library(vip)
library(mapboxapi)
#Package installs -------------------------------------------------------------
load.fun <- function(x) {
x <- as.character(x)
if(isTRUE(x %in% .packages(all.available=TRUE))) {
eval(parse(text=paste("require(", x, ")", sep="")))
print(paste(c(x, " : already installed; requiring"), collapse=''))
} else {
#update.packages()
print(paste(c(x, " : not installed; installing"), collapse=''))
eval(parse(text=paste("install.packages('", x, "')", sep="")))
print(paste(c(x, " : installed and requiring"), collapse=''))
eval(parse(text=paste("require(", x, ")", sep="")))
}
}
########### Required Packages ###########
packages = c("bayesplot", "lme4","RcppEigen",
"tidyverse", "tidyr", "AmesHousing", "broom", "caret", "dials", "doParallel", "e1071", "earth",
"ggrepel", "glmnet", "ipred", "klaR", "kknn", "pROC", "rpart", "randomForest",
"sessioninfo", "tidymodels","ranger", "recipes", "workflows", "themis","xgboost",
"sf", "nngeo", "mapview", "spatialsample", "DiagrammeR", "vetiver", "pins")
for(i in seq_along(packages)){
packge <- as.character(packages[i])
load.fun(packge)
}
session_info()
#########################################
set.seed(717)
theme_set(theme_bw())
"%!in%" <- Negate("%in%")
g <- glimpse
# Import Packages (if you don't have one of them please do install.packages("name of packages"))
library(tidyverse)
library(sf)
library(tidycensus)
library(riem)
library(gridExtra)
library(rgeos)
library(FNN)
library(spdep)
library(caret)
library(ckanr)
library(FNN)
library(grid)
library(ggcorrplot) # plot correlation plot
library(corrplot)
library(corrr) # another way to plot correlation plot
library(kableExtra)
library(jtools) # for regression model plots
library(ggstance) # to support jtools plots
library(ggpubr) # plotting R^2 value on ggplot point scatter
library(broom.mixed) # needed for effects plots
library(knitr)
library(rmarkdown)
library(RSocrata)
library(viridis)
library(ggplot2)
library(stargazer)
library(XML)
library(data.table)
#library(ggpmisc)
library(patchwork)
library(raster)
library(classInt) # for KDE and ML risk class intervals
library(tableHTML)
library(exactextractr)
library(sp)
library(units)
library(lubridate)
library(pscl)
library(cvms)
library(yardstick)
library(plotROC)
library(gganimate)
library(gifski)
library(mapview)
library(ggspatial)
library(RColorBrewer)
library(leaflet)
library(leafsync)
library(leaflet.extras)
library(leafem)
library(shiny)
# Etc
options(scipen=999)
options(tigris_class = "sf")
# Themes
mapTheme_a <- function(base_size = 10, title_size = 12, small_size = 8) {
theme(
text = element_text( color = "black"),
plot.title = element_text(size = title_size, colour = "black", hjust = 0.5),
plot.subtitle=element_text(size = base_size, colour = "black", hjust = 0.5, face="italic"),
plot.caption=element_text(size = small_size, colour = "black", hjust = 0.5),
axis.ticks = element_blank(),
axis.title = element_blank(),
axis.text = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
panel.background = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
strip.text.x = element_text(size = small_size, face="italic"),
strip.text.y = element_text(size = small_size, face="italic"),
strip.background = element_rect(colour="transparent", fill="transparent"),
legend.title = element_text(size = small_size),
legend.text = element_text(size = small_size),
legend.key.size = unit(0.4, "cm"))
}
mapTheme_b <- function(base_size = 8, title_size = 10, small_size = 6) {
theme(
text = element_text( color = "black"),
plot.title = element_text(size = title_size, colour = "black", hjust = 0.5),
plot.subtitle=element_text(size = base_size, colour = "black", hjust = 0.5, face="italic"),
plot.caption=element_text(size = small_size, colour = "black", hjust = 0.5),
axis.ticks = element_blank(),
axis.title = element_blank(),
axis.text = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
panel.background = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
strip.text.x = element_text(size = base_size),
strip.text.y = element_text(size = base_size),
strip.background = element_rect(colour="transparent", fill="transparent"),
legend.title = element_text(size = small_size),
legend.text = element_text(size = small_size),
legend.key.size = unit(0.25, "cm"))
}
plotTheme_a <- function(base_size = 10, title_size = 12, small_size = 8){
theme(axis.text = element_blank(),
axis.ticks = element_blank(),
text = element_text(size = 10),
panel.background = element_rect(fill = greyPalette5[1]),
axis.title.x = element_text(size = small_size),
axis.title.y = element_text(size = small_size),
plot.subtitle = element_text(hjust = 0.5, size = base_size),
plot.title = element_text(hjust = 0.5, size = title_size),
plot.caption=element_text(size = small_size, colour = "black", hjust = 0.5),
strip.background = element_rect(colour="transparent", fill="transparent"))
}
plotTheme_b <- function(base_size = 10, title_size = 12, small_size = 8){
theme(axis.text = element_text(size = small_size),
text = element_text(size = 10),
panel.background = element_rect(fill = greyPalette5[1]),
axis.title.x = element_text(size = small_size),
axis.title.y = element_text(size = small_size),
plot.subtitle = element_text(hjust = 0.5, size = base_size, face="italic"),
plot.title = element_text(hjust = 0.5, size = title_size),
plot.caption=element_text(size = small_size, colour = "black", hjust = 0.5),
strip.background = element_rect(colour="transparent", fill="transparent"),
strip.text.x = element_text(size = small_size, face="italic"),
strip.text.y = element_text(size = small_size, face="italic"))
}
plotTheme_c <- function(base_size = 8, title_size = 10, small_size = 6){
theme(axis.text = element_text(size = small_size),
text = element_text(size = 8),
panel.background = element_rect(fill = greyPalette5[1]),
axis.title.x = element_text(size = small_size),
axis.title.y = element_text(size = small_size),
plot.subtitle = element_text(hjust = 0.5, size = base_size, face="italic"),
plot.title = element_text(hjust = 0.5, size = title_size),
plot.caption=element_text(size = small_size, colour = "black", hjust = 0.5),
strip.background = element_rect(colour="transparent", fill="transparent"),
strip.text.x = element_text(size = small_size, face="italic"),
strip.text.y = element_text(size = small_size, face="italic"))
}
#Functions
q5 <- function(variable) {as.factor(ntile(variable, 5))}
qBr_a <- function(df, variable, rnd) {
if (missing(rnd)) {
as.character(quantile(round(df[[variable]],0),
c(.01,.2,.4,.6,.8), na.rm=T))
} else if (rnd == FALSE | rnd == F) {
as.character(formatC(quantile(df[[variable]],
c(.01,.2,.4,.6,.8), na.rm=T),
digits = 3))
}
}
qBr_b <- function(df, variable, rnd) {
if (missing(rnd)) {
as.character(quantile(round(df[[variable]]*100,0)/100,
c(.01,.2,.4,.6,.8), na.rm=T))
} else if (rnd == FALSE | rnd == F) {
as.character(formatC(quantile(df[[variable]],
c(.01,.2,.4,.6,.8), na.rm=T),
digits = 3))
}
}
nn_function <- function(measureFrom,measureTo,k) {
measureFrom_Matrix <- as.matrix(measureFrom)
measureTo_Matrix <- as.matrix(measureTo)
nn <- get.knnx(measureTo, measureFrom, k)$nn.dist
output <- as.data.frame(nn) %>%
rownames_to_column(var = "thisPoint") %>%
gather(points, point_distance, V1:ncol(.)) %>%
arrange(as.numeric(thisPoint)) %>%
group_by(thisPoint) %>%
summarize(pointDistance = mean(point_distance)) %>%
arrange(as.numeric(thisPoint)) %>%
dplyr::select(-thisPoint) %>%
pull()
return(output)
}
myCrossValidate <- function(dataset, id, dependentVariable, indVariables) {
allPredictions <- data.frame()
cvID_list <- unique(dataset[[id]])
for (i in cvID_list) {
thisFold <- i
cat("This hold out fold is", thisFold, "\n")
fold.train <- filter(dataset, dataset[[id]] != thisFold) %>% as.data.frame() %>%
dplyr::select(id, geometry, indVariables, dependentVariable)
fold.test <- filter(dataset, dataset[[id]] == thisFold) %>% as.data.frame() %>%
dplyr::select(id, geometry, indVariables, dependentVariable)
regression <-
glm(countMVTheft ~ ., family = "poisson",
data = fold.train %>%
dplyr::select(-geometry, -id))
thisPrediction <-
mutate(fold.test, Prediction = predict(regression, fold.test, type = "response"))
allPredictions <-
rbind(allPredictions, thisPrediction)
}
return(st_sf(allPredictions))
}
`%!in%` <- Negate(`%in%`)
source("https://raw.githubusercontent.com/urbanSpatial/Public-Policy-Analytics-Landing/master/functions.r")
#color pallettes
bluePalette5 <- c("#eff3ff","#bdd7e7","#6baed6","#3182bd","#08519c")
blue2Palette5 <- c("#08519c","#3182bd","#6baed6","#bdd7e7","#eff3ff")
orangePalette5 <- c("#FFFBCB","#FFDCA0","#FEBE75","#FE9F4A","#FD801F")
greyPalette5 <- c("#f7f7f7","#cccccc","#969696","#636363","#252525")
greenPalette5 <- c("#edf8e9","#bae4b3","#74c476","#31a354","#006d2c")
purplePalette5 <- c("#f2f0f7","#cbc9e2","#9e9ac8","#756bb1","#54278f")
indego_green <- "#AFF708"
igreenPalette52 <- c("#F2FFA1","#C2ED79","#92DC51","#61CA28","#31B800")
bluePalette10 <- c('white', '#f7fbff','#deebf7','#c6dbef','#9ecae1','#6baed6','#4292c6','#2171b5','#08519c','#08306b', '#081d58') #, '#000000')
blue2Palette10 <- c('#000000', '#081d58', '#08306b', '#08519c', '#2171b5', '#4292c6', '#6baed6', '#9ecae1', '#c6dbef', '#deebf7', '#f7fbff')
#mapview palletes
MV_greenPalette <- colorRampPalette(c("#edf8e9", "#006d2c"))
# ============= Set up analysis grid =====================
# Read in Philadelphia boundary
phl_boundary <- st_read("https://opendata.arcgis.com/datasets/063f5f85ef17468ebfebc1d2498b7daf_0.geojson") %>%
st_transform("ESRI:102729") #%>% # transform so measurements are in feet
#st_union()
# set grid cell size (in feet)
gridSize <- 1000 #1000 feet
# Create grid based on Philadelphia boundary
fishnet <- st_make_grid(phl_boundary,
cellsize = gridSize,
square = TRUE) %>%
.[phl_boundary] %>%
st_sf() %>%
mutate(unique_id = 1:n()) %>%
st_transform("EPSG:4326")
# Create centroids of fishnet for later measuring
fishnet_centroids <- st_centroid(fishnet)
# Plot grid
ggplot() +
annotation_map_tile("cartolight", forcedownload = FALSE, zoom = 12)+
geom_sf(data = fishnet, fill = "transparent", color=bluePalette5[3], size=1)+
labs(title = "Analysis Grid",
subtitle = "size 1000*1000ft")+
mapTheme_a()
2-1-2. Census data
As a brief introduction to Philadelphia, we collected 2020 census data for all block groups. By clicking on the layers tab in the right hand corner of the map below, you can interact with various data including population density, employment rates, and types of transit users. In terms of racial distribution, we can see that there is a higher density of white residents in Center City, in certain East neighborhoods facing the Delaware river and in the North-West of the city moving closer towards Wissahikon Valley Park. Most Black residents live in West Philadelphia and in central northern parts of the city. The majority of Hispanic and Latinx communities reside in North Philadelphia, while Asian residents concentrate in central and southern neighborhoods.
# Call the census to view a list of available variables:
acs_variable_list.2020 <- load_variables(2020, #year
"acs5", #five year ACS estimates
cache = TRUE)
acs_vars <- c("B01001_001E", #ACS total Pop estimate
"B03002_012E", # Estimated total Hispanic or Latino
"B03002_003E", # Estimated total White, not Hispanic or Latino
"B03002_004E", # Estimated total Black, not Hispanic or Latino
"B03002_005E", # Estimated total American Indian and Alaska Native population, not Hispanic or Latino
"B03002_006E", # Estimated total Asian population, not Hispanic or Latino
"B03002_007E", # Estimated total Native Hawaiian and Other Pacific Islander alone, not Hispanic or Latino
"B03002_008E", # Estimated total other race alone, not Hispanic or Latino
"B03002_009E", # Estimated total two or more races, not Hispanic or Latino
"B03002_010E", # Estimated total two or more races including other race, not Hispanic or Latino
"B03002_011E", # Estimated total two or more races excluding other race & 3 or more races, not Hispanic or Latino
"B23025_004E", #Estimated Total population who are employed
"B19013_001E", #Estimated Median Household Income
"B01002_001E", # Median age
"B01001_031E", "B01001_032E", "B01001_033E", "B01001_034E", # Female bachelor age 18 to 24
"B01001_007E", "B01001_008E", "B01001_009E", "B01001_010E", # Male bachelor age 18 to 24
"B08301_001E", # Means of transportation to work total
"B08301_002E", # Means of transportation to work: Car
"B08301_010E", # Means of transportation to work: Public transit (bus, subway, rail)
"B08301_018E", # Means of transportation to work: Bicycle
"B08301_019E", # Means of transportation to work: walked
"B14007_002E" # Total:!!Enrolled in school
)
blockgrps <- get_acs(geography = "block group",
year = 2020,
variables = acs_vars,
geometry = TRUE,
state = 42,
county = 101,
output = "wide") %>%
dplyr::select (GEOID, all_of(acs_vars)) %>%
mutate (total_other_race = B03002_007E + B03002_005E + B03002_008E + B03002_009E + B03002_010E + B03002_011E,
fml_Bch = B01001_031E + B01001_032E + B01001_033E + B01001_034E, # Female bachelor age 18 to 24
ml_Bch = B01001_007E + B01001_008E + B01001_009E + B01001_010E) %>%
dplyr::select(-B01001_031E, -B01001_032E, -B01001_033E, -B01001_034E, -B01001_007E, -B01001_008E, -B01001_009E, -B01001_010E, -B03002_007E, -B03002_005E, -B03002_008E, -B03002_009E, -B03002_010E, -B03002_011E) %>%
rename (total_pop = B01001_001E, #ACS total Pop estimate
total_white = B03002_003E, #Estimated total White only population
total_Black = B03002_004E, # Estimated total Black population
total_Asian = B03002_006E, # Estimated total Asian population
total_HispLat = B03002_012E, # Estimated total Hispanic or Latino
total_employed = B23025_004E, #Total employed
medIncome = B19013_001E, #Median HH income
medAge = B01002_001E, #Median age
total_commutes = B08301_001E, # Means of transportation to work total
commute_car = B08301_002E, # Means of transportation to work: Car
commute_pub = B08301_010E, # Means of transportation to work: public transit
commute_bic = B08301_018E, # Means of transportation to work: Bicycle
commute_wlk = B08301_019E, # Means of transportation to work: Walked
total_school = B14007_002E # Total:!!Enrolled in school
)
blockgrps <- blockgrps %>%
mutate(area = st_area(blockgrps),
pct_white = round((total_white / total_pop) * 100,2),
pct_Black = round((total_Black / total_pop) * 100,2),
pct_Asian = round((total_Asian / total_pop) * 100,2),
pct_HispLat = round((total_HispLat / total_pop) * 100, 2),
pct_other_race = round((total_other_race / total_pop) * 100, 2),
pct_employed = round((total_employed / total_pop) *100,2),
pct_commute_car = round((commute_car/total_commutes) * 100,2),
pct_commute_pub = round((commute_pub/total_commutes) * 100,2),
pct_commute_bic = round((commute_bic/total_commutes) * 100,2),
pct_commute_wlk = round((commute_wlk/total_commutes) * 100,2),
pct_school = round((total_school / total_pop)*100,2),
pop_density = round(total_pop / area * 1000000,2) #km2
) %>%
st_transform("EPSG:4326") %>%
drop_units() %>%
dplyr::select(-total_white, -total_Black, - total_Asian, -total_HispLat, -total_employed, -total_commutes, -commute_car, -commute_pub, -commute_bic, -commute_wlk, -total_school, -total_other_race, -fml_Bch, -ml_Bch)
# Extract centroids of block groups
blockgrps_centroids <- st_centroid(blockgrps) %>%
st_sf() %>%
st_transform(st_crs(blockgrps))
# ==== Plot Census Variables ==========
pal <- colorNumeric(palette = "Blues", domain = c(0,100), na.color = NA)
pal2 <- colorNumeric(palette = "Blues", blockgrps$pop_density, na.color = NA)
m222 <- leaflet(blockgrps) %>%
addProviderTiles(providers$CartoDB.Positron) %>%
addControl(html = "<p style='font-size: 10pt; text-align:center; margin:0px; color:black;'>Figure 2-2-2. Census data</p>", position = "topright") %>%
#pct_school
addPolygons(fillColor = ~pal(pct_school), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_school, "%"), group = "Percentage of Students") %>%
#pct_commute_wlk
addPolygons(fillColor = ~pal(pct_commute_wlk), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_commute_wlk, "%"), group = "Percentage of walk commuters") %>%
#pct_commute_bic
addPolygons(fillColor = ~pal(pct_commute_bic), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_commute_bic, "%"), group = "Percentage of bicycle commuters") %>%
#pct_commute_pub
addPolygons(fillColor = ~pal(pct_commute_pub), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_commute_pub, "%"), group = "Percentage of transit commuters") %>%
#pct_commute_car
addPolygons(fillColor = ~pal(pct_commute_car), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_commute_car, "%"), group = "Percentage of car commuters") %>%
#pct_employed
addPolygons(fillColor = ~pal(pct_employed), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_employed, "%"), group = "Percentage of employed") %>%
#pct_white
addPolygons(fillColor = ~pal(pct_white), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_white, "%"), group = "Percentage of white residents") %>%
#pct_Black
addPolygons(fillColor = ~pal(pct_Black), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_Black, "%"), group = "Percentage of Black Residents") %>%
#pct_Asian
addPolygons(fillColor = ~pal(pct_Asian), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_Asian, "%"), group = "Percentage of Asian Residents") %>%
#pct_HispLat
addPolygons(fillColor = ~pal(pct_HispLat), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_HispLat, "%"), group = "Percentage of Hispanic and Latino Residents") %>%
#pct_other_race
addPolygons(fillColor = ~pal(pct_other_race), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pct_other_race, "%"), group = "Percentage of other race") %>%
#pop_density
addPolygons(fillColor = ~pal2(pop_density), fillOpacity = 0.8, weight = 1, color = "transparent", highlightOptions = highlightOptions(color = "white", weight = 2, bringToFront = TRUE), label = paste0(blockgrps$pop_density, "ppl/km2"), group = "Population density") %>%
addLayersControl(overlayGroups = c("Percentage of Students", "Percentage of walking commuters", "Percentage of bicycle commuters", "Percentage of transit commuters", "Percentage of car commuters", "Percentage of rent spent", "Percentage of employed", "Percentage of white Residents", "Percentage of Black Residents", "Percentage of Asian Residents", "Percentage of Hispanic and Latino Residents", "Percentage of other race", "Population density"), options = layersControlOptions(collapsed = TRUE)) %>%
hideGroup(c("Percentage of walk commuters", "Percentage of bicycle commuters", "Percentage of transit commuters", "Percentage of car commuters", "Percentage of employed", "Percentage of white Residents", "Percentage of Black Residents", "Percentage of Asian Residents", "Percentage of Hispanic and Latino Residents", "Percentage of other race", "Population density")) %>%
addMeasure(position = "bottomleft") %>%
addMouseCoordinates()
#sync(m222, m222)
m222