読者です 読者をやめる 読者になる 読者になる

疫学と医療統計学と遺伝学と時々、大学院生活

疫学を専門とする大学院生の研究に関する備忘録的ページ。

R package 'maptools'と'spsurvey'を利用して、カラーマップを作成する。

R programming visualization

今回は様々な指標をカラーマップを使用して、地図上に見やすく可視化するコードについてまとめる。
(特別なパッケージとしてはmaptoolsとspsurveyが挙げられる。)
今回は事例として、愛知県の男女の平均寿命を地図上にプロットすることを最終目的とする。

概要

多くの場合、データは表として表され、非常に解釈が困難かつ地域データを用いるときには傾向も分かりにくい現状がある。その解決法として、カラースケールを利用して、地図上にプロットすることで地域間の傾向や関係を明らかにできる。

それぞれのパッケージの説明

maptools

ESRI形式の地理的なファイルデータを読み込み、加工するためのパッケージ。
PBSmapping,spatstat, maps, RArcInfo, Stata tmap, WinBUGS, Mondrianなどの地理データの形式に変換するコードもあり。

spsurvey

probability surveyに使用されるデザインや解析のためのアルゴリズムを提供する。

コード

1. 最初に必要なパッケージを読み込む。
(すでにインストール済みであることを前提とする。)

library(maptools)
library(spsurvey)
library(ggplot2)
library(RColorBrewer)

2. まず全国(都道府県レベル)の地図を描いてみる。地図データはあらかじめダウンロードする必要あり。

jpn <- readShapePoly("JPN_adm1.shp")
jpn.dbf <- read.dbf("JPN_adm1.dbf") 
par(cex=0.4)
plot(jpn, xlim=c(128,146), ylim=c(30,45), col="grey", border="white")

f:id:ryosukefujii0320:20160118191340p:plain


3. 市区町村レベルデータを読み込みます。あらかじめダウンロードしておく必要あり。

#Japan_ver72.shpは2014年8月にESRI JAPANによって更新されたshpファイルです。
jpn <- readShapePoly("Japan_ver72.shp")  #シェープファイル(xy座標の組み込まれたデータ)を読み込みます。
jpn.dbf <- read.dbf("Japan_ver72.dbf")  #dbfファイルを読み込みます。

4. 愛知県のデータをJCODE(それぞれの市区町村に割り当てられているコード)を元に選択する。

aichipoly <- jpn[which(jpn$JCODE=="23101"): which(jpn$JCODE=="23563"),]
aichidata <- jpn.dbf[1060:1128,]
aichicoord <- coordinates(aichipoly)

5. 簡単に地図を描きます。土地をグレーで、境界線を白に指定します。

plot(aichipoly, col="gray", border="white")

f:id:ryosukefujii0320:20160118191405p:plain


6. 次に、市区町村の上にそれぞれの名前(aichipoly$CITY_ENG)をプロットします。

plot(aichipoly)
text(aichicord[,1], aichicord[,2], aichipoly$CITY_ENG, cex=0.5, col="black")

f:id:ryosukefujii0320:20160118192027p:plain


7. ここで政府の統計資料(e-stat)から愛知県の統計データをダウンロードしたcsvを読みこむ。*1

aichicsv <- read.csv("StatChiiki_201507291321554.csv", sep=",", header=T, fileEncoding="cp932")
colnames(aichicsv)
## [1] "NAME"      "JCODE"     "total"     "X65over"   "agingrate" "LEmen"    
## [7] "LEwomen"   "hospital"  "doctor"

8. 次に、aichidataとaichicsvを結合する。その際にJCODEを共通の列として選択する。

aichidata <- merge(aichicsv, aichidata, sort=FALSE, by="JCODE")

9. 最後にカラーマップに必要な色分けを指定する。ここではLEmen(平均余命:男性)とLEwomen(平均余命:女性)を使用する。

DDm <- aichidata$LEmen
cutnum <- 8
classesm <-cut(DDm, seq(min(DDm), max(DDm), length=cutnum+1), include.lowest=TRUE)
DDw <- aichidata$LEwomen
classesw <-cut(DDw, seq(min(DDw), max(DDw), length=cutnum+1), include.lowest=TRUE)

10. 愛知県の男性と女性の平均余命についてカラーマップを作成する。

#男性
par(cex=0.9)
cols <- rev(heat.colors(cutnum))
plot(aichipoly, col=cols[classesm], xlab="", ylab="", axes=FALSE)
box()

#女性
par(cex=0.9)
cols <- rev(heat.colors(cutnum))
plot(aichipoly, col=cols[classesw], xlab="", ylab="", axes=FALSE)
box()

男性
f:id:ryosukefujii0320:20160118191435p:plain
女性
f:id:ryosukefujii0320:20160118191442p:plain

まとめ

このように可視化することで結果の解釈を手助けする可能性は大いにあるので、今後の研究に使える見込みがある。

データの出典

愛知県の男女の平均余命(平成22年市区町村別生命表

*1:データ名はそれぞれダウンロードした統計資料のファイル名を指定