空間分析方法與應用 期中考
授課教師：溫在弘 課程助教：杜承軒
考試時間：2018年6月27日中午12:00～6月28日下午6:00
請利用提供的圖資，依照題目中的敘述繪製圖表，並以文字加以說明。
（答案以提供的圖資為準）
※請保留程式碼以利批閱
※以R Markdown編寫的html格式繳交
※作答後請上傳至ceiba作業區之期中考上傳區，逾時不以計分。
※圖資說明：
>>tpe_sqr_bnd.csv 台北市邊界資料
>>tpe_subtown.shp 台北市次分區資料（Sum_CENSUS欄位為人口數資料）
>>school.shp 台北市各級學校點位資料（TYPE欄位為學校種類）
>>>高中：TYPE為「高級中等學校」
>>>國中：TYPE為「國民中學」及「附設國民中學」
>>>國小：TYPE為「國民小學」及「附設國民小學」
library(GISTools);library(rgdal);library(splancs);library(raster);library(dbscan);library(spdep)
rm(list=ls())
setwd("C:/Users/user/Desktop/fin")
#台北市人口邊界資料
tpe=readOGR(dsn = ".", layer = "tpe_subtown", encoding="utf8",verbose = F)
bondary=read.table("tpe_sqr_bnd.csv", header=TRUE, sep=",")
BDP=as.points(bondary[,2], bondary[,3])
#台北市各級學校資料
school=readOGR(dsn = ".", layer = "school", encoding="utf8",verbose = F)
high=subset(school,school@data$TYPE=="高級中等學校")
junior=subset(school,school@data$TYPE %in% c("附設國民中學","國民中學"))
primary=subset(school,school@data$TYPE %in% c("附設國民小學","國民小學"))
school_dt=school@data
high_dt=high@data
junior_dt=junior@data
primary_dt=primary@data
school.pt=as.points(school_dt$XCOORD,school_dt$YCOORD)
high.pt=as.points(high_dt$XCOORD,high_dt$YCOORD)
junior.pt=as.points(junior_dt$XCOORD,junior_dt$YCOORD)
primary.pt=as.points(primary_dt$XCOORD,primary_dt$YCOORD)一、(20%)
「『高中職社區化』是指教育部為均衡高中職發展、整合高中職資源，使國中畢業生能就近升學高中、高職的一種策略性規劃，目的在達成均衡高中職教育品質、學生適性學習以及就近入學，以建構高中職『就學社區』的理想。」
在高中職社區化構想下，想知道台北市內高中與國中之間的分布狀況。請使用Bivariate F function回答「高中附近是否有國中群聚」，即探討國中到最近高中的距離之遠近關係。
參數設定： 請找出3公里以內的F值，以100公尺為間隔。
作答要求： 請畫出Bivariate F function的圖，同時包含隨機情形的包絡線，並以文字加以說明。#F s=seq(0,3000,100) NP.F=Fhat(high.pt,junior.pt, s) nsim=99; U.Fhat=-9999;L.Fhat=99999 for(i in 1:nsim) { NP.CSR=csr(BDP, npts(high.pt)) NP.F.CSR=Fhat(NP.CSR, junior.pt, s) U.Fhat=pmax(NP.F.CSR, U.Fhat) L.Fhat=pmin(NP.F.CSR, L.Fhat) } plot(s, NP.F, type="l", xlab="距離（單位：公尺）", ylab="F值", main="F函數：高中附近是否有國中群聚") lines(s, U.Fhat, col="red", lty=3) lines(s, L.Fhat, col="blue", lty=3)
二、(20%)
配合台北市各村里的人口數或人口密度，使用DBSCAN找出國中的空間群聚關係。
參數設定： 搭配各村里人口資料，來調整DBSCAN設定的參數。
作答要求： 請解釋如何設定參數，繪製出點群聚的地圖，並以文字說明國中的空間群聚關係。#DBSCAN shade=auto.shading(tpe@data$Sum_CENSUS/area(tpe),n=5,cols = paste(brewer.pal(5, "YlGn"),"88",sep="")) colors=c("black",brewer.pal(7, "Set1"),brewer.pal(11, "Paired")) xdbscan=function(pt,u,v){ res=dbscan(pt, eps = u, minPts = v) hullplot(pt, res, asp=1,col =colors,xlim=c(300000,315000), ylim=c(2760000,2790000),main ="",axes =F,cex=NA,xlab="",ylab="") choropleth(tpe,tpe@data$Sum_CENSUS/area(tpe),shade,add=T,border="#BBBBBB") pointmap(pt, col = colors[res$cluster +1] , pch = res$cluster +1 ,cex=1.5,add=T) plot(tpe,col=rgb(0,0,0,0),add=T) } par(mfrow=c(1,2)) xdbscan(junior.pt,1200,3) xdbscan(junior.pt,1300,3)
三、(30%)
加總國中與國小點位至台北市次分區圖層，以每個次分區中，國中小學校之「數量」與「密度」（學校數量／面積），分別繪製出LISA空間自相關地圖。
參數設定： 各次分區的空間鄰近使用距離來判斷(distance-based)，在五公里以內即算是互為鄰近關係。LISA的顯著水準為 α=0.1 (雙尾)。
作答要求： 分別以學校數量和密度，繪製出兩張LISA地圖，並請在地圖上以顏色區分HH、LL、HL、LH及不顯著的地區，以文字加以說明空間群聚關係。並回答你認為使用數量或密度哪個比較適合，以及比較第二題使用點資料作群聚的方法結果有什麼差異。##Q3. xLISA=function(vec,poly,nbw,alpha=0.05,alternative ="two.sided"){ LISA=localmoran(vec,nbw,alternative ="two.sided") diff=vec - mean(vec) z=LISA[,4] quad=vector(mode="numeric",length=nrow(LISA)) quad[diff>0 & z>0]=1 # H-H quad[diff<0 & z>0]=2 # L-L quad[diff>0 & z<0]=3 # H-L quad[diff<0 & z<0]=4 # L-H quad[LISA[, 5]>alpha]=5 #因為在localmoran就設定雙尾 colors=c("red", "blue", "lightpink", "skyblue2", rgb(.95, .95, .95)) plot(poly, border="grey", col=colors[quad], main ="LISA Map") legend("bottomright",legend=c("High-High","Low-Low","High-Low","Low-High","Not Significant"), fill=colors,bty="n",cex=0.7,y.intersp=1,x.intersp=1) } coords=coordinates(tpe) tpe.dnn.nb=dnearneigh(coords, d1=0, d2=5000, row.names=tpe$OBJECTID) tpe.dnn.nb.w=nb2listw(tpe.dnn.nb) sch.count=poly.counts(junior,tpe)+poly.counts(primary,tpe) sch.density=sch.count/area(tpe) xLISA(sch.count,tpe,tpe.dnn.nb.w,0.1)
xLISA(sch.density,tpe,tpe.dnn.nb.w,0.1)
四、(30%)
欲探討台北市各次分區中，人口與國中小學的學校數量的關係。因變數受到空間自相關的影響，故選擇空間迴歸分析，使用Spatial Lag Model作為分析方法：
\(y={\rho}Wy+x{\beta}+{\varepsilon}\)
其中 y 代表各次分區之國中小學的學校數量（單位：間）， x 代表各次分區人口數（單位：每萬人）， W 為空間鄰近矩陣，以Queen 一階鄰近的方式來定義，並經過列標準化後得到的結果。
經過SLM模型計算後得到 ρ=0.4826 ，β=0.2594。請計算spatial multiplier 矩陣後回答下列問題。tpe.nb=poly2nb(tpe) tpe.nb.w=nb2listw(tpe.nb) tpe.nb.wmat <-nb2mat(tpe.nb, style="W") CENSUS=tpe@data$Sum_CENSUS/10000 sar=lagsarlm(sch.count~CENSUS, listw=tpe.nb.w) rho=0.4826 #coef(sar)[1] beta=0.2594 #coef(sar)[3] N=length(CENSUS) eye=diag(N)（一） 照ID順序（即不更動順序下），列出spatial multiplier前5列×前5行的矩陣
sp_multiplier=solve(eye - rho * tpe.nb.wmat) sp_multiplier[1:5,1:5]## 0 1 2 3 4 ## [1,] 1.051821533 0.09066754 0.01786745 0.09471168 0.009978498 ## [2,] 0.126934553 1.04799806 0.12282584 0.12531073 0.022026055 ## [3,] 0.020845356 0.10235487 1.04764594 0.10562754 0.106113732 ## [4,] 0.094711683 0.08950766 0.09053789 1.05428853 0.096766675 ## [5,] 0.008731186 0.01376628 0.07958530 0.08467084 1.060612622（二） 若每個次分區人口都均等地成長10%，請繪製SLM估計每個次分區學校數量增加的趨勢地圖，並列出可能成長前十名的次分區（名字與預測的增加量）。
cvec=CENSUS*0.1 rus.est=sp_multiplier%*% cvec * beta rus.est=data.frame(tpe@data$SEC_NAME,rus.est) rus.est[rev(order(rus.est$rus.est)),][1:10,]## tpe.data.SEC_NAME rus.est ## 62 石牌地區 0.7478363 ## 41 興隆地區 0.6671393 ## 5 三張犁地區 0.6386085 ## 2 東社地區 0.6280708 ## 52 金龍地區 0.6275514 ## 49 東湖地區 0.6220419 ## 42 木柵地區 0.6133745 ## 1 三民地區 0.5924985 ## 60 天母地區 0.5719952 ## 40 景美地區 0.5586238par(mfrow=c(1,1)) tpe$rus.est=rus.est[,2] colR=colorRampPalette(c("white","orange","red"), space = "rgb") spplot(tpe, zcol="rus.est", col.regions=colR(20), main="Spillover Effects")
（三） 假設某發生事件使台北市人口有不規則的改變，在公館地區、學府地區、興隆地區、東門地區的人口各增加一萬人，在大直地區、西湖地區、金龍地區、灣仔地區的人口各減少一萬人，其餘地方保持不變。請繪製SLM估計每個次分區學校數量增減之趨勢地圖（請用視覺化方式呈現出增加或減少的差異），並列出預測增加與減少量最多的各前五名之次分區（名字與預測的增減量）。
cvec=rep(0,N) cvec[c(11,25,29,41)]=1 cvec[c(20,50,52,53)]=-1 rus.est=sp_multiplier%*% cvec * beta rus.est=data.frame(tpe@data$SEC_NAME,rus.est) rus.est[rev(order(rus.est$rus.est)),][1:5,]## tpe.data.SEC_NAME rus.est ## 25 公館地區 0.3321525 ## 11 學府地區 0.3254932 ## 41 興隆地區 0.3175596 ## 29 東門地區 0.2800971 ## 40 景美地區 0.1248943rus.est[rev(order(rus.est$rus.est)),][N:(N-5),]## tpe.data.SEC_NAME rus.est ## 50 西湖地區 -0.35814339 ## 52 金龍地區 -0.33707942 ## 53 灣仔地區 -0.31168098 ## 20 大直地區 -0.29981125 ## 51 紫陽地區 -0.13126086 ## 49 東湖地區 -0.08104363tpe$rus.est=rus.est[,2] colR=colorRampPalette(c("blue","white","red"), space = "rgb") spplot(tpe, zcol="rus.est", col.regions=colR(20), main="Spillover Effects")
