MUSA 508 Assignment 1: TOD - Los Angeles City

1. Introduction

This assignment investigated Transit Oriented Development (TOD) potential in Los Angeles City, CA, through the analysis of space/time indicators wrangled from Census data with respect to LA Metro stations. We selected Population, Median Rent, Percent of Bachelor, and Percent of Poverty of the census tracts within 0.5 mile of each Metro stations in 2009 and 2018 as TOD indicators for comparison purpose. Besides, we also visualized the spatial patterns of rape crime incidents with stations and median rent data respectively.

2. Set up

Load Libraries

library(tidyverse)
library(tidycensus)
library(sf)
library(kableExtra)
options(scipen=999)
options(tigris_class = "sf")

Load Styling options

# Set up map theme function
mapTheme <- function(base_size = 12) {
  theme(
    text = element_text( color = "black"),
    plot.title = element_text(size = 16,colour = "black"),
    plot.subtitle=element_text(face="italic"),
    plot.caption=element_text(hjust=0),
    axis.ticks = element_blank(),
    panel.background = element_blank(),axis.title = element_blank(),
    axis.text = element_blank(),
    axis.title.x = element_blank(),
    axis.title.y = element_blank(),
    panel.grid.minor = element_blank(),
    panel.border = element_rect(colour = "black", fill=NA, size=2),
    strip.text.x = element_text(size = 14))
}
# Set up plot theme function
plotTheme <- function(base_size = 12) {
  theme(
    text = element_text( color = "black"),
    plot.title = element_text(size = 16,colour = "black"),
    plot.subtitle = element_text(face="italic"),
    plot.caption = element_text(hjust=0),
    axis.ticks = element_blank(),
    panel.background = element_blank(),
    panel.grid.major = element_line("grey80", size = 0.1),
    panel.grid.minor = element_blank(),
    panel.border = element_rect(colour = "black", fill=NA, size=2),
    strip.background = element_rect(fill = "grey80", color = "white"),
    strip.text = element_text(size=12),
    axis.title = element_text(size=12),
    axis.text = element_text(size=10),
    plot.background = element_blank(),
    legend.background = element_blank(),
    legend.title = element_text(colour = "black", face = "italic"),
    legend.text = element_text(colour = "black", face = "italic"),
    strip.text.x = element_text(size = 14)
  )
}

Load Quantile Break Functions

qBr <- 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]]), digits = 3),
                 c(.01,.2,.4,.6,.8), na.rm=T)
  }
}
q5 <- function(variable) {as.factor(ntile(variable, 5))}

Load Hexadecimal Color Palette

palette5 <- c("#f0f9e8","#bae4bc","#7bccc4","#43a2ca","#0868ac")

Load Census API Key

census_api_key("61f5998ae3ee7a7ab62b29f55769a3b864515ac4", overwrite = TRUE, install = TRUE)
readRenviron("~/.Renviron")

Load Multiple Rings Buffer Function

