使用 Plot.jl 在日本地图上显示数据

动机

我想用 Julia 的 Plots.jl 显示日本统计数据的地图。我认为“应该使用 Makie.jl”，但我想用熟悉的工具（如 StatsPlots）来实现它。我想显示这样的图表。
Plot.jlで日本地図にデータを表示する

获取地图数据

日本地理空间信息局、国土交通省下载（zip 文件）并将所有图层数据解压缩为 .根据维基百科，地图数据以一对三个文件的形式存在。

.shp — 形状标准：地形信息主体。
.shx — 形状索引标准：用于快速向前和向后搜索地形数据的位置索引。
.dbf — 属性标准：每个形状的垂直表格属性信息。符合 dBASE IV 格式。

如果你有 shp 数据，你可以只显示一个地图，所以如果你有 Shapfile.jl，它是一个处理这种格式的 Libiray，你就可以做到。

引入 Shapefile.jl 并读取地图数据

   using Plots, DataFrames, CSV, Shapefile
   coast_table = Shapefile.Table("geodata-jpn-all_u_2_2\\coastl_jpn.shp")
   table_jpn = Shapefile.Table("geodata-jpn-all_u_2_2\\polbnda_jpn.shp")

Table 函数将 shp 数据读取为表数据。表数据包括区域名称、ID等。形状函数将表格中的多行收集到一个多边形数据中。如果传递给 Plots，收集的多边形数据将按原样显示。例如，上面的coast_shape数据可以显示如下。

    coast_shape = Shapefile.shapes(coast_table)
    plot(coast_shape,color=:gray)

Plot.jlで日本地図にデータを表示する

看一眼表中的数据

表格数据可以轻松转换为 DataFrame。

Julia> df_pola = DataFrame(table_jpn);
5×9 DataFrame
 Row │ f_code  coc     nam        laa            pop      ypc    adm_code  salb    soc    
     │ String  String  String     String         Int64    Int64  String    String  String 
─────┼────────────────────────────────────────────────────────────────────────────────────
   1 │ FA001   JPN     Hokkai Do  Sapporo Shi    1930496   2014  01100     UNK     JPN
   2 │ FA001   JPN     Hokkai Do  Hakodate Shi    274485   2014  01202     UNK     JPN
   3 │ FA001   JPN     Hokkai Do  Otaru Shi       127224   2014  01203     UNK     JPN
   4 │ FA001   JPN     Hokkai Do  Asahikawa Shi   349057   2014  01204     UNK     JPN
   5 │ FA001   JPN     Hokkai Do  Muroran Shi      91276   2014  01205     UNK     JPN

Plots.jl 可以通过从 Table 中提取相应的行数据并根据 Table 的元数据将其传递给 Shape (Shapes if multiple) 来显示。 adm_code 是日本地理调查研究所最重要的元数据。这是城市代码它用于各种政府统计数据的区域分类。总务省页面有一个城市代码列表，所以请参考它。

市政代码由 5 位数字和最后 1 位奇偶校验组成（防止输入错误）。
日本地理空间信息管理局元数据仅使用前五位数字。
是一个数字但不是整数 Integer 而是一串数字 String。

注意到以上几点，可以想到一个函数，它用自治市代码返回相关区域的形状。

从市政代码返回形状的函数

function selectshape(table, tg)
    rs = empty(Shapefile.shapes(table_jpn))
    for row in table_jpn
        if typeof(tg) == Regex &&  match(tg, row.adm_code) !== nothing
            push!(rs,Shapefile.shape(row))
        elseif  typeof(tg) == String && tg ==row.adm_code
            push!(rs,Shapefile.shape(row))
        end
    end
    return rs
end

之所以可以通过正则表达式进行搜索，是因为市町村代码使用前两位数字表示都道府县。您可以使用正则表达式来搜索：

Julia> hokkaido = selectshape(table_jpn, r"^01d*");　#01~は北海道
Julia> plot(hokkaido,alpha=0.3,color=:green)

Plot.jlで日本地図にデータを表示する

像关东这样的广阔区域可以如下进行。

kanto = begin
    kanto_rgx = [r"^08d*",r"^09d*",r"^10d*",r"^11d*",r"^12d*",r"^13d*",r"^14d*"]
    # 08~茨木、09~栃木、10~群馬、11~埼玉、12~千葉、13～東京都, 14~神奈川
    [ selectshape(table_jpn, rgx) for rgx in kanto_rgx]
end
pl = plot(xlims=(138,142),ylims=(33.5,37.5))
for (e,c) in zip(kanto,[:yellow,:red,:green,:blue,:orange,:skyblue,:lightgreen])
    plot!(pl, e,color=c,lw=0)
end
display(pl)

Plot.jlで日本地図にデータを表示する

高级——显示实际统计数据

内阁府网站让我们在 1980 年和 2010 年重现人口统计指数（普通生育率：每 1000 人的出生人数）的映射示例。

首先，从下载的数据中删除不必要的数据并以 CSV 格式读取。此时市政府代码将类型指定为字符串必须。

Julia> birth_da=CSV.read("市区町村別人口1000人当たり年間出生数.csv",
                DataFrame ; types=Dict(:adm_code => String) )
5×6 DataFrame
 Row │ nam            adm_code  y1980    y1990    y2000    y2010   
     │ String         String    Float64  Float64  Float64  Float64 
─────┼─────────────────────────────────────────────────────────────
   1 │ 北海道 札幌市  011002       15.0     10.5      8.5      7.6
   2 │ 北海道 函館市  012025       13.3      9.1      7.5      6.5
   3 │ 北海道 小樽市  012033       10.8      7.1      6.4      5.4
   4 │ 北海道 旭川市  012041       13.9      9.4      8.3      7.4
   5 │ 北海道 室蘭市  012050       13.7      8.1      7.6      6.7

绘制

n,m=size(birth_da)
pl_80=plot()
pl_10=plot()
year80 = "y1980"
year10 = "y2010"

for i = 1:n
    streetnam = birth_da[i,:adm_code][1:5]
    damap = selectshape(table_jpn, streetnam)
    if length(damap) > 0
        plot!(pl_80,damap, 
            lw=0, 
            clims=(1,25),
            fill_z=birth_da[i, Symbol(year80)],
            fillcolor=:inferno)
        plot!(pl_20,damap, 
            lw=0, 
            clims=(1,25),
            fill_z=birth_da[i, Symbol(year10)], 
            fillcolor=:inferno)
    end
end

L=@layout [a{0.5w} b{0.5w}]
plot(pl_80, pl_20,layout=L,size=(1200,800),dpi=120)

您现在可以比较地图。
Plot.jlで日本地図にデータを表示する

未来的任务

我想整齐地布置冲绳和小笠原。
我希望能够使用 GIS 地理信息进行地图绘制。
我希望能够使用 OpenStreetMapX.jl。
我想通过投影法反映正确的投影。

原创声明：本文系作者授权爱码网发表，未经许可，不得转载;

原文地址：https://www.likecs.com/show-308623766.html