graphics | 基础绘图系统（九）——栅格图、点密度图、等高线（填充）图、三维图...

先插播一则消息：刚建的交流群莫名其妙地被判违规，不仅无法加新人，而且群内的消息可能只能自己看到。先看看过几天能不能解封，不行再另想办法。

本篇是“基础绘图系统”系列的第九篇推文，也是倒数第二篇，介绍的几种图形都与地理上的一些概念有关。

rasterImage()

该函数绘制的是栅格图，栅格是指由颜色构成的矩阵。语法结构如下：

rasterImage(image,
            xleft, ybottom, xright, ytop,
            angle = 0, interpolate = TRUE, ...)

• image：栅格对象；

• xleft、ybottom、xright、ytop：绘制栅格的坐标四至。

下面的代码就是生成栅格对象的一个例子：

col <- hcl(seq(0, 360, 24), 80, 80)
mat <- matrix(col, nrow = 4, ncol = 4)
raster <- as.raster(mat)

raster
##      [,1]      [,2]      [,3]      [,4]     
## [1,] "#FFA2BC" "#BBD143" "#00E2E2" "#EBAEFF"
## [2,] "#FFAC91" "#89DA67" "#00DBFF" "#FFA1FF"
## [3,] "#FEB965" "#37E192" "#6ECFFF" "#FF9DE3"
## [4,] "#E1C641" "#00E4BC" "#B9BFFF" "#FFA2BC"

• hcl函数为颜色生成函数，三个参数分别表示颜色的色相（hue）、彩度（chroma）和明度（luminance）；后续会有专门推文介绍R语言中的各类颜色生成方法。

绘制栅格：

par(plt = c(0,1,0,1))
plot(c(100, 250), c(300, 450), type = "n",
     axes = F, ann = F)
rasterImage(raster, 100, 300, 250, 450)

由xleft、ybottom、xright、ytop四个参数控制的绘图区域有无数个点，而栅格对象包含的元素必然是有限的。rasterImage函数在绘制栅格时默认会进行插值（上图）。而设置interpolate = F的效果如下：

par(plt = c(0,1,0,1))
plot(c(100, 250), c(300, 450),
     type = "n", axes = F, ann = F)
rasterImage(raster, 100, 300, 250, 450,
            interpolate = F)

smoothScatter()

该函数绘制的是点密度图。语法结构如下：

smoothScatter(x, y = NULL, nbin = 128, bandwidth,
              colramp = colorRampPalette(c("white", blues9)),
              nrpoints = 100, ret.selection = FALSE,
              pch = ".", cex = 1, col = "black",
              transformation = function(x) x^.25,
              postPlotHook = box,
              xlab = NULL, ylab = NULL, xlim, ylim,
              xaxs = par("xaxs"), yaxs = par("yaxs"), ...)

• x和y参数分别是点的横坐标和纵坐标（与其他散点图函数一样）；

• bandwidth参数为带宽，nbin为“带”的个数，二者设置一个即可，后者默认值为128。

该函数的默认效果：

set.seed(123)
x <- rpois(1000, 10)
y <- rpois(1000, 8)
smoothScatter(x, y)

带宽bandwidth越大，点密度图看起来越平滑。比较下面两幅图：

smoothScatter(x, y, bandwidth = 0.5)

smoothScatter(x, y, bandwidth = 2)

contour()

该函数绘制等高线。语法结构如下：

contour(x = seq(0, 1, length.out = nrow(z)),
        y = seq(0, 1, length.out = ncol(z)),
        z,
        nlevels = 10, levels = pretty(zlim, nlevels),
        labels = NULL,
        xlim = range(x, finite = TRUE),
        ylim = range(y, finite = TRUE),
        zlim = range(z, finite = TRUE),
        labcex = 0.6, drawlabels = TRUE, method = "flattest",
        vfont, axes = TRUE, frame.plot = axes,
        col = par("fg"), lty = par("lty"), lwd = par("lwd"),
        add = FALSE, ...)

• z参数为矩阵，表示“高程”；

• x、y参数可以理解成矩阵z的各行、各列的坐标，其长度分别与它的行数、列数相同，可省略。

• x、y省略时，默认值在0-1之间均匀取值；此时z参数可放置在第一个参数位置；

• 等高线的“高”不一定是真实的高度数据，也可以是各类属性数据。

