【Rで近接行列】エリア間の「近接関係」を調べる

この記事は、「Rではじめる地理空間データの統計解析入門」の第6章を読んで自分なりの理解を整理するために書いたものです。
間違った理解もあるかもしれませんので、正確な情報、詳しい説明は「[Rではじめる地理空間データの統計解析入門](https://www.kspub.co.jp/book/detail/5273036.html)」をお読みください。

# 空間相関
空間データが持つ性質の１つに「空間相関」というのがある。

あるエリアの情報（地価や人口、気温など）が周囲と似た傾向、もしくは正反対の傾向を持っているといった関係性のこと。
例えば、地価の高いエリアの近くのエリアはやはり地価が高くなるなど。（この場合は正の空間相関）
逆に、あるエリアの数値が高いが、その周囲のエリアの数値は低くなるのは負の空間相関（大型スーパーの売り上げは高く、その周囲の小売り商店の売り上げは低い、など）。
空間相関がない場合、エリアごとの値は空間的な規則性を持たず、ランダムな分布となる。



## 近さの定義
あるエリアと「近い」というのを定義する場合、以下のいずれかのパターンを考える。
- a：境界線を共有している（ルーク型）
- b：境界線または点を共有している（クイーン型）
- c：最近傍のk地点（最近隣kゾーン）
- d：一定の距離内


aとbはチェスの駒になぞらえてルーク型、クイーン型と呼ばれるらしい。
グリッド上に並んだエリアで考えると、a（ルーク型）はグリッドに対して、境界線を共有しているグリッドを「隣接している」とする。
b（クイーン型）は境界線を共有しているグリッドと、４隅の頂点を共有しているグリッドを「隣接している」とする。
それぞれイメージすると以下のような画像の状態で、グレーのエリアに対して、ピンク色のエリアが隣接していることになる。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/a76c3c8a-3a9c-4699-8963-a8e92a148a74.png)



# 近接行列
エリア同士が近いかどうかを行列で示したもの。
例えば、札幌市の各区は画像のように位置しており、それぞれのエリアがどのエリアと接しているかを1or0で示すと表のようになる。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/c9eb56ad-8090-4f11-847b-09373acdcdd8.png)

中央区は、北区、東区、白石区、豊平区、南区、西区と接していることになる。



# 近接行列を行基準化する
上で示した近接行列の行の和を1に揃えることを行基準化という。
例えば、中央区は6つのエリアと接しているので、それぞれのエリアの近接行列の数値として1を6等分した0.1666...という値が得られる。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/26d48352-3809-42d0-950c-b7c4c9e40367.png)

こうすることで、各エリアの隣接性の強さが評価できる（ということだと思う）。

# Rで近接行列を求める
これらのことをRで求めてみる。
今回はエリアのデータとして北海道の市区町村別のエリアで求めてみる。
データは国土数値情報の[行政区域データ](https://nlftp.mlit.go.jp/ksj/gml/datalist/KsjTmplt-N03-2024.html)（一部加工）を使用してみます。



## データの読み込み
```R
library(sf); library(spdep)

# ポリゴンデータを読み込む
hokkaido <- read_sf("hokkaido.shp")

# 投影法の定義
st_crs(hokkaido) <-6668

#投影変換
hokkaido <- st_transform(hokkaido,crs=6680)
```

## 近接行列（ルーク型、クイーン型）
ルーク型、クイーン型で近接行列を作成するには以下を実行する。
```R
nb1 <- poly2nb(hokkaido)
```
この時、デフォルトではクイーン型で判定する。ルーク型にするには`queen = FALSE`をつける。

この時nb1にはnb形式という近接情報を格納する形式のオブジェクトが入っている。
`nb1`を実行すると、近接の判定結果を見ることができる
```R
## nb1の結果
Neighbour list object:
Number of regions: 179 
Number of nonzero links: 904 
Percentage nonzero weights: 2.821385 
Average number of links: 5.050279 
2 regions with no links:
3, 52
4 disjoint connected subgraphs
```

- Number of regions：全ポリゴンの数。179ポリゴン存在する
- Number of nonzero links：全エリア間の近接ノードの総数
- Percentage nonzero weights：全接続のうち、実際に近接関係があるリンクが存在する割合、ということらしい。
- Average number of links：１つのエリアあたりの平均接続数
- regions with no links：接続しないエリアの数。今回は２エリア（今回の場合、礼文島と奥尻町）
- disjoint connected subgraphs：エリア全体を地域群に分類した時のグループの数。今回は大きく４つの地域群に分類される

## プロット
```R
coords <-st_coordinates(st_centroid(hokkaido))
plot(st_geometry(hokkaido), col="white", border="grey")
plot(nb1, coords,add=TRUE, col="red",cex=0.01,lwd=1.5)
```

プロットした結果はこのように各エリアと隣接するノードが伸びていることがわかる

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/67e8f215-9ea3-46ab-a94a-3127cd9dc2df.png)


