前言

进行地理空间分析及可视化目前有多种软件，如ArcGIS，QGIS,Geoda等，与传统的图形用户界面（GUI）类型的软件相比，R更强调命令行界面（CUI），在结果的可重复性、外部资源的调用性方面具有较大的优势，如通过leaflet^[1]包中的几行代码就可以绘制出一个全球夜间灯光的交互式地图。当前R语言中用于空间统计分析的package呈快速增加趋势，越来越多的人参与到R中地理空间分析的维护与开发，R语言空间分析在未来快速发展大数据、数据科学领域仍有较大的发展空间。

#全球夜间灯光可视化图
install.packages("leaflet")
library(leaflet)
leaflet() %>%
  addProviderTiles("NASAGIBS.ViirsEarthAtNight2012")

tmap简介

tmap^[2]主要用于空间数据的可视化，tmap中提供了灵活且易于操作的方法来进行地图的创建，与ggplot2有一定相似性，可通过不同图层的叠加形式，不断修饰和完善图形。本文主要以tmap包为例绘制一个基础分级色彩地图，以及气泡图、分面图。

利用tmap进行空间可视化

准备工作。首先准备一个地区的shp底图，其实相关数据。本文演示数据为河南省shp底图，其实为河南省2014-2019年14个地级市GDP（单位：亿元），数据来源于河南省统计局（未纳入直管县）。数据如下图所示。

#     地级市   y2014   y2015   y2016   y2017    y2018    y2019
#1    郑州市 6776.99 7311.52 8113.97 9193.77 10670.14 11589.72
#2    开封市 1492.06 1605.84 1755.10 1887.55  2157.70  2364.14
#3    洛阳市 3284.57 3469.03 3820.11 4290.19  4613.49  5034.85
#4  平顶山市 1637.17 1686.01 1825.14 1994.66  2170.86  2372.64
#5    安阳市 1791.81 1872.35 2029.85 2249.85  2141.87  2229.29
#6    鹤壁市  682.20  715.65  771.79  827.65   921.18   988.69
#7    新乡市 1917.81 1975.03 2166.97 2357.76  2671.60  2918.18
#8    焦作市 1844.31 1926.08 2095.08 2280.10  2501.76  2761.11
#9    濮阳市 1253.61 1328.34 1449.56 1585.47  1442.65  1581.49
#10   许昌市 2087.23 2171.16 2377.71 2632.92  3140.93  3395.68
#11   漯河市  941.16  992.59 1081.93 1165.04  1435.90  1578.44
#12 三门峡市 1240.06 1251.04 1325.86 1447.42  1319.01  1443.82
#13   南阳市 2675.57 2866.82 3114.97 3345.30  3500.56  3814.98
#14   商丘市 1697.64 1812.16 1989.15 2195.55  2659.52  2911.20
#15   信阳市 1757.34 1879.67 2037.80 2194.51  2534.47  2758.47
#16   周口市 1989.75 2089.70 2263.86 2459.70  2939.59  3198.49
#17 驻马店市 1691.30 1807.69 1972.99 2175.04  2485.26  2742.06
#18   济源市  480.46  492.54  538.91  600.12   630.46   686.96

主要思路。分为：导入数据；合并数据；基础绘图；修改图形；参数理解；交互式地图；气泡图；其他部分

1.安装本文所需要的包

#在运行之前如果没有相关的包，需要首先进行安装
install.packages("tmap")
install.packages("rgdal")
install.packages("vioplot")
install.packages("RColorBrewer")

2.设置工作路径、加载包、读取数据

setwd("F:/地理空间分析")#修改为你的工作路径
library(tmap) #用来绘图
library(rgdal) #用来读取shp数据
library(vioplot)#用来绘制小提琴图
library(RColorBrewer)#用来修改颜色

3.读取数据

HN.Areas <- readOGR('.', "河南省_地级市", encoding = "UTF-8")#读取shp文件，"河南省_地级市"为shp文件名称
HN.GDP <- read.csv("HNGDP.csv", header = TRUE)#读取csv文件

4.查看文件并画出shp底图

names(HN.Areas)#查看列名
#[1] "Name"  "Layer"
head(HN.Areas)
#class       : SpatialPolygonsDataFrame 
#features    : 6 
#extent      : 12353941, 12802933, 3776148, 4253525  (xmin, xmax, ymin, #ymax)
#crs         : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +uni#ts=m +no_defs 
#variables   : 2
#names       :   Name, Layer 
#min values  : 济源市,    市 
#max values  : 郑州市,    市 
plot(HN.Areas)#画出shp底图

