空间数据可视化与simple future模型应用

2018-07-25 11:48:59 浏览数 (1)

这是一篇关于关于空间地理信息数据可视化与simple feature 模型应用的笔记小结。

之前关于simple feature地理信息数据模型的分享已经有两篇了,这里会继续分享simple feature模型在构建 Point/MutiPoint、LineString/MutiString、Polygon/MutiPolygons。

Polygon/MutiPolygons的应用其实就是数据地图模型的新拓展,在R语言的ggplot2中使用geom_polygon图层函数制作数据地图,需要使用基于sp包的 SpatialPolygonsDataFrame模型,从中提取所需要的行政区划信息和地理特征信息结合业务数据合并,最终实现可视化需求。

simple feature模型构建了新的基于特征的空间地理信息数据存储格式,详细的介绍及其实现可以参考以下资源:

空间数据可视化笔记——simple features空间对象基础

左手用R右手Python系列12——空间数据可视化与数据地图

基于R语言sf包simple feature案例应用:

代码语言:javascript复制
devtools::install_github("tidyverse/ggplot2")
#如果提示无geom_sf函数则需下载github上开发板的ggplot2
library("ggplot2")
library("magrittr")
library("reshape2")
library("ggthemes")
library('dplyr')
library("sf")

1、素材导入:

代码语言:javascript复制
china_map3 <- st_read("D:/R/mapdata/State/china.geojson",stringsAsFactors=FALSE) 
str(china_map3)
class(china_map3)
class(china_map3$geometry)

业务数据:

代码语言:javascript复制
mydata <- data.frame(
  id = 1:34,
  name = china_map3$name,
  scale = runif(34,100,200) %>% round(),
  Scope = rep(LETTERS[1:5],length = 34)
)

将业务数据与地图数据合并

代码语言:javascript复制
final_namdata2 <- left_join(china_map3,mydata[,c(1,3,4)],by = 'id')
final_namdata2 <- st_transform(final_namdata2, 3857) #转换编码

连续型指标

代码语言:javascript复制
ggplot()  
  geom_sf(data=final_namdata2,aes(fill = scale))  
  coord_sf()  
  scale_fill_distiller(palette ='BrBG' )  
  theme_void()

离散指标

代码语言:javascript复制
ggplot()  
  geom_sf(data=final_namdata2,aes(fill=Scope))  
  coord_sf()  
  scale_fill_brewer(palette ='BrBG' )  
  theme_void()

点图(point):

代码语言:javascript复制
province_city <- data.table::fread("D:/R/rstudy/Province/chinaprovincecity.csv") 

cities <- c("北京","上海","天津","重庆","沈阳","呼和浩特","太原","郑州","西安","兰州","合肥","南京","杭州","长沙","武汉")
map_point <- st_as_sf(province_city,coords = c("jd", "wd")) %>% filter(city %in% cities)
st_crs(map_point) = 4326 #设置投影信息str(map_point)
class(map_point)
class(map_point$geometry)

叠加散点图

代码语言:javascript复制
ggplot()  
  geom_sf(data=final_namdata2, aes(fill = Scope),alpha = .3)  
  geom_sf(data=map_point,aes(colour = class,size=zhibiao))  
  coord_sf()  
  scale_fill_brewer(palette ='BrBG' )  
  scale_colour_wsj()  
  theme_void()

线图(line):

代码语言:javascript复制
mutiline_data <- data.frame(city = rep(NA,28),group = rep(NA,28))
mutiline_data$group <- rep(1:14,each = 2)
for (i in 1:28){
  if (i%%2) {
    mutiline_data$city[i] = '郑州'
  } else {
    mutiline_data$city[i] = cities[cities!='郑州'][i/2]
  }
}

转换为sf特有的LINESTRING格式

代码语言:javascript复制
mutiline_data <- left_join(mutiline_data,province_city[,c("city","jd","wd")],by = "city") %>% 
  st_as_sf(coords = c("jd", "wd")) %>%
  group_by(group) %>% 
  summarize() %>% 
  st_cast("LINESTRING")

str(mutiline_data)
class(mutiline_data)
class(mutiline_data$geometry)

st_crs(mutiline_data) <- 4326 #设置地理信息投影
代码语言:javascript复制
ggplot()  
  geom_sf(data = final_namdata2, aes(fill = Scope),alpha = .3)   
  geom_sf(data = mutiline_data,aes(group = group),colour = "black",size = .8)  
  geom_sf(data = map_point,aes(colour = class,size=zhibiao))  
  coord_sf()  
  scale_fill_brewer(palette ='BrBG' )  
  scale_colour_wsj()  
  theme_void()  

基于Python GeoDataFrame包的simple feature应用实现:


