使用 Stata 绘制三变量填充中国地图

1952 字 约 7 分钟 课程 可视化与绘图地图Stata

在之前「使用 Stata 绘制历年中国各省市区县地图（小地图版本 + 长版）」课程的基础上，我们今天再来学习下如何绘制三变量填充地图，之前的课程链接：

  1. 使用 Stata 绘制历年中国省级行政区划（小地图版本 + 长版）：https://rstata.duanshu.com/#/brief/course/e194d75fe1674f7c8921daec1ece7f7d
  2. 使用 Stata 绘制历年中国市级行政区划（小地图版本 + 长版）：https://rstata.duanshu.com/#/brief/course/0e9d78633ae64fefa684d6aad89633bc
  3. 使用 Stata 绘制历年中国县级行政区划（小地图版本 + 长版）：https://rstata.duanshu.com/#/brief/course/99c2e97b88ea401bbdb8da4511341cc8

另外平台上也有双变量填充地图绘制的教程：

今天我们将会进一步讲解如何使用 Stata 制作三变量填充地图，同样也包含两种版本：

小地图版本

由于三变量填充地图的图例制作较为费事，所以我只提供了 10 阶的。在制作地图之前我们需要简单了解下这种三元图如何阅读：

为了给对应的格子生成特定的颜色，我设计了一个“三变量填充地图颜色生成器”：https://observablehq.com/d/ac2462c27dac3e35

该生成器生成的颜色顺序就是按照图中的 1～55:

三元图看起来很绕，但是抛去具体的值，直接看颜色还是很直观的，例如上图中的蓝色就表示这组数据的 B 含量“满满”，A 和 C 都较少。

所以这个图例的使用方法很简单，下面我们结合具体数据来看下。

各城市人口年龄结构

读取和处理个城市人口年龄结构数据（2020 年第七次人口普查结果）：

use "七普各市人口数据_各年龄组人口比重、有老年人口的户数、户口登记地在外乡镇街道人口.dta", clear
keep 市 市代码 各年龄组人口占总人口比重_0_14岁人口比重_百分比 各年龄组人口占总人口比重_15_64岁人口比重_百分比 各年龄组人口占总人口比重_65岁及以上人口比重_百分比
ren 各年龄组人口占总人口比重_0_14岁人口比重_百分比 年龄0_14岁
ren 各年龄组人口占总人口比重_15_64岁人口比重_百分比 年龄15_64岁
ren 各年龄组人口占总人口比重_65岁及以上人口比重_百分比 年龄65岁及以上
save mydata, replace

把图例的 shp 文件转换成 dta 文件：