multipleRingBuffer <- function(inputPolygon, maxDistance, interval) 
{
  #create a list of distances that we'll iterate through to create each ring
  distances <- seq(0, maxDistance, interval)
  #we'll start with the second value in that list - the first is '0'
  distancesCounter <- 2
  #total number of rings we're going to create
  numberOfRings <- floor(maxDistance / interval)
  #a counter of number of rings
  numberOfRingsCounter <- 1
  #initialize an otuput data frame (that is not an sf)
  allRings <- data.frame()
  
  #while number of rings  counteris less than the specified nubmer of rings
  while (numberOfRingsCounter <= numberOfRings) 
  {
    #if we're interested in a negative buffer and this is the first buffer
    #(ie. not distance = '0' in the distances list)
    if(distances[distancesCounter] < 0 & distancesCounter == 2)
    {
      #buffer the input by the first distance
      buffer1 <- st_buffer(inputPolygon, distances[distancesCounter])
      #different that buffer from the input polygon to get the first ring
      buffer1_ <- st_difference(inputPolygon, buffer1)
      #cast this sf as a polygon geometry type
      thisRing <- st_cast(buffer1_, "POLYGON")
      #take the last column which is 'geometry'
      thisRing <- as.data.frame(thisRing[,ncol(thisRing)])
      #add a new field, 'distance' so we know how far the distance is for a give ring
      thisRing$distance <- distances[distancesCounter]
    }
    
    #otherwise, if this is the second or more ring (and a negative buffer)
    else if(distances[distancesCounter] < 0 & distancesCounter > 2) 
    {
      #buffer by a specific distance
      buffer1 <- st_buffer(inputPolygon, distances[distancesCounter])
      #create the next smallest buffer
      buffer2 <- st_buffer(inputPolygon, distances[distancesCounter-1])
      #This can then be used to difference out a buffer running from 660 to 1320
      #This works because differencing 1320ft by 660ft = a buffer between 660 & 1320.
      #bc the area after 660ft in buffer2 = NA.
      thisRing <- st_difference(buffer2,buffer1)
      #cast as apolygon
      thisRing <- st_cast(thisRing, "POLYGON")
      #get the last field
      thisRing <- as.data.frame(thisRing$geometry)
      #create the distance field
      thisRing$distance <- distances[distancesCounter]
    }
    
    #Otherwise, if its a positive buffer
    else 
    {
      #Create a positive buffer
      buffer1 <- st_buffer(inputPolygon, distances[distancesCounter])
      #create a positive buffer that is one distance smaller. So if its the first buffer
      #distance, buffer1_ will = 0. 
      buffer1_ <- st_buffer(inputPolygon, distances[distancesCounter-1])
      #difference the two buffers
      thisRing <- st_difference(buffer1,buffer1_)
      #cast as a polygon
      thisRing <- st_cast(thisRing, "POLYGON")
      #geometry column as a data frame
      thisRing <- as.data.frame(thisRing[,ncol(thisRing)])
      #add teh distance
      thisRing$distance <- distances[distancesCounter]
    }  
    
    #rbind this ring to the rest of the rings
    allRings <- rbind(allRings, thisRing)
    #iterate the distance counter
    distancesCounter <- distancesCounter + 1
    #iterate the number of rings counter
    numberOfRingsCounter <- numberOfRingsCounter + 1
  }
  
  #convert the allRings data frame to an sf data frame
  allRings <- st_as_sf(allRings)
}

3. Data Wrangling

Wrangling Census Tracts Data

options(tigris_use_cache = TRUE)

city_bound_LA <- st_read('https://opendata.arcgis.com/datasets/09f503229d37414a8e67a7b6ceb9ec43_7.geojson') %>%
  st_transform('ESRI:102241')

tracts09 <- 
  get_acs(geography = "tract", variables = c("B25026_001E","B02001_002E","B15001_050E",
                                             "B15001_009E","B19013_001E","B25058_001E",
                                             "B06012_002E"), 
          year=2009, state=06, county=037, geometry=T, output="wide") %>%
  st_transform('ESRI:102241') %>%
  rename(TotalPop = B25026_001E, 
         Whites = B02001_002E,
         FemaleBachelors = B15001_050E, 
         MaleBachelors = B15001_009E,
         MedHHInc = B19013_001E, 
         MedRent = B25058_001E,
         TotalPoverty = B06012_002E) %>%
  dplyr::select(-NAME, -starts_with("B")) %>%
  mutate(pctWhite = ifelse(TotalPop > 0, Whites / TotalPop,0),
         pctBachelors = ifelse(TotalPop > 0, ((FemaleBachelors + MaleBachelors) / TotalPop),0),
         pctPoverty = ifelse(TotalPop > 0, TotalPoverty / TotalPop, 0),
         year = "2009") %>%
  dplyr::select(-Whites, -FemaleBachelors, -MaleBachelors, -TotalPoverty) 
tracts09 <- tracts09[city_bound_LA,] 

