空間分析  期末考

助教 杜承軒 2020.06.14

圖表解讀

甲、[10%]
※重點:兩模型「殘差」與空間自相關的關係。

圖中畫出OLS最小平方法(圓形與虛線)、SAM自迴歸模型(方形)的殘差Moran's I分布圖。
模型殘差代表未能解釋的部分,顯示出OLS沒能考慮空間的鄰近相關性,因此OLS的殘差之Moran's I有顯著不為0。且隨著延遲距離越大,空間自相關下降至接近0。
而SAM已考慮空間自相關,因此殘差之Moran's I接近0。

乙、[10%]
※重點:Moran散布圖的意義,以及兩條斜線的意涵。

動態Moran散布圖,此處「動態」是指放入不同觀察值會有不同結果。
在Moran散布圖分成四個象限,分別代表自己與鄰居和平均的關係,四個象限分別為:HH、LH、LL、HL。
圖中黑點是槓桿點(leverage points),指的是自變量因子的離群點。
在納入所有點來計算時,Moran's I為0.51 (實線);
但把離群點移除(這裡就是動態的概念)後,計算新的Moran's I為0.61 (虛線)。

實作題

圖資讀取

library(rgdal);library(GISTools);library(spdep);library(gstat)
setwd("D:/1082SA/Fin")
KH=readOGR("KH.shp",encoding='utf8')
KH=spTransform(KH,CRS("+init=epsg:3826"))
for(i in 3:6) KH@data[,i]=as.numeric(as.character(KH@data[,i]))
KH.pt=gCentroid(KH,byid=T)

Q1.罷免傾向度

Q1-1 [5%]
Moran’s I 數值[2%]+Moran 散布圖[3%]

KH$AD=KH$AGREE/KH$DISAGREE #罷免傾向度

KH.Qnb = poly2nb(KH) 
KH.Qw  = nb2listw(KH.Qnb)
KH.QM  = nb2mat(KH.Qnb)

M=moran.test(KH$AD,KH.Qw)
cat("Moran’s I =",M$estimate[1])
Moran’s I = 0.2502881
moran.plot(KH$AD,KH.Qw,labels = KH$TOWNNAME,xlab="罷免傾向度",ylab="鄰居罷免傾向度平均")

Q1-2 [10%]
扣分:沒有使用單尾(右尾)檢定(-2%),只找出正值正相關(HH) (-2%)

AD.lisa = localmoran(KH$AD,KH.Qw)
plot(KH,col='grey',main="LISA 地圖\n罷免傾向度顯著正相關")
plot(KH[AD.lisa[, 5]<.05,],col='red',add=T)
legend("bottomright",legend=c("顯著正相關","不顯著"), fill=c('red','grey'),bty="n",cex=0.7)

Q1-3 [10%]
計算[2%]+列出正確行政區[4*2%]

id = KH$AD<mean(KH$AD) & KH.QM%*%KH$AD>mean(KH.QM%*%KH$AD)
KH$TOWNNAME[id]%>%as.character
[1] "橋頭區" "路竹區" "茄萣區" "彌陀區"

Q2.投票率

Q2-1 [10%]
半變異元圖:x軸[2%]+y軸[2%]+點[3%]+線[3%]
x軸:注意cutoff。y軸:注意投票率要乘上100,結果大部分會落在0~20。
趨勢類似即可,會因為width設定而有擬和上的不同,width越小,y軸範圍會較大一些。

KH$VR=rowSums(KH@data[,4:6])/KH$TOTAL*100 #投票率

vgm=variogram(KH$VR~1,KH.pt,cutoff=15000) #預設width=1000 (cutoff/15) 
fit=fit.variogram(vgm, model = vgm(0,"Exp",9000,20))
plot(vgm,fit,main="投票率半變異元圖")

Q2-2 [5%]
扣分:忘記include.self(-2%),閾值設定錯誤(qnorm(.95)=1.645)(-2%)

KH.Dnb    = dnearneigh(KH.pt@coords,d1=0,d2=fit$range[2])
KH.Dnb.in = include.self(KH.Dnb)
KH.Dw.in  = nb2listw(KH.Dnb.in)

Gi=localG(KH$VR,KH.Dw.in)%>%as.vector

plot(KH,col='grey',main="投票熱區地圖")
plot(KH[Gi>=qnorm(.95),],col='red',add=T)
legend("bottomright",legend=c("投票熱區","不顯著"), fill=c('red','grey'),bty="n",cex=0.7)

Q2-3[10%]
數值計算正確才給分[10%]

vgm2=variogram(KH$VR~1,KH.pt,cutoff=10000,width=2000) 
vgm2$gamma[2]
[1] 6.388584

d=dist(KH.pt@coords)
z=dist(KH$VR)^2
mean(z[d>2000 & d<4000])/2
[1] 6.388584

Q3.鄰居同意票數差值

Q3[10%]
計算(設定knn、計算鄰居差)[4%]+三個行政區名稱和數值[6%]

KH.Knb  = knn2nb(knearneigh(KH.pt@coords, k=5)) 
KH.KB = nb2listw(KH.Knb,style='B')

KH$diff=unlist(lapply(KH.KB$neighbours, function(.) diff(range(KH$AGREE[.]))))
KH@data[order(KH$diff)[1:3],c('TOWNNAME','diff')]
   TOWNNAME diff
26   桃源區 2820
27 那瑪夏區 2820
19   旗山區 8066

Q4.與鄰居票數的差異

Q4-1[10%]
計算過程[4%]+兩個行政區名稱和相差數值[3*2%]

KH.QB = nb2listw(KH.Qnb,style='B')
KH.QM = nb2mat(KH.Qnb,style='B')
KH.QMN = ifelse(KH.QM==1,1,NA)

KH.max=c()
for(i in 1:38)  KH.max[i]=max(abs(KH.QMN[i,]*KH$AGREE-KH$AGREE[i]),na.rm=T)
Max=max(KH.max)
DISTs=KH$TOWNNAME[KH.max==Max]%>%as.character
cat(DISTs,"相差",Max,"票")
鹽埕區 三民區 相差 111204 票

Q4-2[10%]
重複排列計算[5%]+計算p-value3%+結論[2%]
P.S.此處p-value會是隨機的,但應該都會大於0.9

RAND=c()
for(i in 1:1000){
rs=sample(KH$AGREE)
RAND[i]=max(KH.QM*abs(matrix(rs,38,38,T)-matrix(rs,38,38,F)))
}

p.sim=sum(RAND>Max)/1000
p.norm=pnorm(Max,mean(RAND),sd(RAND),lower.tail = F)
cat("p=value 為",p.sim,"或",p.norm,",票數差異沒有顯著性大於隨機情形")
p=value 為 0.929 或 0.9842702 ,票數差異沒有顯著性大於隨機情形