利用R包ggmap进行空间可视化

ggmap 是在R环境里调用地图作用可视化的利器。它的语法结构跟ggplot2非常相似，也使R语言的用户可以迅速上手。 ggmap 结合 ggplot 可以方便快速绘制基于地图的可视化图表。下面的文章里，我将用两个例子 （”三藩市的犯罪记录” 和 “Taxi in Porto“） 演示一下它的使用。案例数据均取自kaggle。

初步了解ggmap

get_map: 是最基本的功能，可以用于从（google/openstreetmap）下载地图。

地图下载参数

location: 如果我们知道数据的准确经纬度，我们可以直接使用。或者我们可以直接用geocode. (如 ‘San Francisco’ 用表示三藩市)

maptype: 可以根据可视化的需要，选择 ‘terrain’, ‘satellite’, ‘roadmap’, ‘hybrid’， 通常“terrain”就可以达到不错的效果。

zoom: 可以选择从 3 (大陆范围) 到 21 (建筑物范围). 因为这里的例子都是针对城市，所以我选择10-14 （城市范围）

ggmap: 用作绘制地图，可以和ggplot叠加

旧金山的犯罪记录

这组数据 包含了旧金山从2013年至2015年的各种犯罪记录有及发生地点。让我以2014年的数据做以下的可视化演示。

了解数据

这里我们只需要使用犯罪的类别(Category)，经度(X)，纬度(Y)

head(data_2014[c('Category','DayOfWeek','PdDistrict','X','Y','Year')])

##            Category DayOfWeek PdDistrict         X        Y Year
## 27585 LARCENY/THEFT    Sunday    MISSION -122.4356 37.76080 2014
## 27586 LARCENY/THEFT    Sunday   NORTHERN -122.4336 37.80004 2014
## 27587 VEHICLE THEFT    Sunday    CENTRAL -122.4128 37.78683 2014
## 27588       ASSAULT    Sunday   NORTHERN -122.4217 37.77364 2014
## 27589       ASSAULT    Sunday   NORTHERN -122.4217 37.77364 2014
## 27590       ASSAULT    Sunday   NORTHERN -122.4217 37.77364 2014

发生最为频繁的是Larceny, other offences, non-criminal, assault and vehicle theft, 然而在这些的类别里面other offences and non criminal并不能给我们带来太多信息。因此，我会在地图上显示Larceny, Assault and Vehicle Theft，这三种犯罪的空间密度。

Crime=data.frame(table(data_2014$Category))
Crime=Crime[order(Crime$Freq,decreasing=T),]
head(Crime)

##              Var1  Freq
## 17  LARCENY/THEFT 18901
## 22 OTHER OFFENSES 10290
## 21   NON-CRIMINAL  9669
## 2         ASSAULT  6254
## 35  VEHICLE THEFT  3564
## 34      VANDALISM  3538

地图下载和绘制

Use ggmap to plot the map and use geom_jitter to plot the crimes in San Francisco:

使用get_map语句下载得到旧金山市区的地图。

library(ggmap)
map=get_map(location='San Fransico',maptype='roadmap',zoom=12)

## Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=San+Fransico&zoom=12&size=640x640&scale=2&maptype=roadmap&language=en-EN&sensor=false[/code] 
使用ggmap语句绘制地图，geom_jitter画化不互相叠加的点，以表示每个犯罪的发生。

plot_crime <- ggmap(map)+geom_jitter(aes(X,Y,colour=Category),
data=data_2014[data_2014$Category %in% Top3,],alpha=0.2)+geom_jitter(aes(X,Y,colour=Category),
data=data0[data0$Category %in% Top3,],alpha=1)+labs(x='Longitude',y='Latitude')
plot_crime




一些发现
在2014年里，Larceny/Theft （绿色）的案件要远远多过别的犯罪事件。大多数的犯罪发生在北海岸。

The assault crimes (红色)主要集中发生在 Tenderloin, 而该地区就是旧金山一个知名的high-crime neighbourhood 。

The Vehicle theft crime (蓝色) 的发生率比起前两种要少很多，没有特别严重的地区。

波尔图的出租车
这组数据 包含了波尔图出租车的位置轨迹信息。波尔图为葡萄牙的第二大城市。

了解数据
head(test[c('TRIP_ID','TAXI_ID','TIMESTAMP','DAY_TYPE','MISSING_DATA')])


##   TRIP_ID  TAXI_ID  TIMESTAMP DAY_TYPE MISSING_DATA
## 1      T1 20000542 1408039037        A        False
## 2      T2 20000108 1408038611        A        False
## 3      T3 20000370 1408038568        A        False
## 4      T4 20000492 1408039090        A        False
## 5      T5 20000621 1408039177        A        False
## 6      T6 20000607 1408037146        A        False


在这里，我将只使用POLYLIN列的数据（出租车的运动轨迹）。这里我们定义POLYLIN的第一个位置为上客点，最后一个位置为下客点。

数据预处理
test$POLYLINE[1]


## [1] "[[-8.585676,41.148522],[-8.585712,41.148639],[-8.585685,41.148855],[-8.58573,41.148927],[-8.585982,41.148963],[-8.586396,41.148954],[-8.586072,41.14872],[-8.586324,41.147847],[-8.586999,41.14746],[-8.586576,41.147154],[-8.584884,41.146623]]"


getCoord=function(x){
x=gsub('[/[]','',x)
x=gsub(']','',x)
x=strsplit(x,',')[[1]]
n=length(x)
lon=as.numeric(x[seq(1,n,2)])
lat=as.numeric(x[seq(2,n,2)])
df=data.frame(lon=lon,lat=lat)
df$status='moving'
df$status[1]='pickup'
df$status[nrow(df)]='dropoff'
return(df)
}


通过一些处理，我们可以从PO
a4ab
LYLINE数据中，提取经纬度，以及出租车的运动状态。

loc=NULL
for (i in 1:length(polyline)){
loc=rbind(loc,getCoord(polyline[i]))
}
head(loc)


##         lon      lat status
## 1 -8.585676 41.14852 pickup
## 2 -8.585712 41.14864 moving
## 3 -8.585685 41.14885 moving
## 4 -8.585730 41.14893 moving
## 5 -8.585982 41.14896 moving
## 6 -8.586396 41.14895 moving


运用同样的方法， 我们用ggmap画出波尔图的城市地图以及出租车的运动轨迹。

library(ggmap)
map=get_map(location=c(lon=median(loc$lon),lat=median(loc$lat)),
maptype='roadmap',zoom=13)


## Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=41.161833,-8.612033&zoom=13&size=640x640&scale=2&maptype=roadmap&language=en-EN&sensor=false[/code] 
plot_taxi=ggmap(map)+labs(x='Longitude',y='Latitude')+
geom_jitter(aes(lon,lat,colour=factor(status)),data=loc[loc$status=='moving',],alpha=0.4)+
geom_point(aes(lon,lat,colour=factor(status)),data=centers,size=10,alpha=0.6)+
geom_jitter(aes(lon,lat,colour=factor(status)),data=loc[loc$status!='moving',],alpha=0.8)
plot_taxi




一些发现
我从未去过葡萄牙，但是我也可以猜到Vitoria是商业中心或者是风景名胜，因为出租车的在那里的运动密度特别大。

我们发现在东北和西北有两个主要的上客和下客地点，推测可能是居住的密集区。如果你注意到图中几个大的点，那是用kmeans聚类的上下客的地点。

the completed R code is available: here, if you are interested to check out.