tracts18 <- 
  get_acs(geography = "tract", variables = c("B25026_001E","B02001_002E","B15001_050E",
                                             "B15001_009E","B19013_001E","B25058_001E",
                                             "B06012_002E"), 
          year=2018, state=06, county=037, geometry=T, output="wide") %>%
  st_transform('ESRI:102241') %>%
  rename(TotalPop = B25026_001E, 
         Whites = B02001_002E,
         FemaleBachelors = B15001_050E, 
         MaleBachelors = B15001_009E,
         MedHHInc = B19013_001E, 
         MedRent = B25058_001E,
         TotalPoverty = B06012_002E) %>%
  dplyr::select(-NAME, -starts_with("B")) %>%
  mutate(pctWhite = ifelse(TotalPop > 0, Whites / TotalPop,0),
         pctBachelors = ifelse(TotalPop > 0, ((FemaleBachelors + MaleBachelors) / TotalPop),0),
         pctPoverty = ifelse(TotalPop > 0, TotalPoverty / TotalPop, 0),
         year = "2018") %>%
  dplyr::select(-Whites, -FemaleBachelors, -MaleBachelors, -TotalPoverty) 

tracts18 <- tracts18[city_bound_LA,] 
allTracts <- rbind(tracts09, tracts18)

Wrangling Transit Open Data

station_2018 <- st_read('Metro_Rail_Lines_Stops.shp')%>%
  select(MetroLine, Station) %>%
  st_transform(st_crs(tracts09))  

station_2009 <- 
  station_2018[!(station_2018$MetroLine=='Expo Line' 
                 | station_2018$Station=='Arcadia Station'
                 | station_2018$Station=='Monrovia Station'
                 | station_2018$Station=='Duarte / City of Hope Station'
                 | station_2018$Station=='Irwindale Station'
                 | station_2018$Station=='Azusa Downtown Station'
                 | station_2018$Station=='APU / Citrus College Station'
                 | station_2018$Station=='Pico Station'),]

metroBuffers_2009 <- 
  rbind(
    st_buffer(station_2009, 804.67) %>%
      mutate(Legend = "Buffer09") %>%
      dplyr::select(Legend),
    st_union(st_buffer(station_2009, 804.67)) %>%
      st_sf() %>%
      mutate(Legend = "Unioned Buffer 2009"))

metroBuffers_2018 <- 
  rbind(
    st_buffer(station_2018, 804.67) %>%
      mutate(Legend = "Buffer18") %>%
      dplyr::select(Legend),
    st_union(st_buffer(station_2018, 804.67)) %>%
      st_sf() %>%
      mutate(Legend = "Unioned Buffer 2018"))

buffer_2009 <- filter(metroBuffers_2009, Legend=="Unioned Buffer 2009")
buffer_2018 <- filter(metroBuffers_2018, Legend=="Unioned Buffer 2018")

Relating Metro Stations and Tracts

# 2009

Tracts.2009 <- 
  rbind(
    st_centroid(tracts09)[buffer_2009,] %>%
      st_drop_geometry() %>%
      left_join(allTracts) %>%
      st_sf() %>%
      mutate(TOD = "TOD"),
    st_centroid(tracts09)[buffer_2009, op = st_disjoint] %>%
      st_drop_geometry() %>%
      left_join(tracts09) %>%
      st_sf() %>%
      mutate(TOD = "Non-TOD"))

Tracts.2009.inf <- 
  rbind(
    st_centroid(tracts09)[buffer_2009,] %>%
      st_drop_geometry() %>%
      left_join(allTracts) %>%
      st_sf() %>%
      mutate(TOD = "TOD"),
    st_centroid(tracts09)[buffer_2009, op = st_disjoint] %>%
      st_drop_geometry() %>%
      left_join(tracts09) %>%
      st_sf() %>%
      mutate(TOD = "Non-TOD"))%>%
  mutate(MedRent.inf = ifelse(year == "2009", MedRent * 1.17, MedRent)) 

# 2018

Tracts.2018 <- 
  rbind(
    st_centroid(tracts18)[buffer_2018,] %>%
      st_drop_geometry() %>%
      left_join(tracts18) %>%
      st_sf() %>%
      mutate(TOD = "TOD"),
    st_centroid(tracts18)[buffer_2018, op = st_disjoint] %>%
      st_drop_geometry() %>%
      left_join(tracts18) %>%
      st_sf() %>%
      mutate(TOD = "Non-TOD"))

