R:创建选定的加拿大省份和美国各州的 map

Question

我正在尝试创建选定的加拿大省/地区和选定的美国州的 map 。到目前为止，最好的 map 似乎是使用 GADM 数据生成的 map :http://www.gadm.org/

但是，我无法在同一张 map 上绘制美国和加拿大，或者仅绘制选定的省/地区和州。例如，我对阿拉斯加、育空地区、西北地区、不列颠哥伦比亚省、艾伯塔省和蒙大拿州等感兴趣。

此外，美国 map 似乎是沿着国际日期变更线分割的。

有人可以帮我吗:

1. 在一张 map 上绘制上述省/地区和州
2. 避免美国沿着国际日期变更线 split
3. 叠加经纬度网格
4. 选择一个特定的投影，可以是多圆锥投影。

也许 spplot 不允许用户指定投影。我在 spplot 帮助页面上没有看到选择投影的选项。我知道如何使用 map 包中的 map 函数选择投影，但这些 map 看起来不太好，而且我也无法使用该函数绘制所需的省/地区和州子集。

我不知道如何开始添加经纬度网格。然而，文件“sp.pdf”的第 3.2 节似乎解决了这个主题。

下面是我到目前为止想出的代码。我已经加载了我偶然发现的每个与 map 相关的包，并注释掉了除省/地区或州边界之外的 GADM 数据。

不幸的是，到目前为止我只成功绘制了加拿大或美国的 map

library(maps)
library(mapproj)
library(mapdata)
library(rgeos)
library(maptools)
library(sp)
library(raster)
library(rgdal)

# can0<-getData('GADM', country="CAN", level=0) # Canada
  can1<-getData('GADM', country="CAN", level=1) # provinces
# can2<-getData('GADM', country="CAN", level=2) # counties

plot(can1)    
spplot(can1, "NAME_1") # colors the provinces and provides
                       # a color-coded legend for them
can1$NAME_1            # returns names of provinces/territories
# us0 <- getData('GADM', country="USA", level=0)
  us1 <- getData('GADM', country="USA", level=1)
# us2 <- getData('GADM', country="USA", level=2)
plot(us1)              # state boundaries split at 
                       # the dateline
us1$NAME_1             # returns names of the states + DC
spplot(us1, "ID_1")
spplot(us1, "NAME_1")  # color codes states and
                       # provides their names
#
# Here attempting unsuccessfully to combine U.S. and Canada on one map.
# Attempts at selecting given states or provinces have been unsuccessful.
#
plot(us1,can1)
us.can1 <- rbind(us1,can1)

感谢您的帮助。到目前为止，我在上述步骤 2 - 4 上还没有取得任何进展。也许是我要求太多了。也许我应该切换到 ArcGIS 并尝试该软件。

我已阅读此 StackOverflow 帖子:

Can R be used for GIS?

编辑

我现在借用了 Bevand 等人的《R 应用空间数据分析》的电子版。 (2008) 并从本书的网站下载(或找到)相关的 R 代码和数据:

http://www.asdar-book.org/

我还在这里发现了一些漂亮的 GIS 相关 R 代码:

https://sites.google.com/site/rodriguezsanchezf/news/usingrasagis

如果我学会了如何实现预期目标，我将在这里发布解决方案。尽管如果我无法实现 R 中的目标，我最终可能会转向 ArcGIS。

Answer 1

要在同一设备上绘制多个 SpatialPolygons 对象，一种方法是先指定要绘制的地理范围，然后使用 plot(..., add=TRUE) 。这将仅将那些感兴趣的点添加到 map 中。

使用投影(例如多圆锥投影)进行绘图需要首先使用 rgdal 包中的 spTransform() 函数，以确保所有图层都位于相同的位置投影。

## Specify a geographic extent for the map
## by defining the top-left and bottom-right geographic coordinates
mapExtent <- rbind(c(-156, 80), c(-68, 19))

## Specify the required projection using a proj4 string
## Use http://www.spatialreference.org/ to find the required string
## Polyconic for North America
newProj <- CRS("+proj=poly +lat_0=0 +lon_0=-100 +x_0=0 
            +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs")

## Project the map extent (first need to specify that it is longlat) 
mapExtentPr <- spTransform(SpatialPoints(mapExtent, 
                  proj4string=CRS("+proj=longlat")),
                  newProj)

## Project other layers
can1Pr <- spTransform(can1, newProj)
us1Pr <- spTransform(us1, newProj) 

## Plot each projected layer, beginning with the projected extent
plot(mapExtentPr, pch=NA)
plot(can1Pr, border="white", col="lightgrey", add=TRUE)
plot(us1Pr, border="white", col="lightgrey", add=TRUE)

使用相同的方法可以轻松地向 map 添加其他功能，例如突出显示感兴趣的管辖区:

## Highlight provinces and states of interest
theseJurisdictions <- c("British Columbia",
                        "Yukon",
                        "Northwest Territories",
                        "Alberta",
                        "Montana",
                        "Alaska")

plot(can1Pr[can1Pr$NAME_1 %in% theseJurisdictions, ], border="white", 
    col="pink", add=TRUE)

plot(us1Pr[us1Pr$NAME_1 %in% theseJurisdictions, ], border="white", 
    col="pink", add=TRUE)

结果如下:

我认为，在使用投影时添加网格线非常复杂，需要另一篇文章。看起来好像是@Mark Miller 在下面添加的!