如何使用 Stata 绘图展示保险机构的地理分布

1295 字 约 5 分钟 课程 金融可视化与绘图Stata

之前给大家分享过一份保险许可证持有机构列表数据，在数据中我根据地址解析了经纬度，今天我们来看一下如何使用 Stata 绘图展示这些保险机构的分布。

首先我们需要准备下面几个数据：

  1. chinaprov2021mini_db.dta
  2. chinaprov2021mini_coord.dta
  3. chinaprov2021mini_line_coord2.dta
  4. chinaprov2021mini_label2.dta

这些地图数据都来自（里面有视频讲解）：

×

还有保险机构数据：

  • 保险机构地理位置数据.dta

这个数据的完整版来源于：

×

这里为了减少文件的大小，仅仅保留了经纬度及其所处省份变量。

首先读取数据，生成 id 变量（用以标志观测值，防止数据在之后的处理中乱掉）：

use "保险机构地理位置数据.dta", clear
gen id = _n
encode 省, gen(prov)
keep id 纬度 经度 prov
drop if missing(经度)

encode 命令可以把“省”变量转换生成因子变量 prov。

由于上面的几个 mini 数据都是 “+proj=aea +lat_0=0 +lon_0=105 +lat_1=25 +lat_2=47 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs” 坐标系的，所以这里需要对经纬度进行坐标系转换：

* ssc install geo2xy
geo2xy 纬度 经度, gen(lat lon) projection(albers, 6378137 298.257223563 25 47 0 105) replace
gen class = "main"
save mainpoint, replace

这样就完成了坐标转换。

不过因为这个地图还带一个南海诸岛的小地图，所以为了更好的绘图效果，我们再提取小地图范围的散点并平移：

* 提取小地图框格的
use mainpoint, clear
keep if inrange(lon, 120000, 1766004.1) & inrange(lat, 320000, 2557786.0)
replace lon = lon * 0.5 + 2100000
replace lat = lat * 0.5 + 1665139
replace class = "smallbox"

* 把两部分数据合并起来
append using mainpoint
save pointdata, replace

下面就可以把这些散点绘制到地图上了。

首先绘制一副空白地图：

use chinaprov2021mini_db.dta, clear
spmap using chinaprov2021mini_coord.dta, id(ID) ///
 ocolor("black" ...) osize(vvthin ...) ///
 line(data(chinaprov2021mini_line_coord2.dta) ///
  select(keep if inlist(group, 1, 2, 3, 5, 6)) ///
  by(group) size(vvthin *1 *0.5 *0.5 *0.5) ///
  pattern(solid ...) ///
  color(white /// 省界颜色
     black /// 国界线颜色
     "0 85 170" /// 海岸线颜色
     black /// 小地图框格颜色
     black /// 比例尺和指北针颜色
     )) ///
 polygon(data(polygon2) fcolor(black) ///
  osize(vvthin)) ///
 label(data(chinaprov2021mini_label2) x(X) y(Y) label(cname) length(20) size(*0.8))

chinaprov2021mini_line_coord2.dta 数据包含了绘图所需的所有线条，其中我把线条分成了其中类型（group 变量）：

  1. 省界
  2. 国界线（含九段线）
  3. 海岸线
  4. 秦岭-淮河线
  5. 小地图框格
  6. 比例尺和指北针
  7. 胡焕庸线

这里不需要秦岭-淮河线和胡焕庸线，所以 select(keep if inlist(group, 1, 2, 3, 5, 6)) 中删除了这两个；size(vvthin *1 *0.5 *0.5 *0.5) 里面的 5 种 size 分别是剩下五种线条的粗细；pattern(solid ...) 表示所有的线条都使用 solid（实线）；color() 中的五种颜色也分别是五种类型线条的颜色。

然后添加 point() 选项：

spmap using chinaprov2021mini_coord.dta, id(ID) ///
 ocolor("black" ...) osize(vvthin ...) ///
 line(data(chinaprov2021mini_line_coord2.dta) ///
  select(keep if inlist(group, 1, 2, 3, 5, 6)) ///
  by(group) size(vvthin *1 *0.5 *0.5 *0.5) ///
  pattern(solid ...) ///
  color(white /// 省界颜色
     black /// 国界线颜色
     "0 85 170" /// 海岸线颜色
     black /// 小地图框格颜色
     black /// 比例尺和指北针颜色
     )) ///
 polygon(data(polygon2) fcolor(black) ///
  osize(vvthin)) ///
 label(data(chinaprov2021mini_label2) x(X) y(Y) label(cname) length(20) size(*0.8)) ///
 point(data(pointdata) by(prov) ///
  fcolor("80 80 255" "206 61 50" "116 155 88" "240 230 133" ///
   "70 105 131" "186 99 56" "93 177 221" "128 34 104" ///
   "107 215 107" "213 149 167" "146 72 34" "131 123 141" ///
   "199 81 39" "213 143 92" "122 101 165" "228 175 105" ///
   "59 27 83" "205 222 183" "97 42 121" "174 31 99" ///
   "231 199 111" "90 101 94" "204 153 0" "153 204 0" ///
   "169 169 169" "204 153 0" "153 204 0" "51 204 0" ///
   "0 204 51" "0 204 153" "0 153 204") ///
  x(lon) y(lat) ///
  size(*0.1 ...) shape(O ...)) ///
 leg(off)

point() 中的 by(prov) 表示分省份绘制，fcolor() 里面的颜色分别表示各省的散点颜色，至少需要 31 种颜色。

颜色的选择可以使用这个应用：https://tidyfriday.cn/colors

最后添加标题、副标题等：

spmap using chinaprov2021mini_coord.dta, id(ID) ///
 ocolor("black" ...) osize(vvthin ...) ///
 line(data(chinaprov2021mini_line_coord2.dta) ///
  select(keep if inlist(group, 1, 2, 3, 5, 6)) ///
  by(group) size(vvthin *1 *0.5 *0.5 *0.5) ///
  pattern(solid ...) ///
  color(white /// 省界颜色
     black /// 国界线颜色
     "0 85 170" /// 海岸线颜色
     black /// 小地图框格颜色
     black /// 比例尺和指北针颜色
     )) ///
 polygon(data(polygon2) fcolor(black) ///
  osize(vvthin)) ///
 label(data(chinaprov2021mini_label2) x(X) y(Y) label(cname) length(20) size(*0.8)) ///
 point(data(pointdata) by(prov) ///
  fcolor("80 80 255" "206 61 50" "116 155 88" "240 230 133" ///
   "70 105 131" "186 99 56" "93 177 221" "128 34 104" ///
   "107 215 107" "213 149 167" "146 72 34" "131 123 141" ///
   "199 81 39" "213 143 92" "122 101 165" "228 175 105" ///
   "59 27 83" "205 222 183" "97 42 121" "174 31 99" ///
   "231 199 111" "90 101 94" "204 153 0" "153 204 0" ///
   "169 169 169" "204 153 0" "153 204 0" "51 204 0" ///
   "0 204 51" "0 204 153" "0 153 204") ///
  x(lon) y(lat) ///
  size(*0.1 ...) shape(O ...)) ///
 leg(off) ///
 ti("保险许可证持有机构地理分布（截止2023-01-31）", color(black)) ///
 subti("绘制：微信公众号 RStata") ///
 graphr(margin(medium)) ///
 caption("数据来源：银保监会，使用高德地图地理编码接口解析经纬度", size(*0.8))

点击这里跳转到 RStata 短书平台获取附件：如何使用 Stata 绘图展示保险机构的地理分布

相关文章

评论