在 R 中使用 ggplot2 在 Choropleth 地图上绘制数据不正确答案

【问题标题】：Incorrect plotting of data on Choropleth map with ggplot2 in R在 R 中使用 ggplot2 在 Choropleth 地图上绘制数据不正确
【发布时间】：2015-03-28 13:08:10
【问题描述】：

我正在尝试使用 ggplot2 制作德国的 Choropleth 地图。我的数据是一个 .csv 文件，有 2 行（RS= 每个德国州包含 1 到 16 的数字，Tariff= 16 个随机正数和负数）。

    RS  Tariff
1   01  -5.25
2   02  7.16
3   03  6.65
4   04  3.10
5   05  3.69
6   06  2.49
7   07  1.89
8   08  3.93
9   09  -5.84
10  10  -2.61
11  11  -0.21
12  12  2.35
13  13  -5.94
14  14  -7.54
15  15  -3.27
16  16  -8.75

我还有一个形状文件Germany shape file。我想要做的是将这个正数和负数绘制到每个州的德国地图上，有 2 种颜色（正数=绿色和负数=红色）。下面是我的代码

library(XLConnect)
library(sp)
library(rgdal)
library(ggplot2)
library(plyr)
library(RColorBrewer)
library(DataCombine)
library(rgeos)
library(maptools)

#### EEG Data Read ####
eeg<-read.csv(file = "data/testdata1.csv", skip = 0, sep = ",", dec=".",       header=TRUE)
colnames(eeg)<-c("RS", "Tariff")
eeg$RS<-   c("01","02","03","04","05","06","07","08","09","10","11","12","13","14","15","16")
eeg$RS<-as.factor(eeg$RS)
eeg$Tariff<- as.numeric(eeg$Tariff)

#### Shape Data Read ####
bundesl<-readOGR("data/new_shape/vg2500_bld.shp", "vg2500_bld")
bundesl@data<- bundesl@data[order(bundesl$RS, na.last=NA),]

### Rearrange shape data for better merging ###
levels(bundesl$GEN)<- c("Schleswig-Holstein", "Mecklenburg-Vorpommern",   "Hamburg", "Bremen", "Niedersachsen", "Sachsen-Anhalt", "Brandenburg", 
                    "Berlin", "Nordrhein-Westfalen", "Hessen","Thüringen","Sachsen", "Rheinland-Pfalz", "Saarland", "Baden-  Württemberg", "Bayern")
bundesl$GEN<- c("Schleswig-Holstein", "Mecklenburg-Vorpommern", "Hamburg", "Bremen", "Niedersachsen", "Sachsen-Anhalt", "Brandenburg", 
            "Berlin", "Nordrhein-Westfalen", "Hessen","Thüringen","Sachsen",  "Rheinland-Pfalz", "Saarland", "Baden-Württemberg", "Bayern")
bundesl$SHAPE_LENG<-  c("1217255.7","1780980.8","175253.8","154971.6","2016496.4","949096.8",
                   "1295460.4","180751.2","1352108","1105092.8","961942.7","979294.3","910650.4",
                   "282910.8","1298891.7","2046039.3")
bundesl$SHAPE_AREA<- c("15857425536","23044684847","760539820","405480872","47716406483","20494982327","29653902483","886480139","34047269991","21092318103","16178531941","18401642456","19834907486","2578541706","35801397076","70550070623")

# #### Shape Data und EEG Data join ####
bundesl@data<-merge(bundesl@data, eeg, by="RS",  all=TRUE)

# #### Shapes Plot ####
bundesl@data$id <- (as.numeric(rownames(bundesl@data))-1)
bundesl.df<-fortify(bundesland)
bundesl.df <- join(bundesl.df, bundesl@data, by="id")


 ggp <- ggplot(data=bundesl.df, aes(x=long, y=lat, group=group))
 ggp <- ggp + geom_polygon(aes(fill=Tariff), col="black") 
 ggp <- ggp + coord_map() 
 ggp <- ggp + scale_fill_continuous(name=expression(Tariff), low = "red", high = "green", space = "Lab", na.value = "white", guide = "colourbar")
ggp <- ggp + theme_minimal()
ggp <- ggp + theme(axis.title=element_blank(), axis.ticks=element_blank(),   axis.text=element_blank())
 ggp

到目前为止，我设法绘制了地图，但数据映射错误。我的意思是像石勒苏益格-荷尔斯泰因州这样关税为正的州应该是绿色的，但应该是红色的，巴伐利亚州应该是红色的，但应该是绿色的。

我的猜测是强化功能有问题。我的数据只有 16 行，但在强化后打印 1000 多行。为什么？？这是导致数据不匹配的原因。我在互联网上进行了所有可能的搜索以找到解决方案。如果有人能回答我为什么会出现这个问题，我将不胜感激。

提前感谢您的帮助！

【问题讨论】：

标签： r ggplot2 fortify choropleth

【解决方案1】：

fortify 将 shapefile 中的多边形放入ggplot 可以绘制的东西中，因此有 1,000 多行。虽然您可以为加固的多边形附加值，但这不是必需的。

所以，你真的不必为合唱团经历所有的麻烦。看看下面的内容。我添加了一些额外的位来显示哪些值映射到哪个RS：

library(rgdal)
library(ggplot2)
library(gridExtra)