## 近接行列（最近隣kゾーン）
ルーク型やクイーン型では、島など陸続きでなっていないエリアは隣接していない判定となる。
島などでエリアがくっついていなくても隣接していると判定するには最近隣kゾーンを使う。
```R
knn <- knearneigh(coords,4)
nb2 <- knn2nb(knn)

#プロット
plot(st_geometry(hokkaido),border="grey")
plot(nb2,coords,add = TRUE, col="red",cex=0.01,lwd=1.5)
```

この場合、最近隣の４つのエリアと隣接しているものとして近接行列を作成する。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/2c2a8834-6793-42be-b439-a06def234c80.png)

```R
knn <- knearneigh(coords,8)
nb3 <- knn2nb(knn)

#プロット
plot(st_geometry(hokkaido),border="grey")
plot(nb2,coords,add = TRUE, col="red",cex=0.01,lwd=1.5)
```
こうすると、最近隣の８つのエリアで近接行列を作成する。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/d8fe8b4c-c5a8-4d65-a5f2-418e4bb9bada.png)


## 近接行列（一定距離内のエリアを隣接しているとする場合）
```R
nb4 <- dnearneigh(coords,d1=0,d2=30000)

#プロット
plot(st_geometry(hokkaido),border="grey")
plot(nb4,coords, add=TRUE,col="red",cex=0.01,lwd=1.5)
```

この場合、最小距離を0km、最大距離を30kmとして近接行列を作成している。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/07470190-a403-484b-ac56-4a4e2bb2fcb7.png)


## 近接行列を行基準化する
```R
w1 <- nb2listw(nb1)
w1mat <- listw2mat(w1)
```
- `nb2listw`：nb形式の近接情報を行基準化する。この時どことも隣接していないエリアがある場合はエラーが出るので、そんな時は`zero.policy = TRUE`をつける。
- `listw2mat`:nb2listwで行基準化したものを表形式に変換する

`w1mat`を確認すると、エリアごとの行基準化された値が整理されている(画像は札幌市に隣接しているエリアのみ表示）。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/0c636b0e-33f7-4bfe-b22b-6ba3bcc42721.png)


## 距離減衰を加味した近接行列

```R
glist <- nbdists(nb1,coords)
glist <- lapply(glist,function(x) 1/x)
w1 <- nb2listw(nb1,glist=glist)
w1mat <- listw2mat(w1)
```
- 1行目で、各エリア間の距離を評価する
- 2行目で、距離減衰関数を評価している
- 3行目で、その評価結果をもとに近接行列を行基準化している

こうすることで、各エリア間の距離に応じて行基準化の数値を調整してくれる。つまり、距離が近いエリア間では遠いエリア間より数値が高くなっている。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/61492/590159a5-5b98-41b4-b233-83987cd9afd0.png)



# まとめ
- 隣接するエリア間である数値（気温や地価など）においt関係性が見られることを空間相関という
- 空間相関を考える際に、どれくらいのエリア間を「近い」とするのかを定義づける必要がある
- 「近い」という定義には、ルーク型やクイーン型、最近隣kゾーン、一定距離内　などがある
- それぞれのエリア間を「近い」か「そうでないか」の判定を行列で示したものを近接行列という
- その近接行列の行の和を1に揃えることを行基準化という
- それぞれの近さの定義における近接行列や行基準化はRで簡単に求めることができる