之前给大家分享过使用 Stata 绘制中国地图的课程和数据资料：
使用 Stata 绘制地图课程汇总：https://mp.weixin.qq.com/s/wCVsicJaq-lXzrQmCRKXgQ
其中强烈建议学习这个课程：使用 Stata 绘制历年中国市级行政区划（小地图版本 + 长版）：https://rstata.duanshu.com/#/brief/course/0e9d78633ae64fefa684d6aad89633bc
前不久也给大家分享了台风路径数据：
在数据介绍中我绘制了这么一幅图：
今天的课程中我们将分享这幅图的绘制方法。
首先是绘制各区县台风过境次数：
clear all
use "1949～2024年各区县每天是否有台风（有台风的天）.dta", clear
gen year = yofd(日期)
collapse (sum) 是否有台风, by(省 省代码 市 市代码 县 县代码 year)
keep if year == 2024
merge 1:1 省 省代码 市 市代码 县 县代码 using chinacounty2021mini_db.dta
replace 是否有台风 = 0 if mi(是否有台风)
drop _m
grmap 是否有台风 using chinacounty2021mini_coord.dta,
id(ID) osize(vvthin ...) ocolor(white ...)
clmethod(unique)
fcolor("242 240 247" "203 201 226"
"158 154 200" "106 81 163")
graphr(margin(medium))
leg(order(`r(legorder)'))
line(data(chinaprov2021mini_line_coord2.dta) by(group)
size(vvthin *1 *0.5 *0.5 *0.5) pattern(solid ...)
select(drop if inlist(group, 4, 7))
color(gs20
black
"0 85 170"
black
black
))
polygon(data(polygon2) fcolor(black)
osize(vvthin))
label(data(chinacounty2021mini_label2) x(X) y(Y) label(cname) length(20) size(*0.8))
ti("2024年各区县台风（热带气旋）过境天数")
subti("数据整理 & 绘制：微信公众号 RStata")
caption("数据来源：中国气象局热带气旋资料中心" "<https://tcdata.typhoon.org.cn/zjljsjj.html>", size(*0.8))
gr export "2024年各区县台风（热带气旋）过境天数.png", replace width(2400)
|
这部分代码学习上面的课程就可以理解了。
由于现有的代码已经使用了 line 数据，也就是 chinaprov2021mini_line_coord2.dta 数据，所以如果还想添加台风路径线条，就需要把台风路径的线条数据和该数据合并。
先把台风路径的 shp 文件转换成 dta 文件：
shp2dta using "台风分日期路径/line.shp", database(tfline_db) coord(tfline_coord) replace
|
这里我们截取 2024 年以来的数据：
use tfline_db, clear
keep if date >= 20240101
keep _ID
joinby _ID using tfline_coord
keep if inrange(_X, -180, 180) | mi(_X)
|
由于使用的地图底图是 aea 坐标系，所以这里我们也需要把台风路径的数据也转换成和底图相同的坐标系：
geo2xy _Y _X, gen(y x) projection(albers, 6378137 298.257223563 25 47 0 105) replace
drop _Y _X
ren (y x) (_Y _X)
*- 保留中国周边的（这个范围就是中国地图主图的坐标范围）
keep if (inrange(_X, -2625586, 2206965) & inrange(_Y, 1836814.1, 5921583)) | mi(_X)
save 2024tfline_coord, replace
|
台风路径的最后两个点可以用来构造箭头：
use 2024tfline_coord, clear
bysort _ID: keep if _n == _N | _n == _N - 1
drop if mi(_X)
bysort _ID: gen order = _n
gather _X _Y
tostring order, replace
replace var = var + order
drop order
spread var value
save 2024arrow, replace
|
箭头数据的格式是这样的：
list in 1/10
*> +-------------------------------------------------------+
*> | _ID _X1 _X2 _Y1 _Y2 |
*> |-------------------------------------------------------|
*> 1. | 18804 2012746.9 2074545.5 1840813.5 1928910.8 |
*> 2. | 18805 2126568.9 2183960.5 2004878.6 2104597 |
*> 3. | 18811 759256.16 757734.94 1885722.1 1907261.7 |
*> 4. | 18812 685930.26 683071.68 2174476 2218075.4 |
*> 5. | 18813 866485.23 953791.25 2540490.9 2592759.1 |
*> |-------------------------------------------------------|
*> 6. | 18814 1193659.1 1227891 2640711.9 2611436.2 |
*> 7. | 18821 2183960.5 . 2104597 . |
*> 8. | 18822 2102173.1 2022566.5 2342652.3 2471601.2 |
*> 9. | 18823 1703833.2 1714694.1 2602877.1 2660714.3 |
*> 10. | 18824 1427029.5 1407064.9 2727280.2 2724425.5 |
*> +-------------------------------------------------------+
|
在原有线条数据 chinaprov2021mini_line_coord2.dta 里面，group 变量用以对线条数据分组，取值为 1-7。合并台风路径后可以把台风路径设置为 group = 8：
use chinaprov2021mini_line_coord2.dta, clear
append using 2024tfline_coord
replace group = 8 if mi(group)
save mylinecoord, replace
|
然后就可以在地图上添加台风路径和方向箭头了：
use "1949～2024年各区县每天是否有台风（有台风的天）.dta", clear
gen year = yofd(日期)
collapse (sum) 是否有台风, by(省 省代码 市 市代码 县 县代码 year)
keep if year == 2024
merge 1:1 省 省代码 市 市代码 县 县代码 using chinacounty2021mini_db.dta
replace 是否有台风 = 0 if mi(是否有台风)
drop _m
gr export "2024年我国周边台风（热带气旋）路径.png", replace width(2400)
|
这样就搞定了这个图。
点击这里跳转到 RStata 短书平台获取附件：使用 Stata 绘制中国地图+台风路径（线条+箭头）
评论