egg <- read.table(text="RS  Tariff
01  -5.25
02  7.16
03  6.65
04  3.10
05  3.69
06  2.49
07  1.89
08  3.93
09  -5.84
10  -2.61
11  -0.21
12  2.35
13  -5.94
14  -7.54
15  -3.27
16  -8.75", header=TRUE, colClasses=c("character", "numeric"))

bundesl <- readOGR("vg2500_geo84/vg2500_bld.shp", "vg2500_bld")
bundesl@data<- bundesl@data[order(bundesl$RS, na.last=NA),]

# good projection for germany but if you intende to draw additional
# lines or points you'll have to project them before plotting so this
# may be more trouble than it's worth and you can just use 
# coord_map("mollweide") or something else that works for you besides mercator

bundesl <- spTransform(bundesl, CRS("+proj=utm +zone=33 +ellps=WGS84 +datum=WGS84 +units=m +no_defs "))

bundesl_map <- fortify(bundesl, region="RS")

# only doing this bit to plot the RS # at the center of each polygon
# totally not necessary for the choropleth

egg <- cbind(egg, data.frame(gCentroid(bundesl, byid=TRUE)))

gg <- ggplot()

# this bit ensures you have the outlines

gg <- gg + geom_map(data=bundesl_map, map=bundesl_map,
                    aes(x=long, y=lat, map_id=id),
                    color="#7f7f7f", size=0.15)

# this bit here does your choropleth

gg <- gg + geom_map(data=egg, map=bundesl_map,
                    aes(fill=Tariff, map_id=RS),
                    color="#7f7f7f", size=0.15)
gg <- gg + geom_text(data=egg, aes(x=x, y=y, label=RS), size=3)
gg <- gg + coord_equal() # we already projected it
gg <- gg + scale_fill_continuous(name=expression(Tariff), 
                                 low="red", high="green", space="Lab", 
                                 na.value="white", guide="colourbar")
gg <- gg + labs(x=NULL, y=NULL)

# decent map theme 

gg <- gg + theme_bw()
gg <- gg + theme(panel.grid=element_blank())
gg <- gg + theme(panel.border=element_blank())
gg <- gg + theme(axis.ticks=element_blank())
gg <- gg + theme(axis.text=element_blank())


gt <- tableGrob(cbind(bundesl@data[,c(2,4)], egg[,2]))

grid.arrange(gg, gt, ncol=2)

08 & 16 中包含 unicode，因此没有转换就无法显示。我也意识到在质心上绘制RS 数字对于柏林和勃兰登堡来说是有问题的，但这是为了给出一个总体思路，而不是完美的。

我强烈建议使用 cut 为值定义 5 或 6 个标准化中断，而不是使用连续刻度。

【讨论】：

感谢您的帮助@hrbmstr，您的回答澄清了问题所在。但我仍然无法解决实际问题。上帝帮助我:( P.S.抱歉这么晚才回复，生活受到了阻碍。再次感谢。
@hrbemstr 你能告诉我你为什么使用 spTransform 吗？
正如 cmets 所说，它将坐标预先投影到国际公认的适当投影上，以供德国使用。你可以摆脱它，要么只使用coord_map()（而不是coord_equal()），这将使ggplot为你使用投影（默认情况下它将是墨卡托），或者你可以使用建议的替代方案cmets中的一个coord_map("mollweide") 。如果没有coord_map 或spTransform + coord_equal，您将拥有一个偏向绘图尺寸的地图。
谢谢。我想我会继续使用 sptransform。

【解决方案2】：

这有点老了，但由于这是关于德国等值线地图的最佳问题，所以我想添加一些我在遵循 @hrbrmstr 的出色答案时学到的东西。

正如您在他绘制的地图中看到的那样，柏林刚刚获得勃兰登堡的颜色。要解决此问题，必须编辑 bundesl_map 中的顺序，以确保柏林 (10) 位于勃兰登堡 (11) 之后。所以地图的完整处理应该是这样的：

library(rgdal)
library(rgeos)
library(maptools)

bundesl <- readOGR("vg2500_geo84/vg2500_bld.shp", "vg2500_bld")
bundesl@data<- bundesl@data[order(bundesl$RS, na.last=NA),]
bundesl_map <- fortify(bundesl, region="RS")
bundesl_map <- rbind(
    bundesl_map[bundesl_map$id != 10, ], 
    bundesl_map[bundesl_map$id == 10, ]
)
saveRDS(bundesl_map, "bundesland")

在最后一步中，我们保存地图以供将来使用 (bundesl_map <- readRDS("budesland"))。 Here is a copy of the file I created with named ids.

绘图也可以更简洁一些，如下所示：

library(magrittr)
library(ggplot2)
library(viridis)

egg %>% ggplot(aes(fill=Tariff, map_id=RS)) +
    geom_map(map=bundesl_map, color="white", size=0.2) + 
    geom_text(aes(x=x, y=y, label=egg$RS), size=2, color="white") +
    coord_map("mercator") +
    expand_limits(x=bundesl_map$long, y=bundesl_map$lat) +
    scale_fill_viridis(begin=0.4, end=0.9, breaks=-8:7, guide=guide_legend(reverse=T)) + 
    theme_map(base_size=8)

其中theme_map 定义为：

theme_map <- function(...) {
    theme_classic(...) %+replace% 
    theme(
        axis.ticks = element_blank(), 
        axis.text = element_blank(),
        axis.title = element_blank(), 
        line=element_blank()
    )
}

这将产生一个像这样的地图：

【讨论】：