R语言实战（2）-制作地区分布图及其动画

跟着xiangyun练习R语言

R
ggplot2
Author

Lee

Published

November 16, 2022

因为在工作工程中，经常要接触到一些地图的数据，并需要对其进行可视化操作，所以 我希望学习ggplot2包制作地图的详细过程。原文为xiangyun的博文。xiangyun主要使用的是baseR的技术，我则尝试使用tidy风格对文章代码进行复现。

1 地图制作-地区分布图

这是最常用也最常见的一个地图的种类。

1.1 数据准备

使用maps包中内置的vote.repub数据集。

library(tidyverse)
library(maps)
library(sf)
data(votes.repub)
str(votes.repub)
 num [1:50, 1:31] NA NA NA NA 18.8 ...
 - attr(*, "dimnames")=List of 2
  ..$ : chr [1:50] "Alabama" "Alaska" "Arizona" "Arkansas" ...
  ..$ : chr [1:31] "1856" "1860" "1864" "1868" ...

数据集为51行31列的矩阵，行代表各州，列代表大选年份。进一步查看votes.repub的维度信息：

attributes(votes.repub)
$dim
[1] 50 31

$dimnames
$dimnames[[1]]
 [1] "Alabama"        "Alaska"         "Arizona"        "Arkansas"      
 [5] "California"     "Colorado"       "Connecticut"    "Delaware"      
 [9] "Florida"        "Georgia"        "Hawaii"         "Idaho"         
[13] "Illinois"       "Indiana"        "Iowa"           "Kansas"        
[17] "Kentucky"       "Louisiana"      "Maine"          "Maryland"      
[21] "Massachusetts"  "Michigan"       "Minnesota"      "Mississippi"   
[25] "Missouri"       "Montana"        "Nebraska"       "Nevada"        
[29] "New Hampshire"  "New Jersey"     "New Mexico"     "New York"      
[33] "North Carolina" "North Dakota"   "Ohio"           "Oklahoma"      
[37] "Oregon"         "Pennsylvania"   "Rhode Island"   "South Carolina"
[41] "South Dakota"   "Tennessee"      "Texas"          "Utah"          
[45] "Vermont"        "Virginia"       "Washington"     "West Virginia" 
[49] "Wisconsin"      "Wyoming"       

$dimnames[[2]]
 [1] "1856" "1860" "1864" "1868" "1872" "1876" "1880" "1884" "1888" "1892"
[11] "1896" "1900" "1904" "1908" "1912" "1916" "1920" "1924" "1928" "1932"
[21] "1936" "1940" "1944" "1948" "1952" "1956" "1960" "1964" "1968" "1972"
[31] "1976"

为方便后续操作，将原数据的矩阵格式更改为数据框格式。

rownames <- rownames(votes.repub)
voted_repub <- as_tibble(votes.repub,
  # 保留原数据的行名
  rownames = pkgconfig::get_config("tibble::rownames", NA)
)

再增加两列state和state_name列，为后续与地图数据集关联，从而将数据映射到提图上。

state列可与 ggplot2::map_data() 生成的地图数据合并，state_name 列可与 tigris 包获取的数据合并。
voted_repub <-
  voted_repub %>%
  mutate(
    state_names = rownames(voted_repub),
    region = str_to_lower(rownames(voted_repub))
  ) %>%
  arrange(state_names)
Note

需要什么格式的州名取决于所能获得的地图数据的形式，一般需要对地图数据有所了解才能知晓。

接下来为方便后续的操作(如ggplot2绘图)，就是数据重塑。

us_voted_repub <-
  voted_repub %>%
  pivot_longer(-c(region, state_names),
    names_to = "year",
    values_to = "repub_prop"
  )
us_voted_repub
# A tibble: 1,550 × 4
   state_names region  year  repub_prop
   <chr>       <chr>   <chr>      <dbl>
 1 Alabama     alabama 1856       NA   
 2 Alabama     alabama 1860       NA   
 3 Alabama     alabama 1864       NA   
 4 Alabama     alabama 1868       51.4 
 5 Alabama     alabama 1872       53.2 
 6 Alabama     alabama 1876       40.0 
 7 Alabama     alabama 1880       37.0 
 8 Alabama     alabama 1884       38.4 
 9 Alabama     alabama 1888       32.3 
10 Alabama     alabama 1892        3.95
# ℹ 1,540 more rows

可以看到，us_vote_repub数据集存在很多缺失值，这是由于美国大选有些历史数据没有记录，我们需要关注从何时起，大选数据基本全覆盖了。

