《经济评论》2024年第3期论文「高铁网络如何促进城市间合作创新——基于高铁网络通达性与合作专利的实证分析」中介绍了这样的一个指标：
论文全文：https://mp.weixin.qq.com/s/jDl_v0-gCHWgoTUaI86BCA
今天我们一起来学习下如何在 R 语言中高效的计算这一指标，也就是城市间合作申请专利数量。
前不久我们分享的「1985～2024年共同申请专利的所有申请人经纬度及其所属的省市区县」数据正好可以用来计算这一指标：
1985～2024年共同申请专利的所有申请人经纬度及其所属的省市区县：https://rstata.duanshu.com/#/brief/course/6ce40daa5ebe4164af66b07f68374c49
该数据是对所有申请人为企业的进行地理编码，由于直接对企业名称地理编码会错误率很高，所以先匹配工商注册信息，使用工商注册地址解析，然后匹配不到工商注册信息的使用企业名称解析获得经纬度，有了经纬度之后就可以根据经纬度判断每个企业所处的省市区县了。对于实在无法解析的企业，尽可能使用其名称中包含的城市信息。
在附件中我选择了 2015～2017 年的数据作为示例：
2015～2017年共同申请专利的所有申请人经纬度及其所属的省市区县.dta
首先读取该数据，去除重复的专利：
library(tidyverse)
haven::read_dta("2015～2017年共同申请专利的所有申请人经纬度及其所属的省市区县.dta") -> df
df %>%
mutate(公开公告号 = str_remove_all(公开公告号, "[A-Z]$")) %>%
distinct(公开公告号, .keep_all = T) %>%
distinct(newipzlid) -> unique_patent
|
这是因为该专利数据中同时包含了专利的引用与授权公告，直接删除授权公告又担心有些授权专利的申请公告没有包含在数据里面：
统计的时候可以先去除公开公告号里面的 A、B、U、S。其中 A 代表发明专利的申请公开，B 代表发明专利的授权公告，U 代表实用新型专利的授权公告，S 代表外观设计专利的授权公告。这里我索性去除了所有结尾的大写字母。
上面获取的是所有互不重复的专利编号，再和原始数据匹配，保留匹配成功的，就得到了所有互不重复的专利申请人信息：
df %>%
inner_join(unique_patent) %>%
mutate(专利类型 = if_else(str_detect(专利类型, "发明"), "发明", 专利类型)) -> df
df
|
#> # A tibble: 787,363 × 16
#> 年份 newipzlid 申请人编号 企业名称 经度 纬度 省 省代码 市 市代码
#> <dbl> <chr> <chr> <chr> <dbl> <dbl> <chr> <dbl> <chr> <dbl>
#> 1 2015 20150000221 申请人1 北京比动商贸有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> 2 2015 20150000221 申请人2 北京比动广告有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> 3 2015 20150000222 申请人1 北京比动商贸有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> 4 2015 20150000222 申请人2 北京比动广告有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> 5 2015 20150000223 申请人1 北京比动商贸有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> 6 2015 20150000223 申请人2 北京比动广告有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> 7 2015 20150000224 申请人1 北京比动商贸有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> 8 2015 20150000224 申请人2 北京比动广告有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> 9 2015 20150000225 申请人1 北京比动商贸有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> 10 2015 20150000225 申请人2 北京比动广告有限公… 116. 39.9 北京市…… 110000 北京市…… 110000
#> # ℹ 787,353 more rows
#> # ℹ 6 more variables: 县 <chr>, 县代码 <dbl>, 公开公告号 <chr>, 专利类型 <chr>,
#> # IPC <chr>, 授权公告日 <date>
|
然后按照论文中介绍的方法进行一些筛选：
然后就可以分年、城市对统计合作申请专利的数量了。不过在此之前我们先了解两个函数：
# combn() 函数：生成所有的城市对组合，通过 sort 确保城市对的一致性（如 "北京 - 上海" 和 "上海 - 北京" 被视为同一对）
combn(c("a", "b", "c"), 2)
|
#> [,1] [,2] [,3]
#> [1,] "a" "a" "b"
#> [2,] "b" "c" "c"
|
combn(sort(c("a", "b", "c")), 2)
|
#> [,1] [,2] [,3]
#> [1,] "a" "a" "b"
#> [2,] "b" "c" "c"
|
# 如果考虑反向情况，可以用 expand.grid
expand.grid(c("a", "b", "c"), c("a", "b", "c"))
|
#> Var1 Var2
#> 1 a a
#> 2 b a
#> 3 c a
#> 4 a b
#> 5 b b
#> 6 c b
#> 7 a c
#> 8 b c
#> 9 c c
|
使用 combn()：
考虑反向情况：
df %>%
group_by(newipzlid, 年份) %>%
summarise(城市 = list(unique(市)), .groups = "drop") %>%
filter(map_int(城市, length) > 1) %>%
mutate(城市对 = map(城市, ~expand.grid(.x, .x))) %>%
unnest(城市对) %>%
rename(城市1 = Var1, 城市2 = Var2) %>%
filter(城市1 != 城市2) %>%
group_by(年份, 城市1, 城市2) %>%
summarise(合作专利数量 = n(), .groups = "drop") %>%
arrange(年份, desc(合作专利数量)) -> df2b
df2b
|
#> # A tibble: 22,944 × 4
#> 年份 城市1 城市2 合作专利数量
#> <dbl> <fct> <fct> <int>
#> 1 2015 北京市 南京市 2531
#> 2 2015 南京市 北京市 2531
#> 3 2015 北京市 天津市 2260
#> 4 2015 天津市 北京市 2260
#> 5 2015 北京市 上海市 2250
#> 6 2015 上海市 北京市 2250
#> 7 2015 北京市 济南市 1852
#> 8 2015 济南市 北京市 1852
#> 9 2015 北京市 苏州市 1682
#> 10 2015 苏州市 北京市 1682
#> # ℹ 22,934 more rows
|
如果想要分专利类型统计合作数量，加入专利类型变量即可：
df3 %>%
spread(专利类型, 合作专利数量, fill = 0) %>%
mutate(总合作申请数量 = 发明 + 外观设计 + 实用新型)
|
#> # A tibble: 22,944 × 7
#> 年份 城市1 城市2 发明 外观设计 实用新型 总合作申请数量
#> <dbl> <fct> <fct> <dbl> <dbl> <dbl> <dbl>
#> 1 2015 北京市 深圳市 714 21 172 907
#> 2 2015 北京市 金华市 136 0 127 263
#> 3 2015 北京市 上海市 1822 32 396 2250
#> 4 2015 北京市 苏州市 901 65 716 1682
#> 5 2015 北京市 盐城市 39 0 43 82
#> 6 2015 北京市 温州市 49 1 38 88
#> 7 2015 北京市 青岛市 630 4 266 900
#> 8 2015 北京市 佛山市 28 14 12 54
#> 9 2015 北京市 潍坊市 294 2 315 611
#> 10 2015 北京市 临沂市 295 11 129 435
#> # ℹ 22,934 more rows
|
这样我们就完成了这个操作～
点击这里跳转到 RStata 短书平台获取附件：名师讲堂｜R 语言：如何根据专利申请数据构造城市间合作申请数量矩阵？
评论