空間分析　第三章

第三章探索地理空間的幾何關係

讀取資料與初步整理

library(sf);library(raster)
setwd("D:/SA_Book/Lab_Data_Ch3")
hospital=st_read(dsn = "Taipei_hospital.shp")
village=st_read(dsn = "Taipei_village.shp")
flood.raster=raster("Taipei_flood200.tif")
windowsFonts(JH = windowsFont("微軟正黑體"))

問題3-1：幾何中心點是否落在淹水區內

依據淹水潛勢資料作評估，當台北市發生兩百年重現期規模的降雨事件時，哪些里會是潛在的淹水受災區？

#st_centroid函數：幾何中心點
centroids=st_centroid(village)
#找出淹水受災區，捨去淹水不到0.5m的網格，網格轉向量並融合
values(flood.raster)[values(flood.raster)<0.5]=NA
flood=st_as_sf(as(flood.raster, "SpatialPolygonsDataFrame")) #網格轉向量
flood=st_union(flood)#融合面資料以利計算
flood=st_transform(flood,st_crs(village)) #投影座標轉換成一致的座標
#篩選出被影響村里（st_covered_by函數：點是否在面之中）
cen_id=which(st_covered_by(centroids,flood,sparse=F))
#繪圖
par(family="JH",mar=c(1,1,1,1),mfrow=c(1,2)) #圖框與字型設定
plot(st_geometry(village))
plot(st_geometry(flood),col = 'blue',add=T,border=NA)
plot(st_geometry(centroids[cen_id,]),col='red',pch=20,add=T)
plot(st_geometry(village))
plot(st_geometry(village[cen_id,]),col='red',add=T)

問題3-2：行政區是否落在淹水區內

承問題3-1，然而對於里長而言，上述的操作型定義可能不夠完善，因為只要該里有部分範圍淹水就需要進行處理。已知發生兩百年重現期規模的降雨事件時，一共需通知哪些台北市的里長，該里為受災區（淹水深度達50公分以上）？

#st_is_within_distance函數：兩圖資之間的距離是否在閾值以內
vill_id=st_is_within_distance(flood,village,dist=0,sparse =F)
#繪圖
par(family="JH",mar=c(1,1,1,1)) #圖框與字型設定
plot(st_geometry(village))
plot(st_geometry(village[vill_id,]),col='red',add=T)
plot(st_geometry(flood),col='blue',add=T)

問題3-3：評估受災影響人數

災害受影響人數的估計是災害風險評估與救災資源準備的重要參考資訊。請問臺北市該次淹水事件的潛在受影響人數為多少？

#村里中心點
sprintf("以村里中心點估計，有%d個村里受災，受災人數為%d人",length(cen_id),sum(village$Total[cen_id]))

## [1] "以村里中心點估計，有78個村里受災，受災人數為453736人"

#村里面空間
sprintf("以村里面空間估計，有%d個村里受災，受災人數為%d人",length(vill_id),sum(village$Total[vill_id]))

## [1] "以村里面空間估計，有456個村里受災，受災人數為2018074人"

#以各里的淹水面積估計受影響人數
#st_intersection函數：兩圖資之間交集
VF=st_intersection(village,flood)
village$flood_area=0
village$flood_area[vill_id]=st_area(VF) #每個村里的淹水面積
village$damaged_pop=round(village$Total*village$flood_area/st_area(village)) #每個村里的受災人數
sprintf("以村里淹水面積估計受災人數為%d人",sum(village$damaged_pop))

## [1] "以村里淹水面積估計受災人數為481266人"

問題3-4：脆弱人口密度分布圖

如何建立受災的脆弱人口密度主題圖？

A65_up.df=st_drop_geometry(village[,25:32])
village$A65_up=rowSums(A65_up.df) #計算各村里高齡（脆弱）人口
village$vulpop=round(village$A65_up*village$flood_area/st_area(village)) #每個村里的受災脆弱人數
village$vden=village$vulpop/st_area(village)*10^6 #計算受災脆弱人口密度

library(sf);library(cartography)
#面量圖繪製
par(family="JH",mar=c(1,1,1,1)) #圖框與字型設定
plot(st_geometry(village),bg = "lightblue1") #繪製底圖
choroLayer(village,var='vden',method="fisher-jenks",nclass=7,col=carto.pal(pal1 = "red.pal", n1 = 7),legend.values.rnd=3,add=T,legend.title.txt ='脆弱人口密度')
layoutLayer(title="台北市200年重現期降雨淹水事件受災脆弱人口密度面量圖",north=T,postitle = "center")

問題3-5：環域＋疊圖分析

假設在災害期間，各家醫院能提供的緊急醫療服務為周圍的1公里範圍，台北市有多少個里完全沒有可近的醫療資源。

#st_buffer函數：建立環域範圍
hospital_1km=st_buffer(hospital,dist=1000) #建立醫院1公里環域範圍
vill.host_id=which(lengths(st_is_within_distance(village,hospital_1km,dist=0))==0) #找出和醫院1公里沒有交集的村里
#繪圖
par(mar=c(1,1,1,1))
plot(st_geometry(village))
plot(st_geometry(village[vill.host_id,]),col='orange',add=T)
plot(st_geometry(hospital_1km),col="#FF004455",pch=20,border=NA,cex=2,add=T)
plot(st_geometry(hospital),col='red',pch=20,cex=1,add=T)

問題3-6：鄰近性分析

若能事先評估醫院災期時的服務人口總數，可有助於提升災害應變時醫療資源的分配效率。假設醫院應服務其所在位置最近的5個里的里民，請問馬偕醫院的服務人口總數有多少？

MMH=subset(hospital,attr_ANNO=="馬偕醫院") #選出馬偕醫院
#st_distance：計算兩圖層之間的最短距離
cen_MMH=st_distance(MMH,centroids)#計算各里中心到馬偕醫院的距離
near5_MMH=order(cen_MMH)[1:5]#依距離排序後挑選離馬偕醫院最近的五個里
sprintf("馬偕醫院最近的五個里的人口數總和為%d人",sum(village$Total[near5_MMH]))

## [1] "馬偕醫院最近的五個里的人口數總和為23206人"

第三章 探索地理空間的幾何關係

第三章探索地理空間的幾何關係