local name = "triscale10"
shp2dta using `name'/`name', database(`name'_db) coordinates(`name'_coord) genid(ID2) gencentroids(centroids) replace

use triscale10_db.dta, clear
drop ID2
order ID
save triscale10_db.dta, replace

图例的标签也需要替换成我们需要的：

import excel using "triscale10/triscale10_label.xlsx", clear first
replace cname = "15～64岁" in 1 // B
replace cname = "0～14岁" in 2  // A
replace cname = "65 岁及以上" in 3 // C
save "triscale10_label", replace

这里要非常注意 B A C 三个的对应关系。

这里我准备绘制的是 市级填充地图 + 省级线条 + 散点图，需要下面的这些数据：

  • chinacity2020mini_coord.dta
  • chinacity2020mini_db.dta
  • chinacity2020mini_label2.dta
  • chinaprov2020mini_line_coord2.dta
  • polygon2.dta

这些在之前的课程中都有提供。

可以看到我准备的这个图例数据包含两部分，多边形和文本标签。由于 Stata 绘制地图的时候只能添加一组多边形数据和一组文本标签数据，所以我们需要把 triscale10 合并到 polygon2.dta 数据上面：

use polygon2, clear
gen class = "指北针和比例尺"
append using triscale10_coord
replace class = "三变量图例" if missing(class)

如果使用的双变量图例阶数很高，可能会出现 _ID 重复的情况，所以这里我们最后给指北针和比例尺的 _ID 调大点：

replace _ID = _ID + 100 if class == "指北针和比例尺"
save polygon2_with_triscale, replace

还需要把 triscale10_label 合并到 chinacity2020mini_label2.dta 文件中：

use chinacity2020mini_label2.dta, clear
append using triscale10_label
save chinacity2020mini_label2_triscale, replace

把数据的每个变量分成 10 组：

use mydata, clear
drop if missing(年龄65岁及以上)
*- B: 年龄15_64岁
*- A: 年龄0_14岁
*- C: 年龄65岁及以上

replace 年龄65岁及以上 = 100 - 年龄0_14岁 - 年龄15_64岁

egen groupB = cut(年龄15_64岁), at(0 10 20 30 40 50 60 70 80 90 100) icodes
replace groupB = groupB + 1
tab groupB

egen groupA = cut(年龄0_14岁), at(0 10 20 30 40 50 60 70 80 90 100) icodes
replace groupA = groupA + 1
tab groupA

egen groupC = cut(年龄65岁及以上), at(0 10 20 30 40 50 60 70 80 90 100) icodes
replace groupC = groupC + 1
tab groupC

gen v = groupB + groupA + groupC
tab v

可以看到由于分组的“损失”，并不是所有城市的 groupA + groupB + groupC 都是 12，因此我们这里需要使点“诈”，可以把自己希望强调的变量增加 1 单位，例如老年人口：

replace groupC = 12 - groupA - groupB
drop v

这种分组的“损失”原理是这样的，例如 9.8、73.92 和 16.28，三者的和是 100，但是分组的时候 9.8 和 73.92 分别因为不足 10 和 不足 80 被归类到第 1 组和第 8 组（都是舍去了），这就导致最终的总和 1 + 8 + 2 < 12 了，而我们的绘图数据要求这三个变量之后应该是 12（例如老年人口最多的城市，三个变量应该是 10 + 1 + 1）。

tostring groupB groupA groupC, replace

*- 这里一定要注意是 B + A + C
gen groupclass = groupB + "-" + groupA + "-" + groupC
drop groupB groupA groupC

codebook groupclass

*- 对 groupclass 进行有序因子化
merge m:1 groupclass using triscale10_db
drop if _m == 2
drop _m
*- 安装 labmask: net install gr0034.pkg, from("http://www.stata-journal.com/software/sj8-2") replace
labmask ID, val(groupclass)
ren ID groupclass2
drop x_centroids - C

和 chinacity2020mini_db.dta 匹配：

merge 1:1 市 市代码 using "chinacity2020mini_db.dta"
replace groupclass2 = -1 if missing(groupclass2)
save "mapdata", replace

绘制地图

然后就可以绘制地图了：

use mapdata, clear
gsort 市代码

local colorlist = `""0 64 255" "30 41 233" "18 96 205" "58 24 210" "49 73 185" "35 123 164" "85 13 187" "78 56 165" "67 101 147" "51 145 131" "112 6 164" "106 44 145" "98 83 129" "86 123 116" "70 165 104" "137 2 143" "135 36 126" "129 70 112" "120 106 101" "108 142 91" "92 181 81" "163 1 123" "163 30 108" "160 60 96" "154 91 87" "146 124 79" "135 159 70" "120 196 60" "189 2 107" "191 26 93" "191 52 83" "189 79 75" "185 108 69" "178 139 61" "168 173 53" "155 210 41" "216 2 94" "220 22 82" "223 44 73" "225 68 67" "224 93 61" "222 121 55" "217 152 47" "210 186 36" "200 224 22" "244 1 86" "251 18 75" "255 36 67" "255 57 61" "255 79 56" "255 104 51" "255 132 44" "255 164 34" "255 199 20" "255 238 0""'
grmap groupclass2 using chinacity2020mini_coord.dta, ///
  id(ID) osize(vvthin ...) ocolor(white ...) ///
  clmethod(custom) clbreaks(-1 0 1(1)55) ///
    fcolor("gs12" `colorlist') ///
    leg(off) ///
  graphr(margin(medium)) ///
  line(data(chinaprov2020mini_line_coord2.dta) by(group) ///
    size(vthin *1 *0.5 *0.5 *0.5) pattern(solid ...) ///
    select(drop if inlist(group, 4, 7)) ///
    color(gs5 /// 省界颜色
        "162 154 196" /// 国界线颜色
        "0 85 170" /// 海岸线颜色
        black /// 小地图框格颜色
        black /// 比例尺和指北针颜色
        )) ///
  polygon(data(polygon2_with_triscale) by(_ID) ///
    osize(vvvthin ...) ocolor(white ...) ///
    fcolor(`colorlist' black ...)) ///
  label(data(chinacity2020mini_label2_triscale) x(X) y(Y) ///
    select(keep if inlist(cname, "N", "1000km") | index(cname, "岁")) ///
    label(cname) length(40 ...) size(*0.7 ...)) ///
  ti("第七次人口普查各城市分年龄组人口比重", size(*0.8)) ///
  subti("数据处理&绘制：微信公众号 RStata", size(*0.8)) ///
  caption("数据来源：2020年分县市人口统计年鉴", size(*0.6))

gr export pic10x10x10.png, replace width(4800)

可以看到地图上蓝色的城市表示这些城市的劳动力人口比重较高（15～64岁），偏黄色的则表示少年人口比重较老年人口比重高，偏红色则相反。

长版

附件中的“长版”文件夹中也提供了长版地图的绘制数据和代码，绘制方法类似，这里就不再赘述了：

点击这里跳转到 RStata 短书平台获取附件：使用 Stata 绘制三变量填充中国地图

相关文章

评论