RでGIS（ディゾルブ）

# ディゾルブとは
あるデータの地物を、特定の属性の同一となる値の地物ごとに結合する処理です。

イメージ図
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/587bd30b-0ee4-bf65-804a-6de49e795f43.png)



# 本編
ディゾルブしたいデータを読み込みます
```
# Shapefileの読み込み
library(sf)
shapefile_data <- st_read("B002005212020DDSWC01101/r2kb01101.shp")
```

今回はe-Statより境界データ（小地域）の北海道札幌市中央区のポリゴンを使用します。

https://www.e-stat.go.jp/gis/statmap-search?type=2

このようなポリゴンになっています。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/86ddca17-22ba-0bcc-22eb-decfe57037d7.png)


属性には`CITY_NAME`に札幌市中央区という同じ文字列が入っているので、それをキーとしてディゾルブしていきます。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/ad0f5a9d-7430-e306-f793-d6ff5cfd8698.png)
```
library(dplyr)

# 属性に基づいてグループ化し、ディゾルブ
dissolved_shapefile <- shapefile_data  %>%
  group_by("CITY_NAME")%>%
  summarise(geometry = st_union(geometry))
```
group_byでディゾルブする際のキーとなるカラムを指定します。
さらにdplyrのsummariseでジオメトリをグループごとにユニオンします。


処理結果は下図の通りで、全ての地物が１まとまりになりました。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/a8a4857c-65a9-b4be-da8d-511bfd2042ef.png)

属性テーブルの情報は削ぎ落とされる
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/e8b0c23a-ae55-9ddf-dd05-23956e308b77.png)



## ジオメトリエラーが起きた時
もし
```
! Loop 0 is not valid: Edge 4 is degenerate (duplicate vertex)
Run `rlang::last_trace()` to see where the error occurred.
```
といったエラーが出た場合、重複した頂点などがあり、ジオメトリエラーが起きているので、ディゾルブを実行する前に以下を実行してみてください。

```
# ジオメトリを修正
shapefile_data <- st_make_valid(shapefile_data)
```

## ディゾルブ時に他の属性の値を合計したりする場合
今回のデータの属性にはJINKOというカラムに各地区の人口が入力されていたので、ディゾルブ時にその人口の数値を合計したいと思います。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/895e92ab-db2f-86ce-c582-0d0de20e0a69.png)


```
# 属性に基づいてグループ化し、ジオメトリをディゾルブし、JINKOのカラムを合計
dissolved_shapefile <- shapefile_data %>%
  group_by("CITY_NAME") %>%
  summarise(
    geometry = st_union(geometry),
    total_JINKO = sum(JINKO, na.rm = TRUE)
  )
```

最初のコードから増えたのは`total_JINKO = sum(JINKO, na.rm = TRUE)`の部分でsummariseする際にJINKOのカラムをsum（合計）して、total_JINKOというカラムに結果を入力するようにしています。

この結果、処理実行後の属性テーブルは以下の通りtotal_JINKOというカラムにJINKOの合計値が入っています。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/ac092876-b4b5-71fc-2cde-285aadd63026.png)
札幌市中央区の人口は248,680人であることがわかります。

summariseの集計処理はsum（合計）の他に以下のような集計が可能です
- 合計 (sum): 指定されたカラムの合計
- 平均 (mean): 指定されたカラムの平均値
- 中央値 (median): 指定されたカラムの中央値
- 最大値 (max): 指定されたカラムの最大値
- 最小値 (min): 指定されたカラムの最小値
- 標準偏差 (sd): 指定されたカラムの標準偏差
- 分散 (var): 指定されたカラムの分散
- 行数のカウント (n): グループごとの行数をカウントする



それぞれ実行する際には以下のように書く。
```
dissolved_shapefile <- shapefile_data %>%
  group_by("CITY_NAME") %>%
  summarise(
    geometry = st_union(geometry),
    total_JINKO = sum(JINKO, na.rm = TRUE),
    average_JINKO = mean(JINKO, na.rm = TRUE),
    median_JINKO = median(JINKO, na.rm = TRUE),
    maximum_JINKO = max(JINKO, na.rm = TRUE),
    minimum_JINKO = min(JINKO, na.rm = TRUE),
    sd_JINKO = sd(JINKO, na.rm = TRUE),
    var_JINKO = var(JINKO, na.rm = TRUE),
    count_JINKO = n()
  )
```
データ数をカウントする`n()`のみ引数は使用しなくて大丈夫です。



## データのエクスポート
ディゾルブした結果をGeoPackageなどにエクスポートしたい場合は、以下のようにする。
```
output_path <- "path_to_save/hogehoge.gpkg"
# GeoPackage形式で保存
st_write(dissolved_shapefile, output_path, driver = "GPKG")
```