4

我正在使用 raster 包中的 getData 函数来检索阿根廷的地图。我想使用 ggplot2 绘制生成的地图,所以我使用 broom 包中的 tidy 函数转换为数据框。这很好用,但我不知道如何保留联邦区的名称,以便我可以在地图上使用它们。

这是我不保留地区名称的原始代码:

# Original code: ##################################
# get the map data from GADM.org and then simplify it
arg_map_1 <- raster::getData(country = "ARG", level = 1, path = "./data/")     %>% 
  # simplify
  rmapshaper::ms_simplify(keep = 0.01) %>% 
  # tidy to a dataframe
  broom::tidy()

# plot the map
library(ggplot2)
ggplot(data=arg_map_1) +
  geom_map(map=arg_map_1, aes(x=long, y=lat, map_id=id, fill=id),
       color="#000000", size=0.25)

下面是从 SPDF 中提取地区名称并将其用作地图 ID 的 hack 代码:

# Code with a hack to keep the district names: ################################
# get the map data from GADM.org and then simplify it
arg_map_1 <- raster::getData(country = "ARG", level = 1, path = "./data/") %>% 
  # simplify
  rmapshaper::ms_simplify(keep = 0.01)  

for(region_looper in seq_along(arg_map_1@data$NAME_1)){
  arg_map_1@polygons[[region_looper]]@ID <- 
    as.character(arg_map_1@data$NAME_1[region_looper]) 
}

# tidy to a dataframe
arg_map_1 <- arg_map_1 %>% 
  broom::tidy()

library(ggplot2)
ggplot(data=arg_map_1) +
  geom_map(map=arg_map_1, aes(x=long, y=lat, map_id=id, fill=id),
           color="#000000", size=0.25)

我一直在想,必须有某种方法可以使用保留名称的 tidy 功能,但对于我的生活,我无法弄清楚。

4

2 回答 2

7

您可以使用joinpackage 中的功能plyr。这是一个通用的解决方案(看起来很长,但实际上很容易):

  1. 加载 shapefile:假设my_shapefile.shp您的工作目录中有一个 shapefile。让我们加载它:

    shape <- readOGR(dsn = "/my_working_directory", layer = "my_shapefile")
    

    请注意,在这个 shapefile 中有一个数据框,可以使用shape@data. 例如,此数据框可能如下所示:

    > head(shape@data)
           code                   region     label
    0 E12000006          East of England E12000006
    1 E12000007                   London E12000007
    2 E12000002               North West E12000002
    3 E12000001               North East E12000001
    4 E12000004            East Midlands E12000004
    5 E12000003 Yorkshire and The Humber E12000003
    
  2. 从 shapefile 创建新数据框:使用broom包来调整 shapefile 数据框:

    new_df <- tidy(shape)
    

这会导致这样的结果:

> head(new_df)
      long      lat order  hole piece group id           
1 547491.0 193549.0     1 FALSE     1   0.1  0 
2 547472.1 193465.5     2 FALSE     1   0.1  0 
3 547458.6 193458.2     3 FALSE     1   0.1  0 
4 547455.6 193456.7     4 FALSE     1   0.1  0 
5 547451.2 193454.3     5 FALSE     1   0.1  0 
6 547447.5 193451.4     6 FALSE     1   0.1  0

不幸的是,tidy()丢失了变量名(在本例中为“区域”)。相反,我们得到了一个新变量“id”,从 0 开始。幸运的是,“id”的顺序与存储在shape@data$region. 让我们用它来恢复名称。

  1. 使用行名创建辅助数据框:让我们使用行名创建一个新的数据框。此外,我们将添加一个“id”变量,与tidy()创建的变量相同:

    # Recover row name 
    temp_df <- data.frame(shape@data$region)
    names(temp_df) <- c("region")
    # Create and append "id"
    temp_df$id <- seq(0,nrow(temp_df)-1)
    
  2. 使用“id”将行名称与新数据框合并:最后,让我们将名称放回新数据框:

    new_df <- join(new_df, temp_df, by="id")
    

而已!您甚至可以使用join命令和“id”索引向新数据框添加更多变量。最终结果将类似于:

> head(new_df)
      long      lat order  hole piece group id            name    var1    var2 
1 547491.0 193549.0     1 FALSE     1   0.1  0 East of England   0.525   0.333   
2 547472.1 193465.5     2 FALSE     1   0.1  0 East of England   0.525   0.333   
3 547458.6 193458.2     3 FALSE     1   0.1  0 East of England   0.525   0.333   
4 547455.6 193456.7     4 FALSE     1   0.1  0 East of England   0.525   0.333   
5 547451.2 193454.3     5 FALSE     1   0.1  0 East of England   0.525   0.333   
6 547447.5 193451.4     6 FALSE     1   0.1  0 East of England   0.525   0.333   
于 2017-05-15T12:25:54.710 回答
2

alistaire 的评论促使我继续推动region=参数。我尝试了许多迭代,并在此线程https://github.com/tidyverse/ggplot2/issues/1447中找到了一些想法。

这是获取地区名称的代码:

# load the magrittr library to get the pipe
library(magrittr)
# load the maptools library to get the rgeos object
library(maptools)

arg_map_1 <- raster::getData(country = "ARG", level = 1, path = "./data/") %>% 
  # simplify
  rmapshaper::ms_simplify(keep = 0.01) %>% 
  # tidy to a dataframe
  broom::tidy(region="NAME_1")

# plot the map
library(ggplot2)
ggplot(data=arg_map_1) +
  geom_map(map=arg_map_1, aes(x=long, y=lat, map_id=id, fill=id),
           color="#000000", size=0.25)

首先,请注意必须加载 maptools 库才能使 tidy 操作正常工作。另外,我想强调的是,从中提取区域信息的变量必须用引号引起来。我一直错误地假设 broom 会以与其他 tidyverse 包(如 dplyr)识别未引用或被反引号包围的列名相同的方式识别变量名称。

于 2016-11-13T19:07:25.987 回答