INPUTしたらOUTPUT!

忘れっぽいんでメモっとく

{flexdashboard}と{leaflet}でGISっぽいダッシュボードを作る

仕事で{flexdashboard} + {leaflet}( + {plotly})でGISっぽいダッシュボードを作ることがあったのでメモしとく。

以下参考にしたサイト


0. モチベーション

地図上にデータを可視化する際に

  • 都道府県単位だと粗い
  • 国勢調査の小地域単位だと細かすぎる
    • しかもキーが紐づかない

といった課題がある。粒度的には国土数値情報で公開されている市区町村ごとのポリゴンが良さそうだが、社内のデータとの組み合わせはどうやるんだっけ?というのが今回の動機。


1. 塗り分けデータのダウンロード

社内のデータは公開できないため、サンプルとしてestatから取得できる住民基本台帳人口移動報告(統計表ID0003157723)の転入超過数を元に塗り分けを行う。 データの取得は@yutannihilation氏謹製の{estatapi}パッケージ*1を使用する。

> library(estatapi)
> library(dplyr)

> appId <- "**********"

> # 統計表ID0003157723のデータを取得
> d <- estat_getStatsData(appId = appId, statsDataId = "0003157723")
> glimpse(d)
Observations: 93,093
Variables: 11
$ cat01_code                     <chr> "0", "0", "0", "0", "0", "0", "0", "0", "0", "0", "0", "0", "0", "0", "0...
$ 性別                           <chr> "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総...
$ cat02_code                     <chr> "000", "000", "000", "000", "000", "000", "000", "000", "000", "000", "0...
$ 年齢                           <chr> "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総数", "総...
$ cat03_code                     <chr> "001", "001", "001", "001", "001", "001", "001", "001", "001", "001", "0...
$ 都市間移動者数・その他         <chr> "他市町村からの転入者数", "他市町村からの転入者数", "他市町村からの転入者数", "他市町村からの転入者数", "他市町村からの転入者数", "他市町村からの...
$ cat04_code                     <chr> "61000", "61000", "61000", "61000", "61000", "61000", "61000", "61000", ...
$ 国籍                           <chr> "日本人移動者", "日本人移動者", "日本人移動者", "日本人移動者", "日本人移動者", "日本人移動者", "日本人移動者", "日本人...
$ area_code                      <chr> "00000", "00001", "00002", "01000", "01001", "01002", "01100", "01101", ...
$ 全国・都道府県・市区町村2016~ <chr> "全国", "全国市部", "全国郡部", "北海道", "北海道市部", "北海道郡部", "札幌市", "中央区", "北区", "東区", "白石区", "豊平区"...
$ value                          <dbl> 4880967, 4530526, 350441, 240769, 204729, 36040, 118336, 20274, 15580, 1...


性別、年齢、国籍などの項目が縦持ちのフォーマットとなっているため総数のレコードのみに絞り込む。ついでに東京のデータにだけにしておく。

> d <- d %>% 
+   filter(cat01_code == "0" & cat02_code == "000" & cat03_code == "003") %>% 
+   filter(grepl('^13', area_code)) %>% 
+   select(area_code, value)
> glimpse(d)
Observations: 71
Variables: 2
$ area_code <chr> "13000", "13001", "13002", "13100", "13101", "13102", "13103", "13104", "13105", "13106", "13...
$ value     <dbl> 74177, 74218, -41, 58207, 881, 4787, 3231, 1417, 2246, 1924, 2735, 2539, 3169, 1128, 4263, 58...
> d %>% head() %>% knitr::kable()


area_code value
13000 74177
13001 74218
13002 -41
13100 58207
13101 881
13102 4787


これで市区町村ごとの転入超過数のデータができた。


2. シェープファイルのダウンロードと読み込み

次に地図に描画するシェープファイルの取り込みを行う。シェープファイル国土数値情報 行政区域データの詳細 からダウンロードできる行政境界を使用する。


windowsで試した際はencodingを指定しないと文字化けとなった。*2


> library(rgdal)

> layers <- ogrListLayers("N03-170101_GML/N03-17_170101.shp")

> # macの場合はencodingの指定はしなくても文字化けしない
> shp <- readOGR("N03-170101_GML/N03-17_170101.shp", layer=layers[1], 
+                stringsAsFactors = FALSE)
OGR data source with driver: ESRI Shapefile 
Source: "N03-170101_GML/N03-17_170101.shp", layer: "N03-17_170101"
with 116024 features
It has 5 fields

> # windowsの場合はencodingを指定しないと文字化けする
> Encoding(layers[1]) <- "UTF-8"
> shp <- readOGR("N03-170101_GML/N03-17_170101.shp", layer=layers[1], 
+                stringsAsFactors = FALSE, encoding = "UTF-8")

