第二章 地圖繪製：空間分布的呈現
讀取資料與初步整理
library(sf);library(cartography)
setwd('D:/SA_Book/Lab_Data_Ch2')
temple=st_read(dsn='Taiwan_temple.shp')
town=st_read(dsn='Taiwan_town.shp')
county=st_read(dsn='Taiwan_county.shp')
windowsFonts(JH=windowsFont("微軟正黑體"))問題2-1：實際點位資料的空間分布
請問台灣媽祖廟宇的空間分布狀況為何？
par(family="JH",mar=c(1,1,1,1)) #圖框與字型設定
#subset函數：從所有廟宇中選出媽祖
mazhou=subset(temple,God=='媽祖')
plot(st_geometry(county),col='lightgreen',bg = "lightblue1") #繪製底圖
plot(st_geometry(mazhou),add=T,cex=0.5,col='red',pch=20) #繪製媽祖廟宇點位
layoutLayer(title="台灣媽祖廟宇分布圖",north=T,postitle = "center")
legend('bottomleft',legend='媽祖廟宇',col='red',pch=20,bty ='n')
問題2-2：空間統計資料的分布－點子圖
如果手邊現有的資料只有各行政區的媽祖廟宇數量，而沒有各廟宇的實際點位，透過點子圖，我們可以從媽祖廟宇的面資料（行政區）表現出其空間分布及疏密。比較台灣中南部地區以兩種行政區尺度單位－鄉鎮市區及縣市統計資料產生的點子圖，兩者的空間分布有何差異？
par(family="JH",mar=c(1,1,1,1)) #圖框與字型設定
#計算各鄉鎮媽祖廟宇數量
#st_contains函數：計算多邊形中有哪些點
town$mazhou_count=lengths(st_contains(town,mazhou))
plot(st_geometry(county),col='lightgreen',bg = "lightblue1") #繪製底圖
dotDensityLayer(town,var='mazhou_count',n=1,cex=0.5,pch=20,col='red',legend.txt='1點代表1間廟',legend.pos = 'bottomleft') #繪製點子圖
layoutLayer(title="台灣媽祖廟宇點子圖",north=T,postitle = "center") 
#計算各縣市媽祖廟宇數量
county$mazhou_count=lengths(st_contains(county,mazhou))
plot(st_geometry(county),col='lightgreen',bg = "lightblue1") #繪製底圖
dotDensityLayer(county,var='mazhou_count',n=1,cex=0.5,pch=20,col='red',legend.txt='1點代表1間廟',legend.pos = 'bottomleft') #繪製點子圖
layoutLayer(title="台灣媽祖廟宇點子圖（縣市）",north=T,postitle = "center") 
問題2-3：空間統計資料的分布－分級符號圖
問題2-1及問題2-2透過實際點位資料或點子圖呈現媽祖廟宇的空間分布及疏密。現在，如果同樣只有縣市統計的媽祖數量資料，應如何繪製地圖，以突顯各縣市間媽祖廟宇數量上的差異？
par(family="JH",mar=c(1,1,1,1)) #圖框與字型設定
plot(st_geometry(county),col='lightgreen',bg = "lightblue1") #繪製底圖
propSymbolsLayer(county,var='mazhou_count',inches=0.2,legend.title.txt ='媽祖廟宇數量') #繪製分級符號圖
layoutLayer(title="台灣媽祖廟宇分級符號圖",north=T,postitle = "center") 
par(mar=c(4,4,2,4),mfrow=c(2,1)) #繪製各縣市媽祖廟宇數量長條圖
barplot(county$mazhou_count,names.arg = county$COUNTYNAME,las=2)
barplot(sort(county$mazhou_count,decreasing = T),names.arg = county$COUNTYNAME[order(county$mazhou_count,decreasing = T)],las=2)
問題2-4：空間統計資料的分布－面量圖（合適的分級方法）
在問題2-2中，利用點子圖點子的疏密來表現媽祖廟宇的分布密度。然而，若想以面為單位來呈現媽祖廟宇的分布密度，則應如何繪製地圖？此外，分級組界應如何劃定較為合適？
par(family="JH",mar=c(1,1,1,1)) #圖框與字型設定
#計算人口密度（st_area函數：計算密度）
town['mazhou_density']=town$mazhou_count/st_area(town)*10^6
plot(st_geometry(county),bg = "lightblue1") #繪製底圖
choroLayer(town,var='mazhou_density',method="quantile",nclass=5,col=carto.pal(pal1 = "red.pal", n1 = 5),legend.values.rnd=3,add=T,legend.title.txt ='媽祖廟宇密度')
layoutLayer(title="台灣媽祖廟宇分布密度面量圖",north=T,postitle = "center") 
#不同的分級方法可以透過method參數來變更
#自然分級法：fisher-jenks、等量分級法：quantile、等距分級法：equal
methods_name=c("自然分級法","等量分級法","等距分級法")
methods=c("fisher-jenks","quantile","equal")
par(mfrow=c(1,3),mar=c(2,0,2,0))
for(i in 1:3){
plot(st_geometry(county),bg = "lightblue1")
choroLayer(town,var='mazhou_density',method=methods[i],nclass=5,col=carto.pal(pal1 = "red.pal", n1 = 5),legend.values.rnd=3,add=T,legend.title.txt ='媽祖廟宇密度')
layoutLayer(title=methods_name[i],north=T,postitle = "center")
}
問題2-5：空間統計資料的分布－面量圖（合適的數值單位）
繪製面量圖時，應以推衍數值（如密度、比例）作為分層設色的數值單位。請問若以不同的數值單位製作面量圖，如媽祖廟宇數量及媽祖廟宇數量等級，可能會有哪些效果？可能衍伸什麼樣的誤導？
par(family="JH",mfrow=c(1,3),mar=c(2,0,2,0)) #圖框與字型設定
# 媽祖廟宇數量
plot(st_geometry(county),bg = "lightblue1") #繪製底圖
choroLayer(town,var='mazhou_count',method="fisher-jenks",nclass=5,col=carto.pal(pal1 = "red.pal", n1 = 5),legend.values.rnd=3,add=T,legend.title.txt ='媽祖廟宇數量')
layoutLayer(title="媽祖廟宇數量面量圖",north=T,postitle = "center")
# 媽祖廟宇數量等級
class=getBreaks(town$mazhou_count,method = "fisher-jenks",n=3) #以自然分級法分成三級做分級範例
town$level=cut(town$mazhou_count,class,include.lowest=T)
levels(town$level)=c("Level1","Level2","Level3")
plot(st_geometry(county),bg = "lightblue1") #繪製底圖
typoLayer(town,var='level',col=carto.pal(pal1 = "red.pal", n1 = 3),add=T,legend.title.txt ='媽祖廟宇數量等級')
layoutLayer(title="媽祖廟宇數量等級面量圖",north=T,postitle = "center")
# 媽祖廟宇數量密度
plot(st_geometry(county),bg = "lightblue1") #繪製底圖
choroLayer(town,var='mazhou_density',method="fisher-jenks",nclass=5,col=carto.pal(pal1 = "red.pal", n1 = 5),legend.values.rnd=3,add=T,legend.title.txt ='媽祖廟宇密度')
layoutLayer(title="媽祖廟宇密度面量圖",north=T,postitle = "center") 