向抗疫勇士们致敬!

向抗疫勇士们致敬!

昨天是清明节,感觉不适合发推文。

其实我昨天是想仿照徐小洋同学(微信搜索公众号“走天涯徐小洋地理数据科学”)一样,绘制一幅图为逝者默哀,向抗疫勇士们致敬。

截至 4 月 4 日,已有 388 位勇士牺牲在了抗击新冠疫情的岗位上。

下面我对这幅图的绘制方法进行讲解。

首先从 清明,地图祭奠抗疫勇士 这篇推文里面爬取这份名单:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
library(rvest)
library(tidyverse)
library(lubridate)
library(sf)
html <- read_html('https://mp.weixin.qq.com/s/pZKM6LjCrjG7b57hKOOq1A')
html %>%
html_table(header = TRUE) -> table

# 删除数据框中的 ?
rmfun <- function(x){
str_remove_all(x, pattern = "?")
}

table[[1]] %>%
as_tibble() %>%
mutate_all(rmfun) %>%
type_convert() %>%
mutate(时间 = ymd(paste0("2020年", 时间))) -> herosdf

统计地理分布,等下绘制地图的时候我们只使用省份的前两个字匹配就好:

1
2
3
4
5
6
7
8
9
10
11
12
13
# 地理分布
herosdf %>%
mutate(省份 = str_sub(所在地, 1, 2)) %>%
mutate(省份 = case_when(
省份 == "乌鲁" ~ "新疆",
省份 == "呼伦" ~ "内蒙",
省份 == "哈尔" ~ "黑龙",
省份 == "广州" ~ "广东",
省份 == "建设" ~ "新疆",
省份 == "齐齐" ~ "黑龙",
T ~ 省份
)) %>%
count(省份) -> herosdf_prov

下面我们使用高德地图接口解析每个勇士牺牲的地点的经纬度,并把这个经纬度转换成 lcc 投影坐标系下的坐标。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
# 经纬度
# 编写一个根据地址解析经纬度的函数
library(jsonlite)

# 高德地图密钥
key <- Sys.getenv("amap.key")
get_loc <- function(address){
fromJSON(paste0('https://restapi.amap.com/v3/geocode/geo?address=', address, '&key=', key)) -> start
start$geocodes$location -> start
return(start)
}

herosdf %>%
count(所在地) -> herosdf_longlat

herosdf_longlat$tude = ""
for(i in 1:nrow(herosdf_longlat)){
try((herosdf_longlat$tude[i] = get_loc(herosdf_longlat$所在地[i])))
}

herosdf_longlat %>%
separate(tude, into = c("lon", "lat"), sep = ",") %>%
type_convert() %>%
dplyr::filter(!is.na(lon)) %>%
st_as_sf(coords = c("lon", "lat")) %>%
st_set_crs(4326) %>%
st_transform(crs = "+proj=lcc +lat_1=30 +lat_2=62 +lat_0=0 +lon_0=105 +x_0=0 +y_0=0 +ellps=krass +units=m +no_defs") -> herosdf_sf
herosdf_sf

cbind(st_drop_geometry(herosdf_sf),
st_coordinates(herosdf_sf)) -> herosdf_coords

最后我们就可以绘制地图 + 柱状图了,最后使用 cowplot 进行拼图:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
# 读取中国省份数据
library(sf)
boundary <- read_sf('chinamap/国界线.shp')
read_sf('chinamap/中国省界.shp') %>%
mutate(省份 = str_sub(NAME, 1, 2)) %>%
left_join(herosdf_prov) %>%
mutate(group = cut(n,
breaks = c(1, 4, 8, 13, 32, 60),
labels = c("1 - 4", "5 - 8", "9 - 13",
"14 - 32", "33 - 60"),
include.lowest = TRUE)) -> provmap

library(ggimage)
ggplot(provmap) +
geom_sf(data = subset(provmap, !is.na(group)),
aes(fill = group), size = 0.1) +
geom_sf(data = boundary, size = 0.2) +
geom_sf_text(aes(label = NAME), family = cnfont) +
geom_image(data = herosdf_coords,
aes(X, Y, image = "candle.png"),
size = 0.02) +
theme_ipsum(base_family = cnfont) +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.text.x = element_blank(),
axis.text.y = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
legend.box = "horizontal") +
scale_fill_brewer(palette = "Greys") +
labs(fill = "各省牺牲人数") -> p1
p1

provmap %>%
subset(!is.na(n)) %>%
mutate(NAME = forcats::fct_reorder(NAME, -n)) %>%
st_drop_geometry() %>%
ggplot(aes(x = NAME, y = n, fill = n, color = n)) +
geom_col() +
coord_flip() +
theme_ipsum(base_family = cnfont) +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position = "none") +
scale_fill_gradientn(colors = brewer.pal(9, "Greys")[4:9]) +
scale_color_gradientn(colors = brewer.pal(9, "Greys")[4:9]) +
scale_y_continuous(breaks = seq(from = 0, to = 60, by = 5)) +
labs(x = "", y = "抗疫牺牲人数") +
theme(panel.background = element_blank(),
plot.background = element_blank()) -> p2
p2

library(cowplot)
ggdraw() +
draw_plot(p1, x = 0.1, y = 0.03) +
draw_plot(p2)

这样就得到了开头绘制的那幅图了。

知识星球附件链接:https://t.zsxq.com/qnmIEYF

#

评论

Your browser is out-of-date!

Update your browser to view this website correctly. Update my browser now

×