下列代码利用泊松随机生成函数生成了矩阵z：

set.seed(123)
z <- rpois(100, 10)
z <- matrix(z, ncol = 10, byrow = T)
contour(z)

volcano数据集记录了某火山周围的海拔高度，示意如下图：

利用该数据集生成等高线图：

contour(volcano)

• nlevels参数可以控制等高线水平的个数，levels参数可以指定需要绘制出的等高线；二者设置一个即可，或使用默认值。

下列的代码来自谢益辉的《现代统计图形》:

library(MSG)
library(KernSmooth)
data(ChinaLifeEdu)
x <- ChinaLifeEdu
head(x)
##   Life.Expectancy High.Edu.NO
## 1           76.10       48001
## 2           74.91       18601
## 3           72.54       40036
## 4           71.65       23197
## 5           69.87       23660
## 6           73.34       44404

est <- bkde2D(x, apply(x, 2, dpik))
contour(est$x1, est$x2, est$fhat,
  nlevels = 15, col = "darkgreen",
  vfont = c("sans serif", "plain"),
  xlab = "预期寿命",
  ylab = "高学历人数"
)
points(x, pch = 20)

• 上图绘制过程中利用了KernSmooth工具包中平滑函数对原始数据进行了处理；

• 散点表示原始数据，每个点代表一个中国大陆省级行政区的数据。

filled.contour

该函数与contour函数类似，不同点是它使用颜色填充来区分不同的等高线水平。语法结构如下：

filled.contour(x = seq(0, 1, length.out = nrow(z)),
               y = seq(0, 1, length.out = ncol(z)),
               z,
               xlim = range(x, finite = TRUE),
               ylim = range(y, finite = TRUE),
               zlim = range(z, finite = TRUE),
               levels = pretty(zlim, nlevels), nlevels = 20,
               color.palette = function(n) hcl.colors(n, "YlOrRd", rev = TRUE),
               col = color.palette(length(levels) - 1),
               plot.title, plot.axes, key.title, key.axes,
               asp = NA, xaxs = "i", yaxs = "i", las = 1,
               axes = TRUE, frame.plot = axes, ...)

下列代码的效果可与contour(volcano)进行比较：

filled.contour(volcano)

persp()

该函数绘制的是三维图，可以理解成将等高线立体起来。语法结构如下：

persp(x = seq(0, 1, length.out = nrow(z)),
      y = seq(0, 1, length.out = ncol(z)),
      z, xlim = range(x), ylim = range(y),
      zlim = range(z, na.rm = TRUE),
      xlab = NULL, ylab = NULL, zlab = NULL,
      main = NULL, sub = NULL,
      theta = 0, phi = 15,
      r = sqrt(3), d = 1,
      scale = TRUE, expand = 1,
      col = "white", border = NULL,
      ltheta = -135, lphi = 0,
      shade = NA, box = TRUE,
      axes = TRUE, nticks = 5,
      ticktype = "simple", ...)

• x、y、z：同contour函数；

• theta、phi：控制三维图的观察方向，使用角度表示；前者为方位角（x-y平面旋转角度），后者为俯视角；

• r：观察点与绘图中心之间的距离；

• d：控制透视强度；大于1减少透视效果，小于或等于1夸大透视效果；

• scale：是否对坐标进行标准化；若为TRUE，则对x、y、z坐标分别进行标准化（此时各方向的比例将改变）；若为FALSE，则进行等比例的缩放；

• expand：取0-1之间的数值以对z进行收缩；

• ltheta、lphi：控制灯光方向，绘图效果如同从该方向对三维图照亮；

• shade：控制阴影生成的参数。

下列代码使用的大部分参数都是默认值：

set.seed(123)
z <- rpois(100, 10)
z <- matrix(z, ncol = 10, byrow = T)
par(plt = c(0,1,0.1,1))
persp(volcano, col = "green3")

下列代码来自官方文档（也是对上图绘图要素的调整）：

z <- 2 * volcano
x <- 10 * (1:nrow(z))
y <- 10 * (1:ncol(z))
par(bg = "slategray", plt = c(0,1,0,1))
persp(x, y, z,
      theta = 135, phi = 30,
      col = "green3", scale = FALSE,
      ltheta = -120, shade = 0.75,
      border = NA, box = FALSE)