查看原文
其他

R语言可视化——关于ggplot所支持的数据地图素材类型

2017-05-31 杜雨 数据小魔方

做了这么多数据地图,是时候该总结一些心得和理念了,今天这篇讨论ggplot2所支持的数据地图素材格式。


library("plyr")

library("dplyr")

library("rgdal")

library("sf")

library(maptools)

library("ggplot2")

library("ggthemes")

library("geojsonio")


options(stringsAsFactors=FALSE,warn=FALSE,encoding="UTF-8")


今天来跟大家分享一下关于ggplot2系统制作数据地图的源数据支持问题,一直觉得这个问题很重要。


其实分享过这么多期的数据地图,我所用的数据地图素材无外乎以下三种:


  • 传统的shp素材;

  • json素材;

  • 地图包内置地图素材。


其中因为第三种素材的格式跟我们导入的shp空间数据集格式一致,这里重点讲解前两种数据源。


虽然从数据存储格式上来讲我们分为shp素材、json素材,但是由于在R语言中使用ggplot2作图,所支持的数据集对象大致又可分为两类,它们都可以由shp、json数据文件转化而来。


  • sp:SpatialPolygonDataFrame

  • sf:Simple feature list column


所以说数据文件格式和空间数据集对象格式的关系可以这么表述:



这两种格式的数据集所描述的信息差不多是一致的。第一种格式(sp)是R语言绘图比较传统的数据格式,它将地理信息数据分割为两大块:描述层和映射层。


在数据存放时,描述层记录各个地理区域的名称、ID、编号、简写、iOS编码,以及其他标识信息和度量变量,描述层是一个dataframe,我们可以用data@data来提取描述层的数据框。


而对应的几何映射层,是每一个行政区域的多边形边界点,这些边界点按照order排序,按照group分组。多边形边界点信息是一个多层嵌套的list结构,但是我们仍然可以通过fortity函数将其转化为数据框。


即sp空间数据对象是一个dataframe(描述层)和polygons(几何映射层)两个对象的组合对象。


而sf对象将这种控件数据格式件进行了更加整齐的布局,使用st_read()导入的空间数据对象完全是一个整齐的数据框,拥有整齐的行列,这些行列中包含着数据描述和几何多边形的边界点信息。其中最大的特点是,它将每一个行政区划所对应的几何边界点封装成了一个list对象的记录,这条记录就像其他普通的文本记录、数值记录一样,被排列在对应行政区划描述的单元格中。


这样做的好处是,我们不必要自己做这种从描述层到几何映射层的对应关系的链接,因为对应关系本身就已经存在。然后如果是第一种sp格式的话,在制作ggplot2地图过程中,我们需要分离描述层和几何映射层,并为两者指定连接的id(主键),如果算上你要将自己的业务数据和描述层数据合并这一动作的话,那么总共我们需要合并两次数据。(倘若描述层均没有对应的id,你需要为其构造虚拟id,这一次合并算上的话,那么就需要三次合并)。


然而在sf对象中我们仅需指定一次合并即可,即描述层和业务指标数据的合并。


接下来通过案例演示来解释以上原理:


通常我们制作一个数据地图的方式如下:


shapefile文件导入:

setwd("D:/R/rstudy/CHN_adm")

china_map<-readShapePoly("CHN_adm1.shp")

Warning message:

use rgdal::readOGR or sf::st_read 


当我使用sp包导入shp数据集时,R提示以上warning,也就意味着这种传统的方式在不久的将来就会被遗弃,而且建议使用rgdal::readOGR和sf::st_read 来导入。(这也是我今天讲这篇内容的意义所在,真不好说sp包哪天就停止服务了,之前那那些写过的代码可能全部都要挂掉!)


rgdal包可以替代sp包导入shp数据(事实上它也支持json数据的导入,似不似很强大),而sf包则是新崛起,支持sf对象格式导入的包,而且ggplot2率先给这个包开了绿灯,直接创建了geom_sf图层(之后会讲到)