us_voted_repub %>%
  group_by(year) %>%
  summarise(n = sum(is.na(repub_prop))) %>%
  ggplot(aes(x = year, y = n, group = 1)) +
  geom_line() +
  scale_x_discrete(breaks = seq(1856, 1976, by = 4)) +
  scale_y_continuous(
    expand = c(0, 1),
    limits = c(0, 30)
  ) +
  theme(
    axis.text.x = element_text(angle = 45)
  )
图 1: 历届大选未收集到得票率数据的州有多少

如 图 1 所示，大选在20世纪初就基本覆盖了全美国。

接下来，以某一年的数据为基础，绘制一张静态地区分布图：

us_voted_repub_1976 <- us_voted_repub %>%
  filter(year == 1976)

# 获取美国洲际地图数据
usa_map <- map_data(map = "state")
usa_map %>%
  slice_head(n = 6)
       long      lat group order  region subregion
1 -87.46201 30.38968     1     1 alabama      <NA>
2 -87.48493 30.37249     1     2 alabama      <NA>
3 -87.52503 30.37249     1     3 alabama      <NA>
4 -87.53076 30.33239     1     4 alabama      <NA>
5 -87.57087 30.32665     1     5 alabama      <NA>
6 -87.58806 30.32665     1     6 alabama      <NA>

usa_map中，每一对经度long和纬度lat对应州边界多边形的一个点，按一定顺序连接起来，就可以形成每个州的行政区划边界。注意，由于usa_map是州级地图数据，不包含各郡的边界，故而 subregion列都为缺失值。

下面我们需要将观测数据合并到usa_map数据中，注意，地图数据在左边，观测数据在右侧，因为右侧数据是不完整的，而绘制地图时最好保证地图数据的完整。

usa_voting_df <- usa_map %>%
  left_join(us_voted_repub_1976, by = "region")
usa_voting_df %>%
  slice_head(n = 6)
       long      lat group order  region subregion state_names year repub_prop
1 -87.46201 30.38968     1     1 alabama      <NA>     Alabama 1976      43.48
2 -87.48493 30.37249     1     2 alabama      <NA>     Alabama 1976      43.48
3 -87.52503 30.37249     1     3 alabama      <NA>     Alabama 1976      43.48
4 -87.53076 30.33239     1     4 alabama      <NA>     Alabama 1976      43.48
5 -87.57087 30.32665     1     5 alabama      <NA>     Alabama 1976      43.48
6 -87.58806 30.32665     1     6 alabama      <NA>     Alabama 1976      43.48

由于地图上的每个区域都有编号顺序，不能随意打乱，原文中使用baseR需要对数据进行排序。而使用tidy风格的代码，order列的顺序没有改变，故不需要进行重新排序的操作。

1.2 数据展示

与普通的图形不同，绘制地图需要注意底图和投影坐标，ggplot2提供了geom_map()和coord_map()两个函数实现对应配置。如 图 2，是根据maps包的北美州级地图数据绘制而成的州级边界地图。

ggplot(
  usa_voting_df,
  aes(x = long, y = lat, group = subregion)
) +
  geom_map(
    map = usa_map, aes(map_id = region),
    color = "grey80", fill = "grey30", linewidth = 0.3
  )
图 2: 北美州级地图

接着需要设置坐标系。maps包内置的 state 数据集，其经纬度坐标来自常见 WGS 84 参考系，即所谓的世界大地测量系统 World Geodetic System，于 1984 制定，全球定位系统 GPS 即采用此坐标参考系。

为了在二维的地图上更加形象的展示地球的球形特点，使用mapproj包，以北极为中心的方位角进行投影：纬度为同心圆弧线，经度则是等距径向线，详细描述可参看?mapproj::mapproject具体如下：

ggplot(usa_voting_df, aes(long, lat, group = subregion)) +
  geom_map(
    map = usa_map, aes(map_id = region),
    fill = "grey30", color = "grey80", linewidth = 0.3
  ) +
  coord_map(
    projection = "orthographic",
    orientation = c(
      latitude = 39,
      longitude = -98, ortation = 0
    )
  )

接着，利用多边形图层geom_polygon()添加得票数据，采用不同颜色对得票多少进行区分，并对图形进行一定美化。

