—
title: “2015 US Fatal Vehicle Accidents Analysis”
author: “Matthew Sun”
date: “2/01/2018”
output:
html_document: default
—

#### Introduction

There are over 5 millions of vehicle accidents reported in the US in 2015, according to data from National Highway Traffic Safety Administration (NHTSA). About 7% of them are fatal accidents. As the new technology of artificial intelligence going to mature in the next one or two decades, self-driving vehicles will become part of a normal life in human society. I am just wondering how this new technology can help reduce the car accidents on the road.

#### Data Overview

There is a dataset on the NHTSA website for 2015 fatal accident: ftp.nhtsa.dot.gov. It reports the fatal incidents in major US cities, counties and states. It also includes the longitude and latitude of the accidents so the data can be visualized on the map. The state, county and city are recorded in FIPS codes. After proper data cleaning and re-arrangement, the following data fields are used in the analysis.

* State
* Population
* County
* City
* Longitude
* Latitude
* Month
* Day
* Hour
* Time of Day
* Day Week
* Day of Week
* Collision Manner
* Road Function
* Interstate
* Junction
* Trafficway
* Vehicles Involved
* Drunken
* Road Departure
* Weather

#### Geographical Analysis

Each fatal accident can be represented on the US geographical map as a blue circle dot, by its longitude and latitude number. By slicing the data according to the combination of chosen state, including all states as an option, the fatal accidents can be visually inspected on US map in each state. The data can also be sliced by a timeframe of continuous number of hours ranging from one to 24. All the vehicle accident volumes on the map are aggregated into annual summation. The volume of the accidents can be differentiated by the intensity of blue color in each dot. By hovering the mouse over each dot, it will show a pop up message with the county name, a total number of vehicle accidents within the chosen timeframe and the state name.

“`{r map, message=FALSE, include=FALSE}
knitr::opts_chunk$set(eval = TRUE, echo=TRUE)
library(dplyr)
library(data.table)
library(geojsonio)
library(leaflet)
library(shiny)
library(ggplot2)
setwd(‘~/NYC_data_science/projects/’)
accidents <- fread(file = './project1/vehicle_accidents.csv')
collision_total <- accidents %>% filter(Hour < 24) %>%
select(Hour, Collision_manner, Road_function)

hour_range <- seq(from = 0, to = 1)
state_range <- c('California', 'New York')
accidents_selection <- accidents %>% filter(State %in% state_range & Hour %in% hour_range) %>%
select(State, State_ab, Longitude, Latitude, County, Hour, N_state, N_nation,
Density_state, Density_nation) %>% group_by(State) %>%
mutate(N_hour_state = length(State), Ratio_state = N_hour_state/N_state) %>%
ungroup()

accidents_temp <- accidents_selection %>% select(State, N_hour_state, N_nation) %>%
distinct() %>% mutate(N_hour_nation = sum(N_hour_state),
Ratio_nation = N_hour_nation/N_nation) %>%
select(-N_nation, -N_hour_state)
accidents_selection <- accidents_selection %>% left_join(., accidents_temp, by = ‘State’) %>%
group_by(State) %>% arrange(desc(N_hour_state))

states <- geojson_read('./project1/us_population_density.json', what = "sp")
bins <- c(0, 10, 20, 50, 100, 200, 500, 1000, Inf)
pal <- colorBin("YlOrRd", domain = states$density, bins = bins)

leaflet(data = states) %>% addProviderTiles(“Esri.WorldStreetMap”) %>%
setView(lng = -93.85, lat = 37.45, zoom = 4) %>%
clearMarkers() %>%
clearShapes() %>%
addPolygons(
fillColor = ~pal(density),
weight = 2,
opacity = 1,
color = “white”,
dashArray = “3”,
fillOpacity = 0.5) %>%
addCircleMarkers(data = accidents_selection, ~Longitude, ~Latitude, radius = 2,
stroke = F, fillOpacity = 1, fillColor = ‘blue’,
label = paste0(‘County: ‘, accidents_selection$County,
‘: ‘, accidents_selection$N_state))
“`

Since the population of each state varies, the pure volume comparison is only one side of the story. The US basemap is at the same time, plotted as choropleth map according to the 2012 population density of each state.

“`{r graph2, echo=FALSE, fig.height=5, fig.width=8, dev=’png’}
qty_vs_density <- accidents %>%
select(State, N_state, Density_state, Density_nation) %>%
group_by(State) %>% distinct() %>% arrange(desc(N_state))

g <- ggplot(qty_vs_density, aes(x=factor(State, levels = unique(State)), y=N_state+Density_state))
g + geom_col(aes(fill=State)) + guides(fill=F) +
theme(axis.text.x = element_text(angle = 50, hjust = 1, size = 6.5)) +
ylab(‘No. of Accidents’) +
xlab(‘US States’) +
ggtitle(‘US Fatal Accidents Distribution – Quantity vs Density’) +
theme(plot.title = element_text(hjust = 0.5)) +
geom_col(aes(x=factor(State, levels = unique(State)),y=Density_state),
fill = ‘black’, alpha = 0.9, size = 0.25) +
geom_hline(aes(yintercept =Density_nation),
linetype = ‘solid’, color = ‘red’, size = 0.5) +
annotate(‘text’, x=’Vermont’, y = 480, label = ‘National Density\nLine (pp3m) – Red’,
size = 2.4, color = ‘red’)
“`