# Combine data

allTracts.group <- 
  rbind(Tracts.2009,Tracts.2018)%>%
  mutate(MedRent.inf = ifelse(year == "2009", MedRent * 1.17, MedRent)) 

4. TOD Indicators Analytics

TOD - Census Variables

Bachelor Per 10K People

allTracts.group.visual <- 
  rbind(Tracts.2009,Tracts.2018)%>%
  mutate(MedRent.inf = ifelse(year == "2009", MedRent * 1.17, MedRent))%>%
  mutate(pctPoverty = pctPoverty * 100) %>%
  mutate(pctBachelors = pctBachelors * 10000)

ggplot(allTracts.group)+
  geom_sf(data = st_union(Tracts.2009))+
  geom_sf(data = st_union(Tracts.2018))+
  geom_sf(aes(fill = q5(pctBachelors))) +
  geom_sf(data = buffer_2018, fill = "transparent", color = "red") +
  geom_sf(data = buffer_2009, fill = "transparent", color = "blue")+
  scale_fill_manual(values = palette5,
                    labels = qBr(allTracts.group.visual, "pctBachelors"),
                    name = "Bachelors Per 10K Population\n(Quintile Breaks)") +
  labs(title = "Bachelors Per 10K Population 2009-2018", 
       subtitle = "1/10000; Blue-area close to metro stations (2009); red-area close to new stations (2018)",
       caption = "TOD Indicators - Bachelor Per 10K | Data: US Cencus Bureau, City of Los Angeles Open Data") +
  facet_wrap(~year)+
  mapTheme() + 
  theme(plot.title = element_text(size=22))

Percentage of Poverty

# Poverty
ggplot(allTracts.group)+
  geom_sf(data = st_union(Tracts.2009))+
  geom_sf(data = st_union(Tracts.2018))+
  geom_sf(aes(fill = q5(pctPoverty))) +
  geom_sf(data = buffer_2018, fill = "transparent", color = "red") +
  geom_sf(data = buffer_2009, fill = "transparent", color = "blue")+
  scale_fill_manual(values = palette5,
                    labels = qBr(allTracts.group.visual, "pctPoverty"),
                    name = "Percent_Poverty\n(Quintile Breaks)") +
  labs(title = "Percent of Poverty 2009-2018", 
       subtitle = "Percentage(%);  Blue-area close to metro stations (2009); red-area close to new stations (2018)",
       caption = "TOD Indicators - Percentage of Poverty | Data: US Cencus Bureau, City of Los Angeles Open Data") +
  facet_wrap(~year)+
  mapTheme() + 
  theme(plot.title = element_text(size=22))

Total Population

# Population
ggplot(allTracts.group)+
  geom_sf(data = st_union(Tracts.2009))+
  geom_sf(data = st_union(Tracts.2018))+
  geom_sf(aes(fill = q5(TotalPop))) +
  geom_sf(data = buffer_2018, fill = "transparent", color = "red") +
  geom_sf(data = buffer_2009, fill = "transparent", color = "blue")+
  scale_fill_manual(values = palette5,
                    labels = qBr(allTracts.group, "TotalPop"),
                    name = "Population\n(Quintile Breaks)") +
  labs(title = "Total Population 2009-2018", 
       subtitle = "Counts of People;  Blue-area close to metro stations (2009); red-area close to new stations (2018)", 
       caption = "TOD Indicators - Total Population | Data: US Cencus Bureau, City of Los Angeles Open Data") +
  facet_wrap(~year)+
  mapTheme() + 
  theme(plot.title = element_text(size=22))

Median Rent

