演習の概要

  • 空間隣接行列、空間的自己相関

注意事項

  • 慶應義塾大学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)

データの読み込み

神奈川県の市区町村境界データを国土数値情報から入手する

kokudosuuchi::getKSJURL("N03", prefCode = 14) %>%
  filter(year==2015) %>%  select(zipFileUrl) %>% 
  pull(1) -> url
tf <- tempfile(fileext = ".zip")
curl::curl_download(url, tf)
unzip(tf, exdir = getwd(), junkpaths = TRUE)
unlink(tf)
rm(tf)

神奈川県の行政境界データを読み込み、更に横浜市を抽出する。その際、同じ行政コードのポリゴンを統合する。

ad <- sf::st_read('N03-15_14_150101.shp')%>%
  dplyr::group_by(N03_007) %>%
  dplyr::summarise(geometry = st_union(geometry)) %>%
  dplyr::ungroup()
ad %>%
  filter(N03_007 %in% c(14101:14118)) -> y_ad

横浜市の行政境界を描画する。

y_ad %>%
  st_geometry() %>%
  plot(border="white", col="grey")

空間隣接行列

ドロネー三角網による空間離接性

y_ad_cen <- sf::st_centroid(y_ad)
y_ad_coord <- sf::st_coordinates(y_ad_cen)
y_ad.tri.nb <- spdep::tri2nb(y_ad_coord)

y_ad %>%
  st_geometry() %>%
  plot(border="white", col="grey")
plot(y_ad.tri.nb, y_ad_coord, add=TRUE)
st_centroid assumes attributes are constant over geometries of xst_centroid does not give correct centroids for longitude/latitude data

最近隣k地点による空間隣接性

k=4の場合

y_ad.knn4.nb <- spdep::knn2nb(knearneigh(y_ad_coord, k=4))
y_ad %>%
  st_geometry() %>%
  plot(border="white", col="grey")
plot(y_ad.knn4.nb, y_ad_coord, add=TRUE)

k=3の場合

y_ad.knn3.nb <- spdep::knn2nb(knearneigh(y_ad_coord, k=3))

y_ad %>%
  st_geometry() %>%
  plot(border="white", col="grey")
plot(y_ad.knn3.nb, y_ad_coord, add=TRUE)

k=3とk=4の違い

diffs <- diffnb(y_ad.knn3.nb, y_ad.knn4.nb)

y_ad %>%
  st_geometry() %>%
  plot(border="white", col="grey")
plot(y_ad.knn3.nb, y_ad_coord, add=TRUE)
plot(diffs, y_ad_coord, col="red", lty=2, lwd=2, add=TRUE)

距離により隣接関係を定義する

y_ad.r.nb <- spdep::dnearneigh(y_ad_coord, 0, 0.06)

y_ad %>%
  st_geometry() %>%
  plot(border="white", col="grey")
plot(y_ad.r.nb, y_ad_coord, add=TRUE)

面オブジェクトの隣接関係で隣接性を定義

y_ad.poly.nb <- spdep::poly2nb(as(y_ad, "Spatial"))

y_ad %>%
  st_geometry() %>%
  plot(border="white", col="grey")
plot(y_ad.poly.nb, y_ad_coord, add=TRUE)

空間的自己相関：事例１

次に、空間的自己相関を示す指標（Moran’s I）を計算します。

ここでは、関東圏の市区町村別住宅地平均地価に関する空間的自己相関を算出する例を示します。

データの準備

関東圏の市区町村境界データに、市区町村別住宅地平均地価データを結合します。

kanto <- sf::st_read('kanto_area.shp')
lph <- readr::read_csv("lph.csv")
kanto_lph <- dplyr::inner_join(kanto, lph, by=c("JCODE"="JCODE"))
ggplot2::ggplot(data = kanto, aes(fill=LPH)) + 
  geom_sf() + 
  scale_fill_distiller(palette="Blues", trans = "reverse") + 
  theme_bw() +
  labs(title = "Land Price of House in Kanto Region")

Moran’s Iを計算する

空間隣接行列の定義

ここでは、ドロネー三角網を使って空間隣接行列を定義する

kanto.tri.nb <- sf::st_coordinates(st_centroid(kanto)) %>%
                    spdep::tri2nb()

Moran’s Iの算出

moran.test()関数を使ってMoran’s Iを計算します。

空間重み付け行列のスタイルをstyle=で指定する。ここでは、style="W"とした場合の例を示します。

moran.test(kanto$LPH,nb2listw(kanto.tri.nb, style="W"))

空間的自己相関：事例２

事例２として、都道府県別死亡者数データの例を示します

データの準備

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

pref.pnt <- sf::st_read('pref_gov.shp')
jpn.pref <- sf::st_read('jpn_pref.shp')
jpn.pref <- readr::read_csv("pref_gov.txt")

Moran’s Iを計算する

空間隣接行列の定義

ここでは、ドロネー三角網を使って空間隣接行列を定義する