The bar chart visualizes the fatal vehicle accidents volume against the accident density, which is normalized by the population of each state. Since the density number is extremely low, I used parts per 3 million for visualization comparison. The density of each state is plotted as black color, against a red line, which is national average density number.

From the bar chart, only South Carolina and Mississippi present both relative high volume and high density. The remaining majority of high volume states either align with or lower than the national average, after normalized with the population.

#### Cause Analysis

These three charts demonstrate the three major reasons that have caused the fatal accidents.

“`{r graph4, echo=FALSE, fig.height=4, fig.width=7, dev=’png’}
collision_rdeparture <- accidents %>% group_by(Collision_manner, Road_departure) %>%
summarise(Num_collisions = n()) %>% group_by(Collision_manner) %>%
mutate(Total = sum(Num_collisions)) %>%
arrange(Total, Num_collisions)

g <- ggplot(collision_rdeparture, aes(x = factor(Collision_manner,
levels = unique(Collision_manner)), y = Num_collisions))
g + geom_bar(aes(fill = Road_departure), stat = ‘identity’, position = ‘dodge’) +
theme(axis.text.x = element_text(angle = 50, hjust = 1, size = 8.5)) +
ylab(‘No. of Accidents’) +
xlab(‘Collision Manner’) +
ggtitle(‘Impact of Road Departure on Fatal Accidents’) +
theme(plot.title = element_text(hjust = 0.5, size = 11)) +
guides(fill = guide_legend(‘Road Departure’)) +
theme(legend.title = element_text(size = 8),
legend.position = ‘bottom’,
legend.key.size = unit(0.2,’cm’),
legend.text = element_text(size = 7))
“`

Road Departure away from the road is the major contributor towards the stationary collision and Head-on collision.

“`{r graph5, echo=FALSE, fig.height=4, fig.width=7, dev=’png’}

collision_drunken <- accidents %>% group_by(Collision_manner, Drunken) %>%
summarise(Num_collisions = n()) %>% group_by(Collision_manner) %>%
mutate(Total = sum(Num_collisions)) %>%
arrange(Total, Num_collisions)

g <- ggplot(collision_drunken, aes(x = factor(Collision_manner,
levels = unique(Collision_manner)), y = Num_collisions))
g + geom_bar(aes(fill = Drunken), stat = ‘identity’, position = ‘dodge’) +
theme(axis.text.x = element_text(angle = 50, hjust = 1, size = 8.5)) +
ylab(‘No. of Accidents’) +
xlab(‘Collision Manner’) +
ggtitle(‘Impact of Drunken on Fatal Accidents’) +
theme(plot.title = element_text(hjust = 0.5, size = 11)) +
guides(fill = guide_legend(‘Drunken – BAC’)) +
theme(legend.title = element_text(size = 8),
legend.position = ‘bottom’,
legend.key.size = unit(0.2,’cm’),
legend.text = element_text(size = 7))
“`

Drunken represented only 25% of stationary collision.

“`{r graph6, echo=FALSE, fig.height=4, fig.width=7, dev=’png’}
collision_weather <- accidents %>% group_by(Weather, Collision_manner) %>%
summarise(Num_collisions = n()) %>% group_by(Weather) %>%
mutate(Total = sum(Num_collisions)) %>%
arrange(Total, Num_collisions)

top_5_weather <- accidents %>% select(Weather, Collision_manner) %>%
group_by(Weather) %>% summarise(Num_collisions = n()) %>%
distinct() %>% top_n(5, wt = Num_collisions)

collision_weather <- collision_weather %>% filter(Weather %in% top_5_weather$Weather)

g <- ggplot(collision_weather, aes(x = factor(Weather,
levels = unique(Weather)), y = Num_collisions))
g + geom_bar(aes(fill = Collision_manner), stat = ‘identity’, position = ‘dodge’) +
theme(axis.text.x = element_text(angle = 50, hjust = 1, size = 8.5)) +
ylab(‘No. of Accidents’) +
xlab(‘Weather’) +
ggtitle(‘Impact of Weather on Fatal Accidents – Top 5’) +
theme(plot.title = element_text(hjust = 0.5, size = 11)) +
guides(fill = guide_legend(‘Collision Manner’)) +
theme(legend.title = element_text(size = 8),
legend.position = ‘bottom’,
legend.key.size = unit(0.2,’cm’),
legend.key.height = unit(0.25,’cm’),
legend.text = element_text(size = 7))
“`

