8

我打算在瑞士做一个调查。将询问 NPA。

NPA（邮政编码）包含 4 个数字。

  • 例如 1227 是 Carouge 的 NPA（瑞士日内瓦州的一部分）。
  • 例如 1784 年是 Courtepin 的 NPA（瑞士弗里堡州的一部分）。
  • 等等

我想知道如何在地图上表示所有观察（大约 1500 个）。我正在考虑使用 ggplot，因为我将它用于其他图表（我认为 ggplot 是“美丽的”）。但是，我对任何其他建议持开放态度。

以下是一些假数据： http: //pastebin.com/HsuQnLP3

瑞士地图的输出应该有点像美国地图（来源：http ://www.openintro.org ）

在此处输入图像描述

更新：

我试图创建一些代码：

library(sp)
test <-  url("https://dl.dropboxusercontent.com/u/6421260/CHE_adm3.RData")
print(load(test))
close(test)

gadm$NAME_3
gadm$TYPE_3

但似乎http://gadm.org/没有提供公社的 NPA ......

新更新：

我找到了（感谢@yrochat）一个带有 NPA 的 shapefile：http: //www.cadastre.ch/internet/cadastre/fr/home/products/plz/data.html

我的 ZIP 文件名为：Shape LV03

然后我试过了

library("maptools")
swissmap <- readShapeLines("C:/Users/yourName/YourPath/PLZO_SHP_LV03/PLZO_PLZ.shp")
plot(swissmap)
data <- data.frame(swissmap)
data$PLZ #the row who gives the NPA

由于我在 shapefile 上有 PLZ，我如何在地图上为我的观察着色？我提供了一些关于数据http://pastebin.com/HsuQnLP3的假数据

在此处输入图像描述

谢谢

4

1 回答 1

9

好的，有了 shapefile，我们可以很容易地绘制东西。

work.dir <- "directory_name_no_trailing slash"

# open the shapefile
require(rgdal)
require(rgeos)
require(ggplot2)
ch <- readOGR(work.dir, layer = "PLZO_PLZ")

# convert to data frame for plotting with ggplot - takes a while
ch.df <- fortify(ch)

# generate fake data and add to data frame
ch.df$count <- round(runif(nrow(ch.df), 0, 100), 0)

# plot with ggplot
ggplot(ch.df, aes(x = long, y = lat, group = group, fill = count)) +
    geom_polygon(colour = "black", size = 0.3, aes(group = group)) +
    theme()

# or you could use base R plot
ch@data$count <- round(runif(nrow(ch@data), 0, 100), 0)
plot(ch, col = ch@data$count)

就我个人而言，我发现ggplot使用起来比默认输出要容易得多plot，而且默认输出看起来要好得多。

截屏

并ggplot使用简单的数据框，这使得子集化变得容易。

# plot just a subset of NPAs using ggplot
my.sub <- ch.df[ch.df$id %in% c(4,6), ]
ggplot(my.sub, aes(x = long, y = lat, group = group, fill = count)) +
    geom_polygon(colour = "black", size = 0.3, aes(group = group)) +
    theme()

结果：

截图2

于 2013-08-06T23:37:03.417 回答