使用 ggmap 绘制形状文件：当形状文件大于 ggmap 时进行裁剪答案

【问题标题】：Plotting shape files with ggmap: clipping when shape file is larger than ggmap使用 ggmap 绘制形状文件：当形状文件大于 ggmap 时进行裁剪
【发布时间】：2015-04-23 13:48:02
【问题描述】：

当我尝试将 ggmap 与 shape 文件结合时，我遇到了裁剪问题。 Kahle and Wickham (2013: 158) 中的示例运行良好，因为来自 ggmap 的光栅图像覆盖了整个形状文件。下面是一个示例，说明当我尝试在覆盖较小区域的 ggmap 图上绘制美国各州的形状文件时会发生什么。 ggmap 显示了纽约市，我想将其与美国各州的边界重叠（仅作为示例）。生成的地图没有任何意义。问题是形状文件被剪裁了，而 ggplot 连接了未剪裁的点。下面是代码。形状文件来自here。我只是在这里展示最后一个情节。

我该如何解决这个问题？

path <- "PATH TO SHAPEFILE"
library("ggmap")
library("rgdal")

# shapefile
states <- readOGR(dsn = path, layer = "states")
states_df <- fortify(states)
# plot shapefile
plot(states, lwd = 0.1)
ggplot(states_df, aes(long, lat, group = group)) +
    geom_polygon(colour = "black", fill = NA, size = 0.1)


# combine ggmap with shapefile
map <- get_map("new york city", zoom = 10, source = "stamen")
ggmap(map, extent = "device")

ggmap(map, extent = "device") +
    geom_polygon(aes(long, lat, group=group), data = states_df, colour = "red", fill = NA, size = 1)

Kahle、David 和 Hadley Wickham。 2013.“Ggmap：使用 ggplot2 进行空间可视化”。 R 杂志 5(1):144-61。

【问题讨论】：

这里的问题与这个问题有关吗：stackoverflow.com/questions/13469566/…?

标签： r ggplot2 ggmap

【解决方案1】：

这是我的尝试。我经常使用 GADM shapefile，您可以使用 raster 包直接导入。我对 NY、NJ 和 CT 的形状文件进行了子集化。您最终可能不必这样做，但减少数据量可能会更好。当我绘制地图时，ggplot 自动删除了留在 ggmap 图像的 bbox 之外的数据点。因此，我不必做任何额外的工作。我不确定您使用的是哪个 shapefile。但是，GADM 的数据似乎适用于 ggmap 图像。希望这对您有所帮助。

library(raster)
library(rgdal)
library(rgeos)
library(ggplot2)

### Get data (shapefile)
us <- getData("GADM", country = "US", level = 1)

### Select NY and NJ
states <- subset(us, NAME_1 %in% c("New York", "New Jersey", "Connecticut"))

### SPDF to DF
map <- fortify(states)

## Get a map
mymap <- get_map("new york city", zoom = 10, source = "stamen")


ggmap(mymap) +
geom_map(data = map, map = map, aes(x = long, y = lat, map_id = id, group = group))

如果你只想要线条，下面就是你所追求的。

ggmap(mymap) +
geom_path(data = map, aes(x = long, y = lat, group = group))

【讨论】：

不错的结果，使用 geom_path 处理多边形会遇到的一个问题是，如果形状文件中的顶点数不够高，那么在地图的边缘。
这是否解决了截断问题ggmap(mymap) + geom_map(data = map, map = map, aes(x = long, y = lat, map_id = id, group = group), color = "black", fill=NA, size = 0.5)？我在边界线的中上端看到了问题，但geom_map 和fill=NA 似乎不是这种情况。
@user2503795 看来您的代码运行良好。我为你感到高兴。 :)
@jazzurro：谢谢！我认为根据您的解决方案对我的问题最简单的答案是：使用geom_map 而不是geom_polygon 并调整参数（将map = [fortified dataframe] 和map_id = id 添加为aes 调用的美学）。
@user2503795 我很高兴听到您根据此答案找到了解决方案。感谢您在这里分享您自己的想法。这将对将来的其他一些用户有所帮助。 :)

【解决方案2】：

我会查看这个答案，当您放大时，您所期望的 ggmap 似乎无法以理想的方式处理多边形，即不在绘图上的项目被截断，导致形状文件的“有趣”结果.

Polygons nicely cropping ggplot2/ggmap at different zoom levels

# transform for good measure
states <- spTransform(states,CRS("+datum=WGS84 +proj=longlat") )

# combine ggmap with shapefile
states_df <- fortify(states)

# get your map
map <-get_map("new york city", zoom = 10, source = "stamen")

a <- ggmap(map, # this is where we get our raster
       base_layer=ggplot(aes(x=long, y=lat), data=states_df), # this defines the region where things are plotted
       extent = "normal",  # this won't work with device, you need normal (see examples in ggmap documentation)
       maprange=FALSE
       ) +
coord_map( # use map's bounding box to setup the 'viewport' we want to see
  projection="mercator",
  xlim= c(attr(map, "bb")$ll.lon, attr(map, "bb")$ur.lon),
  ylim=c(attr(map, "bb")$ll.lat, attr(map, "bb")$ur.lat)
) +
geom_polygon( # plot the polygon
  aes(x=long, y=lat,group=group), data =states_df, color = "red", fill=NA, size = 1)

print(a)

有输出：

作为旁注，您可能想使用美国人口普查数据查看州地图，它们似乎比 ESRI 数据集质量更高。

ftp://ftp2.census.gov/geo/pvs/tiger2010st/tl_2010_us_state10.zip

最后一点，两极附近的 ggmap 存在问题，因此我还将按您感兴趣的州对您的数据进行子集化。

【讨论】：