全国基準地域メッシュデータの利用

村上貴志 氏が公開している「全国基準地域メッシュデータ」を利用する手順を示す。

1. 分割ファイルからの ZIP ファイル mesh05-jgd-01-shp.zip の生成

   #!/bin/bash
   cd /tmp
   rm -f mesh05-jgd-01-shp.zip
   cat mesh05-jgd-01-shp.0 mesh05-jgd-01-shp.1 mesh05-jgd-01-shp.2 > mesh05-jgd-01-shp.zip
   

2. 解凍

   #!/bin/bash
   cd /tmp
   for i in mesh05-jgd-*-shp.zip; do
     echo $i
     unzip $i
   done
   

3. R の shapefiles パッケージを利用した読み込み

   #!/bin/bash
   cd /tmp
   script=/tmp/a.$$.r
   cat > "$script" <<'RCOMMAND'
   library(shapefiles)
   s <- read.shapefile("mesh05-jgd-47")
   RCOMMAND
   Rscript "$script"
   

シェープファイルの読み込みと地図表示に用いる R パッケージ

R では shapefiles, maptools, sp などのパッケージが利用できる。用途は次のとおりである。

• maptools：シェープファイルを読み込むため
• spsurvey：dbf ファイル読み込みのため
• ggplot2：地図を表示するため
• RColorBrewer：地図の表示に適した色を使いたい場合

注記：maptools, rgdal, rgeos の各パッケージは、2023年10月16日に CRAN から撤去（アーカイブ）された。現在はこれらに代わり sf（および terra）の利用が推奨される。sf ではシェープファイルの読み込みに st_read を用い、地図表示には ggplot2 の geom_sf を用いる。以下では従来の maptools を用いた例を残すが、新規に作成する場合は sf を用いる。

参考：https://sudillap.hatenablog.com/entry/2013/03/26/210202

GADM（Global Administrative Areas）のシェープファイルによる地図表示

GADM が配布する行政区域データを用いて地図を表示する。GADM のファイルは次の階層に対応する。

• JPN_adm0：日本列島の全体地図
• JPN_adm1：都道府県地図
• JPN_adm2：市区町村地図

注記：ここで用いる JPN_adm*.RData（sp 形式）は旧バージョンの GADM（2.x／3.x）で配布されていた形式である。現在の GADM（バージョン 4.1 以降）は GeoPackage・シェープファイル・R（sf）用の rds 形式で配布されており、標準の配布形式は GeoPackage である。旧バージョンのデータは GADM の旧版ページから入手できる。

require(maptools)   # 注記参照：現在は sf を推奨
load("JPN_adm1.RData")
plot(gadm)

require(RColorBrewer)
col <- sample(1:8, size=47, replace=TRUE)   # 県ごとの色
plot(gadm, col=brewer.pal(8,"Accent")[col])

上で読み込んだ属性データ（dbf 由来）の NAME_1 カラムに県名が入っている。

require(ggplot2)
require(maptools)

jpn <- readShapePoly("JPN_adm1.shp")  # シェープファイルはカレントディレクトリにあるとする
map <- fortify(jpn)  # 日本全体
#map <- fortify(jpn[jpn$NAME_1=="Tokyo",])  # 東京都だけ表示したい場合

xlim <- c(128, 146) # 経度
ylim <- c(30, 46)   # 緯度

col <- sample(1:8, size=47, replace=TRUE)  # 県ごとの色

ggplot() + geom_polygon(aes(long,lat,group=group,fill=as.character(col[as.integer(id)+1])),color="black",data=map) +
  coord_fixed(ratio=1) +
  geom_point(aes(x=139.7036, y=35.69389),color="blue",size=10) +
  geom_text(aes(x=139.7036, y=35.69389),label="新宿",hjust=-0.2,color="red",size=10) +
  guides(col=FALSE, fill=FALSE) + xlim(xlim) + ylim(ylim) +
  scale_fill_brewer(type="qual",palette="Accent")


plot(jpn[jpn$NAME_1=="Tokyo",])  # 東京都のデータだけを抽出。JPN_adm2.shp を読み込んでいるとする。
points(139.7036, 35.69389, lwd=20, col="red")  # 新宿区役所の位置に打点
text(139.7036, 35.69389, "新宿", col="blue", adj = c(-0.3,0.5), cex=2)  # 文字を表示