1. Preprocess data

自政府資料開放平臺：https://data.gov.tw/dataset/12197 下載109年度交通事故資料合計共3個檔案（109年A1交通事故資料.csv、109年A2交通事故資料(109年1月-6月).csv、109年A2交通事故資料(109年7月-12月).csv。 預先將檔案從中篩選出「新竹市」的交通事故、合併三個檔案並留下4個變數，分別為交通事故類型A1/A2(type)、新竹市東區、北區、香山區三個區(dist)、經度(longitude)、緯度(latitude)，並另存成hsinchu_109_accident.csv。

2. Install packages

## Loading required package: rstudioapi
## Loading required package: ggplot2
## Loading required package: ggmap
## Google's Terms of Service: https://cloud.google.com/maps-platform/terms/.
## Please cite ggmap if you use it! See citation("ggmap") for details.
## Loading required package: factoextra
## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBa
## Loading required package: fpc
## Loading required package: dbscan
## 
## Attaching package: 'dbscan'
## The following object is masked from 'package:fpc':
## 
##     dbscan

3. Read data

data<- read.csv("hsinchu_109_accident.csv", header = TRUE)
head(data)
##   type   dist longitude latitude
## 1   A1 香山區  120.9372 24.80055
## 2   A1   北區  120.9692 24.80585
## 3   A1   北區  120.9754 24.82175
## 4   A1 香山區  120.9344 24.79041
## 5   A1   北區  120.9555 24.80477
## 6   A1 香山區  120.9332 24.79846
df<-data[,3:4]
head(df)
##   longitude latitude
## 1  120.9372 24.80055
## 2  120.9692 24.80585
## 3  120.9754 24.82175
## 4  120.9344 24.79041
## 5  120.9555 24.80477
## 6  120.9332 24.79846

4. Cluster

利用dbscan做分群，先訂初始值eps=0.015，換算距離大約等於1.5公里(x=40057km地球圓周*0.015度/360度)，MinPts=5 最少以5個為一群。 利用兩種套件做圖：fpc& factoextra

#利用fpc套件畫圖
db=fpc::dbscan(df, eps=0.015, MinPts=5)
plot(db, df, main= "DBSCAN", frame= FALSE)

#利用factoextra套件畫圖
fviz_cluster(db, df, stand=FALSE, geom="point")

print(db)
## dbscan Pts=1866 MinPts=5 eps=0.015
##        0    1
## border 2    1
## seed   0 1863
## total  2 1864

兩者所呈現的圖皆沒有分群效果，原因可能為台灣太小又只選新竹市，經緯度範圍不夠大。

利用kNN嘗試找出最佳的距離。

dbscan::kNNdistplot(df, k=5)
abline(h=0.005, lty=2)

由圖可大略看出0.005是較佳的數字，因此下圖再將eps從0.015改為0.005，MinPts一樣維持5，再產出新的圖。

db=fpc::dbscan(df, eps=0.005, MinPts=5)
fviz_cluster(db, df, stand=FALSE, geom="point")

5. Make a Map (Hsinchu)

#Get Hsinchu City's Map
map<- ggmap(get_googlemap(center = c(lon=120.9667, lat=24.8), 
                        zoom = 13, scale = 2,
                        maptype = "terrain",
                        color = "color",
                        extent = "device"))
## Source : https://maps.googleapis.com/maps/api/staticmap?center=24.8,120.9667&zoom=13&size=640x640&scale=2&maptype=terrain&key=xxx-R0

6. Create a map with all accident locations plotted

#Create color variable for graphing
col1 = "#b3cde0" 
map+
  geom_point(
    aes(x = longitude, y = latitude, color= col1),
    data= data, 
    size= 0.5)+
  theme(legend.position = "bottom")
## Warning: Removed 82 rows containing missing values (geom_point).

#Spatical Heatmap
map+
  stat_density2d(
    aes(x = longitude, y = latitude, fill = ..level.., alpha= 0.05), 
    bins= 10, 
    data= data, 
    geom= "polygon")
## Warning: Removed 82 rows containing non-finite values (stat_density2d).

7. Analysis and Conclusion

此次作業需分析政府公開資料集中的交通事故資料，找出事故熱區。 資料集選擇109年新竹市的A1(有死亡)、A2(僅受傷無死亡)之資料分析該年度新竹市的交通事故熱區，利用dbscan分群，調整eps、Minpts參數，最後試出eps=0.005、Minpts=5為較滿意之參數。 由於事故熱區資料多為經緯度資訊，若以地圖方式呈現會更加清楚，因此額外利用ggplot2做出spatical heatmap、point graph 以利對照分析。依上圖所示，易發生事故熱區主要集中在中間大塊的區域，對照地圖應為新竹市區火車站附近，再者為光復路上老爺飯店附近。