china_map<-readOGR("CHN_adm1.shp",stringsAsFactors=FALSE)

OGR data source with driver: ESRI Shapefile 

Source: "CHN_adm1.shp", layer: "CHN_adm1"

with 32 features

It has 9 fields


同样我用rgdal包导入该素材,不仅没有警示,而且还给出详细的数据格式描述信息。


此时按照旧方法,我需要分别提取出描述层的dataframe和几何映射层的数据框。


mydata<-china_map@data

mymapdata<-fortify(china_map)


这样通过提取和转化之后,描述层是一个行政区域描述信息的数据框,几何映射层转化为一个包含经纬度指标,group分组标量,order排序变量以及ID的数据框。


如果我们需要制作填色地图,那么我们需要先将自己的业务指标和描述层数据进行整理和并,并将合并后的数据与几何映射层的数据框进行合并。


这里我虚拟一个指标。

mydata$zhibiao<-runif(32,10,20)

将其与几何映射层进行合并:

mynewmapdata<-merge(mymapdata[,c(-4,-5)],mydata[,c("ID_1","NAME_1","zhibiao")],by.x="id",by.y="ID_1")


此时即可进行多边形的填色映射。

ggplot()+

geom_polygon(data=mynewmapdata,aes(x=long,y=lat,group=group,fill=zhibiao),col="grey95")+

scale_fill_gradient(low="white",high="steelblue") +

coord_map("polyconic") +

theme_map()



如果是要添加点映射的话,只需再增加一个geom_point()的图层。


这样的步骤看起来确实很繁杂,针对这种sp格式的地图数据,如果你不想做复杂的合并整理构成,ggplot2几何图层对象中有一个简化版的函数——geom_map()。


这个函数虽然使用了两个分离的数据框:描述层,几何映射层,但是它可以允许你通过指定两者之间的关键词(主键)来完成这种合并对接过程。而不需要做复杂的拼接合并。(几何映射层的关键词必须命名为region,而描述层的关键字命名无要求,但是需显声明指定)


ggplot(mydata,aes(map_id=ID_1))+

geom_map(aes(fill=zhibiao),map=mymapdata %>% rename(region=id))+

scale_fill_gradient(low="white",high="steelblue") +

expand_limits(x=mynewmapdata$long,y=mynewmapdata$lat)+

coord_map("polyconic")+

theme_map()



以上代码免去了繁杂的合并转化过程,节省了大量代码,是数据格式制作地图的极佳替代方案。


接下来谈论下sf格式,这种格式的数据既可以来源于json格式数据,也可以来源于shp格式数据,非常自由。


china_map<-st_read("CHN_adm1.shp",stringsAsFactors=FALSE,quiet=TRUE)


使用这种格式数据进行地图制作,仅需合并一下业务数据(这里我就直接生成了):


china_map$zhibiao<-runif(32,10,20)

china_map <- sf::st_transform(china_map, "+init=epsg:4326")

ggplot(china_map)+

geom_sf(aes(fill=zhibiao),col="grey95")+

scale_fill_gradient(low="white",high="steelblue") +

coord_sf()+

theme_map()



关于json素材的导入,我在之前讲leaflet数据地图素材时已经讲过这几种情况,这里不再赘述。


R语言中可以制作数据地图的包有很多(我长用到的):


#需要自己准备数据地图素材:

map()

ggplot(china_map1)+geom_line()  #只能做轮廓图

ggplot(china_map1)+geom_ploygon()

ggplot(china_map1)+geom_sf()

ggplot(china_map1)+gg_map()


可以调用在线地图库素材:

ggmap()

leaflet()

REmap()


关于地理信息可视化涉及到的范围实在是太广了,我也只能洞悉到这个地步,以后的路还有很长,会陆续跟大家分享。


欢迎关注魔方学院QQ群





您可能也对以下帖子感兴趣

文章有问题?点此查看未经处理的缓存