環境ガバナンスのデータサイエンス/空間モデリング特論　R演習5

演習の概要

リスク分析：粗率、相対危険度、ポアソン確率、相対リスクのベイズ推定

注意事項

慶應義塾大学SFCで開講している「環境ガバナンスのデータサイエンス」「空間モデリング特論」の授業履修者向けの演習用ページです。
必ずしも全てのバージョンのRやOSで動作確認を行っていません。この演習用ページを作成している段階では、R3.5.3を使っています。
Rの更新などにより、Rコードなどを予告無しに変更する場合があります。

データ分析の準備

データのダウンロード

演習で用いるデータはここからダウンロードしてください。
ここでは./直下にgisdataフォルダを作成し、setwd("gisdata")とディレクトリを指定しています

パッケージのインストール

install.packages()でパッケージをインストールし、library()でパッケージを呼び出す
ここでは以下のパッケージを用います

install.packages("spdep")
install.packages("sf")
install.packages("spatstat")
install.packages("tidyverse")
install.packages("kokudosuuchi")

library(spdep)
library(sf)
library(spatstat)
library(tidyverse)
library(kokudosuuchi)

以下の方法でパイプ %>% の優先順位を高める

needs::prioritize(magrittr)

データの読み込み

都道府県境界データ（ポリゴン）、都道府県庁所在地データ（ポイント）及び都道府県別属性データを読み込みます。

pref_pnt <- sf::st_read('pref_gov.shp')
jpn_pref <- sf::st_read('jpn_pref.shp')
hd06 <- readr::read_csv("hd06.csv")

都道府県境界ポリゴンデータに属性データhd06を結合します。

jpn_pref_hd06 <- dplyr::inner_join(jpn_pref, hd06, by=c("PREF_CODE"="PREF_CODE"))

リスク分析

粗率

確率地図を作成する

spdepのprobmap()で確率地図を作成する。粗率は以下のように計算する。

jpn_hd06_pm <- spdep::probmap(jpn_pref_hd06$HD06, jpn_pref_hd06$POPJ06/100)
summary(jpn_hd06_pm)

      raw            expCount        relRisk            pmap       
 Min.   : 95.15   Min.   :  822   Min.   : 69.45   Min.   :0.0000  
 1st Qu.:137.58   1st Qu.: 1581   1st Qu.:100.42   1st Qu.:0.5122  
 Median :152.10   Median : 2381   Median :111.02   Median :1.0000  
 Mean   :151.09   Mean   : 3677   Mean   :110.28   Mean   :0.7424  
 3rd Qu.:169.73   3rd Qu.: 3729   3rd Qu.:123.89   3rd Qu.:1.0000  
 Max.   :199.93   Max.   :16995   Max.   :145.93   Max.   :1.0000

主題図を作成する

まず、都道府県別心疾患死亡者数を地図上に可視化する

ggplot2::ggplot(data = jpn_pref_hd06, aes(fill=HD06)) + 
  geom_sf() + 
  scale_fill_distiller(palette="Reds", trans = "reverse") + 
  theme_bw() +
  labs(title = "Number of Death by Heart Disease")

粗率（ここでは、心疾患死亡者数/人口数x100）の計算結果はjpn_hd06_pmのrawに格納されている

ggplot2::ggplot(data = jpn_pref_hd06, aes(fill=jpn_hd06_pm$raw)) + 
  geom_sf() + 
  scale_fill_distiller(palette="Reds", trans = "reverse") + 
  theme_bw() +
  labs(title = "Raw Rate of Death by Heart Disease")

環境ガバナンスのデータサイエンス/空間モデリング特論 R演習5