使用美国县级数据创建 Choropleth 地图
Creating a Choropleth map with US county level data
我正在尝试使用 R 生成县级 COVID-19 感染数据的等值线图。我是 R 的新手,所以....
我用 ggmap 做了一些相当基本的东西来绘制空间数据,但从来没有像这样的东西。通常我只需要将兴趣点覆盖在地图上,因此我可以使用 geom_point 及其 lat/lon。在这种情况下,我需要构建底层地图,然后填充区域,而我在 ggplot 世界中努力做到这一点。
我遵循了一些我发现的在线示例:
library(ggplot2)
library(broom)
library(geojsonio)
#get a county level map geoJSON file
counties <- geojson_read("https://eric.clst.org/assets/wiki/uploads/Stuff/gz_2010_us_050_00_500k.json", what = "sp")
#filter our alaska and Hawaii
lower48 <- counties[(counties@data$STATE != "02" & counties@data$STATE != "15") ,]
#turn it into a dataframe for ggmap
new_counties <- tidy(lower48)
# Plot it
print(ggplot() +
geom_polygon(data = new_counties, aes( x = long, y = lat, group = group), fill="#69b3a2", color="white") +
theme_void() +
coord_map())
产生这个情节的是:
到目前为止一切顺利。但是我的 new_counties 数据框现在看起来像这样:
head(new_counties)
# A tibble: 6 x 7
long lat order hole piece group id
<dbl> <dbl> <int> <lgl> <chr> <chr> <chr>
1 -85.4 33.9 1 FALSE 1 0.1 0
2 -85.4 33.9 2 FALSE 1 0.1 0
3 -85.4 33.9 3 FALSE 1 0.1 0
4 -85.4 33.9 4 FALSE 1 0.1 0
5 -85.4 33.9 5 FALSE 1 0.1 0
6 -85.4 33.8 6 FALSE 1 0.1 0
所以我丢失了任何可以与县级感染数据联系起来的东西。
我的数据有每个县的 5 位 FIPS 代码。前两位数字是州,后三位是县。我的 geoJSON 文件有更详细的 FIPS 代码。我尝试只抓取前 5 个并创建我自己的数据元素,我可以映射回
library(ggplot2)
library(broom)
library(geojsonio)
#get a county level map geoJSON file
counties <- geojson_read("https://eric.clst.org/assets/wiki/uploads/Stuff/gz_2010_us_050_00_500k.json", what = "sp")
#filter our alaska and Hawaii
lower48 <- counties[(counties@data$STATE != "02" & counties@data$STATE != "15") ,]
#add my own FIPS code
lower48@data$myFIPS <- substr(as.character(lower48@data$GEO_ID),1,5)
#turn it into a dataframe for ggmap
new_counties <- tidy(lower48, region = "myFIPS")
# Plot it
print(ggplot() +
geom_polygon(data = new_counties, aes( x = long, y = lat, group = group), fill="#69b3a2", color="white") +
theme_void() +
coord_map())
但这会产生这个情节
而且我不得不说我对 broom::tidy 还不够熟悉,无法确切知道原因。
我还注意到,在我键入此内容时,我需要过滤掉波多黎各!
如果有人能给我指点有用的方向....我不拘泥于当前的方法,但我想坚持使用 ggplot2 或 ggmap。我的老板最终希望我覆盖某些功能。最终目标是按照示例 here 制作一个显示随时间变化的数据的动画地图,但我显然离那个还有很长的路要走。
有很多方法可以做到这一点,但概念基本相同:查找包含国家级 FIPS 代码的地图并将它们用于 link 数据源,也包含相同的 FIPS 代码作为绘图的变量(这里是每天的 covid-19 病例数)。
#devtools::install_github("UrbanInstitute/urbnmapr")
library(urbnmapr) # For map
library(ggplot2) # For map
library(dplyr) # For summarizing
library(tidyr) # For reshaping
library(stringr) # For padding leading zeros
# Get COVID cases, available from:
url <- "https://static.usafacts.org/public/data/covid-19/covid_confirmed_usafacts.csv
?_ga=2.162130428.136323622.1585096338-408005114.1585096338"
COV <- read.csv(url, stringsAsFactors = FALSE)
names(COV)[1] <- "countyFIPS" # Fix the name of first column. Why!?
数据以宽格式存储,每个县的每日病例分布在各列中。这需要在与地图合并之前收集。日期需要转换为正确的日期。 FIPS 代码以整数形式存储,因此需要将其转换为具有前导 0 的字符,以便与地图数据合并。我使用 urbnmap 包作为地图数据。
Covid <- pivot_longer(COV, cols=starts_with("X"),
values_to="cases",
names_to=c("X","date_infected"),
names_sep="X") %>%
mutate(date_infected = as.Date(date_infected, format="%m.%d.%Y"),
countyFIPS = str_pad(as.character(countyFIPS), 5, pad="0"))
# Obtain map data for counties (to link with covid data) and states (for showing borders)
states_sf <- get_urbn_map(map = "states", sf = TRUE)
counties_sf <- get_urbn_map(map = "counties", sf = TRUE)
# Merge county map with total cases of cov
counties_cov <- inner_join(counties_sf, group_by(Covid, countyFIPS) %>%
summarise(cases=sum(cases)), by=c("county_fips"="countyFIPS"))
counties_cov %>%
ggplot() +
geom_sf(mapping = aes(fill = cases), color = NA) +
geom_sf(data = states_sf, fill = NA, color = "black", size = 0.25) +
coord_sf(datum = NA) +
scale_fill_gradient(name = "Cases", trans = "log", low='pink', high='navyblue',
na.value="white", breaks=c(1, max(counties_cov$cases))) +
theme_bw() + theme(legend.position="bottom", panel.border = element_blank())
对于动画,您可以使用 gganimate 包并在白天进行过渡。这些命令与上面类似,只是不应汇总 covid 数据。
library(gganimate)
counties_cov <- inner_join(counties_sf, Covid, by=c("county_fips"="countyFIPS"))
p <- ggplot(counties_cov) + ... # as above
p <- p + transition_time(date_infected) +
labs(title = 'Date: {frame_time}')
animate(p, end_pause=30)
我正在尝试使用 R 生成县级 COVID-19 感染数据的等值线图。我是 R 的新手,所以....
我用 ggmap 做了一些相当基本的东西来绘制空间数据,但从来没有像这样的东西。通常我只需要将兴趣点覆盖在地图上,因此我可以使用 geom_point 及其 lat/lon。在这种情况下,我需要构建底层地图,然后填充区域,而我在 ggplot 世界中努力做到这一点。
我遵循了一些我发现的在线示例:
library(ggplot2)
library(broom)
library(geojsonio)
#get a county level map geoJSON file
counties <- geojson_read("https://eric.clst.org/assets/wiki/uploads/Stuff/gz_2010_us_050_00_500k.json", what = "sp")
#filter our alaska and Hawaii
lower48 <- counties[(counties@data$STATE != "02" & counties@data$STATE != "15") ,]
#turn it into a dataframe for ggmap
new_counties <- tidy(lower48)
# Plot it
print(ggplot() +
geom_polygon(data = new_counties, aes( x = long, y = lat, group = group), fill="#69b3a2", color="white") +
theme_void() +
coord_map())
产生这个情节的是:
到目前为止一切顺利。但是我的 new_counties 数据框现在看起来像这样:
head(new_counties)
# A tibble: 6 x 7
long lat order hole piece group id
<dbl> <dbl> <int> <lgl> <chr> <chr> <chr>
1 -85.4 33.9 1 FALSE 1 0.1 0
2 -85.4 33.9 2 FALSE 1 0.1 0
3 -85.4 33.9 3 FALSE 1 0.1 0
4 -85.4 33.9 4 FALSE 1 0.1 0
5 -85.4 33.9 5 FALSE 1 0.1 0
6 -85.4 33.8 6 FALSE 1 0.1 0
所以我丢失了任何可以与县级感染数据联系起来的东西。
我的数据有每个县的 5 位 FIPS 代码。前两位数字是州,后三位是县。我的 geoJSON 文件有更详细的 FIPS 代码。我尝试只抓取前 5 个并创建我自己的数据元素,我可以映射回
library(ggplot2)
library(broom)
library(geojsonio)
#get a county level map geoJSON file
counties <- geojson_read("https://eric.clst.org/assets/wiki/uploads/Stuff/gz_2010_us_050_00_500k.json", what = "sp")
#filter our alaska and Hawaii
lower48 <- counties[(counties@data$STATE != "02" & counties@data$STATE != "15") ,]
#add my own FIPS code
lower48@data$myFIPS <- substr(as.character(lower48@data$GEO_ID),1,5)
#turn it into a dataframe for ggmap
new_counties <- tidy(lower48, region = "myFIPS")
# Plot it
print(ggplot() +
geom_polygon(data = new_counties, aes( x = long, y = lat, group = group), fill="#69b3a2", color="white") +
theme_void() +
coord_map())
但这会产生这个情节
而且我不得不说我对 broom::tidy 还不够熟悉,无法确切知道原因。 我还注意到,在我键入此内容时,我需要过滤掉波多黎各!
如果有人能给我指点有用的方向....我不拘泥于当前的方法,但我想坚持使用 ggplot2 或 ggmap。我的老板最终希望我覆盖某些功能。最终目标是按照示例 here 制作一个显示随时间变化的数据的动画地图,但我显然离那个还有很长的路要走。
有很多方法可以做到这一点,但概念基本相同:查找包含国家级 FIPS 代码的地图并将它们用于 link 数据源,也包含相同的 FIPS 代码作为绘图的变量(这里是每天的 covid-19 病例数)。
#devtools::install_github("UrbanInstitute/urbnmapr")
library(urbnmapr) # For map
library(ggplot2) # For map
library(dplyr) # For summarizing
library(tidyr) # For reshaping
library(stringr) # For padding leading zeros
# Get COVID cases, available from:
url <- "https://static.usafacts.org/public/data/covid-19/covid_confirmed_usafacts.csv
?_ga=2.162130428.136323622.1585096338-408005114.1585096338"
COV <- read.csv(url, stringsAsFactors = FALSE)
names(COV)[1] <- "countyFIPS" # Fix the name of first column. Why!?
数据以宽格式存储,每个县的每日病例分布在各列中。这需要在与地图合并之前收集。日期需要转换为正确的日期。 FIPS 代码以整数形式存储,因此需要将其转换为具有前导 0 的字符,以便与地图数据合并。我使用 urbnmap 包作为地图数据。
Covid <- pivot_longer(COV, cols=starts_with("X"),
values_to="cases",
names_to=c("X","date_infected"),
names_sep="X") %>%
mutate(date_infected = as.Date(date_infected, format="%m.%d.%Y"),
countyFIPS = str_pad(as.character(countyFIPS), 5, pad="0"))
# Obtain map data for counties (to link with covid data) and states (for showing borders)
states_sf <- get_urbn_map(map = "states", sf = TRUE)
counties_sf <- get_urbn_map(map = "counties", sf = TRUE)
# Merge county map with total cases of cov
counties_cov <- inner_join(counties_sf, group_by(Covid, countyFIPS) %>%
summarise(cases=sum(cases)), by=c("county_fips"="countyFIPS"))
counties_cov %>%
ggplot() +
geom_sf(mapping = aes(fill = cases), color = NA) +
geom_sf(data = states_sf, fill = NA, color = "black", size = 0.25) +
coord_sf(datum = NA) +
scale_fill_gradient(name = "Cases", trans = "log", low='pink', high='navyblue',
na.value="white", breaks=c(1, max(counties_cov$cases))) +
theme_bw() + theme(legend.position="bottom", panel.border = element_blank())
对于动画,您可以使用 gganimate 包并在白天进行过渡。这些命令与上面类似,只是不应汇总 covid 数据。
library(gganimate)
counties_cov <- inner_join(counties_sf, Covid, by=c("county_fips"="countyFIPS"))
p <- ggplot(counties_cov) + ... # as above
p <- p + transition_time(date_infected) +
labs(title = 'Date: {frame_time}')
animate(p, end_pause=30)