空間分析　 實習八

助教 杜承軒 2021.05.16

利用KDE方法，分別呈現台灣南部地區（嘉南高屏)「媽祖」與「觀音菩薩」寺廟密度地圖

＊資料前處理

library(sf);library(tmap);library(SpatialKDE)
setwd("D:/1092SA/Data/")
Tw=st_read('Taiwan_county.shp',options="ENCODING=Big5")
SOUTH=Tw[Tw$COUNTY%in%c("嘉義縣","嘉義市","台南市","高雄市","屏東縣"),]

temple=st_read('Tempcycle_twd97.shp')
temple=temple[SOUTH,]

temple=temple[!is.na(temple$主祭神祇),]
Mazu=temple[temple$主祭神祇=="媽祖",]
Guan=temple[temple$主祭神祇=="觀音菩薩",]

SOUTH.lyr=tm_shape(SOUTH)+tm_borders()+tm_layout(frame = F)
Mazu.lyr=qtm(Mazu,symbols.col='blue',symbols.size=.03,symbols.border.lwd=0)
Guan.lyr=qtm(Guan,symbols.col='red',symbols.size=.03,symbols.border.lwd=0)

KDE(grid)

grid=create_grid_rectangular(SOUTH,cell_size=2000)
Mazu.kde=kde(Mazu,5000,grid=grid)
Guan.kde=kde(Guan,5000,grid=grid)

Mazu.map=tm_shape(Mazu.kde)+tm_polygons("kde_value",palette="Blues",border.alpha=0,title="Mazu")+SOUTH.lyr+Mazu.lyr
Guan.map=tm_shape(Guan.kde)+tm_polygons("kde_value",palette="Reds",border.alpha=0,title="Guanyin")+SOUTH.lyr+Guan.lyr
tmap_arrange(Mazu.map,Guan.map)

KDE(raster)

raster=create_raster(SOUTH,2000)
Mazu.KDE=kde(Mazu,5000,grid=raster)
Guan.KDE=kde(Guan,5000,grid=raster)

Mazu.rmap=tm_shape(Mazu.KDE) + tm_raster(palette="Blues",title="Mazu")+SOUTH.lyr+Mazu.lyr
Guan.rmap=tm_shape(Guan.KDE) + tm_raster(palette="Reds",title="Guanyin")+SOUTH.lyr+Guan.lyr
tmap_arrange(Mazu.rmap,Guan.rmap)

Dual KDE(grid)

Dual.kde=Mazu.kde
Dual.kde$kde_value=Mazu.kde$kde_value-Guan.kde$kde_value
tm_shape(Dual.kde)+tm_polygons("kde_value",palette="RdBu",border.alpha=0,title="Dual")+SOUTH.lyr+Mazu.lyr+Guan.lyr

Dual KDE(raster)

Dual.KDE=Mazu.KDE-Guan.KDE
tm_shape(Dual.KDE) + tm_raster(palette="RdBu",title="Dual")+SOUTH.lyr+Mazu.lyr+Guan.lyr