名师讲堂｜使用 Stata 测算城市数字经济发展的工具变量

今天给大家分享如何使用 Stata 基于上市公司数据和行政区划边界，测算城市数字经济发展的地理工具变量（IV）的方法。工具变量定义参考杨本建、唐金汶《数字经济与区域产业布局》（2026 年第 3 期），其采用「各城市到杭州的球面距离 × 全国与杭州上市公司数字化转型程度之比」作为城市数字经济发展的地理工具变量，以缓解内生性。

基于早期杭州在数字经济发展方面的引领地位，计算各城市到杭州的球面距离。本文将这一距离乘以全国上市公司与位于杭州的上市公司数字化转型程度之比来反映杭州数字经济辐射作用随时间的变化，并以这一变量作为城市数字经济发展的工具变量（地理工具变量）。其中，上市公司数字化转型程度参考吴非等（2021）的方法计算而得。

处理流程概览：

步骤	脚本位置	任务
Step 0	第 0 节	shp2dta 读取 2021行政区划/市.shp，用 gencentroid() 生成各市域质心坐标
Step 1	第 1 节	读入上市公司数字化转型关键词总词频（公司-年度）
Step 2	第 2 节	读入注册/办公地址文件，提取「注册地址所在市」「办公地址所在市」
Step 3	第 3 节	计算年度全国均值、杭州（两口径）均值，并求「全国/杭州」之比
Step 4	第 4 节	geodist 计算各城市到杭州的 WGS84 椭球球面距离
Step 5	第 5 节	交叉匹配构造注册地址 / 办公地址两套城市-年度面板与工具变量
Step 6 + (B)	第 6 节及标签段	保存年度比值与质心文件，并为全部 dta 添加中文变量标签

数据来源说明：测算依赖三类数据——上市公司数字化转型关键词总词频（吴非 2021，公司-年度）、上市公司注册与办公地址（含所在地级市，公司-年度）、2021 年行政区划市域边界（shapefile，用于生成质心）。其中上市公司数据与地址数据来自 RStata 平台分享的「2001～2024 年各上市公司数字转型关键词总词频」与「2000～2024 年上市公司注册地址与办公地址」数据集；行政区划来自「2021 行政区划」数据。

上市公司数字化转型关键词总词频（吴非 2021，公司-年度）数据来自之前分享的：

2001～2024 年各上市公司数字转型关键词总词频及TF-IDF指标计算结果（吴非2021版本，基于管理层讨论分析文本）: https://rstata.duanshu.com/#/brief/course/91b932a04b774f2785e4159ed56a3cd4

上市公司数据与地址数据来自之前分享的：

2000～2024 年上市公司注册地址与办公地址（含经纬度及其所处的省市区县）：https://rstata.duanshu.com/#/brief/course/746d4f595eba41a2a116b5a6edc3ac0a

---

一、指标来源与计算原理

1.1 核心指标：地理工具变量

变量	含义
dist_hz_km	各城市到杭州的球面距离（km，Stata geodist 命令，WGS84 椭球）
ratio_natl_hz	全国 / 杭州 上市公司数字化转型程度之比（随年份变化、对所有城市相同）
iv_product	地理工具变量 = dist_hz_km × ratio_natl_hz

两种“公司所在地”口径：为便于稳健性检验，本结果同时按注册地址所在市（注册地址_市）与办公地址所在市（办公地址_市）各算一套。两者仅在「哪些城市纳入、以及杭州均值取注册地还是办公地在杭州的公司」上不同；全国均值与球面距离口径一致。

1.2 数字化转型程度的测算原理

“上市公司数字化转型程度”由 吴非等（2021） 的方法判定：从上市公司年报的「管理层讨论与分析（MD&A）」文本中，提取数字转型相关关键词并统计其总词频。本讲义取公司-年度「总词频」的年度均值作为该群体的数字化转型程度：

• 全国：年份 t 全部上市公司的总词频均值（两口径相同）；
• 杭州（注册口径）：年份 t 注册地位于杭州市的上市公司总词频均值；
• 杭州（办公口径）：年份 t 办公地位于杭州市的上市公司总词频均值。

