第 4 章 latticeMap

4.1 安装

该包尚未在CRAN发布,需要通过github安装。

remotes::install_github("rpkgs/latticeMap")

4.2 layer_*

latticeMap layer_*系列函数与ggplot中的geom_*类函数使用方法类似。目前已经完成的函数包括layer_title, layer_latFreq, layer_barchart, layer_contourf, layer_statistic

latticeMap绘图主要针对栅格数据,以SpatialGridDataFrameSpatialPixelsDataFrame 这两种数据为主,详情可以翻阅R包sp

以干旱指数(\(Aridity = P/PET\))为例,下文将会演示如何利用latticeMap定制一幅美观的地图。

首先加载包:

library(rgdal)
library(sf)
library(latticeMap)
library(rcolors)
  1. 读取数据
g <- readGDAL("data-raw/AridityIndex_MSWEP-prcp_div_GLEAM-Ep_1980-2020.tif") # Aridity index
## data-raw/AridityIndex_MSWEP-prcp_div_GLEAM-Ep_1980-2020.tif has GDAL driver GTiff 
## and has 180 rows and 360 columns
# continent boundary
file_continent <- "data-raw/shp/continent.shp"
get_continent <- function() {
    read_sf(file_continent) %>% as_Spatial()
}
sp_cont <- list("sp.lines", get_continent(), lwd = 0.5, first = F)
sp_layout = list(sp_cont)
  1. 设置颜色和brks
brks = c(-Inf, 0, 0.05, 0.2, 0.5, 0.65, 1, 2, 5, 10, Inf)
# brks2 = c(0.05, 0.2, 0.5, 0.65, 1, 2, 5)
# sub- humid (0.65 > AI ≥ 0.5), semi- arid (0.5 > AI ≥ 0.2),
# arid (0.2 > AI ≥ 0.05) and hyper- arid (AI < 0.05) regions
nbrk = length(brks) - 1
n_more = 2
cols <- get_color(rcolors$amwg256, nbrk+n_more) %>% .[-(1:n_more)] %>% rev()
p <- sp_plot(g, brks = brks, colors = cols,
             xlim = c(-180, 240),
             ylim = c(-60, 90),
             aspect = 0.5,
             sp.layout = sp_layout,
             key.num2factor = TRUE)
p

## layer_title添加标题

p + layer_title(x = 0, y = 1, labels = c("(a) AI")) 

4.3 layer_latFreq添加随纬度波动图

p + layer_title(x = 0, y = 1, labels = c("(a) AI")) + 
    layer_latFreq(bbox = c(185, 240, -60, 90), zlim = c(-2, 2),
                  unit = "native", is_spatial =TRUE)

4.4 layer_barchart添加频率分布直方图

p + layer_title(x = 0, y = 1, labels = c("(a) AI")) + 
    layer_latFreq(bbox = c(185, 240, -60, 90), zlim = c(-2, 2),
                  unit = "native", is_spatial =TRUE) +
    layer_barchart(x = 0.01, y = 0.05, width = 0.22)

4.5 layer_contourf添加contour

brks2 = c(0.05, 0.2, 0.65, 1)
p + layer_title(x = 0, y = 1, labels = c("(a) AI")) + 
    layer_latFreq(bbox = c(185, 240, -60, 90), zlim = c(-2, 2),
                  unit = "native", is_spatial =TRUE) +
    layer_barchart(x = 0.01, y = 0.05, width = 0.22) +
    layer_contourf(brks = brks2)

4.6 layer_statistic添加均值±方差统计

p + layer_title(x = 0, y = 1, labels = c("(a) AI")) + 
    layer_latFreq(bbox = c(185, 240, -60, 90), zlim = c(-2, 2),
                  unit = "native", is_spatial =TRUE) +
    layer_barchart(x = 0.01, y = 0.05, width = 0.22) +
    layer_contourf(brks = brks2) +
    layer_statistic(x = 0.4, y = 0.1, cex = 1.4)

4.7 layer_grob (deprecated)

今天为大家介绍,如何采用panel.annotation制作下图格式的专业地图,抛砖引玉,以飨读者。

