Rを使って郵便番号界の地図を作る(統計局のShape+日本郵便のCSV)
郵便番号がどんな感じで割り当てられているか、地図上で確認したかったんです。
GPSで取った地点情報がプライバシーに配慮云々で郵便番号になっていた、なんてことがありまして、その粒度がどんなものか確認したかったってのが発端です。
検索してみても無償のものは意外に見つからないもので、国際航業さんが売っていたりします↓
PAREA-Zip郵便番号界 | 国際航業株式会社
ノースリーブな(=無い袖は振れない)わたくしとしては、自分で作るしかないわけで。
郵便番号のリストは↓ここにありました。
郵便番号データダウンロード - 日本郵便
住所文字列と郵便番号が入っていればどれでもいいんですが、とりあえず、「住所の郵便番号(CSV形式)」の「読み仮名データの促音・拗音を小書きで表記するもの」というのをダウンロードしました。
次は地図です。こういうときに使う境界データとしては、統計局のe-StatでShapeファイルが公開されているので、これが便利。入手の方法は下記を参照。
e-Stat(統計局)で公開されているShapeファイルを、Rで表示する - Rプログラミングの小ネタ
今回は、我らがサラリーマンの聖地 新橋を有する港区を対象にしてみました。
とりあえず地図をプロットするところまでやってみましょう。
library(maptools) # GIS系のパッケージ
# 座標系の情報を持つオブジェクトを作る
# longlatは緯度経度指定、WGS84は世界測地系であることを表す
pj <- CRS("+proj=longlat +datum=WGS84")
# シェープファイルの読み込み
map <- readShapePoly("h22ka13103.shp", proj4string=pj)
# とりあえずプロット
plot(map)
この地図上に追加で郵便番号をプロットすればよさそうですね。
読み込んだシェープファイルにはポリゴンの境界情報だけでなく、それぞれのポリゴンの属性の情報も含まれています。
こうやれば確認できます
head(map@data)
AREA (略) KEN_NAME SITYO_NAME GST_NAME (略) MOJI
0 729237.30 東京都 <NA> 港区 元赤坂2丁目
1 148246.50 東京都 <NA> 港区 北青山1丁目
2 89514.52 東京都 <NA> 港区 元赤坂1丁目
3 93129.48 東京都 <NA> 港区 赤坂3丁目
4 151487.30 東京都 <NA> 港区 赤坂4丁目
5 199895.30 東京都 <NA> 港区 北青山2丁目
「MOJI」っていう列に住所が入っているようです。これだけにアクセスしたい場合は、「map$MOJI」と打てばOK。厳密には「map@data$MOJI」なんでしょうけど、省略できるようです。
これを郵便局のCSVファイルとマッチさせればよさそう。
# 郵便番号データの読み込み
post_data <- read.csv("13TOKYO.CSV", head=F)
post_data <- post_data[post_data$V8=="港区", ] # 港区だけ取り出し
# ○丁目の部分をカットして、MOJI2の列に格納
map$MOJI2 <- sub("[0-9].*", "", map$MOJI)
# MOJI2と郵便番号データの文字列がマッチしたら、
# 郵便番号を取り出して設定する
for ( i in 1:nrow(map@data) ) {
# 複数マッチしたら、安直に1つ目を採用
matched_index <- grep(map$MOJI2[i], post_data$V9)[1]
map$POST[i] <- post_data$V3[matched_index]
}
plot(map)
text( coordinates(map), as.character(map$POST), cex=0.5)
一応、最低限のところまではできました。でも、文字が見にくいとか、地名も欲しいとか、色塗ったほうがキャッチーだよね(?)とか、ということで付け足してみたのが以下。
# ハイフンを入れる
post <- paste ( substr(map$POST, 1, 3), substr(map$POST, 4, 7), sep="-" )
# 地名&郵便番号の文字列を作る
labs <- paste(map$MOJI, post, sep="\n")
# 塗りわけのパレット
library(RColorBrewer) # カラーパレット
cols <- brewer.pal(12, "Set3")
cols <- rep(cols, 10) # 12色じゃ足りないので繰り返し
pdf("郵便番号界マップ(港区).pdf", family="Japan1GothicBBB", width=16.54, height=11.69)
# 郵便番号を因子型にすると、番号ごとに塗りわけできる
plot(map, col=cols[factor(map$POST)], border="gray", lwd=0.1)
text( coordinates(map), labs, cex=0.3)
dev.off()
上記のサンプルコードで作ったのが↓これです。
郵便番号界マップ(港区)(PDFファイル)
↓PDFファイルのスクリーンショット
なかなかいい感じになったんじゃないでしょうか。郵便番号ごとに塗り分けるために、因子型にするなんてのは、ちょっとしたtipsです。
が、ツッコミが入る前に白状しておかねばなりますまい。この郵便番号地図、正しくないところがあります。
住所のマッチングが結構甘いです。grepでマッチしたものが複数ある場合、単純に1つ目のものを採用しています。
例えば、↓こういうパターンに対してはうまく機能しています。
"1070052", "赤坂(次のビルを除く)"
"1076090", "赤坂赤坂アークヒルズ・アーク森ビル(地階・階層不明)"
"1076001", "赤坂赤坂アークヒルズ・アーク森ビル(1階)"
・・・
が、↓こんなやつだとうまくいきませんね。
"1050014", "芝(1~3丁目)"
"1080014", "芝(4、5丁目)"
芝1丁目、芝2丁目、・・・、芝5丁目を全部別の行にしてくれていたとしたら、まだプログラム側でなんとかしようなんて気も起きますが、「1~3」とか「4、5」とか、機械的に処理することを全く考慮していないデータなんですもの。今後の日本郵便のがんばりに期待したいところですね。(上から目線)
ということで、こんな間違いだらけのマップじゃ使えないよってかたは、国際航業さんから買ってください。
GPSで取った地点情報がプライバシーに配慮云々で郵便番号になっていた、なんてことがありまして、その粒度がどんなものか確認したかったってのが発端です。
検索してみても無償のものは意外に見つからないもので、国際航業さんが売っていたりします↓
PAREA-Zip郵便番号界 | 国際航業株式会社
ノースリーブな(=無い袖は振れない)わたくしとしては、自分で作るしかないわけで。
郵便番号のリストは↓ここにありました。
郵便番号データダウンロード - 日本郵便
住所文字列と郵便番号が入っていればどれでもいいんですが、とりあえず、「住所の郵便番号(CSV形式)」の「読み仮名データの促音・拗音を小書きで表記するもの」というのをダウンロードしました。
次は地図です。こういうときに使う境界データとしては、統計局のe-StatでShapeファイルが公開されているので、これが便利。入手の方法は下記を参照。
e-Stat(統計局)で公開されているShapeファイルを、Rで表示する - Rプログラミングの小ネタ
今回は、我らがサラリーマンの聖地 新橋を有する港区を対象にしてみました。
とりあえず地図をプロットするところまでやってみましょう。
library(maptools) # GIS系のパッケージ
# 座標系の情報を持つオブジェクトを作る
# longlatは緯度経度指定、WGS84は世界測地系であることを表す
pj <- CRS("+proj=longlat +datum=WGS84")
# シェープファイルの読み込み
map <- readShapePoly("h22ka13103.shp", proj4string=pj)
# とりあえずプロット
plot(map)
shapeをプロットしただけのもの |
この地図上に追加で郵便番号をプロットすればよさそうですね。
読み込んだシェープファイルにはポリゴンの境界情報だけでなく、それぞれのポリゴンの属性の情報も含まれています。
こうやれば確認できます
head(map@data)
AREA (略) KEN_NAME SITYO_NAME GST_NAME (略) MOJI
0 729237.30 東京都 <NA> 港区 元赤坂2丁目
1 148246.50 東京都 <NA> 港区 北青山1丁目
2 89514.52 東京都 <NA> 港区 元赤坂1丁目
3 93129.48 東京都 <NA> 港区 赤坂3丁目
4 151487.30 東京都 <NA> 港区 赤坂4丁目
5 199895.30 東京都 <NA> 港区 北青山2丁目
「MOJI」っていう列に住所が入っているようです。これだけにアクセスしたい場合は、「map$MOJI」と打てばOK。厳密には「map@data$MOJI」なんでしょうけど、省略できるようです。
これを郵便局のCSVファイルとマッチさせればよさそう。
# 郵便番号データの読み込み
post_data <- read.csv("13TOKYO.CSV", head=F)
post_data <- post_data[post_data$V8=="港区", ] # 港区だけ取り出し
# ○丁目の部分をカットして、MOJI2の列に格納
map$MOJI2 <- sub("[0-9].*", "", map$MOJI)
# MOJI2と郵便番号データの文字列がマッチしたら、
# 郵便番号を取り出して設定する
for ( i in 1:nrow(map@data) ) {
# 複数マッチしたら、安直に1つ目を採用
matched_index <- grep(map$MOJI2[i], post_data$V9)[1]
map$POST[i] <- post_data$V3[matched_index]
}
plot(map)
text( coordinates(map), as.character(map$POST), cex=0.5)
郵便番号をtext関数で追加プロットしたもの |
一応、最低限のところまではできました。でも、文字が見にくいとか、地名も欲しいとか、色塗ったほうがキャッチーだよね(?)とか、ということで付け足してみたのが以下。
# ハイフンを入れる
post <- paste ( substr(map$POST, 1, 3), substr(map$POST, 4, 7), sep="-" )
# 地名&郵便番号の文字列を作る
labs <- paste(map$MOJI, post, sep="\n")
# 塗りわけのパレット
library(RColorBrewer) # カラーパレット
cols <- brewer.pal(12, "Set3")
cols <- rep(cols, 10) # 12色じゃ足りないので繰り返し
pdf("郵便番号界マップ(港区).pdf", family="Japan1GothicBBB", width=16.54, height=11.69)
# 郵便番号を因子型にすると、番号ごとに塗りわけできる
plot(map, col=cols[factor(map$POST)], border="gray", lwd=0.1)
text( coordinates(map), labs, cex=0.3)
dev.off()
上記のサンプルコードで作ったのが↓これです。
郵便番号界マップ(港区)(PDFファイル)
↓PDFファイルのスクリーンショット
郵便番号に、色と地名を追加したもの |
なかなかいい感じになったんじゃないでしょうか。郵便番号ごとに塗り分けるために、因子型にするなんてのは、ちょっとしたtipsです。
が、ツッコミが入る前に白状しておかねばなりますまい。この郵便番号地図、正しくないところがあります。
住所のマッチングが結構甘いです。grepでマッチしたものが複数ある場合、単純に1つ目のものを採用しています。
例えば、↓こういうパターンに対してはうまく機能しています。
"1070052", "赤坂(次のビルを除く)"
"1076090", "赤坂赤坂アークヒルズ・アーク森ビル(地階・階層不明)"
"1076001", "赤坂赤坂アークヒルズ・アーク森ビル(1階)"
・・・
が、↓こんなやつだとうまくいきませんね。
"1050014", "芝(1~3丁目)"
"1080014", "芝(4、5丁目)"
芝1丁目、芝2丁目、・・・、芝5丁目を全部別の行にしてくれていたとしたら、まだプログラム側でなんとかしようなんて気も起きますが、「1~3」とか「4、5」とか、機械的に処理することを全く考慮していないデータなんですもの。今後の日本郵便のがんばりに期待したいところですね。(上から目線)
ということで、こんな間違いだらけのマップじゃ使えないよってかたは、国際航業さんから買ってください。
コメント
コメントを投稿