> class(shp)
[1] "SpatialPolygonsDataFrame"
attr(,"package")
[1] "sp"


読み込んだシェープファイルはSpatialPolygonsDataFrameという形式で、以下のとおり5つのスロットをもっている。

> glimpse(shp)
Formal class 'SpatialPolygonsDataFrame' [package "sp"] with 5 slots
  ..@ data       :'data.frame':    116024 obs. of  5 variables:
  .. ..$ N03_001: chr [1:116024] "北海道" "北海道" "北海道" "北海道" ...
  .. ..$ N03_002: chr [1:116024] "石狩振興局" "石狩振興局" "石狩振興局" "石狩振興局" ...
  .. ..$ N03_003: chr [1:116024] "札幌市" "札幌市" "札幌市" "札幌市" ...
  .. ..$ N03_004: chr [1:116024] "中央区" "北区" "東区" "白石区" ...
  .. ..$ N03_007: chr [1:116024] "01101" "01102" "01103" "01104" ...
  ..@ polygons   :List of 116024
  .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
  .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
(中略)
  .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
  .. .. [list output truncated]
  ..@ plotOrder  : int [1:116024] 47228 9163 963 8466 7982 8685 30565 29963 13819 15864 ...
  ..@ bbox       : num [1:2, 1:2] 122.9 20.4 154 45.6
  .. ..- attr(*, "dimnames")=List of 2
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot


dataスロットのdata.frameに塗り分けで使用する項目を追加するとleaflet::addPolygons()でコロプレス図が描ける。

> # 東京のみに限定
> shp.sub <- subset(shp, shp@data$N03_007 %in% d$area_code)

> # マージで項目を追加
> shp.sub@data <- merge(shp.sub@data, d, by.x="N03_007", by.y="area_code")
> glimpse(shp.sub)
Formal class 'SpatialPolygonsDataFrame' [package "sp"] with 5 slots
  ..@ data       :'data.frame':    6086 obs. of  6 variables:
  .. ..$ N03_007: chr [1:6086] "13101" "13102" "13103" "13103" ...
  .. ..$ N03_001: chr [1:6086] "東京都" "東京都" "東京都" "東京都" ...
  .. ..$ N03_002: chr [1:6086] NA NA NA NA ...
  .. ..$ N03_003: chr [1:6086] "千代田区" "中央区" "港区" "港区" ...
  .. ..$ N03_004: chr [1:6086] NA NA NA NA ...
  .. ..$ value  : num [1:6086] 881 4787 3231 3231 3231 ...
  ..@ polygons   :List of 6086
  .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
  .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
(以下略)


3. leafletでプロット

あとはleafletでプロットするだけ。

> library(leaflet)

> # 色塗りのカラーパレット
> pal <- colorNumeric("Spectral", domain=shp.sub@data$value, reverse=TRUE)

> # マウスオーバー時の表示内容
> labels <- sprintf(
+   "<strong>%s</strong><br/>%g",
+   paste0(shp.sub@data$N03_001,
+          ifelse(is.na(shp.sub@data$N03_003), "", shp.sub@data$N03_003),
+          ifelse(is.na(shp.sub@data$N03_004), "", shp.sub@data$N03_004)),
+   shp.sub@data$value
+ ) %>% lapply(htmltools::HTML)

> # 地図にプロット
> shp.sub %>% 
+   leaflet() %>% 
+   setView(lat=35.71239, lng=139.4089, zoom=10) %>% 
+   addProviderTiles(providers$CartoDB.Positron) %>% 
+   addPolygons(fillOpacity = 0.5,
+               weight=1,
+               fillColor = ~pal(value),
+               label = labels,
+               labelOptions = labelOptions(
+                 style = list("font-weight" = "normal", padding = "3px 8px"),
+                 textsize = "15px",
+                 direction = "auto"),
+               highlight = highlightOptions(
+                 weight = 5,
+                 color = "#666",
+                 dashArray = "",
+                 fillOpacity = 0.7,
+                 bringToFront = TRUE)) %>% 
+   addLegend("bottomright", pal = pal, values = ~value,
+             title = "転入超過数"
+   )

f:id:tak95:20170614120935p:plain

東京都内では世田谷区に人が集まっているようですね。



以下のように{flexdashboard}と組み合わせることでGISっぽいダッシュボードを構築することができた。

{flexdashboard}と{leaflet}で簡易GIS

f:id:tak95:20170614170537p:plain


  • ダッシュボードを作るのは{flexdashboad} + {plotly} + {leaflet}が簡単でオススメ
    • supersetも良いけど地図への描画にMapboxのAPI Keyが必要となるのがネック