该图使用的是R包latticeMap中的levelplot2绘制。

作为lattice的用户,想必你对如下形式的panel函数非常熟悉:

panel <- function(x, y, z, subscripts, ..., sp.layout) {
  ...
}

lattice中以panel函数控制着每个面板中的绘图内容。通过修改panel函数,我们可以批量、便捷地做出专业、复杂的地图。

latticeMap中设计了panel.annotation函数,便捷地为子图添加插件。panel.annotation主要是将grob对象绘制到viewpoint的指定位置。调用方法如下:

panel.annotation(
  grob,
  bbox = c(xmin, xmax, ymin, ymax),
  unit = "npc",
  xscale = FALSE,
  yscale = FALSE,
  clip = "on",
  ...
)

主要参数:

  • grob: grob object
  • bbox: The region to plot grob, [xmin, xmax, ymin, ymax].
  • unit : A character vector specifying the units for the corresponding numeric values.
  • xscale: A numeric vector of length two indicating the minimum and maximum on the x-scale. The limits may not be identical.
  • yscale: A numeric vector of length two indicating the minimum and maximum on the y-scale. The limits may not be identical.
  • clip : One of “on,” “inherit,” or “off,” indicating whether to clip to the extent of this viewport, inherit the clipping region from the parent viewport, or turn clipping off altogether. For back-compatibility, a logical value of TRUE corresponds to “on” and FALSE corresponds to “inherit.”

现在已完成的基于panel.annotation的函数,也是上图中用到的panel:

  • panel.barchartFreq: 左下角的直方图
  • panel.horizontalFreq: 右侧的随纬度波动的子图。

借助R包ggplotify,R语言任何绘图内容均可转化为grob对象,因此panel.annotation的功能绝不仅限于此,借助panel.annotation可以很容易扩展其他绘图功能饼状图、箱线图、散点图等。用户脑洞多大,就能做出多漂亮的图。

panel.barchartFreq为例,该函数的编写思路如下:

  1. 先将左下角的直方图通过lattice::panel.barchart绘制出来
  2. 将绘图内容转为grob对象
  3. 采用panel.annotation将grob对象绘制在指定位置
  4. 将上述过程封装成panel函数,以便批量制图。

闲话莫提,直接上代码:

经过前人的努力,绘制一幅专业的地图仅需20行代码左右,远比ArcGIS方便,值得你去学习。 latticeMap处于频繁维护期,尚未发布稳定版。想要尝鲜的用户,需要关注微信号,回复latticeMap。

4.8 趋势及显著性水平

此处使用的是GLEAM ET数据。

首先计算线性趋势和pvalue。

# Aridity index
g <- readGDAL("data-raw/AridityIndex_MSWEP-prcp_div_GLEAM-Ep_1980-2020.tif")
## data-raw/AridityIndex_MSWEP-prcp_div_GLEAM-Ep_1980-2020.tif has GDAL driver GTiff 
## and has 180 rows and 360 columns
# continent boundary
file_continent <- "data-raw/shp/continent.shp"
get_continent <- function() {
    read_sf(file_continent) %>% as_Spatial()
}
sp_cont <- list("sp.lines", get_continent(), lwd = 0.5, first = F)
sp_layout = list(sp_cont)
  1. 设置颜色和brks
brks = c(-Inf, 0, 0.05, 0.2, 0.5, 0.65, 1, 2, 5, 10, Inf)
# brks2 = c(0.05, 0.2, 0.5, 0.65, 1, 2, 5)
# sub- humid (0.65 > AI ≥ 0.5), semi- arid (0.5 > AI ≥ 0.2),
# arid (0.2 > AI ≥ 0.05) and hyper- arid (AI < 0.05) regions
nbrk = length(brks) - 1
n_more = 2
cols <- get_color(rcolors$amwg256, nbrk+n_more) %>% .[-(1:n_more)] %>% rev()
p <- sp_plot(g, brks = brks, colors = cols,
             xlim = c(-180, 240),
             ylim = c(-60, 90),
             aspect = 0.5,
             sp.layout = sp_layout,
             key.num2factor = TRUE)
p