注意：原始数字化文件仅含总词频 ≥ 1 的公司-年度（无 0 值记录），且 99.8% 为 A 股；比值在全国与杭州间采用一致口径，时间维度上的变化有效。

1.3 球面距离的测算原理

geodist 命令可以基于 WGS84 椭球直接计算两点间大圆距离，比球面近似更精确，且命令简洁：

*- 安装：ssc install geodist
geodist centroid_lat centroid_lon `hz_lat' `hz_lon', generate(dist_hz_km) in(km)

其中 centroid_lat / centroid_lon 为各城市市域质心经纬度，hz_lat / hz_lon 为杭州市质心经纬度。

1.4 数据来源与口径

数据	角色
2001～2024年各上市公司数字转型关键词总词频（吴非2021，基于管理层讨论与分析部分文本）.dta	公司-年度数字化转型程度（总词频）
2000～2024年上市公司注册地址与办公地址（含经纬度、所处的省市区县及搬迁距离）.dta	公司-年度所在地级市（注册地址_市 / 办公地址_市）
2021行政区划/市.shp	各地级市行政边界，用于求市域质心坐标

---

二、数据准备工作

2.1 Step 0：生成城市质心坐标

关键命令说明：

• shp2dta …, gencentroid(xvar yvar)：读取 shapefile 的同时直接输出每个多边形的质心经纬度（面积加权质心，自动处理含岛屿/飞地的多部件多边形，洞的有向面积为负自动抵消）；前者为经度（x），后者为纬度（y），与 EPSG:4326 对应。
• gencentroid(centroid_lon centroid_lat)：命名质心变量，便于后续 geodist 调用。

*- ==== 0. 生成城市质心坐标（来自 2021 行政区划 市.shp）====
*- 此处代码需下载讲义材料查看～

说明：gencentroid() 内部完成面积加权质心计算；rename 市 city / rename 省 province 将属性名改为英文以便后续合并；生成的 城市质心坐标_2021行政区划.dta 同时被第 4、5 节复用。

2.2 Step 1–2：读取数字化转型文件与地址文件

主程序先将两张公司-年度表读入内存并存为 tempfile，供后续合并。统一把中文变量名改为英文（year / stkcd / digi / city_reg / city_off），规避后续合并与计算中的字符变量名限制。

*- ==== 1. 读取数字化转型文件（公司-年度）====
use "2001～2024年各上市公司数字转型关键词总词频（吴非2021，基于管理层讨论与分析部分文本）.dta", clear
rename 年份 year
rename 股票代码 stkcd
rename 总词频 digi
keep year stkcd digi
tempfile dt
save `dt'

*- ==== 2. 读取地址文件，保留注册/办公所在地级市 ====
use "2000～2024年上市公司注册地址与办公地址（含经纬度、所处的省市区县及搬迁距离）.dta", clear
rename 年份 year
rename 股票代码 stkcd
rename 注册地址_市 city_reg
rename 办公地址_市 city_off
keep year stkcd city_reg city_off
tempfile addr
save `addr'

说明：keep 仅保留后续需要的变量，减小内存占用；tempfile 临时保存，避免反复读盘。两表通过 stkcd + year 唯一对应。

---

三、主程序：年度测算与面板构造

以下代码对应 计算地理工具变量.do 第 3–5 节。

3.1 年度数字化转型程度：全国均值与杭州均值（两口径）

*- ==== 3. 年度“数字化转型程度”：全国均值(全体) 与 杭州均值(两种口径) ====
*- 此处代码需下载讲义材料查看～

说明：collapse (mean) 按年求均值；注册/办公两口径通过 keep if city_? == "杭州市" 分别筛选；merge 1:1 year 把三组年度指标横向合并；drop if missing(...) 自动剔除杭州双口径都无样本的年份（2006 年原始数据无样本），保证后续面板年份一致。

3.2 城市到杭州的距离（geodist）

