library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.2 ✔ tibble 3.3.0
## ✔ lubridate 1.9.4 ✔ tidyr 1.3.1
## ✔ purrr 1.0.4
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(sf)
## Linking to GEOS 3.13.0, GDAL 3.8.5, PROJ 9.5.1; sf_use_s2() is TRUE
popl_file <- "data/monthly_mdp_mesh1km.csv"
mesh_file <- "data/17_ishikawa_1km.geojson"石川県広域データ連携基盤推進事業 公式note の記事「QGISでCSVとGeoJSONを結合する方法をハンズオンで解説！」をみて、Rで同様のことをやってみようと思ったので、やってみました。
データの準備
元記事と同様に、G空間情報センターから、全国の人流オープンデータ（1kmメッシュ、市区町村単位発地別）（国土交通省）と、都道府県別1kmメッシュ（ネオプランニングラボ株式会社）の、ともに石川県分をダウンロードしました。
人流オープンデータは、「1kmメッシュ別の滞在人口データ」から2021年10月分を使用することにして、そのCSVファイル（17/2021/10/monthly_mdp_mesh1km.csv）をdataフォルダに入れておきます。1kmメッシュデータのGeoJSONファイル（17_ishikawa_1km.geojson）も同様にdataフォルダに入れておきました。
tidyverseパッケージとsfパッケージも読み込んでおきました。
データ読み込み
人口データのCSVファイルはread_csv関数で読み込みます。col_types引数でデータ型を指定しています。コードは因子型、数値（整数）は整数型としました。
1kmメッシュのGeoJSONファイルはsfパッケージのread_sf関数で読み込みます。読み込んだデータのcode列を因子型に変換しています。
popl_data <- read_csv(popl_file, col_types = "fffiiiii")
mesh_data <- read_sf(mesh_file) |>
dplyr::mutate(code = factor(code))データ確認
読み込んだデータを確認します。
人口データの方です。変数（列）の詳細は「オープンデータ（人流データ）定義書」に説明があります。
head(popl_data)
## # A tibble: 6 × 8
## mesh1kmid prefcode citycode year month dayflag timezone population
## <fct> <fct> <fct> <int> <int> <int> <int> <int>
## 1 54365567 17 17201 2021 10 0 0 13
## 2 54365567 17 17201 2021 10 1 0 17
## 3 54365567 17 17201 2021 10 2 0 15
## 4 54365578 17 17201 2021 10 0 0 85
## 5 54365578 17 17201 2021 10 0 1 99
## 6 54365578 17 17201 2021 10 0 2 82メッシュデータです。
head(mesh_data)
## Simple feature collection with 6 features and 1 field
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: 136.7125 ymin: 36.35 xmax: 136.75 ymax: 36.375
## Geodetic CRS: WGS 84
## # A tibble: 6 × 2
## code geometry
## <fct> <POLYGON [°]>
## 1 54364529 ((136.7375 36.35, 136.75 36.35, 136.75 36.35833, 136.7375 36.35833, …
## 2 54364538 ((136.725 36.35833, 136.7375 36.35833, 136.7375 36.36667, 136.725 36…
## 3 54364539 ((136.7375 36.35833, 136.75 36.35833, 136.75 36.36667, 136.7375 36.3…
## 4 54364547 ((136.7125 36.36667, 136.725 36.36667, 136.725 36.375, 136.7125 36.3…
## 5 54364548 ((136.725 36.36667, 136.7375 36.36667, 136.7375 36.375, 136.725 36.3…
## 6 54364549 ((136.7375 36.36667, 136.75 36.36667, 136.75 36.375, 136.7375 36.375…データの結合
left_join関数を使って、メッシュデータ（mesh_data）に人口データ（popl_data）を結合させます。キーとするのはメッシュコードで、これはmesh_dataではcode、popl_dataではmesh1kmidという列名なので、join_by関数でこの関係を指定します。
df_join <- mesh_data |>
dplyr::left_join(popl_data, by = join_by(code == mesh1kmid))結果の確認と修正
結果を確認します。1kmメッシュごとの人口の要約を表示します。
summary(df_join$population)
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 10.0 37.0 191.0 795.6 727.5 28329.0 2609人口が10人未満だったところは欠損値になっています。 ここはそのままにしておきました。
結合したデータの出力
結合したデータから、平日（dayflag == 0）の終日（timezone == 2）の値を抽出して、GeoJSONファイルに出力します。
write_sf関数の場合は、ファイルが存在した場合には置換するのがデフォルトの動作とヘルプにはあるのですが、警告が出たりするので、出力ファイルが存在する場合には、あらかじめ削除するようにしました。
output_file <- "outputs/population-2021-10.geojson"
if (file.exists(output_file)) {
file.remove(output_file)
}
## [1] TRUE
df_join |>
dplyr::filter(dayflag == 0, timezone == 2) |>
write_sf(output_file)QGISでの表示
出力したGeoJSONファイルをQGISで読み込んで表示させてみました。2021年10月、石川県の1kmメッシュごとの平日・終日の滞在人口です。背景は地理院タイル淡色です。