names(HN.GDP)#查看csv列名
#[1] "地级市" "y2014"  "y2015"  "y2016"  "y2017"  "y2018"  "y2019" 
head(HN.GDP)
#    地级市   y2014   y2015   y2016   y2017    y2018    y2019
#1   郑州市 6776.99 7311.52 8113.97 9193.77 10670.14 11589.72
#2   开封市 1492.06 1605.84 1755.10 1887.55  2157.70  2364.14
#3   洛阳市 3284.57 3469.03 3820.11 4290.19  4613.49  5034.85

5.对数据进行初步探索

summary(HN.GDP)#利用summary函数对数据进行初步探索
#    地级市              y2014            y2015            y2016            y2017            y2018             y2019      
# Length:18          Min.   : 480.5   Min.   : 492.5   Min.   : 538.9   Min.   : 600.1   Min.   :  630.5   Min.   :  687  
# Class :character   1st Qu.:1313.2   1st Qu.:1397.7   1st Qu.:1525.9   1st Qu.:1661.0   1st Qu.: 1617.5   1st Qu.: 1743  
# Mode  :character   Median :1727.5   Median :1842.3   Median :2009.5   Median :2195.0   Median : 2493.5   Median : 2750  
                    #Mean   :1957.8   Mean   :2069.6   Mean   :2262.8   Mean   :2493.5   Mean   : 2774.3   Mean   : 3021  
                    #3rd Qu.:1971.8   3rd Qu.:2061.0   3rd Qu.:2239.6   3rd Qu.:2434.2   3rd Qu.: 2872.6   3rd Qu.: 3128  
                    #Max.   :6777.0   Max.   :7311.5   Max.   :8114.0   Max.   :9193.8   Max.   :10670.1   Max.   :11590  

也可以利用基础图形对变量进行可视化，观察变量的分布情况，如绘制变量的箱线图与小提琴图，可进一步参考R语言绘图|多变量箱线图（boxplot）与小提琴图（vioplot）的绘制，关于RColorBrewer的使用可参考R语言绘图|如何调用RColorBrewer包对图形颜色进行修改

#绘制变量箱线图
boxplot(HN.GDP$y2014, HN.GDP$y2015, HN.GDP$y2016, HN.GDP$y2017, HN.GDP$y2018, HN.GDP$y2019,
        names = c("2014年", "2015年", "2016年", "2017年", "2018年", "2019年"),
        xlab = "年份", ylab = "GDP(亿元)",
        col = brewer.pal(6,"Blues"))
#绘制变量条形图
vioplot(HN.GDP$y2014, HN.GDP$y2015, HN.GDP$y2016, HN.GDP$y2017, HN.GDP$y2018, HN.GDP$y2019,
        names = c("2014年", "2015年", "2016年", "2017年", "2018年", "2019年"),
        xlab = "年份", ylab = "GDP(亿元)",
        col = brewer.pal(6,"YlGnBu"))

5.利用tmap进行分级色彩地图绘制

5.1首先，需要将csv文件与shp文件进行链接

map.data <- merge(HN.Areas, HN.GDP, by.x = "Name" ,by.y = "地级市")#合并数据，其中x,y分别为两个合并文件所依据的字段
names(map.data)#查看合并后的数据
#[1] "Name"  "Layer" "y2014" "y2015" "y2016" "y2017" "y2018" "y2019"
head(map.data)#查看合并后的数据
#class       : SpatialPolygonsDataFrame 
#features    : 6 
#extent      : 12353941, 12802933, 3776148, 4253525  (xmin, xmax, ymin, ymax)
#crs         : +proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
#variables   : 8
#names       :   Name, Layer,   y2014,   y2015,   y2016,   y2017,    y2018,    y2019 
#min values  : 济源市,    市,  480.46,  492.54,  538.91,  600.12,   630.46,   686.96 
#max values  : 郑州市,    市, 6776.99, 7311.52, 8113.97, 9193.77, 10670.14, 11589.72 

5.2其次，使用qtm进行简单绘图，如下图所示。

qtm(map.data, fill = "y2014")#其中qtm第一个参数为合并为的数据，第二个参数为填充数据

可以发现该图缺少必要的元素，如指北针、比例尺、标题等。因此下面对图形进行修饰，对图例、颜色、分段、比例尺、指北针等进行调整。