ggplot(allTracts.group)+
  geom_sf(data = st_union(Tracts.2009))+
  geom_sf(data = st_union(Tracts.2018))+
  geom_sf(aes(fill = q5(MedRent))) +
  geom_sf(data = buffer_2018, fill = "transparent", color = "red") +
  geom_sf(data = buffer_2009, fill = "transparent", color = "blue")+
  scale_fill_manual(values = palette5,
                    labels = qBr(allTracts.group, "MedRent"),
                    name = "Rent\n(Quintile Breaks)") +
  labs(title = "Median Rent 2009-2018", 
       subtitle = "Real Dollars,  Blue-area close to metro stations (2009); red-area close to new stations (2018)",
       caption = "TOD Indicators - Median Rent | Data: US Cencus Bureau, City of Los Angeles Open Data") +
  facet_wrap(~year)+
  mapTheme() + 
  theme(plot.title = element_text(size=22))

TOD Indicators Summary Table

allTracts.Summary <- 
  st_drop_geometry(allTracts.group) %>%
  group_by(year, TOD) %>%
  summarize(Population = mean(TotalPop, na.rm = T),
            Median.Rent = mean(MedRent, na.rm = T),
            Percent.Poverty = mean(pctPoverty*100, na.rm = T),
            Bachelors.Per10K = mean(pctBachelors*10000, na.rm = T))

kable(allTracts.Summary) %>%
  kable_styling() %>%
  footnote(general_title = "\n",
           general = "TOD Indicators Summary Table: Los Angeles, 2009 versus 2018 | Data: US Cencus Bureau, City of Los Angeles Open Data")

TOD Indicator Differences (Time & Space)

# Indicator differences across time and space
allTracts.Summary %>%
  gather(Variable, Value, -year, -TOD) %>%
  ggplot(aes(year, Value, fill = TOD)) +
  geom_bar(stat = "identity", position = "dodge") +
  facet_wrap(~Variable, scales = "free", ncol=5) +
  scale_fill_manual(values = c("#bae4bc", "#0868ac")) +
  labs(title = "TOD Indicator differences across time and space",
       caption = "TOD Indicator Differences 2009 vs. 2018 | Data: US Cencus Bureau, City of Los Angeles Open Data") +
  plotTheme() + theme(legend.position="bottom")

Bachelors Per 10K People. Compared with 2009 data, bachelors per 10K people in 2018 appears to have a massive growth in the TOD zones and a slower growth in the Non-TOD zones, resulting in a higher number of bachelors per 10K people in the TOD zones than that in the Non-TOD zones, which reverses the 2009 situations that more bachelors were living in the Non-TOD zones.

Median Rent. Median rent of housing increases considerably after 2009 in both Non-TOD and TOD zones.

Percentage of Poverty. Both Non-TOD and TOD zones observes a decrease in poverty percentage from the 2018 data, but the absolute number is still higher in the TOD zones.

Total Population. Generally the total population has a decreasing tendency in LA City, but still more population resides in the Non-TOD zones in 2018.

Graduated Symbol Maps

# Select
TOD09 <- Tracts.2009.inf[(Tracts.2009.inf$TOD =="TOD"),]
TOD18 <- Tracts.2018[(Tracts.2018$TOD =="TOD"),]

buffer09 <- 
  st_transform(station_2009, st_crs(tracts09)) %>%
  st_buffer(804.67) %>%
  dplyr::select(Station)

buffer18 <- 
  st_transform(station_2018, st_crs(tracts09)) %>%
  st_buffer(804.67) %>%
  dplyr::select(Station)

pop_09<- 
  st_join(Tracts.2009.inf, buffer09) %>%
  filter(!is.na(Station))

pop_18<- 
  st_join(Tracts.2018, buffer18) %>%
  filter(!is.na(Station))

popByStation_09 <-
  pop_09 %>% 
  group_by(Station) %>%
  summarize(sumPop = sum(TotalPop)) %>%
  st_drop_geometry()

popByStation_18 <-
  pop_18 %>% 
  group_by(Station) %>%
  summarize(sumPop = sum(TotalPop)) %>%
  st_drop_geometry()

pop_sta_09<-
  left_join(station_2009, popByStation_09, by = c('Station' = 'Station'))%>%
  mutate(Year='2009')