pref.tri.nb <- sf::st_coordinates(st_centroid(pref.pnt)) %>%
                    spdep::tri2nb()

Moran’s Iの算出

moran.test()関数を使ってMoran’s Iを計算します。

ここでは、diabetesに関するMoran’s Iを算出します。

moran.test(pref.pnt$diabetes, nb2listw (pref.tri.nb,style="W"))

Local Moran’s Iの算出

Local Moran’s Iを計算し、局地的な空間的自己相関を可視化しよう。

LMI1 <- localmoran(pref.pnt$diabetes, nb2listw(pref.tri.nb, style="W"))
pref.lm <- data.frame(cbind(LMI1[,1],
           (pref.pnt$diabetes-mean(pref.pnt$diabetes))/sd(pref.pnt$diabetes)),
           row.names=pref.pnt$KENCODE)
colnames(pref.lm) <- c("Ii","standardized diabetes")
pref.lm
# 図をプロットする
plot(pref.lm,xlab="Local Moran's I", ylab="死亡率（標準化済み）", pch=20,
cex=1, cex.axis=1, cex.lab=1)
lines(x=c(-1,2), y=c(0,0), col="grey", lwd=1)
lines(x=c(0,0), y=c(-3,4), col="grey", lwd=1)
# text(pref.lm,rownames(pref.lm), adj=1.2, cex=0.8)
text(1.55, 3.3, "徳島県", cex=0.7)
text(1.16, 1.63, "香川県", cex=0.7)
text(-0.32, 2.15, "大分県", cex=0.7)
text(-0.8, -1.45, "愛知県", cex=0.7)
text(-0.72, 1.8, "富山県", cex=0.7)
text(0.82, -0.57, "東京都", cex=0.7)
text(0.75, -1.91, "神奈川県", cex=0.7)