tm_shape(map.data) + tm_fill("y2014",#填充变量
                         palette = "YlOrRd", #设置颜色
                         alpha = 0.8,#设置alpha数值
                         title = "GDP of 2014",#图例标题
                         style = "jenks", n = 6, #分类方法
                         legend.hist = FALSE) + #直方图是否显示
  tm_borders(col = "black", lwd = 1, alpha = 1) + #边界颜色、宽度、以及alpha
  tm_compass(type = "arrow", #指北针类型
             text.size = 1.2,#指北针大小
             position=c("right", "top")) + #指北针位置
  tm_scale_bar(position=c("0.3", "0.001"),text.size = 1) + #比例尺位置、大小
  tm_text("Name", size = 1.2) + #显示图层名称及字体大小
  tm_layout(title = "河南省2014年GDP", #图标题
            title.position=c("center", "top"),#图标题位置
            title.fontfamily = "serif",#字体
            legend.text.size = 1, #图例文本大小
            legend.text.fontfamily = "serif",#图例字体
            legend.title.size = 1.5, #图例标题大小
            legend.title.fontfamily = "serif",#图例标题字体
            legend.position = c("left", "bottom"), #图例的位置
            frame = FALSE)
#相关参数解读：
#tm_shape函数：图层文件
#tm_fill函数：对图层进行填充
#tm_borders函数：图层边界函数
#tm_compass函数：指北针图层函数
#tm_scale_bar函数：比例尺图层函数
#tm_text函数：图层名称
#tm_layout函数：布局函数

5.3 此外，也可以对palette和style进行调整，更改分类的颜色及分类的标准，如使用cat分类方法

#palette = "YlOrRd", 
#style = "jenks", n = 6
style = "cat"
palette = "YlGnBu"

5.4分面地图绘制。分面地图绘制与上文部分地图绘制形式基本相同，需要添加一个用于分面的图层，即使用tm_facets函数，其中制定根据Name变量进行分面，nrow表示分为6行。

tm_shape(map.data) + tm_fill("y2014",
                         palette = "YlOrRd", 
                         alpha = 0.8,
                         title = "GDP of 2019",
                         style = "jenks", n = 6,
                         legend.hist = FALSE) + 
  tm_borders(col = "black", lwd = 1, alpha = 1) + 
  tm_facets(by = "Name", nrow = 6,scale.factor = 10) + 
  tm_layout(legend.text.size = 1,        
            legend.title.size = 1, 
            legend.position = c("0.1", "0.4"), 
            frame = FALSE)
#tm_facets()中，by分面所依据的变量，brow分面的行数，scale.factor对分面的标题大小进行设置

此外，也可以参照上文部分对分类方式进行修改

#palette = "YlOrRd", 
#style = "jenks", n = 6
style = "cat"
palette = "YlGnBu"

5.5 交互式地图

library(leaflet)
tmap_mode("view")#Interactive maps。运行该行代码后，再进行绘图则为交互式地图
tmap_mode("plot")# tmap mode set to plotting。运行该行代码则返回之前的地图形式
#交互式地图示例如下，可以自由放大与缩小，同时选择不同的图层

气泡图.函数为tm_bubbles(),函数内可以对size，colour，alpha等进行设置。

tm_shape(map.data) + 
  tm_fill(col = "#FFFFCC") + 
  tm_bubbles(size = "y2014", scale = 3, alpha = 0.7,
             border.col = "black", border.alpha = 1) + 
  tm_borders(col = "black", lwd = 1, alpha = 1) + 
  tm_compass(type = "arrow", 
             text.size = 1.4,
             position=c("right", "top")) + 
  tm_scale_bar(position=c("0.6", "0.001"),text.size = 1) + 
  tm_text("Name", size = 1.2) + 
  tm_layout(title = "河南省2014年GDP", 
            title.position=c("center", "top"),
            title.fontfamily = "serif",
            legend.text.size = 0.8, 
            legend.text.fontfamily = "serif",
            legend.title.size = 1.5, 
            legend.title.fontfamily = "serif",
            legend.position = c("left", "bottom"),
            frame = FALSE)
##tm_bubbles为气泡图图层，其中size表示气泡大小的衡量变量，scale表示气泡实际大小，border调整气泡的边界。

本文主要参考tmap官方文档。此外，tmap内置相关参数设置多达上百条，能够对所绘制的地图进行修饰与完善。后续会继续更新tmap相关绘图参数调整以及绘图示例展示。如需要本文的shp地图数据及统计数据请后台回复【20211227】获得。

参考资料