Surprisingly, the highest number of vehicle accidents, about two thirds of them are under clear weather condition, and one of the involving vehicles was in stationary position. Definitely weather is not the major reason why accidents happened.

#### Trend Analysis

The collision manner time distribution by the hour clearly shows the top three collision manner, which is stationary, angle and head-on. I will investigate these three categories in more details.

“`{r graph3, echo=FALSE, fig.height=4, fig.width=7, dev=’png’}
collision_qty <- collision_total %>%
group_by(Hour,Collision_manner) %>% summarise(Num_collisions =n())

g <- ggplot(collision_qty, aes(x=Hour, y=Num_collisions))
g + geom_area(aes(fill = Collision_manner)) +
ylab(‘No. of Accidents’) +
xlab(‘Hour’) +
ggtitle(‘Collision Manner Time Distribution by Hour’) +
theme(plot.title = element_text(hjust = 0.5, size = 13)) +
guides(fill = guide_legend(‘Collision Manner’)) +
theme(legend.title = element_text(size = 9),
legend.position = ‘bottom’,
legend.key.size = unit(0.3,’cm’),
legend.key.height = unit(0.35,’cm’),
legend.text = element_text(size = 8))
“`

The next three area maps plot how each collision manner volume varies with the time during the day.

“`{r graph7, echo=FALSE, fig.height=4, fig.width=7, dev=’png’}
collision_station <- collision_total %>%
filter(Collision_manner==’Stationary collision’) %>%
group_by(Hour, Road_function) %>%
summarise(Num_collisions =n())

g <- ggplot(collision_station, aes(x=Hour, y=Num_collisions))
g + geom_area(aes(fill = Road_function)) +
ylab(‘No. of Accidents’) +
xlab(‘Hour’) +
ggtitle(‘Trend for Collision Manner – Stationary Collision’) +
theme(plot.title = element_text(hjust = 0.5, size = 11)) +
guides(fill = guide_legend(‘Road Function’)) +
theme(legend.title = element_text(size = 8),
legend.position = ‘bottom’,
legend.key.size = unit(0.2,’cm’),
legend.key.height = unit(0.25,’cm’),
legend.text = element_text(size = 7))
“`

For stationary collision, the peak happened around 8:00 to 9:00pm. The second peak happened at midnight.

“`{r graph8, echo=FALSE, fig.height=4, fig.width=7, dev=’png’}
collision_angle <- collision_total %>%
filter(Collision_manner==’Angle’) %>%
group_by(Hour, Road_function) %>%
summarise(Num_collisions =n())

g <- ggplot(collision_angle, aes(x=Hour, y=Num_collisions))
g + geom_area(aes(fill = Road_function)) +
ylab(‘No. of Accidents’) +
xlab(‘Hour’) +
ggtitle(‘Trend for Collision Manner – Angle’) +
theme(plot.title = element_text(hjust = 0.5, size = 11)) +
guides(fill = guide_legend(‘Road Function’)) +
theme(legend.title = element_text(size = 8),
legend.position = ‘bottom’,
legend.key.size = unit(0.2,’cm’),
legend.key.height = unit(0.25,’cm’),
legend.text = element_text(size = 7))
“`

For angle collision manner, the highest peak happened between 3:00 ~ 6:00pm evening rush hour. The second peak between 10:00am to 12:00pm.

“`{r graph9, echo=FALSE, fig.height=4, fig.width=7, dev=’png’}
collision_headon <- collision_total %>%
filter(Collision_manner==’Head-On’) %>%
group_by(Hour, Road_function) %>%
summarise(Num_collisions =n())

g <- ggplot(collision_headon, aes(x=Hour, y=Num_collisions))
g + geom_area(aes(fill = Road_function)) +
ylab(‘No. of Accidents’) +
xlab(‘Hour’) +
ggtitle(‘Trend for Collision Manner – Head-On’) +
theme(plot.title = element_text(hjust = 0.5, size = 11)) +
guides(fill = guide_legend(‘Road Function’)) +
theme(legend.title = element_text(size = 8),
legend.position = ‘bottom’,
legend.key.size = unit(0.2,’cm’),
legend.key.height = unit(0.25,’cm’),
legend.text = element_text(size = 7))
“`

The head-on collision manner happened during morning and evening rush hours, namely between 6:00 to 8:00am, and between 2:00pm to 6:00pm.

#### Conclusion

The surprising findings from the analysis are that two thirds of the total fatal vehicle accidents happened in stationary manner under clear weather condition. During morning and evening rush hours, when the driving conditions are much more complicated, human beings actually did a much better job navigating through the busy traffic. This pattern fits well with self-driving technology, which can handle simple, repetitive tasks much more competently and effectively than human beings. I would expect the number of vehicle fatal accidents dramatically reduced with the advent of self-driving technology within the next one or two decades.