Polygon:

代码语言:javascript复制
from shapely.geometry 
import Point,LineString
import geopandas
from matplotlib.colors 
import LinearSegmentedColormap
from matplotlib import pyplot as plt
import pandas as pd
import numpy  as np
import string
代码语言:javascript复制
China_map = geopandas.GeoDataFrame.from_file("D:/R/mapdata/State/china.geojson", encoding = 'gb18030')
China_map = China_map.to_crs({'init': 'epsg:4326'})

业务数据:

代码语言:javascript复制
mydata = pd.DataFrame({
  "id":range(1,35),
  "name":China_map["name"],
  "scale":np.random.randint(100,200,34),
  "Scope":list(string.ascii_uppercase[:5])*6   list(string.ascii_uppercase[:4])
},columns = ["id","name","scale","Scope"])

将业务数据与地图数据合并

代码语言:javascript复制
final_namdata = China_map.set_index('id').join(mydata.loc[:,["id","scale","Scope"]].set_index('id'))
final_namdata.plot(column='scale',cmap=plt.get_cmap("BrBG"),figsize=(20,10), legend = True);
plt.gca().xaxis.set_major_locator(plt.NullLocator())
plt.gca().yaxis.set_major_locator(plt.NullLocator())
plt.legend(labels = ['a', 'b'], loc = 'best')

Points:

代码语言:javascript复制
province_city = pd.read_csv("D:/R/rstudy/Province/chinaprovincecity.csv",encoding = "gbk") 
cities = ["北京","上海","天津","重庆","沈阳","呼和浩特","太原","西安","兰州","合肥","南京","杭州","长沙","武汉"]
gd_Point = geopandas.GeoDataFrame(
    data = {
    "division":province_city["city"],
    "class":province_city["class"],
    "scale":province_city["zhibiao"]
    },
    geometry=[Point(i,j) for i,j in zip(province_city["jd"],province_city["wd"])]
)
gd_Point.crs = {'init': 'epsg:4326'}
代码语言:javascript复制
base1 = final_namdata.plot(column='scale',cmap=plt.get_cmap("BrBG"),figsize=(20, 10),legend = True)
gd_Point.plot(ax = base1,column ="scale", marker='o')
plt.scatter(province_city['jd'],province_city['wd'],s = province_city['zhibiao'],c=province_city['zhibiao'],alpha=0.6)

LineStrings:

代码语言:javascript复制
line_data = pd.DataFrame({
    "start_city":["郑州"]*14,
    "end_city":cities,
    "s_jd":list(province_city.loc[province_city.city == "郑州",'jd'])*14,
    "s_wd":list(province_city.loc[province_city.city == "郑州",'wd'])*14,
    "jd":province_city.loc[province_city.city.isin(cities),'jd'],
    "wd":province_city.loc[province_city.city.isin(cities),'wd']
},columns = ["start_city","end_city","s_jd","s_wd","jd","wd"])
代码语言:javascript复制
line_data = geopandas.GeoDataFrame(
    data = line_data.iloc[:,:2],
    geometry=[LineString([tuple(line_data.iloc[j,2:4]),tuple(line_data.iloc[j,4:6])]) for j in range(14)]
)
line_data.crs = {'init': 'epsg:4326'}
base1 = final_namdata.plot(column='scale',cmap=plt.get_cmap("BrBG"),figsize=(20, 10),legend = True)
line_data.plot(ax = base1,figsize=(20,10))

Point Line Polygon:

代码语言:javascript复制
base1 = final_namdata.plot(column='scale',cmap=plt.get_cmap("BrBG"),figsize=(20, 10),legend = True)
base2 = line_data.plot(ax = base1,color = 'red')
gd_Point.plot(ax = base2,column ="scale", marker='o',markersize = 50)
plt.gca().xaxis.set_major_locator(plt.NullLocator())
plt.gca().yaxis.set_major_locator(plt.NullLocator())
plt.legend(labels = ['a', 'b'], loc = 'best')

simple feature结构是空间数据结构模型的新型标准,它简洁易懂,便于存储,和诸多开源工具都有api结构,具备良好的扩展性和兼容性,实乃空间可视化的福音,本篇文章仅仅就其中基础应用部分做了案例分享,更为高阶的内容留待以后继续深入探索!

参考资料:

R——sf:

https://r-spatial.github.io/sf/

https://r-spatial.github.io/sf/articles/sf1.html

https://r-spatial.github.io/sf/articles/sf2.html

https://r-spatial.github.io/sf/articles/sf3.html

http://ggplot2.tidyverse.org/reference/ggsf.html

Python——geopandas:

http://geopandas.org/index.html

http://toblerity.org/shapely/manual.html#linestrings

0 人点赞