该课程已经讲完，感兴趣的小伙伴可以点击文末的阅读原文跳转到 RStata 平台观看视频讲解～
「财政研究」2019 年第 4 期中有篇论文<经济分权、地方政府竞争与城市全要素生产率> 提到了衡量政府竞争程度的指标。以往的文献通常采用经济增长率、外商直接投资和宏观税负，背后的逻辑是指标数值越大，代表地方政府竞争越激烈。但是三种方式都无法揭示本地政府受相邻地区地方政府决策的影响，因此该论文提出可以用地区经济增长率与空间权重矩阵的乘积 Wgdpg 衡量地方政府竞争。今天我们就一起来看下如何在 Stata 中通过矩阵运算的方法来计算这一指标。
地区经济增长率通常使用地区生产总值的增长率衡量，这一指标在之前分享的城市统计年鉴面板数据里面有，这里提取了出来：
1999～2022 年中国城市统计年鉴面板数据：https://rstata.duanshu.com/#/brief/course/1dcad8873f394c31a4f082bd74e916b2
1998～2021年各地级市地区生产总值面板数据.dta
另外因为 2018 年的城市统计年鉴（上面的数据是 2017 年的）没有提供全市的地区生产总值，所以我又从 EPS 获取了 2017 年的数据进行补充。
首先读取“2019年各市是否接壤矩阵（对角线为0）.xlsx”文件并处理，该文件来源于平台上之前分享的「历年各省市区县邻接矩阵、质心距离矩阵、质心距离倒数矩阵数据」。
历年各省市区县邻接矩阵、质心距离矩阵、质心距离倒数矩阵数据：https://rstata.duanshu.com/#/brief/course/eb331d940dbf4907b75d75b4914da67d
import excel using "2019年各市是否接壤矩阵（对角线为0）.xlsx", clear
foreach i of varlist _all {
cap replace `i' = subinstr(`i', "-", "_", .) if _n == 1
}
nrow 1
destring, replace
replace 城市 = subinstr(城市, "-", "_", .)
save "2019年各市是否接壤矩阵（对角线为0）", replace
因为 Stata 中的变量名不支持带“-”号的，所以这里都替换成了“_”。
因为城市统计年鉴里面没有提供 2017 年的全市地区生产总值，所以再从 EPS 获取一份 2017 年的处理下：
*- 2017 年的数据
use 2017, clear
ren 城市名称 市
ren 城市代码 市代码
ren 国内生产总值万元 地区生产总值_当年价格_亿元_全市
replace 地区生产总值_当年价格_亿元_全市 = 地区生产总值_当年价格_亿元_全市 / 10000
drop if missing(地区生产总值_当年价格_亿元_全市)
drop 市代码
replace 市 = "北京市" if 市 == "北京"
save 2017b, replace
这份数据里面有市代码变量，但是担心可能不是和我们这份城市统计年鉴数据一个版本的，所以匹配下试试：
shp2dta using "2019行政区划/市.shp", database(citycode2019) coord(coord) replace
cap erase coord.dta
use citycode2019, clear
drop *类型 _ID
merge 1:1 市 using 2017b
*- 有三个市没有匹配成功就不要了
keep if _m == 3
drop _m
save "2017c", replace
然后就可以把 2017 年的数据合并到总数据里面，再补充下缺失值：
*- 补充缺失的增长率
use "1998～2021年各地级市地区生产总值面板数据", clear
drop if 年份 == 2017
append using 2017c
gen 城市 = 省 + "_" + 市
encode 城市, gen(city)
xtset city 年份
replace 地区生产总值增长率_百分比_全市 = ((地区生产总值_当年价格_亿元_全市 - L.地区生产总值_当年价格_亿元_全市) / 地区生产总值_当年价格_亿元_全市) * 100 if missing(地区生产总值增长率_百分比_全市)
drop city
save "1998～2021年各地级市地区生产总值面板数据(填补后)", replace
我们再以 2021 年的数据为例看一下如何计算这一指标。
保留 2021 年的数据和矩阵匹配：
use "1998～2021年各地级市地区生产总值面板数据(填补后)", clear
keep if 年份 == 2021
drop if missing(地区生产总值增长率_百分比_全市)
merge 1:1 城市 using "2019年各市是否接壤矩阵（对角线为0）.dta"
keep if _m == 3
drop _m
为了避免矩阵的顺序错乱，我们首先锁定数据的行顺序：
*- 首先锁定行的顺序
gen id = _n
这样如果后面有操作导致行序变了，只需要 gsort id 就可以修正了（不过这个在下面没用到）。
然后根据城市顺序重新排列下列序：
*- 循环获取城市列表
local vlist = ""
forval i = 1/`=_N' {
local vlist `vlist' `=城市[`i']'
}
keep 年份 省代码 省 市代码 市 地区生产总值_当年价格_亿元_全市 地区生产总值增长率_百分比_全市 城市 `vlist'
order 年份 省代码 省 市代码 市 地区生产总值_当年价格_亿元_全市 地区生产总值增长率_百分比_全市 城市 `vlist'
下面就可以使用 mkmat 命令把变量值转换成矩阵了：
mkmat 地区生产总值增长率_百分比_全市, mat(m1)
mkmat `vlist', mat(m2)
注意这里也用到了 vlist，试过 local 变量。local 变量的定义和使用要一并运行，所以这里的代码要和上面的一段代码同时运行才可以。
这样的话，m1 就是个 nxn 的矩阵，m2 是 nx1 的矩阵，两者相乘就得到了 nx1 的矩阵：
mat Wgdpg = m2 * m1
keep 省代码 省 市代码 市
svmat Wgdpg
list in 1/10
*> +------------------------------------------------------+
*> | 省代码 省 市代码 市 Wgdpg1 |
*> |------------------------------------------------------|
*> 1. | 310000 上海市 310000 上海市 26.2 |
*> 2. | 530000 云南省 530900 临沧市 12.8 |
*> 3. | 530000 云南省 530700 丽江市 8.3 |
*> 4. | 530000 云南省 530500 保山市 7.4 |
*> 5. | 530000 云南省 530100 昆明市 31.16 |
*> |------------------------------------------------------|
*> 6. | 530000 云南省 530600 昭通市 39.87 |
*> 7. | 530000 云南省 530800 普洱市 16.46 |
*> 8. | 530000 云南省 530300 曲靖市 38.6 |
*> 9. | 530000 云南省 530400 玉溪市 10.4 |
*> 10. | 150000 内蒙古自治区 150900 乌兰察布市 38.55 |
*> +------------------------------------------------------+
svmat 可以把矩阵再转换成变量。
然后就可以循环所有年份的了：
*- 循环所有年份
cap mkdir "结果数据"
use "1998～2021年各地级市地区生产总值面板数据(填补后)", clear
forval i = 1998/2021 {
preserve
qui {
keep if 年份 == `i'
drop if missing(地区生产总值增长率_百分比_全市)
merge 1:1 城市 using "2019年各市是否接壤矩阵（对角线为0）.dta"
keep if _m == 3
drop _m
local vlist = ""
forval j = 1/`=_N' {
local vlist `vlist' `=城市[`j']'
}
keep 年份 省代码 省 市代码 市 地区生产总值_当年价格_亿元_全市 地区生产总值增长率_百分比_全市 城市 `vlist'
order 年份 省代码 省 市代码 市 地区生产总值_当年价格_亿元_全市 地区生产总值增长率_百分比_全市 城市 `vlist'
*- 构建矩阵
mkmat 地区生产总值增长率_百分比_全市, mat(m1)
mkmat `vlist', mat(m2)
mat Wgdpg = m2 * m1
keep 省代码 省 市代码 市
svmat Wgdpg
}
save "结果数据/`i'", replace
restore
}
最后合并所有年份的就得到了总的计算结果：
use 结果数据/1998.dta, clear
gen 年份 = 1998
forval i = 1999/2021 {
append using 结果数据/`i'.dta
replace 年份 = `i' if missing(年份)
}
order 年份
ren Wgdpg1 Wgdpg
list in 1/10
*> +------------------------------------------------------------+
*> | 年份 省代码 省 市代码 市 Wgdpg |
*> |------------------------------------------------------------|
*> 1. | 1998 310000 上海市 310000 上海市 34 |
*> 2. | 1998 530000 云南省 530100 昆明市 7.2 |
*> 3. | 1998 530000 云南省 530300 曲靖市 19.1 |
*> 4. | 1998 150000 内蒙古自治区 150300 乌海市 8.9 |
*> 5. | 1998 150000 内蒙古自治区 150200 包头市 13.1 |
*> |------------------------------------------------------------|
*> 6. | 1998 150000 内蒙古自治区 150100 呼和浩特市 21 |
*> 7. | 1998 150000 内蒙古自治区 150400 赤峰市 13.1 |
*> 8. | 1998 110000 北京市 110000 北京市 56.5 |
*> 9. | 1998 220000 吉林省 220200 吉林市 68.2 |
*> 10. | 1998 220000 吉林省 220300 四平市 67.4 |
*> +------------------------------------------------------------+
save "Wgdpg计算结果", replace
另外如果是区县的可能使用 Mata 的矩阵运算会更快速点，首先以 2021 年的为例：
clear all
use "1998～2021年各地级市地区生产总值面板数据(填补后)", clear
keep if 年份 == 2021
drop if missing(地区生产总值增长率_百分比_全市)
merge 1:1 城市 using "2019年各市是否接壤矩阵（对角线为0）.dta"
keep if _m == 3
drop _m
*- 循环获取城市列表
local vlist = ""
forval i = 1/`=_N' {
local vlist `vlist' `=城市[`i']'
}
keep 年份 省代码 省 市代码 市 地区生产总值_当年价格_亿元_全市 地区生产总值增长率_百分比_全市 城市 `vlist'
order 年份 省代码 省 市代码 市 地区生产总值_当年价格_亿元_全市 地区生产总值增长率_百分比_全市 城市 `vlist'
mata:
m1 = st_data(., "`vlist'")
m2 = st_data(., "地区生产总值增长率_百分比_全市")
Wgdpg = m1 * m2
st_matrix("Wgdpg", Wgdpg)
end
keep 省代码 省 市代码 市
svmat Wgdpg
循环所有年份可以采用下面的代码：
*- 由于 mata 程序里面有 end，所以不能直接用在 循环里面，首先定义一个 mata 函数：
clear all
mata:
void matfun() {
vlist = st_local("vlist")
m1 = st_data(., vlist)
m2 = st_data(., "地区生产总值增长率_百分比_全市")
Wgdpg = m1 * m2
st_matrix("Wgdpg", Wgdpg)
}
end
cap mkdir "结果数据2"
use "1998～2021年各地级市地区生产总值面板数据(填补后)", clear
forval i = 1998/2021 {
preserve
qui {
keep if 年份 == `i'
drop if missing(地区生产总值增长率_百分比_全市)
merge 1:1 城市 using "2019年各市是否接壤矩阵（对角线为0）.dta"
keep if _m == 3
drop _m
local vlist = ""
forval j = 1/`=_N' {
local vlist `vlist' `=城市[`j']'
}
keep 年份 省代码 省 市代码 市 地区生产总值_当年价格_亿元_全市 地区生产总值增长率_百分比_全市 城市 `vlist'
order 年份 省代码 省 市代码 市 地区生产总值_当年价格_亿元_全市 地区生产总值增长率_百分比_全市 城市 `vlist'
mata: matfun()
keep 省代码 省 市代码 市
svmat Wgdpg
}
save "结果数据2/`i'", replace
restore
}
合并代码就和上面的一样了。
课程信息
该课程已经讲完，感兴趣的小伙伴可以点击文末的阅读原文跳转到 RStata 平台观看视频讲解～
直播地址：腾讯会议(需要购买 RStata 名师讲堂会员参加) 讲义材料：需要购买 RStata 名师讲堂会员，详情可阅读：一起来学习 R 语言和 Stata 啦！学习过程中遇到的问题也可以随时提问！
更多关于 RStata 会员的更多信息可添加微信号 r_stata 咨询：
附件下载（点击文末的阅读原文即可跳转）：
https://rstata.duanshu.com/#/brief/course/2cfc2c6307bc401b861d0b435e286fff