*- ==== 4. 城市到杭州的距离（geodist, WGS84 椭球, km）====
use "城市质心坐标_2021行政区划.dta", clear
capture drop dist_hz_km             // 用 geodist 重新计算，避免变量重名
sum centroid_lat if city == "杭州市"
local hz_lat = r(mean)
sum centroid_lon if city == "杭州市"
local hz_lon = r(mean)
geodist centroid_lat centroid_lon `hz_lat' `hz_lon', generate(dist_hz_km) in(km)
keep city province centroid_lat centroid_lon dist_hz_km
tempfile cent
save `cent'

说明：先取杭州市质心经纬度存入宏 hz_lat / hz_lon；geodist 以该点为基准，对每一行（城市）计算到杭州的椭球距离（km）并写入 dist_hz_km；capture drop 防止重跑时变量重名报错。

3.3 构造两种口径的城市-年度面板

*- ==== 5. 分别构造两种口径的城市-年度面板并生成工具变量 ====
*- 此处代码需下载讲义材料查看～

说明：cross usingyearly`` 把「城市 × 距离」与「年份 × 比值」做笛卡尔积，得到城市-年度面板；gen iv_product = dist_hz_km * ratio_natl_hz_* 即地理工具变量；两口径仅比值变量不同（注册用 ratio_natl_hz_reg，办公用 ratio_natl_hz_off），其余完全一致。随后第 6 节保存年度比值与质心文件。

---

四、变量标签与结果可视化

4.1 添加中文变量标签

合并脚本在生成所有 dta 后，逐一 use 并重开加标签再 save ... replace，无需额外脚本。

*- 注册地址面板
use "$out/城市数字经济发展地理工具变量_注册地址_2001-2024.dta", clear
label variable province      "省份"
label variable city          "城市(注册地址所在市)"
label variable year          "年份"
label variable dist_hz_km    "各城市到杭州的球面距离(km)"
label variable ratio_natl_hz "全国/杭州(注册地址)上市公司数字化转型程度之比"
label variable iv_product    "地理工具变量(距离×比值)"
label data "城市数字经济发展地理工具变量_注册地址口径(2001-2024,2006年缺)"
save "$out/城市数字经济发展地理工具变量_注册地址_2001-2024.dta", replace

*- 办公地址面板（city 标签改为“城市(办公地址所在市)”）、年度比值、质心文件同理……
di as result "全部 dta 已生成并添加中文变量标签。"

说明：label variable 为中文含义，label data 为数据集标签；变量名仍为英文（province/city/year/dist_hz_km/ratio_natl_hz/iv_product），便于跨软件使用，describe 即可看到中文标签。

4.2 结果可视化

*- —— 图1：各年“全国/杭州”数字化转型程度之比（两口径）趋势 ——
*- 此处代码需下载讲义材料查看～

*- —— 图2：2024 年代表性城市地理工具变量（IV）柱状图 ——
use "城市数字经济发展地理工具变量_结果/城市数字经济发展地理工具变量_注册地址_2001-2024.dta", clear
keep if year == 2024
keep if inlist(city, "杭州市","上海市","深圳市","北京市","广州市") | inlist(city, "成都市","武汉市","西安市","重庆市","乌鲁木齐市")
sort iv_product
graph bar iv_product, over(city, sort(iv_product)) ///
    ti("2024 年代表性城市地理工具变量(IV = 距离×比值)") ///
    ytitle("IV") blabel(bar, format(%9.0f))
graph export "输出/2024年代表性城市IV.png", replace width(4800) height(3000)

说明：图 1 展示两口径比值均在 0.6–0.8 区间波动（杭州长期领先，故比值 < 1，随全国追赶比值上行）；图 2 用柱状图对比代表性城市的 IV，其排序与到杭州距离排序一致（同一年份比值对所有城市相同）。下文图 3 进一步以城市分布图直观呈现 IV 的空间梯度。

4.3 结果展示

下图由上述脚本生成：

点击这里跳转到 RStata 短书平台获取附件：名师讲堂｜使用 Stata 测算城市数字经济发展的工具变量