pop_sta_18<-
  left_join(station_2018, popByStation_18, by = c('Station' = 'Station'))%>%
  mutate(Year='2018')

pop_sta.group <- 
  rbind(pop_sta_09,pop_sta_18)

ggplot() + 
  geom_sf(data = allTracts, fill="black")+
  geom_sf(data = st_centroid(pop_sta.group), aes(size = sumPop), shape = 21, 
          fill = "lightblue", alpha = 1, show.legend = "point") + 
  scale_size_continuous(range = c(0.5, 3))+
  facet_wrap(~Year)+mapTheme()+
  labs(title = "Population Within 0.5 Mile of Each Transit 2009-2018", subtitle = "People", caption = "Data: US Cencus Bureau, City of Los Angeles Open Data")

# MedRent
rentByStation_09 <-
  pop_09 %>% 
  group_by(Station) %>%
  summarize(meanRent = mean(MedRent.inf)) %>%
  st_drop_geometry()

rentByStation_18 <-
  pop_18 %>% 
  group_by(Station) %>%
  summarize(meanRent = mean(MedRent)) %>%
  st_drop_geometry()

rent_sta_09<-
  left_join(station_2009, rentByStation_09, by = c('Station' = 'Station'))%>%
  mutate(Year='2009')

rent_sta_18<-
  left_join(station_2018, rentByStation_18, by = c('Station' = 'Station'))%>%
  mutate(Year='2018')

rent_sta.group <- 
  rbind(rent_sta_09,rent_sta_18)

ggplot() + 
  geom_sf(data = allTracts, fill="black")+
  geom_sf(data = st_centroid(rent_sta.group), aes(size = meanRent), shape = 21, 
          fill = "red", alpha = 1, show.legend = "point") + 
  scale_size_continuous(range = c(0.5, 3))+
  facet_wrap(~Year)+mapTheme()+
  labs(title = "Mean Median rent Within 0.5 Mile of Each Transit 2009-2018", subtitle = "Dollars", caption = "Data: US Cencus Bureau, City of Los Angeles Open Data")

Mean Rent & Distance to Metro Stations

# Multirings_2009
multirings_2009 <-
  st_join(st_centroid(dplyr::select(Tracts.2009.inf, GEOID, year)), 
          multipleRingBuffer(st_union(station_2009), 16898.07, 804.67)) %>%
  st_drop_geometry() %>%
  left_join(dplyr::select(Tracts.2009.inf, GEOID, MedRent.inf, year), 
            by=c("GEOID"="GEOID", "year"="year")) %>%
  st_sf() %>%
  mutate(distance = distance / 1609.34)

multi09.Summary <- 
  st_drop_geometry(multirings_2009) %>%
  group_by(year, distance) %>%
  summarize(Rent = mean(MedRent.inf, na.rm = T))

# Multirings_2018
multirings_2018 <-
  st_join(st_centroid(dplyr::select(tracts18, GEOID, year)), 
          multipleRingBuffer(st_union(station_2009), 16898.07, 804.67)) %>%
  st_drop_geometry() %>%
  left_join(dplyr::select(tracts18, GEOID, MedRent, year), 
            by=c("GEOID"="GEOID", "year"="year")) %>%
  st_sf() %>%
  mutate(distance = distance / 1609.34)

multi18.Summary <- 
  st_drop_geometry(multirings_2018) %>%
  group_by(year, distance) %>%
  summarize(Rent = mean(MedRent, na.rm = T))

# plot the geom_line
ggplot()+
  geom_line(data = multi09.Summary,aes(x = distance,y = Rent,colour = "2009"),size=2)+
  geom_point(data = multi09.Summary,aes(x = distance,y = Rent,colour = "2009"),size=4)+
  geom_line(data = multi18.Summary,aes(x = distance,y = Rent,colour = "2018"),size=2)+
  geom_point(data = multi18.Summary,aes(x = distance,y = Rent,colour = "2018"),size=4)+
  scale_colour_manual("",values = c("2009" = "#7bccc4","2018"="#43a2ca"))+
  xlab("Distance")+ylab("US Dollars")+
  plotTheme() + theme(legend.position="bottom")+
  ggtitle("Median Rent as a function of distance to LA Metro stations")+
  labs(caption = "Data: US Cencus Bureau, City of Los Angeles Open Data")

