上次课程主页：https://rstata.duanshu.com/#/brief/course/582c244c7a8b4beeb8823699c88c6cc9
继续上次的内容，今天我们再来学习如何基于真实专利数据计算数实融合水平。
仿照示例数据，我们这里也需要准备三份数据：
- 专利的 IPC 分类号；
- IPC 与国民经济行业的对照；
- 数字经济与实体经济的行业分类。
Stata 如何读取超大的数据文件
首先我们准备专利的 IPC 分类号数据。专利的原始数据可以从之前分享的这里下载：
1985～2024 年专利申请与授权数据（版本 3，含申请人所处的省市区县）：https://rstata.duanshu.com/#/brief/course/2397451274c546d3a36e156ffc865988
原始数据非常大，我们仅仅保留需要的：
cap mkdir "patent_small1"
cap mkdir "patent_small2"
forval y = 1985/2024 {
use /Volumes/rstata2t/newIP专利数据分年/`y'.dta, clear
drop if 专利类型 == "外观设计"
gen value = 公开公告号
replace value = subinstr(value, "A", "", .)
replace value = subinstr(value, "B", "", .)
replace value = subinstr(value, "U", "", .)
replace value = subinstr(value, "S", "", .)
duplicates drop value, force
keep newipzlid IPC 省 省代码 市 市代码 县 县代码
drop if mi(县代码)
save patent_small1/`y', replace
keep newipzlid IPC
save patent_small2/`y', replace
}
|
不过由于部分年份的专利数据非常大，所以大家的电脑可能无法直接运行上面的循环。不过也可以通过每次读取部分行来解决：
Stata 读取数据的大小受到电脑内存的限制，因此对于超大文件可以通过下面的方式读取：
use 命令有种用法是：
use [varlist] [if] [in] using filename [, clear nolabel]
|
例如每次读 100 条：use in 1/100 using 文件名.dta, clear
不过读取最后一部分的时候需要知道文件的总行数。
获取文件的总行数可以使用 describe 命令，该命令有种用法是：
describe [varlist] using filename [, file_options]
|
例如获取 patent_small2/2012.dta 数据的总行数：
desc using patent_small2/2012.dta
*> Contains data 数据处理：微信公众号
*> RStata
*> Observations: 978,031 11 Mar 2025 01:03
*> Variables: 2
*> ---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
*> Variable Storage Display Value
*> name type format label Variable label
*> ---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
*> newipzlid str11 %11s
*> IPC str788 %-9s
*> ---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
*> Sorted by:
|
patent_small2/2012.dta 文件就是按照上面的循环处理得到的结果：
use patent_small2/2012, clear
list in 1/10
|
准备专利 IPC 数据
按照参考文献的方法，删除只有一个 IPC 号的专利，然后每个专利的分类号只使用前 10 个：
drop if !index(IPC, ";")
keep newipzlid IPC
replace IPC = subinstr(IPC, " ", "", .)
egen ipc = ipc_sub00(IPC), parse(;) choose(10)
drop IPC
save ipcdata, replace
|
保留专利分类号的前 10 个。使用传统的方法是先 split，保留 IPC1-IPC10，再 unite 起来。不过效率很低，这里我使用的是自编的一个 egen 函数：ipc_sub00()。该函数可以用于提取字符串经过指定富豪拆分后的 1-n 个子字符串。源代码如下：
*! 微信公众号 RStata
*! 用于统计字符串使用特定字符拆分后，第 1-n 个子字符串
*! 使用示例：egen IPC_main = ipc_sub00(IPC), parse(;) choose(3)
program define _gipc_sub00
version 6.0
gettoken type 0 : 0
gettoken g 0 : 0
gettoken eqs 0 : 0
syntax varlist(max=1 string), parse(string) choose(integer)
if "`parse'" == "" {
local parse = ";"
}
qui mata: ipc_sub00("`varlist'", "`parse'", `choose')
rename newufvar `g'
end
mata:
void ipc_sub00(string x, string parse, numeric choose) {
v = st_sdata(., x)
ipccountres = J(rows(v), 1, "")
for(z=1;z<=rows(v);z++){
a = ustrsplit(v[z], parse)
if (choose <= cols(a)) {
ipccountres[z] = invtokens(a[1..choose], ";")
}
if (choose > cols(a)) {
ipccountres[z] = v[z]
}
}
stata("cap drop newufvar")
st_addvar("strL", "newufvar")
st_sstore(., "newufvar", ipccountres)
}
end
|
使用前需要把 _gipc_sub00.ado 文件放置到工作目录下。
IPC 与 国民经济行业小类对照
这里我们使用的是「国际专利分类与国民经济行业分类参照关系表（2018）」，该文件下载自国家知识产权局官网。首先我们需要借助外部工具把 pdf 文件转换成 xlsx 文件，我使用的是 Adobe Acrobat。
处理整齐：
import excel using "国际专利分类与国民经济行业分类参照关系表（2018）.xlsx", clear
nrow 1
ren C 国际专利分类号
gen 行业门类代码 = 国民经济行业代码 if ustrregexm(国民经济行业代码, "[A-Z]")
order 行业门类代码
carryforward 行业门类代码, replace
replace 国民经济行业代码 = "0" + 国民经济行业代码 if strlen(国民经济行业代码) == 3 & !ustrregexm(国民经济行业代码, "^0")
keep if strlen(国民经济行业代码) == 4 | mi(国民经济行业代码)
carryforward 国民经济行业代码 国民经济行业名称, replace
drop if mi(国际专利分类号)
foreach i of varlist _all {
gen temp = ustrregexs(0) if ustrregexm(`i', "\n")
replace `i' = subinstr(`i', temp, "", .)
drop temp
gen temp = ustrregexs(0) if ustrregexm(`i', "\r")
replace `i' = subinstr(`i', temp, "", .)
drop temp
cap format `i' %10s
}
replace 国民经济行业代码 = 行业门类代码 + 国民经济行业代码
label data "数据处理：微信公众号 RStata"
save 国民经济分类与IPC分类号对照表, replace
list 国民经济行业代码 国际专利分类号 in 1/10
|
uniq_ipc.dta 数据是历年所有专利涉及的 IPC 分类号。例如 2012 年的可以这么处理得到：
use ipcdata, clear
keep ipc
duplicates drop _all, force
mata: ustrsplit("a;b", ";")
mata:
import delimited using "uniq_ipc.csv", clear
ren v1 uniq_ipc
duplicates drop _all, force
list in 1/10
|
传统的 split + gather 方法也是可以的，不过效率会低非常多。
下面我们需要给每个 IPC 对应上行业小类：
use 国民经济分类与IPC分类号对照表.dta, clear
gen dzid = _n
order dzid
ren 国际专利分类号 IPC
replace IPC = subinstr(IPC, "*", "", .)
replace IPC = subinstr(IPC, " ", "", .)
gen len = strlen(IPC)
tab len
cap mkdir "对照表"
forval i = 3/11 {
use uniq_ipc, clear
cap mkdir "行业匹配结果"
forval i = 3/11 {
use 行业匹配结果/3, clear
forval i = 4/11 {
keep uniq_ipc 国民经济行业代码
duplicates drop _all, force
codebook 国民经济行业代码
save 专利数据与行业小类代码简易对照表, replace
list in 1/10
|
可以看到这个数据的格式其实就类似上次课说的：IPC与产业对照表_示例2.dta
数字经济产业与实体产业
按照论文的介绍，实体产业包括下面四类：
use 国民经济分类与IPC分类号对照表.dta, clear
gen 实体产业 = ""
gen 行业门类 = substr(国民经济行业代码, 1, 1)
order 行业门类
replace 实体产业 = "制造业" if 行业门类 == "C"
replace 实体产业 = "农业" if 行业门类 == "A"
replace 实体产业 = "建筑业及其他工业" if inlist(行业门类, "B", "D", "E")
replace 实体产业 = "服务业" if inlist(行业门类, "I", "O")
keep 国民经济行业代码 实体产业
duplicates drop _all, force
save 实体产业分类, replace
|
这里其实论文里面说的有些含糊，例如农业到底是特指 A01 还是 A，也就是农林牧渔？
数字经济产业分类使用的是「数字经济及其核心产业统计分类（2021）」，该文件也是来源于国家知识产权局。
由于 Stata 不支持 docx 文件的解析，所以这里我们还是使用 R 语言读取处理下：
library(tidyverse)
library(docxtractr)
doc <- read_docx("数字经济及其核心产业统计分类（2021）.docx")
doc %>%
docx_extract_all_tbls() %>%
.[[1]] %>%
slice(-1, -2) %>%
set_names("大类", "中类", "小类", "名称", "说明", "国民经济行业代码及名称") %>%
type_convert() %>%
mutate(国民经济行业代码 = str_extract(国民经济行业代码及名称, "\\d{4}")) %>%
haven::write_dta("数字经济及其核心产业统计分类.dta")
|
然后再用 Stata 处理：
use 数字经济及其核心产业统计分类.dta, clear
drop 国民经济行业代码
keep if substr(小类, 1, 2) != "05" & !mi(小类)
keep 国民经济行业代码及名称 小类
gen v1 = ustrregexs(0) if ustrregexm(国民经济行业代码及名称, "d{4}")
replace 国民经济行业代码及名称 = subinstr(国民经济行业代码及名称, v1, "", .)
gen v2 = ustrregexs(0) if ustrregexm(国民经济行业代码及名称, "d{4}")
replace 国民经济行业代码及名称 = subinstr(国民经济行业代码及名称, v2, "", .)
gen v3 = ustrregexs(0) if ustrregexm(国民经济行业代码及名称, "d{4}")
replace 国民经济行业代码及名称 = subinstr(国民经济行业代码及名称, v3, "", .)
gen v4 = ustrregexs(0) if ustrregexm(国民经济行业代码及名称, "d{4}")
replace 国民经济行业代码及名称 = subinstr(国民经济行业代码及名称, v4, "", .)
gen 数字经济产业 = ""
replace 数字经济产业 = "数字产品制造业" if substr(小类, 1, 2) == "01"
replace 数字经济产业 = "数字产品服务业" if substr(小类, 1, 2) == "02"
replace 数字经济产业 = "数字技术应用业" if substr(小类, 1, 2) == "03"
replace 数字经济产业 = "数字要素驱动业" if substr(小类, 1, 2) == "04"
keep v* 数字经济产业
gen id = _n
gather v*
drop if mi(value)
drop id var
ren value industry2
duplicates drop _all, force
save 数字经济核心产业, replace
list in 1/10
|
然后就可以区分数字与实体了：
use 国民经济分类与IPC分类号对照表.dta, clear
keep 国民经济行业代码
duplicates drop 国民经济行业代码, force
gen industry2 = substr(国民经济行业代码, 2, .)
merge 1:1 industry2 using 数字经济核心产业.dta
tab _m
drop if _m == 2
drop _m
merge m:1 国民经济行业代码 using 实体产业分类.dta
tab _m
drop industry2 _m
keep 国民经济行业代码 数字经济产业 实体产业
gen class = (!mi(数字经济产业))
duplicates drop _all, force
collapse (sum) class, by(国民经济行业代码)
gen 类别 = "数字经济产业" if class > 0
replace 类别 = "实体产业" if class == 0
drop class
tab 类别
save 产业数实分类, replace
list in 1/10
|
这里感觉论文里面说的也比较含糊，所以我采取的方案是如果该小类同时属于数字和实体，那么该小类就归类为数字经济行业。
计算 IPC 融合矩阵
使用上次课编写的 matacode1.do 就可以轻松完成这部操作了：
clear all
use ipcdata, clear
contract ipc
qui do matacode1.do
import delimited using "res1.csv", clear
ren (v1 v2 v3) (ipc1 ipc2 value)
destring, replace
compress
collapse (sum) value, by(ipc1 ipc2)
drop if ipc1 == ipc2
save ipcmat, replace
list in 1/10
|
计算行业共现矩阵
分别使用 ipc 匹配各自的行业即可：
use ipcmat, clear
ren ipc1 uniq_ipc
joinby uniq_ipc using 专利数据与行业小类代码简易对照表.dta, unmatched(master)
drop _m uniq_ipc
ren 国民经济行业代码 industry1
ren ipc2 uniq_ipc
joinby uniq_ipc using 专利数据与行业小类代码简易对照表.dta, unmatched(master)
drop _m uniq_ipc
ren 国民经济行业代码 industry2
drop if mi(value)
collapse (sum) value, by(industry1 industry2)
replace value = 0 if industry1 == industry2
drop if mi(industry1) | mi(industry2)
save industrymat, replace
|
计算数实融合
然后就可以计算数实融合水平了：
use industrymat, clear
ren industry1 国民经济行业代码
recast str5 国民经济行业代码
merge m:1 国民经济行业代码 using 产业数实分类
keep if _m == 3
drop _m
ren 国民经济行业代码 industry1
ren 类别 class1
ren industry2 国民经济行业代码
recast str5 国民经济行业代码
merge m:1 国民经济行业代码 using 产业数实分类
keep if _m == 3
drop _m
ren 国民经济行业代码 industry2
ren 类别 class2
drop if value == 0
bysort class1: egen freq1 = nvals(industry1)
bysort class2: egen freq2 = nvals(industry2)
collapse (sum) value (first) freq1 (first) freq2, by(class1 class2)
gen RH = value / (freq1 * freq2)
gen class = "数实融合" if class1 != class2
replace class = "数数融合" if class1 == "数字经济产业" & class2 == "数字经济产业"
replace class = "实实融合" if class1 == "实体产业" & class2 == "实体产业"
keep class RH
duplicates drop _all, force
list
|
由此我们便计算得到了 2012 年的数实融合水平，不过这个结果比参考文献的大挺多，但是仔细检查了很多遍也没发现计算过程的问题，也有可能是原作者的一些处理细节我我的不一样。
循环各年的就可以计算得到各年的结果了。不过如果想要计算分城市的话，靠在 Stata 里面循环就很费时间了，下次课我们再讲解如何计算分城市和产业的～
点击这里跳转到 RStata 短书平台获取附件：名师讲堂｜使用 Stata 测算数实融合水平（二）： 基于真实专利数据
评论