library(maps)
ggplot(usa_voting_df, aes(long, lat, group = subregion)) +
  geom_map(
    map = usa_map, aes(map_id = region),
    fill = "grey30", color = "grey80", linewidth = 0.3
  ) +
  coord_map(
    projection = "orthographic",
    orientation = c(
      latitude = 39,
      longitude = -98, ortation = 0
    )
  ) +
  geom_polygon(
    aes(group = group, fill = repub_prop / 100),
    color = "black"
  ) +
  scale_fill_gradient2(
    low = "blue", mid = "white", high = "red",
    midpoint = 0.5, labels = scales::label_percent()
  ) +
  theme_minimal() +
  labs(
    title = "1976 年共和党在各州的得票率",
    x = "", y = "", fill = "得票率"
  ) +
  theme(
    plot.title = element_text(
      size = 16, face = "bold", color = "red3"
    ),
    legend.position = "bottom"
  )
图 3: 1976 年共和党在各州的得票率

图 3 中的地图已经初具雏形，但还有几个缺陷：

  • maps包的state数据缺少阿拉斯加和夏威夷两个州。
  • 图中缺少南北纬，东西经的方向指示。

我们还需进一步对图形进行优化：

  • 使用tigris包下载美国的完整地图。(由于网络原因无法顺利下载，这里直接使用现成的rds数据文件)
  • 使用sf包对地图进行绘制。
us_state_map <- read_rds("d:/Myblog/datas/us_state_map.rds")
us_state_map
Simple feature collection with 52 features and 9 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -3112200 ymin: -1697728 xmax: 2258154 ymax: 1558935
Projected CRS: USA_Contiguous_Albers_Equal_Area_Conic
First 10 features:
   STATEFP  STATENS    AFFGEOID GEOID STUSPS         NAME LSAD        ALAND
1       53 01779804 0400000US53    53     WA   Washington   00 172117954267
2       46 01785534 0400000US46    46     SD South Dakota   00 196346195316
3       39 01085497 0400000US39    39     OH         Ohio   00 105823701267
4       01 01779775 0400000US01    01     AL      Alabama   00 131174192284
5       05 00068085 0400000US05    05     AR     Arkansas   00 134776580080
6       35 00897535 0400000US35    35     NM   New Mexico   00 314197253999
7       48 01779801 0400000US48    48     TX        Texas   00 676668210823
8       06 01779778 0400000US06    06     CA   California   00 403660088482
9       21 01779786 0400000US21    21     KY     Kentucky   00 102282218059
10      13 01705317 0400000US13    13     GA      Georgia   00 149484604230
        AWATER                       geometry
1  12549727444 MULTIPOLYGON (((-2000238 15...
2   3383460688 MULTIPOLYGON (((-633753 865...
3  10274036690 MULTIPOLYGON (((1081987 544...
4   4593183334 MULTIPOLYGON (((708460.5 -5...
5   2956395922 MULTIPOLYGON (((122656.3 -1...
6    727781442 MULTIPOLYGON (((-1226428 -5...
7  18991880422 MULTIPOLYGON (((-998043.8 -...
8  20305454540 MULTIPOLYGON (((-2066285 -2...
9   2372611005 MULTIPOLYGON (((571924.5 -8...
10  4420380296 MULTIPOLYGON (((939223.1 -2...

同样的，我们将地图数据与观测数据进行连接。注意，此处使用的是NAME。

library(maps)
library(sf)
us_voted_repub_1976 <-
  us_voted_repub_1976 %>%
  rename(NAME = "state_names")

usa_voting_map_1976 <- us_state_map %>%
  left_join(us_voted_repub_1976, by = "NAME")

使用sf包，可以大大简化梯度的绘制过程。

library(sf)
ggplot() +
  geom_sf(
    data = usa_voting_map_1976,
    aes(fill = repub_prop / 100)
  )
图 4: 1976年共和党在各州的得票率

如 图 4 所示，所有的州都显示了，且把阿拉斯加、夏威夷和波多黎各适当移动到了全美地图的下方，方面查看。

那么下一步，我们就要将具体的得票率映射在地图上。

ggplot() +
  geom_sf(
    data = usa_voting_map_1976,
    aes(fill = repub_prop / 100)
  ) +
  scale_fill_gradient2(
    low = "blue", mid = "white", high = "red",
    midpoint = 0.5, labels = scales::label_percent()
  ) +
  theme_minimal() +
  labs(
    title = "1976 年共和党在各州的得票率",
    x = "", y = "", fill = "得票率"
  ) +
  theme(
    plot.title = element_text(size = 16, face = "bold", color = "red3"),
    legend.position = "bottom"
  )
图 5: 1976年共和党在各州的得票率-美化