5. Crime Analytics and TOD

Wrangle Crime Open Data

# ---- Read Crime Open Data -----
crime09 <- st_read('crime_2009.shp')%>%
  select(Date_Rptd, Crm_Cd_Des) %>%
  st_transform(st_crs(tracts09)) 

crime18 <- st_read('crime_2018.shp')%>%
  select(Date_Rptd, Crm_Cd_Des) %>%
  st_transform(st_crs(tracts09))

crime18_b <-crime18[(crime18$Crm_Cd_Des =="RAPE, FORCIBLE"),]%>%
  mutate(year='2018')

crime09_b <-crime09[(crime09$Crm_Cd_Des =="RAPE, FORCIBLE"),]%>%
  mutate(year='2009')

crime.group <- rbind(crime18_b,crime09_b)
buffer.group <- rbind(buffer09, buffer18)

Visualize Crime Data with Metro Stations and Median Rent

# ---- Plot Crime Data -----

ggplot()+
  geom_sf(data = allTracts, fill="black")+
  geom_sf(data = crime.group,  color = "red",alpha=0.5)+
  geom_sf(data = buffer.group, fill = "transparent", color = "white")+
  scale_fill_manual(values = palette5,labels = qBr(Tracts.2018, "MedRent"),name = "Rent\n(Quintile Breaks)") +
  labs(title = "Crime and Metro Station Distribution 2009-2018", subtitle = "Forcible Rape") +
  facet_wrap(~year)+
  mapTheme() +
  labs(caption = "Data: US Cencus Bureau, City of Los Angeles Open Data") +
  theme(plot.title = element_text(size=22))

ggplot(allTracts.group)+
  geom_sf(data = st_union(tracts09))+
  geom_sf(aes(fill = q5(MedRent))) +
  geom_sf(data = crime.group,  color = "red",alpha=0.7)+
  scale_fill_manual(values = palette5,labels = qBr(allTracts.group, "MedRent"),name = "Rent\n(Quintile Breaks)") +
  labs(title = "Crime and Median Rent Distribution 2009-2018", subtitle = "Forcible Rape") +
  facet_wrap(~year)+
  labs(caption = "Data: US Cencus Bureau, City of Los Angeles Open Data") +
  mapTheme() +
  theme(plot.title = element_text(size=22))

6. Discussion

In terms of census data, TOD zones seem to attract more bachelors, or it might the case that a convenient location to transit stations better fits people living in single apartments compared with suburban homes demanding driving. Ratios of poverty population are higher in the TOD zones, but in both zones the percentages are decreasing. From 2009 to 2018, the tendencies in population and rent did not show much difference between TOD and Non-TOD zones. The LA metro system is sparse and it extends into many directions in the city, but still far from "moving passengers wherever they would like to go".

TOD in LA has not shown many advantages over the non-TOD parts of the city, but our data is limited in itself due to the expansion of metro lines in the city area. Also, the LA metro lines grows into the more inland (Pasadena, Arcadia, etc.) and the coastal part (eg., Santa Monica, CA) of the LA county, and so we may not rely on the data of LA city as a whole due to their spatial and social relations.

We did not observe a lot of spatial ariations regarding crime incidents, forcible rape especially, between 2009 and 2018, but the number of incidents shows a plummet. Both crime incidents and population aggregate in the regions with easy access to transits, or we can say that population density might be postively related to transit convenience, so did crime incidents. In terms of median rent, low-rent-regions were often associated with clustering of forcible rape crime incidents. Similar patterns might be seen in other metropolitan centers with high populations.

In generally, TOD possibly add on the city's potential in attracting new residents for Los Angeles, CA. It can also have other restraints like social conditions and human resource policies. It might be interesting to look at the patterns of bachelors in the census data. But we're not certain about who would pay the bill for TOD yet.