# 空間自己回帰モデル

#

# :目的変数,:説明変数,:空間重み行列

# :回帰係数,:空間ラグ(自己相関)パラメータ

# コロンバス市の犯罪データの例(空間従属性の検定)

# oldcolにCOL.OLD(データフレーム)とCOl.nb(nbクラス)が格納

# CRIME : 1,000世帯辺りの犯罪率

# HOVAL : 家の評価額(単位は$1,000)

# INC : 世帯収入

library(spdep)

data(oldcol)

names(COL.OLD)

# 重回帰モデル

col.ols <- lm(CRIME~INC+HOVAL,data=COL.OLD)

summary(col.ols)

# spdepにあるlm.LMtests()を使って,残差に関する空間従属性の検定を行う

# LMlag : の検定

# LMerr : の検定

W <- nb2listw(COL.nb,style="W") # 区域ごとに標準化した空間重み行列

lm.LMtests(col.ols,listw=W,test=c("LMlag","LMerr"))

# 空間ARモデル

# 空間自己回帰モデルで,とした場合,つまり,

#

# 隣接地域の値でを説明する1次の自己回帰モデル

# Rでは,spdepパッケージのlagsarlm()を使う

# コロンバス市の犯罪データの例(空間ARモデル)

col.ar <- lagsarlm(CRIME~1,data=COL.OLD,listw=W,method="eigen")

summary(col.ar)

# 空間同時自己回帰(SAR)モデル

# 空間自己回帰モデルで,とした場合,つまり,

#

# Rでは,spdepパッケージのlagsarlm()を使う

# コロンバス市の犯罪データの例(SARモデル)

col.sar <- lagsarlm(CRIME~INC+HOVAL,data=COL.OLD,listw=W,method="eigen")

summary(col.sar)

# 空間誤差モデル(SEM; spatial error model)

# 空間自己回帰モデルで,とした場合,つまり,

#

# 重回帰モデルの誤差項に空間自己相関がある形となる

# Rでは,spdepパッケージのerrorsarlm()を使う

# コロンバス市の犯罪データの例(SEM)

col.sem <- errorsarlm(CRIME~INC+HOVAL,data=COL.OLD,listw=W,method="eigen")

summary(col.sem)

# 空間Durbinモデル(SDM)

#

# 目的変数と説明変数の空間自己相関を混合したモデル(SAR混合モデルとも呼ばれる)

# spdepパッケージのlagsarlm(...,type="mixed")を使う

# コロンバス市の犯罪データの例(SDM)

col.sdm <- lagsarlm(CRIME~INC+HOVAL,data=COL.OLD,listw=W,method="eigen",type="mixed")

summary(col.sdm)

# anova()によるsarlmクラスの比較

anova(col.ar,col.sar,col.sem,col.sdm,col.ols)

# 条件付き自己回帰(CAR)モデル

# 空間誤差が条件付き分布に基づいている

# Rでは,spdepパッケージのspautolm(...,family="CAR")を使う

# spautolm(...,family="SAR")とすればSARモデルの推定もできる

# コロンバス市の犯罪データの例(CARモデル)

col.car <- spautolm(CRIME~INC+HOVAL,data=COL.OLD,family="CAR",listw=W)

summary(col.car)

# 地理的加重回帰(GWR)モデル

# グローバル検定

# spgwrパッケージに実装されるグローバル検定

# LMZ.F1GWR.test()

# LMZ2.F1GWR.test()

# LMZ3.F1GWR.test() : 説明変数ごとの有意性を調べる

# BFC99.gwr.test()

# BFC02.gwr.test()

# コロンバス市の犯罪データの例(GWRモデル)

library(spgwr)

data(georgia)

f <- formula(PctBach~TotPop90+PctRural+PctEld+PctFB+PctPov+PctBlack)

g.gauss <- gwr.sel(f,gweight=gwr.Gauss,data=gSRDF) # バンド幅の選択

res <- gwr(f,gweight=gwr.Gauss,data=gSRDF,bandwidth=g.gauss) # GWRの実行

str(res,max.level=1) # resに含まれるオブジェクトの表示(第一階層まで)

class(res$SDF) # 回帰係数はSpatialPolygonsDataFrameクラスの属性値として保存されている

# 黒人人口割合(PctBlack)の回帰係数を地図上に視覚化

brks <- c(-0.025,0,0.01,0.025,0.075)

cols <- grey(5:2/6)

plot(res$SDF,col=cols[findInterval(res$SDF$PctBlack,brks,all.inside=T)])

legend("bottomleft",leglabs(brks),fill=cols,bty="n")

# GWRのグローバル検定

res <- gwr(f,gweight=gwr.Gauss,data=gSRDF,bandwidth=g.gauss,hatmatrix=T)

LMZ.F1GWR.test(res)

LMZ.F2GWR.test(res)

LMZ.F3GWR.test(res)

BFC99.gwr.test(res)

BFC02.gwr.test(res)