将空间多边形转换为地图对象

Question

我正在尝试将 R 中的空间多边形对象转换为地图对象。我已经设法使用以前回答的问题来做到这一点，最具体地说是 (In R, how can I convert SpatialPolygons* to map objects)，但我遇到了与此线程中的人相同的问题：我得到奇怪的多边形链接，导致地图无法使用。

我在下面包含了完整的复制代码。

感谢您的帮助。

library(sp)
library(spdep)
library(maps)
library(rgdal)

## load a file from GADM (you just have to specify the countries "special part" of the file name, like "ARG" for Argentina. Optionally you can specify which level you want to have
loadGADM <- function (fileName, level = 0, ...) {
load(url(paste("http://gadm.org/data/rda/", fileName, "_adm", level, ".RData", sep     = "")))
gadm
}

## the maps objects get a prefix (like "ARG_" for Argentina)
changeGADMPrefix <- function (GADM, prefix) {
GADM <- spChFIDs(GADM, paste(prefix, row.names(GADM), sep = "_"))
GADM
}

## load file and change prefix
loadChangePrefix <- function (fileName, level = 0, ...) {
theFile <- loadGADM(fileName, level)
theFile <- changeGADMPrefix(theFile, fileName)
theFile
}

## this function creates a SpatialPolygonsDataFrame that contains all maps you specify in "fileNames".
## E.g.: 
## spdf <- getCountries(c("ARG","BOL","CHL"))
## plot(spdf) # should draw a map with Brasil, Argentina and Chile on it.
getCountries <- function (fileNames, level = 0, ...) {
polygon <- sapply(fileNames, loadChangePrefix, level)
polyMap <- do.call("rbind", polygon)
polyMap
}



indiamap <- getCountries("IND",level=1)

xx <- indiamap
require(reshape)
xx@data$id <- rownames(xx@data)

# Convert to dataframe
xx.df <- as.data.frame(xx)

#Fortfy automagic
require(ggplot2)
xx.fort <- fortify(xx, region="id")

# Join operation - one row per coordinate vector
require(plyr)
xx <- join(xx.fort, xx.df,by="id")

# Split by ID because we need to add NA at end of each set of polygon coordinates to 'break' the line
xxSp <- split(xx, xx$id)

# Need to insert NA at end of each polygon shape to cut off that shape
xxL <- do.call( rbind , (lapply( xxSp , function(x) { j <- x[ nrow(x) , ] ; j[1:2] <- c(NA,NA); rbind( x , j ) })) )


# Create list object with same structure as map object
xxMap <- list( x = xxL$long , y = xxL$lat , range = c( range(xxL$long) , range(xxL$lat) ) , names = as.character(unique( xxL$NAME ) ) )

# Define as a map class object
attr(xxMap , "class") <- "map"

# Plot!!
map( xxMap )

Answer 1

问题是每个ids（省）都由几个多边形组成，例如群岛。 SpatialPolygonsDataFrame 有信息表明这些单独的多边形是相同的id，但不应通过线连接。通过转换为map，该信息将丢失。在您的解决方法中，您只在ids 之间添加NAs，而不是在所有单独的多边形之间。

诀窍是为所有单独的片段生成一个唯一的 id，因此可以在每个小片段之后添加一个 NA。由join 创建的xx 包含一个名为piece 的字段。似乎该字段对于单个 id 中的每个多边形都有一个单独的编号。因此，例如带有id == 4 的多边形由 3 个部分组成，例如3个岛屿。然后可以通过以下方式生成所有多边形的唯一 ID：

xx$myid <- xx$ID_1*1000 + as.numeric(xx$piece)

每个省的ID * 1000加上每件的ID。印度的多边形内的最大块数为 212。完整的代码如下：

indiamap <- getCountries("IND",level=1)

xx <- indiamap
require(reshape)
xx@data$id <- rownames(xx@data)

# Convert to dataframe
xx.df <- as.data.frame(xx)
# xx.df$myid = xx$ID1*100 + xx$piece

#Fortfy automagic
require(ggplot2)
xx.fort <- fortify(xx, region="id")

# Join operation - one row per coordinate vector
require(plyr)
xx <- join(xx.fort, xx.df,by="id")
unique(xx$piece)

# generate unique id for each polygon piece
xx$myid <- xx$ID_1*1000 + as.numeric(xx$piece)

# Split by myid because we need to add NA at end of each set of polygon coordinates to 'break' the line
xxSp <- split(xx, xx$myid)

# Need to insert NA at end of each polygon shape to cut off that shape
# xxL <- do.call( rbind , (lapply( xxSp , function(x) { j <- x[ nrow(x) , ] ; j[1:2] <- c(NA,NA); rbind( x , j ) })) )
xxL <- do.call( rbind , (lapply( xxSp , function(x) { j <- x[nrow(x) , ] ; j[1:2] <- c(NA,NA); rbind( x , j ) })) )


# Create list object with same structure as map object
xxMap <- list( x = xxL$long , y = xxL$lat , range = c( range(xxL$long) , range(xxL$lat) ) , names = as.character(unique( xxL$NAME_1 ) ) )

# Define as a map class object
attr(xxMap , "class") <- "map"

map( xxMap , fill=T, col=2)