2016.09.23
地理データをRで可視化:コロプレス図
コロプレス図はエリア別の集計データ元に地図を色塗り分けして表現する手法です。地理データをRで可視化する基本の表現方法であるコロプレス図について各種パッケージの使用法、見え方を比較します。
テーマ
地理データをRで可視化する方法について考えます。
中でも基本の表現方法であるコロプレス図について各種パッケージの使用法、見え方を比較します。
0 コロプレス図について
コロプレス図は、エリア別の集計データ元に地図を色塗り分けして表現する手法です。集計データに色を割り当てるため、“どのように色を割り当てるか”の方針を決める必要があります。
具体的には、以下のような方針が挙げられます。
・ 連続値に対して割り当て
・ 階級区分に対して割り当て
・ 等分法
- 等間隔法
- 最適区分法
- 各種クラスター分析による階級 ...and more
連続値として扱う場合は、集計データにそのまま連続色を割り当てますが階級区分を用いる場合は、集計データをある基準をもってグループに分割し、分割したグループに対して色を割り当てます。
たとえば、身近な例として等間隔法を挙げると
各地域別チョコレート消費額が最小値1000円、最大値7999円の集計データに対して "1000円ごとに"分割して"1000~1999円","2000~2999円",・・・,"7000~7999"の7階級に分割し, その階級それぞれに色を割り当てる
となります。
1 準備
1.1 データ
描画のために、今回は札幌市中央区の小地域データを使用します。下記e-statのページから、対象となる地域・調査年度・形式などを選択します。
今回はデータをシェープファイルとしてダウンロードしました。
e-stat(統計GIS) 国勢調査/平成22年度/男女別人口総数及び世帯総数/世界測地系
ダウンロードしたシェープファイルを読み込みます。
今回のようなデータは緯度経度情報を元にエリアがポリゴン化されたもの(Polygon)+人口・世帯などの集計情報(Data)の構成となっており、R言語ではSpatialPoygonDataframeと呼ばれる形式で取り扱うことができます。
file<-"C:/Users/kudo_w/SkyDrive/Document/1_manaty/A002005212010DDSWC01101/h22ka01101.shp"
#データの準備
pj<-CRS("+proj=longlat +datum=WGS84") #緯度経度なので+projにlonglatを指定
Chuo<- maptools::readShapePoly(file,proj4string = pj)
#SpatialPolygonDataframeとして読み込み
読み込まれたデータを確認します。
SpatialPoygonDataframeは「data.frame」と「Polygon」を合わせたリスト構造です。
データ部は@data、polygon部は@polygonでアクセスすることができます。
#データ部の下から10レコードを表示 Chuo@data %>% tail(.) %>% kable()
| AREA | PERIMETER | H22KA01_ | H22KA01_ID | KEN | CITY | KEN_NAME | SITYO_NAME | GST_NAME | CSS_NAME | HCODE | KIHON1 | DUMMY1 | KIHON2 | KEYCODE1 | KEYCODE2 | AREA_MAX_F | KIGO_D | N_KEN | N_CITY | N_C1 | KIGO_E | KIGO_I | TATE | DIR | HIGHT | JIKAKU | NMOJI | MOJI | SEQ_NO2 | KSUM | CSUM | JINKO | SETAI | X_CODE | Y_CODE | KCODE1 | KEY_CODE | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 984 | 25773.480 | 793.756 | 17045 | 17044 | 01 | 101 | 北海道 | 石狩振興局 | 札幌市 | 中央区 | 8101 | 6929 | - | 12 | 101692912 | 101692912 | M | NA | NA | NA | 0 | NA | NA | 0 | 0 | 50 | 10 | 10 | 南二十九条西12丁目 | 17044 | 2 | 2 | 4 | 2 | 141.3443 | 43.02097 | 6929-12 | 01101692912 |
| 985 | 26182.490 | 647.058 | 17052 | 17051 | 01 | 101 | 北海道 | 石狩振興局 | 札幌市 | 中央区 | 8101 | 6929 | - | 10 | 101692910 | 101692910 | M | NA | NA | NA | 0 | NA | NA | 0 | 0 | 50 | 10 | 10 | 南二十九条西10丁目 | 17051 | 8 | 8 | 295 | 120 | 141.3471 | 43.02086 | 6929-10 | 01101692910 |
| 986 | 17949.160 | 545.108 | 17057 | 17056 | 01 | 101 | 北海道 | 石狩振興局 | 札幌市 | 中央区 | 8101 | 6929 | - | 09 | 101692909 | 101692909 | M | NA | NA | NA | 0 | NA | NA | 0 | 0 | 50 | 10 | 9 | 南二十九条西9丁目 | 17056 | 2 | 2 | 77 | 26 | 141.3488 | 43.02091 | 6929-09 | 01101692909 |
| 987 | 9581.004 | 514.771 | 17108 | 17107 | 01 | 101 | 北海道 | 石狩振興局 | 札幌市 | 中央区 | 8101 | 6930 | - | 09 | 101693009 | 101693009 | M | NA | NA | NA | 0 | NA | NA | 0 | 0 | 50 | 10 | 8 | 南三十条西9丁目 | 17107 | 1 | 1 | 20 | 12 | 141.3486 | 43.01963 | 6930-09 | 01101693009 |
| 988 | 38896.330 | 830.011 | 17110 | 17109 | 01 | 101 | 北海道 | 石狩振興局 | 札幌市 | 中央区 | 8101 | 6930 | - | 10 | 101693010 | 101693010 | M | NA | NA | NA | 0 | NA | NA | 0 | 0 | 50 | 10 | 9 | 南三十条西10丁目 | 17109 | 3 | 3 | 69 | 35 | 141.3473 | 43.01900 | 6930-10 | 01101693010 |
| 989 | 31386.130 | 845.620 | 17131 | 17130 | 01 | 101 | 北海道 | 石狩振興局 | 札幌市 | 中央区 | 8101 | 6930 | - | 11 | 101693011 | 101693011 | M | NA | NA | NA | 0 | NA | NA | 0 | 0 | 50 | 10 | 9 | 南三十条西11丁目 | 17130 | 3 | 3 | 53 | 36 | 141.3455 | 43.01863 | 6930-11 | 01101693011 |
#ポリゴン部一つ目を表示 Chuo@polygons[1]
[[1]]
An object of class "Polygons"
Slot "Polygons":
[[1]]
An object of class "Polygon"
Slot "labpt":
[1] 141.32838 43.08521
Slot "area":
[1] 3.785358e-06
Slot "hole":
[1] FALSE
Slot "ringDir":
[1] 1
Slot "coords":
[,1] [,2]
[1,] 141.3281 43.08649
[2,] 141.3285 43.08619
[3,] 141.3289 43.08591
[4,] 141.3293 43.08563
[5,] 141.3299 43.08519
[6,] 141.3299 43.08519
[7,] 141.3299 43.08505
[8,] 141.3299 43.08463
[9,] 141.3291 43.08455
[10,] 141.3286 43.08450
[11,] 141.3281 43.08445
[12,] 141.3277 43.08440
[13,] 141.3272 43.08434
[14,] 141.3271 43.08480
[15,] 141.3271 43.08511
[16,] 141.3271 43.08517
[17,] 141.3271 43.08527
[18,] 141.3274 43.08561
[19,] 141.3274 43.08563
[20,] 141.3278 43.08615
[21,] 141.3279 43.08623
[22,] 141.3281 43.08649
Slot "plotOrder":
[1] 1
Slot "labpt":
[1] 141.32838 43.08521
Slot "ID":
[1] "0"
Slot "area":
[1] 3.785358e-06
元データをそのまま展開するとデータ部に含まれる変数が24と多く、表形式で把握するには見にくいので、今回使用する以下のものに整理してデータ部に入れ直します。
・ MOJI = 各小地域名称
・ lng = 小地域代表点(X)
・ lat = 小地域代表点(Y)
・ KEY_CODE = ユニークキー
・ JINKO = 各小地域内人口
・ SETAI = 各小地域内世帯数
#lng,latに代表点を格納
Chuo.d<- data.frame(MOJI=Chuo@data$MOJI,
X_CODE=Chuo@data$X_CODE,
Y_CODE=Chuo@data$Y_CODE,
KEY_CODE=Chuo@data$KEY_CODE,
JINKO=Chuo@data$JINKO,
SETAI=Chuo@data$SETAI)
Chuo@data <- Chuo.d
Chuo@data %>% tail(.) %>% kable(.,format="markdown")
| MOJI | X_CODE | Y_CODE | KEY_CODE | JINKO | SETAI | |
|---|---|---|---|---|---|---|
| 985 | 南二十九条西12丁目 | 141.3443 | 43.02097 | 01101692912 | 4 | 2 |
| 986 | 南二十九条西10丁目 | 141.3471 | 43.02086 | 01101692910 | 295 | 120 |
| 987 | 南二十九条西9丁目 | 141.3488 | 43.02091 | 01101692909 | 77 | 26 |
| 988 | 南三十条西9丁目 | 141.3486 | 43.01963 | 01101693009 | 20 | 12 |
| 989 | 南三十条西10丁目 | 141.3473 | 43.01900 | 01101693010 | 69 | 35 |
| 990 | 南三十条西11丁目 | 141.3455 | 43.01863 | 01101693011 | 53 | 36 |
見通しがよくなったのでこれを描画していきます。
2 描画
2.0 描画概要
「1 準備」で用意したSpatialPoygonDataframeを基に、以下4つのパッケージを使用してコロプレス図を描画していきます。
手元のRにインストールされていない場合は、install.packages関数でインストール可能です。
・ ggplot2パッケージ
・ leafletパッケージ
・ plotパッケージ
・ choroplethrパッケージ
2.1 ggplot2パッケージによるコロプレス図
2.1.1 連続値
まずggplot2パッケージで描画していきます。
現在データはSpatialPoygonDataframeとなっているので、前処理としてfortify()を使用してdata.frameへの変換と整形を行いました。
#フォントファミリーの指定
#メイリオを"MEI"というfamilyとして定義
windowsFonts(MEI=windowsFont("Meiryo"))
#Polygon部の情報をdata.frameに変換
map_df<-ggplot2::fortify(Chuo)
Regions defined for each Polygons
map_df %>% head(.,10) %>% knitr::kable(.,format="markdown")
| long | lat | order | hole | piece | id | group |
|---|---|---|---|---|---|---|
| 141.3281 | 43.08649 | 1 | FALSE | 1 | 0 | 0.1 |
| 141.3285 | 43.08619 | 2 | FALSE | 1 | 0 | 0.1 |
| 141.3289 | 43.08591 | 3 | FALSE | 1 | 0 | 0.1 |
| 141.3293 | 43.08563 | 4 | FALSE | 1 | 0 | 0.1 |
| 141.3299 | 43.08519 | 5 | FALSE | 1 | 0 | 0.1 |
| 141.3299 | 43.08519 | 6 | FALSE | 1 | 0 | 0.1 |
| 141.3299 | 43.08505 | 7 | FALSE | 1 | 0 | 0.1 |
| 141.3299 | 43.08463 | 8 | FALSE | 1 | 0 | 0.1 |
| 141.3291 | 43.08455 | 9 | FALSE | 1 | 0 | 0.1 |
| 141.3286 | 43.08450 | 10 | FALSE | 1 | 0 | 0.1 |
#PolygonとData部の紐づけ map_join<-data.frame(id=as.numeric(row.names(Chuo@data))-1,Chuo@data) map_join$id<-as.character(map_join$id) map_df<- dplyr::left_join(map_df,map_join,by="id") #JOIN map_df %>% head(.,10) %>% knitr::kable(.,format="markdown")
| long | lat | order | hole | piece | id | group | MOJI | X_CODE | Y_CODE | KEY_CODE | JINKO | SETAI |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 141.3281 | 43.08649 | 1 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
| 141.3285 | 43.08619 | 2 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
| 141.3289 | 43.08591 | 3 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
| 141.3293 | 43.08563 | 4 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
| 141.3299 | 43.08519 | 5 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
| 141.3299 | 43.08519 | 6 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
| 141.3299 | 43.08505 | 7 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
| 141.3299 | 43.08463 | 8 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
| 141.3291 | 43.08455 | 9 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
| 141.3286 | 43.08450 | 10 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 |
#連続値
ggplot(map_df,aes(x=long,y=lat,group=group,fill=JINKO))+
geom_polygon()+
scale_fill_continuous(low="#fee0d2",high = "Red")+
ggtitle("札幌市中央区小地域別人口")+
theme(plot.title=element_text(size = rel(2),family = "MEI"))+
coord_map("polyconic")
2.1.2 離散値:等間隔法
次に離散化での表現です。
classInterval()のstyle引数で離散化の方法を指定することができます。指定できる離散化方法には、たとえば以下のようなものがあります。
・ styleなし > 等分法
・ style="equal" > 等間隔法
・ style="kmeans" > Kmeans法
・ style="hclust" > 階層クラスター法
ここでは等間隔法を用いて人口を8つのクラスに分割しました。
#等分法により8つのクラスに離散化,等間隔法
colcode<- map_df$JINKO %>%
classIntervals(.,n=8,style="equal") %>% findCols(.)
#離散化したクラス情報を保持
colname<-map_df$JINKO %>% classIntervals(.,n=8,style="equal")
#クラス値をデータに格納
map_df$class<-colcode %>% as.factor(.)
#凡例用のラベル作成
classlabel<-data.frame(row.names = 1:8,label=1:8)
for(i in 1:8){
classlabel[i,1]<-paste(round(colname$brks[i],0),round(colname$brks[i+1],0),sep="-")
}
#描画
ggplot(map_df,aes(x=long,y=lat,group=group,fill=class))+
geom_polygon()+
scale_fill_manual(values = brewer.pal(8,"Reds"),
name="人口(人)",
label=classlabel$label)+
guides(fill = guide_legend(reverse = TRUE))+
ggtitle("札幌市中央区小地域別人口")+
theme(plot.title=element_text(size = rel(2),family = "MEI"))+
coord_map("polyconic")
連続値の時と同じデータでの描画ですが、離散化すると色の差がはっきりとし「データに大きな差がある」ように見えます。
2.2 leafletによるコロプレス図
2.2.1 連続値
leafletパッケージでも同様の描画を行います。leaflet()にはカラーパレットを作成するための関数があるので今回は連続値、離散値共にそれを利用しました。
・ colorNumeric()・・・連続値のカラーパレットを作成
#カラーパレットの生成(連続値)
pal <- colorNumeric(
palette = "Reds",
domain = Chuo@data$JINKO
)
#コロプレス図だけ描くver
leaflet::leaflet(Chuo) %>%
addProviderTiles(provider = "OpenStreetMap.BlackAndWhite") %>%
addPolygons(stroke = FALSE,
fillOpacity = 0.9,
color = ~pal(Chuo.d$JINKO)) %>%
addLegend(position='bottomright', pal=pal, values=~JINKO,title="人口(人)")
2.2.2 離散化:等分法
今度は離散化用のパレットを用いて描画します。colorQuantile()で等分法による離散化ができます。
・ colorQuantile()・・・離散のカラーパレットを作成
#カラーパレットの生成(離散値、等分法)
pal_q<-leaflet::colorQuantile(
palette = "Reds",
domain = Chuo@data$JINKO,
n=8
)
leaflet::leaflet(Chuo) %>%
addProviderTiles(provider = "OpenStreetMap.BlackAndWhite") %>%
addPolygons(stroke = FALSE,
fillOpacity = 0.9,
color = ~pal_q(Chuo.d$JINKO)) %>%
addLegend(position='bottomright', pal=pal_q, values=~JINKO,title="人口(%)",labels=~classlabel$label)
2.3 plotによるコロプレス図
2.3.1 連続値
基本となるplot()でもコロプレス図が作図できます。今まで見てきた描画に比べて前処理や引数がかなり少なく済むのでとりあえず描きたい、という場合に便利です。
#連続値のカラーパレットを生成
pal_sp<-smoothPalette(as.numeric(Chuo@data$JINKO),pal="Reds")
#描画
plot(Chuo,col=pal_sp)
legend("bottomright",legend=c(max(Chuo@data$JINKO),min(Chuo@data$JINKO)),
fill = c(min(pal_sp),max(pal_sp)))
2.3.2 離散値:k-means法
離散化でも描いてみます。
クラスタリング一種であるk-means法での離散化を行いました。 classIntervals()のstyle引数を指定することで、任意の方法での離散化が可能です。
#k-means法で離散化しカラーパレット作成
colcode_km<- Chuo@data$JINKO %>%
classIntervals(.,n=8,style="kmeans") %>% findColours(.,pal=brewer.pal(8,"Reds"))
#離散化したクラス情報を保持
class_km<-Chuo@data$JINKO %>% classIntervals(.,n=8,style="kmeans")
plot(Chuo,col=colcode_km)
legend("bottomright",legend=str_sub(gsub(",","~",names(attr(colcode_km,"table"))),start = 2,end = -2),
fill = attr(colcode_km,"palette"))
結果としては等間隔よりも等分法に近い図になりました。
2.4 choroplethrによるコロプレス図
2.4.1 連続値
最後にchoroplethrパッケージによる作図です。
今回のように自分独自の地図を描くにはchoroplethrを継承したR6クラス、そしてそれを利用した関数を作成します。choroplethrではggplot2同様data.frameを引数にとるので、若干の前処理が必要です。
今回の人口のような対象となるデータは
・ region = 地域ごとのユニーク キー
・ value = 人口
となるように整形する必要があります。
#region列を追加 map_dfc <- map_df %>% mutate(region=as.character(id)) kable(head(map_dfc,5),format="markdown")
| long | lat | order | hole | piece | id | group | MOJI | X_CODE | Y_CODE | KEY_CODE | JINKO | SETAI | class | region |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 141.3281 | 43.08649 | 1 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 | 1 | 0 |
| 141.3285 | 43.08619 | 2 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 | 1 | 0 |
| 141.3289 | 43.08591 | 3 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 | 1 | 0 |
| 141.3293 | 43.08563 | 4 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 | 1 | 0 |
| 141.3299 | 43.08519 | 5 | FALSE | 1 | 0 | 0.1 | 北二十二条西 | 141.3284 | 43.08521 | 011017922 | 29 | 1 | 1 | 0 |
#R6クラスをつくる
ChuoChoro <- R6Class(
"ChuoChoropleth",
inherit = choroplethr:::Choropleth,
public = list(
initialize = function(user.df)
{
super$initialize(map_dfc, user.df)
}
)
)
#表現するデータを整形
Chuo.dc<- data.frame(district=map_join$MOJI,
region=as.character(map_join$id),
value=as.numeric(map_join$JINKO))
kable(head(Chuo.dc,5),format="markdown")
| district | region | value |
|---|---|---|
| 北二十二条西 | 0 | 29 |
| 北二十一条西 | 1 | 518 |
| 北二十条西 | 2 | 206 |
| 北十六条西16丁目 | 3 | 0 |
| 北十八条西 | 4 | 136 |
#描画のための関数作成
SappoRo_choropleth<-function (df, title = "", legend = "", num_colors = 8, zoom = NULL)
{
c = ChuoChoro$new(df)
c$title = title
c$legend = legend
c$set_num_colors(num_colors)
c$set_zoom(zoom)
c$render()
}
#連続値で描画(num_colorsに1を指定)
p<-SappoRo_choropleth(Chuo.dc,title="札幌市中央区人口",legend="",num_colors = 1)
p+scale_fill_continuous(low="#fff5f0",high = "Red")
2.4.2 離散化:等分法
先ほど作った関数の引数を変えることで8クラスに分けます。(等分法)
#離散値(num_colorsに8を指定) p2<-SappoRo_choropleth(Chuo.dc,title="札幌市中央区人口",legend="",num_colors = 8) p2+scale_fill_brewer(palette = "Reds")
3 まとめ
このように、同じデータ、同じコロプレス図でも用いるパッケージによって様々な特徴があり、離散化の有無、方法でも見え方はかなり変わるので用途に合わせての選択が必要になってきます。
たまにはこのように描き並べてみるのも楽しいのでおすすめです。
使用パッケージ
No package
1 ggplot2
2 leaflet
3 dplyr
4 knitr
5 DiagrammeR
6 maptools
7 sp
8 rgdal
9 classInt
10 RColorBrewer
11 tagcloud
12 stringr
13 choroplethr
14 R6
参考文献
・ 地理空間データ分析 (Rで学ぶデータサイエンス 7)
・ Rグラフィックスクックブック
・ R言語徹底解説
・ Leaflet for R(https://rstudio.github.io/leaflet/)