R 语言课程|相邻城市共同边界附近的工企与共同边界的距离计算 & 绘图展示

之前给大家分享过三份相邻省份、相邻城市、相邻区县附近工企到共同边界距离的数据。

×

×

×

今天就让我们以城市为例学习下这些数据是如何计算的,也就是如何提取相邻城市的共同边界并计算这两个城市的工企距离共同边界的距离。

加载所需的 R 包

library(tidyverse)
library(sf)

读取工企地理位置数据

工企地理位置数据是以 dta 格式存储的,可以使用 haven 包读取:

haven::read_dta("2012~2013工企地理位置.dta") -> df

剔除没有经纬度的:

df %>%
mutate(经度 = as.numeric(经度),
纬度 = as.numeric(纬度)) %>%
dplyr::filter(!is.na(经度), !is.na(纬度)) -> df

df

#> # A tibble: 668,812 × 17
#> group 年份 gqid 县代码 县 省代码 省 市代码 市 是否和省界接触
#> <dbl> <dbl> <dbl> <dbl> <chr> <dbl> <chr> <dbl> <chr> <dbl>
#> 1 129634 2012 2012116033 621025 正宁… 620000 甘肃… 621000 庆阳… 1
#> 2 129642 2012 2012026811 330784 永康… 330000 浙江… 330700 金华… 0
#> 3 129643 2012 2012173727 230621 肇州… 230000 黑龙… 230600 大庆… 0
#> 4 129666 2012 2012028024 350305 秀屿… 350000 福建… 350300 莆田… 0
#> 5 129668 2012 2012071855 330213 奉化… 330000 浙江… 330200 宁波… 0
#> 6 129683 2012 2012175526 360722 信丰… 360000 江西… 360700 赣州… 1
#> 7 129719 2012 2012313057 430104 岳麓… 430000 湖南… 430100 长沙… 0
#> 8 129737 2012 2012039949 230903 桃山… 230000 黑龙… 230900 七台… 0
#> 9 129746 2012 2012203887 340403 田家… 340000 安徽… 340400 淮南… 0
#> 10 129747 2012 2012034316 370404 峄城… 370000 山东… 370400 枣庄… 1
#> # … with 668,802 more rows, and 7 more variables: 是否和海岸线接触 <dbl>,
#> # 是否和陆地国界线接触 <dbl>, 是否为内部县 <dbl>, 与秦岭淮河线的距离 <dbl>,
#> # 北方或南方 <chr>, 经度 <dbl>, 纬度 <dbl>

后面我们会使用 sf 包进行地理距离的计算,所以我们需要把 df 数据转换成 sf 对象:

# df 转 sf 对象
df %>%
st_as_sf(coords = c("经度", "纬度"), crs = 4326) -> dfsf
dfsf

找到所有的相邻城市对

# 读取区市
read_sf("2019行政区划/市.shp") %>%
st_transform(4326) %>%
dplyr::filter(!is.na(市)) %>%
st_make_valid() -> city

# 获取所有的相邻市对
st_intersects(st_make_valid(city),
st_make_valid(city),
sparse = F) -> touchmat
touchmat %>%
as_tibble() %>%
mutate_all(as.numeric) %>%
mutate(市代码 = city$市代码) %>%
select(市代码, everything()) %>%
set_names(c("市代码", city$市代码)) %>%
gather(2:372, key = "key", value = "value") %>%
dplyr::filter(市代码 != key & value == 1) -> pairdf

pairdf

#> # A tibble: 1,946 × 3
#> 市代码 key value
#> <dbl> <chr> <dbl>
#> 1 120000 110000 1
#> 2 130600 110000 1
#> 3 130700 110000 1
#> 4 130800 110000 1
#> 5 131000 110000 1
#> 6 110000 120000 1
#> 7 130200 120000 1
#> 8 130800 120000 1
#> 9 130900 120000 1
#> 10 131000 120000 1
#> # … with 1,936 more rows

# 去除反向重复的(比较大小,仅仅保留前面的大于后面的)
pairdf %>%
mutate(key = as.numeric(key)) %>%
dplyr::filter(市代码 > key) %>%
rename(相邻市 = key) %>%
select(-value) -> pairdf

pairdf %>%
writexl::write_xlsx("2019年中国相邻市对.xlsx")

提取每个相邻相对的工企业,计算它们距离相邻市界的距离

使用 Stata

虽然计算过程是使用 R 语言完成的,但是得到的数据可以使用 Stata 继续处理,另外还可以使用 Stata 绘图:

循环所有的城市对

上面的代码中仅仅演示了第一个城市对的处理过程,只要循环所有的城市对就可以得到所有城市对的数据了,大家可以自己试试。

点击这里跳转到 RStata 短书平台获取附件:R 语言课程|相邻城市共同边界附近的工企与共同边界的距离计算 & 绘图展示

评论