LS0tCnRpdGxlOiAi56m66ZaT44Oi44OH44Oq44Oz44KwUua8lOe/kjQiCiNhdXRob3I6ICJUb21veXVraSBGdXJ1dGFuaSIKI2RhdGU6ICIyMDE5LzMvMjkiCm91dHB1dDoKICBodG1sX25vdGVib29rOiAKICAgIHRvYzogdHJ1ZQogICAgdG9jX2Zsb2F0OiB0cnVlCiAgaHRtbF9kb2N1bWVudDogZGVmYXVsdAotLS0KCgpgYGB7ciBzZXR1cCwgaW5jbHVkZT1GQUxTRX0Ka25pdHI6Om9wdHNfY2h1bmskc2V0KGVjaG8gPSBUUlVFKQprbml0cjo6b3B0c19jaHVuayRzZXQoY2FjaGUgPSBUUlVFKQprbml0cjo6b3B0c19jaHVuayRzZXQoZXJyb3IgPSBUUlVFKQpgYGAKCiMjIOa8lOe/kuOBruamguimgQotIOepuumWk+mao+aOpeihjOWIl+OAgeepuumWk+eahOiHquW3seebuOmWogoKIyMjIOazqOaEj+S6i+mghQotIOaFtuaHiee+qeWhvuWkp+WtplNGQ+OBp+mWi+ism+OBl+OBpuOBhOOCi+OAjOepuumWk+ODouODh+ODquODs+OCsOOAjeOBruaOiOalreWxpeS/ruiAheWQkeOBkeOBrua8lOe/kueUqOODmuODvOOCuOOBp+OBmeOAggotIOW/heOBmuOBl+OCguWFqOOBpuOBruODkOODvOOCuOODp+ODs+OBrlLjgoRPU+OBp+WLleS9nOeiuuiqjeOCkuihjOOBo+OBpuOBhOOBvuOBm+OCk+OAguOBk+OBrua8lOe/kueUqOODmuODvOOCuOOCkuS9nOaIkOOBl+OBpuOBhOOCi+autemajuOBp+OBr+OAgVIzLjUuM+OCkuS9v+OBo+OBpuOBhOOBvuOBmeOAggotIFLjga7mm7TmlrDjgarjganjgavjgojjgorjgIFS44Kz44O844OJ44Gq44Gp44KS5LqI5ZGK54Sh44GX44Gr5aSJ5pu044GZ44KL5aC05ZCI44GM44GC44KK44G+44GZ44CCIAoKIyMg44OH44O844K/5YiG5p6Q44Gu5rqW5YKZCiMjIyDjg4fjg7zjgr/jga7jg4Djgqbjg7Pjg63jg7zjg4kKLSDmvJTnv5LjgafnlKjjgYTjgovjg4fjg7zjgr/jga9b44GT44GTXShodHRwOi8vd2ViLnNmYy5rZWlvLmFjLmpwL35tYXVuei9hc2FrdXJhX3NwL2FzYWt1cmFfc3BfZGF0YV8xMTAxLnppcCnjgYvjgonjg4Djgqbjg7Pjg63jg7zjg4njgZfjgabjgY/jgaDjgZXjgYTjgIIKLSDjgZPjgZPjgafjga9gLi9g55u05LiL44GrYGdpc2RhdGFg44OV44Kp44Or44OA44KS5L2c5oiQ44GX44CBYHNldHdkKCJnaXNkYXRhIilg44Go44OH44Kj44Os44Kv44OI44Oq44KS5oyH5a6a44GX44Gm44GE44G+44GZCgojIyMg44OR44OD44Kx44O844K444Gu44Kk44Oz44K544OI44O844OrCi0gYGluc3RhbGwucGFja2FnZXMoKWDjgafjg5Hjg4PjgrHjg7zjgrjjgpLjgqTjg7Pjgrnjg4jjg7zjg6vjgZfjgIFgbGlicmFyeSgpYOOBp+ODkeODg+OCseODvOOCuOOCkuWRvOOBs+WHuuOBmQotIOOBk+OBk+OBp+OBr+S7peS4i+OBruODkeODg+OCseODvOOCuOOCkueUqOOBhOOBvuOBmQoKYGBge3IgaW5zdGFsbCBwYWNrYWdlcywgZXZhbD1GQUxTRX0KaW5zdGFsbC5wYWNrYWdlcygic3BkZXAiKQppbnN0YWxsLnBhY2thZ2VzKCJzZiIpCmluc3RhbGwucGFja2FnZXMoInNwYXRzdGF0IikKaW5zdGFsbC5wYWNrYWdlcygidGlkeXZlcnNlIikKaW5zdGFsbC5wYWNrYWdlcygia29rdWRvc3V1Y2hpIikKYGBgCgpgYGB7ciBsaWJyYXJ5LCBldmFsPUZBTFNFfQpsaWJyYXJ5KHNwZGVwKQpsaWJyYXJ5KHNmKQpsaWJyYXJ5KHNwYXRzdGF0KQpsaWJyYXJ5KHRpZHl2ZXJzZSkKbGlicmFyeShrb2t1ZG9zdXVjaGkpCmBgYAoK5Lul5LiL44Gu5pa55rOV44Gn44OR44Kk44OXICU+JSDjga7lhKrlhYjpoIbkvY3jgpLpq5jjgoHjgosKCmBgYHtyIHBpcGVfcHJpb3JpdGl6ZSwgZXZhbD1GQUxTRX0KbmVlZHM6OnByaW9yaXRpemUobWFncml0dHIpCmBgYAoKIyMjIOODh+ODvOOCv+OBruiqreOBv+i+vOOBvwrnpZ7lpYjlt53nnIzjga7luILljLrnlLrmnZHlooPnlYzjg4fjg7zjgr/jgpLlm73lnJ/mlbDlgKTmg4XloLHjgYvjgonlhaXmiYvjgZnjgosKYGBge3IgeV9hZF9yZWFkMSwgZXZhbD1GQUxTRX0Ka29rdWRvc3V1Y2hpOjpnZXRLU0pVUkwoIk4wMyIsIHByZWZDb2RlID0gMTQpICU+JQogIGZpbHRlcih5ZWFyPT0yMDE1KSAlPiUgIHNlbGVjdCh6aXBGaWxlVXJsKSAlPiUgCiAgcHVsbCgxKSAtPiB1cmwKdGYgPC0gdGVtcGZpbGUoZmlsZWV4dCA9ICIuemlwIikKY3VybDo6Y3VybF9kb3dubG9hZCh1cmwsIHRmKQp1bnppcCh0ZiwgZXhkaXIgPSBnZXR3ZCgpLCBqdW5rcGF0aHMgPSBUUlVFKQp1bmxpbmsodGYpCnJtKHRmKQpgYGAKCuelnuWliOW3neecjOOBruihjOaUv+Wig+eVjOODh+ODvOOCv+OCkuiqreOBv+i+vOOBv+OAgeabtOOBq+aoqua1nOW4guOCkuaKveWHuuOBmeOCi+OAguOBneOBrumam+OAgeWQjOOBmOihjOaUv+OCs+ODvOODieOBruODneODquOCtOODs+OCkue1seWQiOOBmeOCi+OAggpgYGB7ciB5X2FkX3JlYWQyLTEsIGV2YWw9RkFMU0V9CmFkIDwtIHNmOjpzdF9yZWFkKCdOMDMtMTVfMTRfMTUwMTAxLnNocCcpJT4lCiAgZHBseXI6Omdyb3VwX2J5KE4wM18wMDcpICU+JQogIGRwbHlyOjpzdW1tYXJpc2UoZ2VvbWV0cnkgPSBzdF91bmlvbihnZW9tZXRyeSkpICU+JQogIGRwbHlyOjp1bmdyb3VwKCkKYWQgJT4lCiAgZmlsdGVyKE4wM18wMDcgJWluJSBjKDE0MTAxOjE0MTE4KSkgLT4geV9hZApgYGAKCuaoqua1nOW4guOBruihjOaUv+Wig+eVjOOCkuaPj+eUu+OBmeOCi+OAggpgYGB7ciB5X2FkX3Bsb3QxLCBldmFsPUZBTFNFfQp5X2FkICU+JQogIHN0X2dlb21ldHJ5KCkgJT4lCiAgcGxvdChib3JkZXI9IndoaXRlIiwgY29sPSJncmV5IikKYGBgCgpgYGB7ciB5X2FkX3Bsb3QyLCBlY2hvPUZBTFNFfQp5X2FkICU+JQogIHN0X2dlb21ldHJ5KCkgJT4lCiAgcGxvdChib3JkZXI9IndoaXRlIiwgY29sPSJncmV5IikKYGBgCgojIyDnqbrplpPpmqPmjqXooYzliJcKIyMjIOODieODreODjeODvOS4ieinkue2suOBq+OCiOOCi+epuumWk+mbouaOpeaApwpgYGB7ciBkcmF1bmV5MS0xLCBldmFsPUZBTFNFfQp5X2FkX2NlbiA8LSBzZjo6c3RfY2VudHJvaWQoeV9hZCkKeV9hZF9jb29yZCA8LSBzZjo6c3RfY29vcmRpbmF0ZXMoeV9hZF9jZW4pCnlfYWQudHJpLm5iIDwtIHNwZGVwOjp0cmkybmIoeV9hZF9jb29yZCkKCnlfYWQgJT4lCiAgc3RfZ2VvbWV0cnkoKSAlPiUKICBwbG90KGJvcmRlcj0id2hpdGUiLCBjb2w9ImdyZXkiKQpwbG90KHlfYWQudHJpLm5iLCB5X2FkX2Nvb3JkLCBhZGQ9VFJVRSkKYGBgCgpgYGB7ciBkcmF1bmV5MS0yLCBlY2hvPUZBTFNFLCBtZXNzYWdlPUZBTFNFLCB3YXJuaW5nPUZBTFNFfQp5X2FkX2NlbiA8LSBzZjo6c3RfY2VudHJvaWQoeV9hZCkKeV9hZF9jb29yZCA8LSBzZjo6c3RfY29vcmRpbmF0ZXMoeV9hZF9jZW4pCnlfYWQudHJpLm5iIDwtIHNwZGVwOjp0cmkybmIoeV9hZF9jb29yZCkKCnlfYWQgJT4lCiAgc3RfZ2VvbWV0cnkoKSAlPiUKICBwbG90KGJvcmRlcj0id2hpdGUiLCBjb2w9ImdyZXkiKQpwbG90KHlfYWQudHJpLm5iLCB5X2FkX2Nvb3JkLCBhZGQ9VFJVRSkKYGBgCgojIyMg5pyA6L+R6Zqja+WcsOeCueOBq+OCiOOCi+epuumWk+mao+aOpeaApwprPTTjga7loLTlkIgKYGBge3Iga25uMS0xLCBldmFsPUZBTFNFfQp5X2FkLmtubjQubmIgPC0gc3BkZXA6OmtubjJuYihrbmVhcm5laWdoKHlfYWRfY29vcmQsIGs9NCkpCnlfYWQgJT4lCiAgc3RfZ2VvbWV0cnkoKSAlPiUKICBwbG90KGJvcmRlcj0id2hpdGUiLCBjb2w9ImdyZXkiKQpwbG90KHlfYWQua25uNC5uYiwgeV9hZF9jb29yZCwgYWRkPVRSVUUpCmBgYAoKYGBge3Iga25uMS0yLCBlY2hvPUZBTFNFfQp5X2FkLmtubjQubmIgPC0gc3BkZXA6OmtubjJuYihrbmVhcm5laWdoKHlfYWRfY29vcmQsIGs9NCkpCnlfYWQgJT4lCiAgc3RfZ2VvbWV0cnkoKSAlPiUKICBwbG90KGJvcmRlcj0id2hpdGUiLCBjb2w9ImdyZXkiKQpwbG90KHlfYWQua25uNC5uYiwgeV9hZF9jb29yZCwgYWRkPVRSVUUpCmBgYAoKaz0z44Gu5aC05ZCICmBgYHtyIGtubjItMSwgZXZhbD1GQUxTRX0KeV9hZC5rbm4zLm5iIDwtIHNwZGVwOjprbm4ybmIoa25lYXJuZWlnaCh5X2FkX2Nvb3JkLCBrPTMpKQoKeV9hZCAlPiUKICBzdF9nZW9tZXRyeSgpICU+JQogIHBsb3QoYm9yZGVyPSJ3aGl0ZSIsIGNvbD0iZ3JleSIpCnBsb3QoeV9hZC5rbm4zLm5iLCB5X2FkX2Nvb3JkLCBhZGQ9VFJVRSkKYGBgCgpgYGB7ciBrbm4yLTIsIGVjaG89RkFMU0UsIG1lc3NhZ2U9RkFMU0UsIHdhcm5pbmc9RkFMU0V9CnlfYWQua25uMy5uYiA8LSBzcGRlcDo6a25uMm5iKGtuZWFybmVpZ2goeV9hZF9jb29yZCwgaz0zKSkKCnlfYWQgJT4lCiAgc3RfZ2VvbWV0cnkoKSAlPiUKICBwbG90KGJvcmRlcj0id2hpdGUiLCBjb2w9ImdyZXkiKQpwbG90KHlfYWQua25uMy5uYiwgeV9hZF9jb29yZCwgYWRkPVRSVUUpCmBgYAoKaz0z44Goaz0044Gu6YGV44GECmBgYHtyIGtubjMtMSwgZXZhbD1GQUxTRX0KZGlmZnMgPC0gZGlmZm5iKHlfYWQua25uMy5uYiwgeV9hZC5rbm40Lm5iKQoKeV9hZCAlPiUKICBzdF9nZW9tZXRyeSgpICU+JQogIHBsb3QoYm9yZGVyPSJ3aGl0ZSIsIGNvbD0iZ3JleSIpCnBsb3QoeV9hZC5rbm4zLm5iLCB5X2FkX2Nvb3JkLCBhZGQ9VFJVRSkKcGxvdChkaWZmcywgeV9hZF9jb29yZCwgY29sPSJyZWQiLCBsdHk9MiwgbHdkPTIsIGFkZD1UUlVFKQpgYGAKCmBgYHtyIGtubjMtMiwgZWNobz1GQUxTRSwgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRX0KZGlmZnMgPC0gZGlmZm5iKHlfYWQua25uMy5uYiwgeV9hZC5rbm40Lm5iKQoKeV9hZCAlPiUKICBzdF9nZW9tZXRyeSgpICU+JQogIHBsb3QoYm9yZGVyPSJ3aGl0ZSIsIGNvbD0iZ3JleSIpCnBsb3QoeV9hZC5rbm4zLm5iLCB5X2FkX2Nvb3JkLCBhZGQ9VFJVRSkKcGxvdChkaWZmcywgeV9hZF9jb29yZCwgY29sPSJyZWQiLCBsdHk9MiwgbHdkPTIsIGFkZD1UUlVFKQpgYGAKCiMjIyDot53pm6LjgavjgojjgorpmqPmjqXplqLkv4LjgpLlrprnvqnjgZnjgosKYGBge3IgZG5lYXJuZWlnaDEtMSwgZXZhbD1GQUxTRX0KeV9hZC5yLm5iIDwtIHNwZGVwOjpkbmVhcm5laWdoKHlfYWRfY29vcmQsIDAsIDAuMDYpCgp5X2FkICU+JQogIHN0X2dlb21ldHJ5KCkgJT4lCiAgcGxvdChib3JkZXI9IndoaXRlIiwgY29sPSJncmV5IikKcGxvdCh5X2FkLnIubmIsIHlfYWRfY29vcmQsIGFkZD1UUlVFKQpgYGAKCmBgYHtyIGRuZWFybmVpZ2gxLTIsIGVjaG89RkFMU0UsIG1lc3NhZ2U9RkFMU0UsIHdhcm5pbmc9RkFMU0V9CnlfYWQuci5uYiA8LSBzcGRlcDo6ZG5lYXJuZWlnaCh5X2FkX2Nvb3JkLCAwLCAwLjA2KQoKeV9hZCAlPiUKICBzdF9nZW9tZXRyeSgpICU+JQogIHBsb3QoYm9yZGVyPSJ3aGl0ZSIsIGNvbD0iZ3JleSIpCnBsb3QoeV9hZC5yLm5iLCB5X2FkX2Nvb3JkLCBhZGQ9VFJVRSkKYGBgCgojIyMg6Z2i44Kq44OW44K444Kn44Kv44OI44Gu6Zqj5o6l6Zai5L+C44Gn6Zqj5o6l5oCn44KS5a6a576pCmBgYHtyIHBvbHluZWlnaDEtMSwgZXZhbD1GQUxTRX0KeV9hZC5wb2x5Lm5iIDwtIHNwZGVwOjpwb2x5Mm5iKGFzKHlfYWQsICJTcGF0aWFsIikpCgp5X2FkICU+JQogIHN0X2dlb21ldHJ5KCkgJT4lCiAgcGxvdChib3JkZXI9IndoaXRlIiwgY29sPSJncmV5IikKcGxvdCh5X2FkLnBvbHkubmIsIHlfYWRfY29vcmQsIGFkZD1UUlVFKQpgYGAKCmBgYHtyIHBvbHluZWlnaDEtMiwgZWNobz1GQUxTRSwgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRX0KeV9hZC5wb2x5Lm5iIDwtIHNwZGVwOjpwb2x5Mm5iKGFzKHlfYWQsICJTcGF0aWFsIikpCgp5X2FkICU+JQogIHN0X2dlb21ldHJ5KCkgJT4lCiAgcGxvdChib3JkZXI9IndoaXRlIiwgY29sPSJncmV5IikKcGxvdCh5X2FkLnBvbHkubmIsIHlfYWRfY29vcmQsIGFkZD1UUlVFKQpgYGAKCiMjIOepuumWk+eahOiHquW3seebuOmWou+8muS6i+S+i++8kQrmrKHjgavjgIHnqbrplpPnmoToh6rlt7Hnm7jplqLjgpLnpLrjgZnmjIfmqJnvvIhNb3JhbidzIEnvvInjgpLoqIjnrpfjgZfjgb7jgZnjgIIKCuOBk+OBk+OBp+OBr+OAgemWouadseWcj+OBruW4guWMuueUuuadkeWIpeS9j+WuheWcsOW5s+Wdh+WcsOS+oeOBq+mWouOBmeOCi+epuumWk+eahOiHquW3seebuOmWouOCkueul+WHuuOBmeOCi+S+i+OCkuekuuOBl+OBvuOBmeOAggoKIyMjIOODh+ODvOOCv+OBrua6luWCmQrplqLmnbHlnI/jga7luILljLrnlLrmnZHlooPnlYzjg4fjg7zjgr/jgavjgIHluILljLrnlLrmnZHliKXkvY/lroXlnLDlubPlnYflnLDkvqHjg4fjg7zjgr/jgpLntZDlkIjjgZfjgb7jgZnjgIIKYGBge3Iga2FudG8sIGV2YWw9RkFMU0V9CmthbnRvIDwtIHNmOjpzdF9yZWFkKCdrYW50b19hcmVhLnNocCcpCmxwaCA8LSByZWFkcjo6cmVhZF9jc3YoImxwaC5jc3YiKQprYW50b19scGggPC0gZHBseXI6OmlubmVyX2pvaW4oa2FudG8sIGxwaCwgYnk9YygiSkNPREUiPSJKQ09ERSIpKQpgYGAKCmBgYHtyIGthbnRvX2xwaF9wbG90MS0xLCBldmFsPUZBTFNFfQpnZ3Bsb3QyOjpnZ3Bsb3QoZGF0YSA9IGthbnRvLCBhZXMoZmlsbD1MUEgpKSArIAogIGdlb21fc2YoKSArIAogIHNjYWxlX2ZpbGxfZGlzdGlsbGVyKHBhbGV0dGU9IkJsdWVzIiwgdHJhbnMgPSAicmV2ZXJzZSIpICsgCiAgdGhlbWVfYncoKSArCiAgbGFicyh0aXRsZSA9ICJMYW5kIFByaWNlIG9mIEhvdXNlIGluIEthbnRvIFJlZ2lvbiIpCmBgYAoKYGBge3Iga2FudG9fbHBoX3Bsb3QxLTIsIGVjaG89RkFMU0UsIG1lc3NhZ2U9RkFMU0UsIHdhcm5pbmc9RkFMU0V9CmdncGxvdDI6OmdncGxvdChkYXRhID0ga2FudG8sIGFlcyhmaWxsPUxQSCkpICsgCiAgZ2VvbV9zZigpICsgCiAgc2NhbGVfZmlsbF9kaXN0aWxsZXIocGFsZXR0ZT0iQmx1ZXMiLCB0cmFucyA9ICJyZXZlcnNlIikgKyAKICB0aGVtZV9idygpICsKICBsYWJzKHRpdGxlID0gIkxhbmQgUHJpY2Ugb2YgSG91c2UgaW4gS2FudG8gUmVnaW9uIikKYGBgCgojIyMgTW9yYW4ncyBJ44KS6KiI566X44GZ44KLCiMjIyMg56m66ZaT6Zqj5o6l6KGM5YiX44Gu5a6a576pCuOBk+OBk+OBp+OBr+OAgeODieODreODjeODvOS4ieinkue2suOCkuS9v+OBo+OBpuepuumWk+mao+aOpeihjOWIl+OCkuWumue+qeOBmeOCiwpgYGB7ciBrYW50b190cmkybmIsIGV2YWw9RkFMU0V9CmthbnRvLnRyaS5uYiA8LSBzZjo6c3RfY29vcmRpbmF0ZXMoc3RfY2VudHJvaWQoa2FudG8pKSAlPiUKICAgICAgICAgICAgICAgICAgICBzcGRlcDo6dHJpMm5iKCkKYGBgCgojIyMjIE1vcmFuJ3MgSeOBrueul+WHugpgbW9yYW4udGVzdCgpYOmWouaVsOOCkuS9v+OBo+OBpk1vcmFuJ3MgSeOCkuioiOeul+OBl+OBvuOBmeOAggoK56m66ZaT6YeN44G/5LuY44GR6KGM5YiX44Gu44K544K/44Kk44Or44KSYHN0eWxlPWDjgafmjIflrprjgZnjgovjgILjgZPjgZPjgafjga/jgIFgc3R5bGU9IlciYOOBqOOBl+OBn+WgtOWQiOOBruS+i+OCkuekuuOBl+OBvuOBmeOAggpgYGB7ciBrYW50b19tb3JhbjEtMSwgZXZhbD1GQUxTRX0KbW9yYW4udGVzdChrYW50byRMUEgsbmIybGlzdHcoa2FudG8udHJpLm5iLCBzdHlsZT0iVyIpKQpgYGAKCmBgYHtyIGthbnRvX21vcmFuMS0yLCBlY2hvPUZBTFNFfQptb3Jhbi50ZXN0KGthbnRvJExQSCxuYjJsaXN0dyhrYW50by50cmkubmIsIHN0eWxlPSJXIikpCmBgYAoKIyMg56m66ZaT55qE6Ieq5bex55u46Zai77ya5LqL5L6L77ySCuS6i+S+i++8kuOBqOOBl+OBpuOAgemDvemBk+W6nOecjOWIpeatu+S6oeiAheaVsOODh+ODvOOCv+OBruS+i+OCkuekuuOBl+OBvuOBmQoKIyMjIOODh+ODvOOCv+OBrua6luWCmQrpg73pgZPlupznnIzlooPnlYzjg4fjg7zjgr/vvIjjg53jg6rjgrTjg7PvvInjgIHpg73pgZPlupznnIzluoHmiYDlnKjlnLDjg4fjg7zjgr/vvIjjg53jgqTjg7Pjg4jvvInlj4rjgbPpg73pgZPlupznnIzliKXlsZ7mgKfjg4fjg7zjgr/jgpLoqq3jgb/ovrzjgb/jgb7jgZnjgIIKYGBge3IganBuX2RhdGEsIGV2YWw9RkFMU0V9CnByZWYucG50IDwtIHNmOjpzdF9yZWFkKCdwcmVmX2dvdi5zaHAnKQpqcG4ucHJlZiA8LSBzZjo6c3RfcmVhZCgnanBuX3ByZWYuc2hwJykKanBuLnByZWYgPC0gcmVhZHI6OnJlYWRfY3N2KCJwcmVmX2dvdi50eHQiKQpgYGAKCiMjIyBNb3JhbidzIEnjgpLoqIjnrpfjgZnjgosKIyMjIyDnqbrplpPpmqPmjqXooYzliJfjga7lrprnvqkK44GT44GT44Gn44Gv44CB44OJ44Ot44ON44O85LiJ6KeS57ay44KS5L2/44Gj44Gm56m66ZaT6Zqj5o6l6KGM5YiX44KS5a6a576p44GZ44KLCmBgYHtyIHByZWZfdHJpMm5iLCBldmFsPUZBTFNFfQpwcmVmLnRyaS5uYiA8LSBzZjo6c3RfY29vcmRpbmF0ZXMoc3RfY2VudHJvaWQocHJlZi5wbnQpKSAlPiUKICAgICAgICAgICAgICAgICAgICBzcGRlcDo6dHJpMm5iKCkKYGBgCgojIyMjIE1vcmFuJ3MgSeOBrueul+WHugpgbW9yYW4udGVzdCgpYOmWouaVsOOCkuS9v+OBo+OBpk1vcmFuJ3MgSeOCkuioiOeul+OBl+OBvuOBmeOAggoK44GT44GT44Gn44Gv44CBYGRpYWJldGVzYOOBq+mWouOBmeOCi01vcmFuJ3MgSeOCkueul+WHuuOBl+OBvuOBmeOAggpgYGB7ciBwcmVmX21vcmFuMS0xLCBldmFsPUZBTFNFfQptb3Jhbi50ZXN0KHByZWYucG50JGRpYWJldGVzLCBuYjJsaXN0dyAocHJlZi50cmkubmIsc3R5bGU9IlciKSkKYGBgCgpgYGB7ciBwcmVmX21vcmFuMSwgZWNobz1GQUxTRX0KbW9yYW4udGVzdChwcmVmLnBudCRkaWFiZXRlcywgbmIybGlzdHcgKHByZWYudHJpLm5iLHN0eWxlPSJXIikpCmBgYAoKIyMjIyBMb2NhbCBNb3JhbidzIEnjga7nrpflh7oKTG9jYWwgTW9yYW4ncyBJ44KS6KiI566X44GX44CB5bGA5Zyw55qE44Gq56m66ZaT55qE6Ieq5bex55u46Zai44KS5Y+v6KaW5YyW44GX44KI44GG44CCCgpgYGB7ciBwcmVmX0xNSTEsIGV2YWw9RkFMU0V9CkxNSTEgPC0gbG9jYWxtb3JhbihwcmVmLnBudCRkaWFiZXRlcywgbmIybGlzdHcocHJlZi50cmkubmIsIHN0eWxlPSJXIikpCnByZWYubG0gPC0gZGF0YS5mcmFtZShjYmluZChMTUkxWywxXSwKICAgICAgICAgICAocHJlZi5wbnQkZGlhYmV0ZXMtbWVhbihwcmVmLnBudCRkaWFiZXRlcykpL3NkKHByZWYucG50JGRpYWJldGVzKSksCiAgICAgICAgICAgcm93Lm5hbWVzPXByZWYucG50JEtFTkNPREUpCmNvbG5hbWVzKHByZWYubG0pIDwtIGMoIklpIiwic3RhbmRhcmRpemVkIGRpYWJldGVzIikKcHJlZi5sbQojIOWbs+OCkuODl+ODreODg+ODiOOBmeOCiwpwbG90KHByZWYubG0seGxhYj0iTG9jYWwgTW9yYW4ncyBJIiwgeWxhYj0i5q275Lqh546H77yI5qiZ5rqW5YyW5riI44G/77yJIiwgcGNoPTIwLApjZXg9MSwgY2V4LmF4aXM9MSwgY2V4LmxhYj0xKQpsaW5lcyh4PWMoLTEsMiksIHk9YygwLDApLCBjb2w9ImdyZXkiLCBsd2Q9MSkKbGluZXMoeD1jKDAsMCksIHk9YygtMyw0KSwgY29sPSJncmV5IiwgbHdkPTEpCiMgdGV4dChwcmVmLmxtLHJvd25hbWVzKHByZWYubG0pLCBhZGo9MS4yLCBjZXg9MC44KQp0ZXh0KDEuNTUsIDMuMywgIuW+s+WztuecjCIsIGNleD0wLjcpCnRleHQoMS4xNiwgMS42MywgIummmeW3neecjCIsIGNleD0wLjcpCnRleHQoLTAuMzIsIDIuMTUsICLlpKfliIbnnIwiLCBjZXg9MC43KQp0ZXh0KC0wLjgsIC0xLjQ1LCAi5oSb55+l55yMIiwgY2V4PTAuNykKdGV4dCgtMC43MiwgMS44LCAi5a+M5bGx55yMIiwgY2V4PTAuNykKdGV4dCgwLjgyLCAtMC41NywgIuadseS6rOmDvSIsIGNleD0wLjcpCnRleHQoMC43NSwgLTEuOTEsICLnpZ7lpYjlt53nnIwiLCBjZXg9MC43KQpgYGAKCmBgYHtyIHByZWZfTE1JMiwgZWNobz1GQUxTRSwgbWVzc2FnZT1GQUxTRSwgd2FybmluZz1GQUxTRX0KTE1JMSA8LSBsb2NhbG1vcmFuKHByZWYucG50JGRpYWJldGVzLCBuYjJsaXN0dyhwcmVmLnRyaS5uYiwgc3R5bGU9IlciKSkKcHJlZi5sbSA8LSBkYXRhLmZyYW1lKGNiaW5kKExNSTFbLDFdLAogICAgICAgICAgIChwcmVmLnBudCRkaWFiZXRlcy1tZWFuKHByZWYucG50JGRpYWJldGVzKSkvc2QocHJlZi5wbnQkZGlhYmV0ZXMpKSwKICAgICAgICAgICByb3cubmFtZXM9cHJlZi5wbnQkS0VOQ09ERSkKY29sbmFtZXMocHJlZi5sbSkgPC0gYygiSWkiLCJzdGFuZGFyZGl6ZWQgZGlhYmV0ZXMiKQpwcmVmLmxtCiMg5Zuz44KS44OX44Ot44OD44OI44GZ44KLCnBsb3QocHJlZi5sbSx4bGFiPSJMb2NhbCBNb3JhbidzIEkiLCB5bGFiPSLmrbvkuqHnjofvvIjmqJnmupbljJbmuIjjgb/vvIkiLCBwY2g9MjAsCmNleD0xLCBjZXguYXhpcz0xLCBjZXgubGFiPTEpCmxpbmVzKHg9YygtMSwyKSwgeT1jKDAsMCksIGNvbD0iZ3JleSIsIGx3ZD0xKQpsaW5lcyh4PWMoMCwwKSwgeT1jKC0zLDQpLCBjb2w9ImdyZXkiLCBsd2Q9MSkKIyB0ZXh0KHByZWYubG0scm93bmFtZXMocHJlZi5sbSksIGFkaj0xLjIsIGNleD0wLjgpCnRleHQoMS41NSwgMy4zLCAi5b6z5bO255yMIiwgY2V4PTAuNykKdGV4dCgxLjE2LCAxLjYzLCAi6aaZ5bed55yMIiwgY2V4PTAuNykKdGV4dCgtMC4zMiwgMi4xNSwgIuWkp+WIhuecjCIsIGNleD0wLjcpCnRleHQoLTAuOCwgLTEuNDUsICLmhJvnn6XnnIwiLCBjZXg9MC43KQp0ZXh0KC0wLjcyLCAxLjgsICLlr4zlsbHnnIwiLCBjZXg9MC43KQp0ZXh0KDAuODIsIC0wLjU3LCAi5p2x5Lqs6YO9IiwgY2V4PTAuNykKdGV4dCgwLjc1LCAtMS45MSwgIuelnuWliOW3neecjCIsIGNleD0wLjcpCmBgYAoK