由于借助 AI 工具学习编程已经变得非常容易了,因此之后的课程就不再默认进行视频讲解了,如果特别需要视频讲解也可以联系李老师预约讲解~讲义材料学习过程中遇到的问题也可以及时与李老师联系。
购买 RStata 名师讲堂会员即可参加该课程啦(之前的和未来的都可以参加)!
价格:2800/年 或者 4800/长期
购买会员可以从这里下单:https://rstata.duanshu.com/#/card/list/
名师讲堂会员权益:
- 参加平台上的其他 R 语言和 Stata 的课程;
- 以会员折扣价购买我们分享的数据资料(10 元/份);
* 如果发票可添加小编微信 r_stata2 (RStata 李老师)开具。如需数据资料,购买后可添加小编微信免费领取数据折扣卡。
更多关于 RStata 会员的更多信息可添加微信号 r_stata2 咨询:

课程主页(点击文末的阅读原文即可跳转):https://rstata.duanshu.com/#/brief/course/5315098919dc49a6a0134fc8ffe0beb5
今天给大家分享使用 Python 计算专利颠覆性创新指数及筛选颠覆性专利的方法。
不过需要注意,直接对每个专利计算和之前年份(即使是 3 年内的)所有专利的相似度是非常非常耗时的,建议还是限定在某个领域内,例如只计算该专利与自己所在的部、大类、小类的。
使用 reticulate 创建与管理 Python 虚拟环境
在 R 中通过 reticulate 包来调用 Python,最好的实践是为项目创建一个专属的 Python 虚拟环境,将所需依赖隔离到独立空间,避免与系统 Python(如 Anaconda)发生版本冲突。
重要说明(避免"已初始化"报错):reticulate 在 R 会话中只能绑定一次 Python——一旦某个 {python} 代码块运行,Python 解释器就被锁定,之后再调用 use_virtualenv() 会报错:
ERROR: The requested version of Python cannot be used, as another version has already been initialized.因此,虚拟环境的激活必须在所有 {python} 代码块之前完成。本文档的解决方案是在 setup chunk 中通过 Sys.setenv(RETICULATE_PYTHON = ...) 提前锁定 Python 路径,这是 reticulate 选取 Python 的最高优先级入口。
安装 reticulate(仅首次)
# 设置 CRAN 镜像(knit 时 R 处于非交互模式,不会自动选择镜像)options(repos = c(CRAN = "https://mirrors.tuna.tsinghua.edu.cn/CRAN/"))# 仅在尚未安装时才安装,避免每次 knit 都重装if (!requireNamespace("reticulate", quietly = TRUE)) { install.packages("reticulate") message("reticulate 安装完成!") message("reticulate 已安装,版本:", packageVersion("reticulate"))虚拟环境初始化原理(已在 setup chunk 中完成)
本文档的 setup chunk(隐藏运行)包含如下逻辑:
.venv_python <- virtualenv_python(.venv_name)if(!file.exists(.venv_python)){ virtualenv_create(.venv_name) .venv_python <- virtualenv_python(.venv_name)# 通过环境变量抢先锁定 Python(优先级最高,早于任何 {python} chunk)Sys.setenv(RETICULATE_PYTHON = .venv_python)use_virtualenv(.venv_name, required =TRUE)这样做的关键在于:knitr 在处理第一个 {python} chunk 时,reticulate 已经通过 RETICULATE_PYTHON 环境变量知道要使用 .venv,不会再去碰 Anaconda。
在虚拟环境中安装 Python 包(仅首次)
py_pkgs <- c("numpy", "pandas", "scipy", "jieba")installed <- py_list_packages(".venv")$packageneed_install <- setdiff(py_pkgs, installed)if (length(need_install) > 0) { virtualenv_install(".venv", packages = need_install) message("已安装缺失的包:", paste(need_install, collapse = ", ")) message("所有 Python 包已就绪,无需安装")验证激活状态
# 验证当前绑定的 Python 路径(应指向 .venv 目录)
1. 指标来源与计算方法
参照冉征、刘修岩、陈露(2025)发表于《中国工业经济》2025 年第 10 期的论文《技术集群结构与颠覆式创新——兼论关联性"陷阱"的突破路径》,本文基于 Kelly et al.(2021)设计的文本相似度法,通过专利文本的前后对比测度每项专利的颠覆性指标。
1.1 核心思想
颠覆式创新聚焦突破原有技术模式,形成全新的生产和价值创造模式,即一方面突破已有技术的限制,另一方面对后续的创新活动产生深远影响(Kerr, 2010)。基于这一内涵,一项专利如果:
则说明该项专利突破了已有技术的框架,形成了新的技术模式,具有更强的颠覆性。
1.2 计算步骤
1.3 两种计算方案
1.4 颠覆性专利筛选
考虑到 radical 指数是连续指标且可以跨时期对比,按照论文方法设置前 5% 作为颠覆式创新的识别区间,即全样本中颠覆性指数位于前 5% 的专利定义为颠覆式创新专利。
2. 数据读取与预处理
2.1 导入 Python 包
from collections import Counter2.2 读取专利数据
首先读取专利数据样本文件夹中的所有 CSV 文件。专利数据包含如下变量:newipzlid(专利 ID,前 4 位为年份)、标题、摘要、公开公告号、申请号。
csv_files = sorted([f for f in os.listdir(data_dir) if f.endswith(".csv")])print("\n===== 读取专利数据 =====") [pd.read_csv(os.path.join(data_dir, f)) for f in csv_files],print(f"原始数据总行数: {len(df_raw)}")2.3 提取年份
从 newipzlid 的前 4 位提取专利申请年份:
df["年份"] = df["newipzlid"].astype(str).str[:4].astype(int)print(f"年份范围: {df['年份'].min()} - {df['年份'].max()}")2.4 专利去重
专利数据中可能存在重复专利,需要按年份和申请号进行去重。先去重公开公告号(去掉末尾的类型标识字母),再去重申请号:
print("\n===== 专利去重 =====")print(f"去重前: {len(df)} 行")# 按 年份 + 公开公告号 去重(去掉末尾字母)df["公开公告号_clean"] = df["公开公告号"].str.replace(r"[A-Z]$", "", regex=True)df = df.drop_duplicates(subset=["年份", "公开公告号_clean"], keep="first")df = df.drop(columns=["公开公告号_clean"])# 按 年份 + 申请号 去重(保留第一次出现的记录)df = df.drop_duplicates(subset=["年份", "申请号"], keep="first")print(f"去重后: {len(df)} 行")关键函数对照:
| | |
|---|
str_replace(x, "[A-Z]$", "") | str.replace(r"[A-Z]$", "", regex=True) | |
distinct(年份, 列, .keep_all=TRUE) | drop_duplicates(subset=[...], keep="first") | |
select(-列) | drop(columns=[...]) | |
2.5 文本清洗
合并专利的标题和摘要作为完整文本,去除标点、数字、英文及多余空格:
2.6 中文分词
使用 jieba 进行中文分词,加载停用词表和用户词典,过滤长度小于 2 的词和非中文字符:
分词逻辑对照:
| | |
|---|
worker(stop_word="stopwords.txt", user="dictionary.txt") | jieba.load_userdict() | |
segment(txt, seg) | jieba.cut(txt) | |
nchar(tokens) >= 2 | len(t) >= 2 | |
!str_detect(tokens, "[^\\u4e00-\\u9fff]") | re.match(r"^[\u4e00-\u9fff]+$", t) | |
3. 辅助函数:TF-IDF 矩阵与颠覆性指数
3.1 构建时变 TF-IDF 矩阵
TF-IDF 计算对照:
| | |
|---|
bind_tf_idf(word, newipzlid, n) | | |
cast_sparse(newipzlid, word, tf_idf) | sparse.csr_matrix(...) | |
dtm / row_norms | sparse.diags(1/norms) @ dtm | |
3.2 计算颠覆性指数
对每一年,确定前后窗口内的专利索引,计算后向相似度(BPS)和前向相似度(FPS),并同时输出两种方案的 radical 指数。其中求和法是论文中使用的:
相似度计算对照:
| | |
|---|
dtm[target_idx, ] %*% t(dtm[backward_idx, ]) | dtm[target_idx, :] @ dtm[backward_idx, :].T | |
rowSums(sim_mat) | np.asarray(sim.sum(axis=1)).flatten() | |
ifelse(BPS > 0, FPS/BPS, NA) | np.where(BPS > 0, FPS/BPS, np.nan) | |
4. 逐年计算颠覆性指数
target_years = list(range(2003, 2008))print("\n===== 计算颠覆性创新指数 =====")print(f"目标年份: {', '.join(map(str, target_years))}")print(f"计算窗口: tau = 3 年")print(f"\n--- 正在计算 {yr} 年 ---") result = compute_radical(df, yr, tau=3)print(f" 专利数: {len(result)}")print(f" 有效 radical_sum 指数: {result['radical_sum'].notna().sum()}")print(f" 有效 radical_avg 指数: {result['radical_avg'].notna().sum()}") all_results.append(result)#> 有效 radical_sum 指数: 980#> 有效 radical_avg 指数: 980#> 有效 radical_sum 指数: 974#> 有效 radical_avg 指数: 974#> 有效 radical_sum 指数: 972#> 有效 radical_avg 指数: 972#> 有效 radical_sum 指数: 979#> 有效 radical_avg 指数: 979#> 有效 radical_sum 指数: 984#> 有效 radical_avg 指数: 984
5. 汇总结果与筛选颠覆性专利
5.1 标记颠覆性专利
在全样本期内,将颠覆性指数位于前 5% 的专利定义为颠覆式创新,分别对两种方案计算阈值并标记:
5.2 合并原始专利信息
# 注意:统一 newipzlid 类型为字符型,避免 merge 时类型不匹配df["newipzlid"] = df["newipzlid"].astype(str)df_radical["newipzlid"] = df_radical["newipzlid"].astype(str)df_final = df[["newipzlid", "年份", "标题", "摘要", "公开公告号", "申请号"]].merge( df_radical, on=["newipzlid", "年份"], how="left"
6. 结果展示
6.1 总体描述
print("\n===== 结果摘要 =====")print(f"目标年份: {min(target_years)} - {max(target_years)}")print(f"总专利数(去重后,含窗口期): {len(df)}")print(f"计算了 radical 指数的专利数: {df_final['radical_avg'].notna().sum()}")6.2 各年颠覆性专利分布(方案一:求和法)
6.3 各年颠覆性专利分布(方案二:均值法)
print("\n--- 各年颠覆性专利分布(方案二:均值法)---")summary_avg = df_final[df_final["年份"].isin(target_years)].groupby("年份").agg( 总专利数=("newipzlid", "count"), 有效radical=("radical_avg", lambda x: x.notna().sum()), 颠覆性专利数=("颠覆性_avg", "sum"), radical均值=("radical_avg", "mean"), radical中位数=("radical_avg", "median"),summary_avg["颠覆性占比"] = (summary_avg["颠覆性专利数"] / summary_avg["总专利数"] * 100).round(2)summary_avg["radical均值"] = summary_avg["radical均值"].round(4)summary_avg["radical中位数"] = summary_avg["radical中位数"].round(4)print(summary_avg.to_string(index=False))6.4 两种方案一致性比较
print("\n--- 两种方案一致性比较 ---")df_target = df_final[df_final["年份"].isin(target_years)]both = ((df_target["颠覆性_sum"] == 1) & (df_target["颠覆性_avg"] == 1)).sum()only_sum = ((df_target["颠覆性_sum"] == 1) & (df_target["颠覆性_avg"] == 0)).sum()only_avg = ((df_target["颠覆性_sum"] == 0) & (df_target["颠覆性_avg"] == 1)).sum()print(f" 两种方案均标记为颠覆性: {both}")print(f" 仅求和法标记: {only_sum}")print(f" 仅均值法标记: {only_avg}")print(f" 一致率: {both / (both + only_sum + only_avg) * 100:.1f}%")
7. 保存结果
7.1 保存为 CSV 文件
"newipzlid", "年份", "标题", "摘要", "公开公告号", "申请号","BPS_sum", "FPS_sum", "radical_sum", "颠覆性_sum","BPS_avg", "FPS_avg", "radical_avg", "颠覆性_avg",csv_file = "patent_disruptive_index.csv"df_output.to_csv(csv_file, index=False, encoding="utf-8-sig")print(f"结果已保存到: {csv_file}")7.2 保存为 Stata .dta 文件
由于 Stata 变量名不支持中文字符,需将中文列名映射为英文变量名,并设置变量标签。最终通过 R 的 haven 包写入 .dta 文件,以保证与 R 版本的完全兼容:
# Stata 变量名不支持中文,需映射为英文变量名"newipzlid": "newipzlid","radical_sum": "radical_s","radical_avg": "radical_a",df_stata = df_output.rename(columns=rename_map).copy()for col in ["bps_sum", "fps_sum", "radical_s", "bps_avg", "fps_avg", "radical_a"]: df_stata[col] = pd.to_numeric(df_stata[col], errors="coerce")for col in ["disrup_s", "disrup_a"]: df_stata[col] = df_stata[col].fillna(0).astype(int)df_stata.to_pickle("df_stata_tmp.pkl")print("数据准备完成,即将通过 R 写入 .dta 文件") radical_s = "颠覆性指数(求和法)", radical_a = "颠覆性指数(均值法)",for (v in names(variable_labels)) { attr(df_stata[[v]], "label") <- variable_labels[[v]]write_dta(df_stata, "patent_disruptive_index.dta", label = "数据处理:微信公众号 RStata")cat("Stata 文件已保存至: patent_disruptive_index.dta\n")df_target <- df_stata[df_stata$year %in% 2003:2007, ]write_dta(df_target, "patent_disruptive_index_2003_2007.dta", label = "数据处理:微信公众号 RStata")cat("2003-2007 年子集已保存至: patent_disruptive_index_2003_2007.dta\n")file.remove("df_stata_tmp.pkl")如何参加课程?
购买 RStata 名师讲堂会员即可参加该课程啦(之前的和未来的都可以参加)!
价格:2800/年 或者 4800/长期
购买会员可以从这里下单:https://rstata.duanshu.com/#/card/list/
名师讲堂会员权益:
- 参加平台上的其他 R 语言和 Stata 的课程;
- 以会员折扣价购买我们分享的数据资料(10 元/份);
* 如果发票可添加小编微信 r_stata2 (RStata 李老师)开具。如需数据资料,购买后可添加小编微信免费领取数据折扣卡。
更多关于 RStata 会员的更多信息可添加微信号 r_stata2 咨询:

课程主页(点击文末的阅读原文即可跳转):https://rstata.duanshu.com/#/brief/course/5315098919dc49a6a0134fc8ffe0beb5






