8

我有一个显示澳大利亚偏远地区的 shapefile，从澳大利亚统计局获得：

http://www.abs.gov.au/AUSSTATS/abs@.nsf/DetailsPage/1270.0.55.005July%202011?OpenDocument

同一个 URL 是 PDF“ASGS Remoteness Structure Edition 2011 PDF Maps”——我正在尝试从这个 PDF 文档中复制第一张地图。

我已阅读 shapefile 并将颜色信息添加到data插槽：

ra <- readShapeSpatial("RA_2011_AUST", delete_null_obj = TRUE)
ra@data$COLOUR <- "#FFFFFF"
ra@data$COLOUR[(as.numeric(as.character(ra@data$RA_CODE11)) %% 10) == 0] <- "#006837"
ra@data$COLOUR[(as.numeric(as.character(ra@data$RA_CODE11)) %% 10) == 1] <- "#31A354"
ra@data$COLOUR[(as.numeric(as.character(ra@data$RA_CODE11)) %% 10) == 2] <- "#78C679"
ra@data$COLOUR[(as.numeric(as.character(ra@data$RA_CODE11)) %% 10) == 3] <- "#C2E699"
ra@data$COLOUR[(as.numeric(as.character(ra@data$RA_CODE11)) %% 10) == 4] <- "#FFFFCC"

我唯一要做的就是绘制地图！这就是我卡住的地方...

ra@polygons是一个包含 35 个多边形的列表，每个多边形都有一个插槽ID，它是数据帧的索引ra@data。所以我所要做的就是告诉plot()在ra@data$COLOUR[ID]. 嗯，不完全是。35 个多边形（“多边形”类）中的每一个都有自己的多边形列表（“多边形”类）；总共有 6902 个多边形！！！

我的理解plot()是，我必须以与绘制多边形相同的顺序向它传递一个颜色向量。因此，我相信我必须创建一个长度为 6902 的向量，每个元素都保存相关多边形的颜色值。到目前为止我做得怎么样？

如果按顺序绘制多边形，那将很容易，但事实并非如此。35 个多边形中的每一个都有plotOrder一个整数向量槽，因此颜色向量可能必须按每个向量中的值排序。

在这一点上，这一切似乎有点太复杂了。我在这里完全偏离轨道了吗？

谢谢你的建议！

4

1 回答 1

13

你已经完成了所有的工作！

plot(ra, col=ra@data$COLOUR)

甚至，正如@Spacedman 建议的那样：

plot(ra, col=ra$COLOUR)

就是这样！

在此处输入图像描述

如果你想摆脱多边形边界：

plot(ra, col=ra$COLOUR, border=NA)

在此处输入图像描述

编辑：试图解释这种行为：

根据?SpatialPointsDataFrame：

具有默认 ID 匹配的 SpatialPolygonsDataFrame 根据多边形 ID 槽检查数据框行名称。然后它们必须彼此一致，并且是唯一的（没有 Polygons 对象可以共享 ID）；如果需要匹配多边形 ID，数据框行将重新排序。

这意味着多边形是根据它们的 ID 排序的，因此按照 slot 中数据帧的行的顺序排列@data。

现在，如果您查看函数plot.SpatialPolygons（使用getAnywhere(plot.SpatialPolygons)），在某些时候会出现这些行：

...
    polys <- slot(x, "polygons")
    pO <- slot(x, "plotOrder")
    if (!is.null(density)) {
        if (missing(col)) 
            col <- par("fg")
        if (length(col) != n) 
            col <- rep(col, n, n)
        if (length(density) != n) 
            density <- rep(density, n, n)
        if (length(angle) != n) 
            angle <- rep(angle, n, n)
        for (j in pO) .polygonRingHoles(polys[[j]], border = border[j], 
            xpd = xpd, density = density[j], angle = angle[j], 
            col = col[j], pbg = pbg, lty = lty, ...)
    }
...

输入到的向量与槽col具有相同的顺序polygons，因此与 ID 相同。plotOrder用于以相同的方式索引所有这些。

于 2013-02-07T08:05:15.247 回答