附件中的文献《生产性服务业集聚何以赋能技术扩散》中提到了城市的生产性服务业专业化集聚指标和多样化集聚指标，今天的课程中我们将会讲解如何使用 Stata 语言根据城市统计年鉴数据来测算这两个指标。
本版本特点：不使用 Mata 矩阵编程语言，全部采用纯 Stata 命令和 forvalues 循环实现，更易于理解和调试。
城市统计年鉴数据在这里：
1999～2025 年中国城市统计年鉴面板数据整理结果: https://rstata.duanshu.com/#/brief/course/c6b0aaf3bcba494faadbacacb012dcef
这里选择的是 2017～2019 年的样本作为演示使用。也就是附件中的 城市统计年鉴样本数据.dta。
一、指标来源与计算公式
1.1 专业化集聚指数 SP（公式4）
1.2 多样化集聚指数 DV（公式5）
二、数据与变量说明
2.1 数据来源
城市统计年鉴数据在这里：
1999～2025 年中国城市统计年鉴面板数据整理结果: https://rstata.duanshu.com/#/brief/course/c6b0aaf3bcba494faadbacacb012dcef
这里选择的是 2017～2019 年的样本作为演示使用。也就是附件中的 城市统计年鉴样本数据.dta。
由于计算所需的变量主要集中在 2003～2019 年，所以完整数据仅可计算得到 2003～2019 年的。
2.2 生产性服务业行业界定
参考 Ke et al.（2014）和韩峰和阳立高（2020），生产性服务业包括 7个细分行业：
| 代码 |
行业名称 |
说明 |
| s1 |
交通运输、仓储和邮政业 |
全市口径 |
| s2 |
金融业 |
全市口径 |
| s3 |
科学研究和技术服务业 |
2003-2016含地质勘探，2017-2019不含 |
| s4 |
租赁和商务服务业 |
全市口径 |
| s5 |
信息传输、计算机服务和软件业 |
全市口径 |
| s6 |
批发和零售业 |
全市口径 |
| s7 |
水利、环境和公共设施管理业 |
全市口径 |
注意：科学研究和技术服务业在2003-2016年统计口径含”地质勘探业”，2017-2019年不含。处理时用 coalesce() 合并两套口径。
2.3 就业数据说明
- 数据来自《中国城市统计年鉴》城镇单位从业人员（全市口径）
- 单位：人
- 总就业 = 城镇单位从业人员期末人数（全市）
三、Stata 代码实现
3.1 读取数据与变量重命名
cd "/Users/ac/Desktop/使用 Stata 测算各城市生产性服务业集聚水平（无 Mata 代码版本）"
use "城市统计年鉴样本数据.dta", clear
rename 年份 year
rename 省 prov
rename 省代码 prov_code
rename 市 city
rename 市代码 cityid
rename 城镇单位从业人员期末人数_人_全市 emp_total
rename 第三产业_交通运输仓储和邮政业_人_全市 emp_s1
rename 第三产业_金融业_人_全市 emp_s2
rename 第三产业_科学研究和技术服务和地质勘探业_人_全市 emp_s3a
rename 第三产业_科学研究和技术服务业_人_全市 emp_s3b
rename 第三产业_租赁和商务服务业_人_全市 emp_s4
rename 第三产业_信息传输计算机服务和软件业_人_全市 emp_s5
rename 第三产业_批发和零售业_人_全市 emp_s6
rename 第三产业_水利环境和公共设施管理业_人_全市 emp_s7
gen emp_s3 = emp_s3a
replace emp_s3 = emp_s3b if missing(emp_s3) & !missing(emp_s3b)
drop emp_s3a emp_s3b
|
3.2 数据清洗与筛选
*- 此处代码需下载讲义材料查看～
disp "样本量（城市×年份）: " _N
summ year, meanonly
disp "年份范围: " r(min) " ~ " r(max)
|
3.3 异常值处理
summ emp_total if city == "嘉峪关市" & year == 2013, meanonly
local jyg_raw = r(mean)
summ emp_total if city == "嘉峪关市" & inlist(year, 2012, 2014), meanonly
local jyg_fix = r(mean)
disp "嘉峪关市 2013 年 emp_total 修复："
disp " 原始值: " `jyg_raw'
disp " 替换值（2012/2014均值）: " `jyg_fix'
replace emp_total = `jyg_fix' if city == "嘉峪关市" & year == 2013
disp "(鹤壁市2013年不在2017-2019样本，跳过剔除）"
local n_before = _N
drop if city == "鹤壁市" & year == 2013
local n_after = _N
disp " 剔除前: `n_before' 条 → 剔除后: `n_after' 条"
disp "异常值处理后样本量: " _N
summ year, meanonly
disp "年份范围: " r(min) " ~ " r(max)
save "cleaned_data.dta", replace
disp "已保存：cleaned_data.dta"
|
3.4 计算全国层面基准（按年份）
use "cleaned_data.dta", clear
save "national_benchmark.dta", replace
disp "已保存：national_benchmark.dta"
summ year, meanonly
disp "全国层面基准计算完成，年份数: " _N
|
3.5 测算专业化集聚指数 SP（公式4）
*- ------------------------------------------------------------
*- 四、测算专业化集聚指数 SP（公式4）
*-
*- SP_i = ∑_s(E_is/E_i) / ∑_s(E_s/E)
*- = [城市生产性服务业7行业就业 / 城市总就业] /
*- [全国生产性服务业7行业就业 / 全国总就业]
*-
*- 经济含义：
*- SP > 1：城市生产性服务业就业比重 > 全国平均水平（相对专业化集聚）
*- SP < 1：城市生产性服务业就业比重 < 全国平均水平
*- 数据处理：微信公众号 RStata
*- ------------------------------------------------------------
*- 此处代码需下载讲义材料查看～
|
3.6 测算多样化集聚指数 DV（公式5）- 纯 Stata 实现
*- ------------------------------------------------------------
*- 五、测算多样化集聚指数 DV（公式5）- 纯 Stata 实现
*-
*- 论文公式(5)的精确结构：
*-
*- DV_i = ∑_s (E_{i,s}/E_i) × [ 1/∑_{s'≠s}^n [ E_{i,s'} / (E_i - E_{i,s}) ]^2 ]
*- ÷ [ 1/∑_{s'≠s}^n [ E_{s'} / (E - E_s) ]^2 ]
*-
*- 关键：公式中的份额是"其余行业在剩余就业中的份额"，而非"占总就业份额"
*- - 城市层面：E_{i,s'} / (E_i - E_{i,s}) = 其余行业就业 / (城市总就业 - 行业s就业)
*- - 全国层面：E_{s'} / (E - E_s) = 其余行业就业 / (全国总就业 - 行业s就业)
*-
*- 本版本使用纯 Stata 命令实现，不使用 Mata
*- 方法：使用 forvalues 循环遍历7个行业，为每个行业计算 DV_s，然后加总
*- 数据处理：微信公众号 RStata
*- ------------------------------------------------------------
*- 此处代码需下载讲义材料查看～
*- 检查 DV 计算结果
summ DV, detail
*- 保存 DV 结果
keep year prov prov_code city cityid DV
save "dv_results.dta", replace
disp "已保存：dv_results.dta"
|
纯 Stata 循环实现详解
本版本采用纯 Stata 命令实现 DV 计算，核心思路是使用 forvalues 循环遍历 7。个行业，为每个行业计算其对 DV 的贡献值 DVs，然后加总得到最终 DV。
实现步骤：
- 外层循环（
forvalues s = 1/7）：遍历7个生产性服务业行业
- 内层循环（
foreach s2 of local others）：遍历除行业s外的其他6个行业
- 计算 DV_s：
- 清理临时变量：每个行业循环结束后，删除该行业相关的临时变量，避免变量名冲突
与 Mata 版本的区别：
| 特性 |
Mata 版本 |
纯 Stata 版本 |
| 编程语言 |
Mata（Stata 矩阵语言） |
纯 Stata 命令 |
| 循环方式 |
Mata for 循环 |
Stata forvalues 循环 |
| 数据访问 |
st_view() 读取到 Mata 矩阵 |
直接操作 Stata 变量 |
| 代码复杂度 |
较高，需要理解 Mata 语法 |
较低，纯 Stata 语法 |
| 调试难度 |
较高，Mata 错误难以调试 |
较低，可使用 disp 调试 |
| 执行速度 |
较快（矩阵运算） |
较慢（逐观测循环） |
推荐场景：
- 学习阶段：建议使用纯 Stata 版本，易于理解和调试
- 大数据量：建议使用 Mata 版本，执行速度更快
3.7 缩尾处理与结果保存
*- ------------------------------------------------------------
*- 六、合并 SP 和 DV，并进行缩尾处理（Winsorize）
*- 数据处理：微信公众号 RStata
*- ------------------------------------------------------------
*- 此处代码需下载讲义材料查看～
*- ------------------------------------------------------------
*- 七、保存为 dta 文件
*- 数据处理：微信公众号 RStata
*- ------------------------------------------------------------
*- 重命名变量以符合输出格式
rename year 年份
rename prov 省
rename prov_code 省代码
rename city 市
rename cityid 市代码
rename SP 生产性服务业专业化集聚SP
*- 保存最终结果
label data "数据处理：微信公众号 RStata"
save "2017-2019年各城市生产性服务业集聚水平.dta", replace
disp "已保存：2017-2019年各城市生产性服务业集聚水平.dta"
|
点击这里跳转到 RStata 短书平台获取附件：名师讲堂｜使用 Stata 测算各城市生产性服务业集聚水平（无 Mata 代码版本）
评论