ねこすたっと

ねこの気持ちと統計について悩む筆者の備忘録的ページ。

空間データを可視化する(1)(sfパッケージ) [R]

はじめに

例えば、人口密度などの特性によって都道府県を塗り分けた図や、医療機関の位置を地図上にプロットした図など、地理空間情報を含んだデータをグラフにしたものを見る機会は多いですよね。 最近はExcelでも簡単に作成できるみたいです( マップ グラフを作成する Excel - Microsoft サポート

Rを使う方が痒い所に手が届くと思うので、今回はRのsfパッケージを使って、神戸市の校区境界と学校所在地をプロットしてみます。

GISデータを準備する

GISとは Geographic Information System の略で、地理情報システムと訳されます。いくつかの機関がGISデータを無償で提供してくれています。 地図と聞いて真っ先に浮かぶ国土地理院もGISデータを無償提供していますが、今回は提供しているデータの種類が豊富な国土数値情報(提供:国土交通省 国土政策局)を使います。

nlftp.mlit.go.jp

今回は以下のzipファイルをダウンロードして使ってみました。

  • 兵庫県の中学校区(2021年):A32-21_28_GML.zip
  • 兵庫県の小学校区(2021年):A27-21_28_GML.zip
  • 兵庫県の学校所在地(2021年)P29-21_28_GML.zip

データファイルの概要

GISデータは、位置や区域を表す形状データと、その点・区域の特徴を示す属性データから構成されます。形状データは以下の4種類があります。

  • 点データ:建物の位置
  • 線データ:河川や道路の形を表す線分データ
  • 面データ:行政や校区の境界を表す多角形データ(ポリゴン*1データともいう)
  • メッシュデータ:エリアを格子状に区切ったデータ(ラスタデータ*2ともいう)

ダウンロードされた圧縮ファイルを解凍して得られるフォルダの中には、いくつかファイルが入っています。必須のファイルは以下の3つです。

  • .shp:主に形状データを格納しているファイル
  • .shx:データの検索効率を向上するためのインデックス情報を格納するファイル
  • .dbf:属性データを格納しているファイル

データの読み込みと確認

まず必要なパッケージを読み込みます。

library(sf)
library(tidyverse)

次に、read_sf( )を使って.shpファイルを読み込みます。他の2つの必須ファイル(.shx, .dbf)も同じ場所に置いておきます。

MS_poly <- read_sf("A32-21_28.shp") # 中学校区境界
ES_poly <- read_sf("A27-21_28.shp") # 小学校区境界
S_point <- read_sf("P29-21_28.shp") # 学校所在地

glimpse( )で中学校区データ(A32-21_28.shp)の構造の概要を確認します。

> MS_poly %>% glimpse()
Rows: 287
Columns: 6
$ A32_001  <chr> "28100", "28100", "28100", "28100", "28100", "28100", "28100", "28100"$ A32_002  <chr> "神戸市", "神戸市", "神戸市", "神戸市", "神戸市", "神戸市", "神戸市",$ A32_003  <chr> "C128210000514", "C128210000373", "C128210000104", "C128210000337", "C…
$ A32_004  <chr> "伊川谷中学校", "井吹台中学校", "烏帽子中学校", "雲雀丘中学校", "塩屋…
$ A32_005  <chr> "神戸市西区伊川谷町上脇字鬼神山1005-2", "神戸市西区井吹台西町2-3", "神…
$ geometry <MULTIPOLYGON [°]> MULTIPOLYGON (((135.0291 34..., MULTIPOLYGON (((135.0296 …

ここでgeometryは形状データで、"MULTIPOLYGON"とあるので多角形データとわかります。 他の5つは属性データです。どの変数がどんな属性を表しているかは、国土数値情報のダウンロードサイトに書いてあります。

  • A32_001:中学校が存在する行政区域コード*3
  • A32_002:設置主体
  • A32_003:学校に設定された固有の学校コード
  • A32_004:学校の名称
  • A32_005:学校の所在地住所

小学校区データ(A27-21_28.shp)に保存されている属性も同じ内容ですが、変数名が異なります(例:A27_001)。

学校所在地データ(P29-21_28.shp)についても見てみます。

> S_point %>% glimpse()
Rows: 2,622
Columns: 8
$ P29_001  <chr> "28203", "28228", "28207", "28100", "28100", "28100", "28100", "28100"$ P29_002  <chr> "A128110000011", "A128110000020", "A128210000019", "A128210000028", "A…
$ P29_003  <int> 16011, 16011, 16011, 16011, 16011, 16011, 16011, 16011, 16011, 16011, …
$ P29_004  <chr> "神戸大学附属幼稚園", "兵庫教育大学附属幼稚園", "伊丹幼稚園ありおか分…
$ P29_005  <chr> "明石市山下町3-4", "加東市山国2013-4", "伊丹市伊丹7-1-30", "神戸市東灘…
$ P29_006  <int> 1, 1, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3…
$ P29_007  <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
$ geometry <POINT [°]> POINT (134.9986 34.64958), POINT (134.9808 34.91444), POINT (135…

今度は形状データが"POINT"、つまり点データになっていますね。属性は以下のとおりです。

  • P29_001:行政区域コード
  • P29_002:学校コード
  • P29_003:学校分類(16001=小学校、16002=中学校など)
  • P29_004:学校の名称
  • P29_005:学校の所在地
  • P29_006:管理者コード(1=国、2=都道府県、3=市区町村、4=民間、0=その他)
  • P29_007:休校区分(0=調査なし、1=開校中、2=休校中)

データの一部を抽出する

今回は神戸市だけプロットするので、兵庫県全体の中から該当する部分を抽出してきます。行政区域コードや設置主体で条件を設定します。

MS_poly_kobe <- MS_poly %>% 
  filter(A32_001 %in% c("28100", "28101", "28102", "28105", "28106", "28107", "28109", "28110", "28111"))

ES_poly_kobe <- ES_poly %>% 
  filter(A27_001 %in% c("28100", "28101", "28102", "28105", "28106", "28107", "28109", "28110", "28111"))

各データの説明ページにあるリンクから行政区域コード一覧が収められたExcelファイルを入手できます。

ここから、28100(神戸市), 28101(神戸市東灘区), ... , 28111(神戸市西区)を選んで指定しました。

学校所在地データについては、以下のコードで神戸市内の公立中学校のみ抽出しました。

MS_point_kobe <- S_point %>% 
  filter(P29_001 %in% c("28100", "28101", "28102", "28105", "28106", "28107", "28109", "28110", "28111")) %>% 
  filter(P29_003 %in% c(16002,16003, 16014)) %>% 
  filter(P29_006 == 3) %>%
  filter(P29_007 == 1)

指定した条件は以下のとおりです。

  • 行政区域コード(P29_001):前述のとおり
  • 学校分類コード(P29_003):16002=中学校, 16003=中等教育学校*4, 16014=義務教育学校*5
  • 管理者コード(P29_006):3=市区町村
  • 休校区分(P29_007):1=開校中

属性データを付け加える

提供されている属性以外の情報を加えて、地図を塗り分けたいということがあると思います。ここではランダムなカテゴリー変数(rnd_cat)とランダムな連続変数(rnd_cont)を追加しておいて、後で塗り分けに使いたいと思います。

set.seed(1234)

ES_poly_kobe <- ES_poly_kobe %>% 
  rowwise() %>% 
  mutate(rnd_cat = sample(c("A","B","C"),prob=c(0.2,0.3,0.5),size=1,replace=TRUE)) %>% 
  mutate(rnd_cont = rnorm(n=1,mean=0,sd=5))

MS_point_kobe <- MS_point_kobe %>% 
  rowwise() %>% 
  mutate(rnd_cat = sample(c("A","B","C"),prob=c(0.2,0.3,0.5),size=1,replace=TRUE)) %>% 
  mutate(rnd_cont = rnorm(n=1,mean=0,sd=5))

おわりに

  • 中学校区数79、中学校数84と乖離がありました。設置主体で条件付けても同じでした。一部は分校によるものでしたが、それ以外にも統廃合などの影響がありそうです。
  • 次回は実際にプロットしてみます。

参考資料

  • GISデータに関する基礎的な解説がわかりやすいです。

www.esrij.com

*1:「多角形」という意味

*2:線・面のように点を繋いで形状を表したものをベクタデータという

*3:都道府県コードと市区町村コードを組み合わせて定義されている

*4:中高一貫なので多分私立のみだが念のため

*5:小中一貫