• 沒有找到結果。

33-02 使用R 軟體 繪製統計地圖的介紹與應用(1)

N/A
N/A
Protected

Academic year: 2021

Share "33-02 使用R 軟體 繪製統計地圖的介紹與應用(1)"

Copied!
27
0
0

加載中.... (立即查看全文)

全文

(1)

1

使用 R 軟體 繪製統計地圖的介紹與應用(1)

蔡靜雯 副統計分析師 地區性的統計資料結合地圖方式呈現,可以更容易觀察到不同區域間的變化, 和是否有地埋位置的影響。本期的 eNews 內容,分成 2 大部分,第 1 部分介紹 地圖的相關資料以及 shp 地圖檔案格式。第 2 部分則是使用台灣地圖檔案和戶 政各地的統計資料,使用 R 軟體實際操作,繪製台灣統計地圖。

第1部分 地圖資料介紹

1-1、坐標系統 目前台灣常見的地圖座標系統有以下3種:

• TWD67 - 採用國際上Geodetic Reference System 1967 橢球體, 作為臺灣地區 參考橢球體,以南投埔里虎子山為大地基準。

• TWD97 - 採用國際上Geodetic Reference System 1980 橢球體, 作為臺灣地區 參考橢球體,以八個衛星追蹤站為大地基準,1997年完成以 GPS 重 新計算坐標基準。 • WGS84 - 美國國防部於西元1984採用新的地球標準物理模型,以地球的質量 中心為中心點,加上了分布在全世界各地的1500個地理座標參考點, 為GPS而制定的世界座標系統,目前為國際通用的經緯度坐標。 TWD67 是與該地區最密合的「區域性」大地基準,只適用於台灣地區的座標系統, 與 WGS84 的坐標基準相差約1公里,TWD97 的坐標基準則與 WGS84 的坐標基 準差異不大,相差約幾公分至數十公分。目前線上的地圖系統幾乎都使用WGS84的 座標系統(Google Map, UrMap, Google Earth)。

(2)

2 1-2、地圖資料格式

ESRI Shapefile (shp),或簡稱shapefile,是美國環境系統研究所公司(ESRI)開發的 空間資料開放格式。一個 Shapefile 檔案通常包含多個檔案,同組地圖資料其主檔 名應該相同,如下圖範例所示,相同主檔名代表為同一組資料, gadm36_TWN_0、 gadm36_TWN_1、gadm36_TWN_2 分別為三組地圖資料。 Shapefile 檔案中有三個必須檔案,其副檔名和儲存的資料內容,分別為 - .shp:圖形格式,用以儲存資料的圖形, - .shx:圖形索引格式,記錄圖形的位置, - .dbf:屬性資料格式,儲存每個圖形的屬性資料。  如果要讀取 shapefile 檔,以上三個檔案都必須同時具備,否則會無法讀取 成功。 有時可能還會有其它非必須的檔案: - .prj:投影格式,用於儲存地理座標系統與投影資訊., - .sbn:空間索引檔,可以加快開啟shapefile的速度。.

(3)

3

第2部分 繪製統計地圖

2-1、台灣地圖資料下載

以下兩個地方的網址都可以下載台灣的地圖資料

1. GADM maps and data https://gadm.org/download_country_v3.html

GADM是全球行政區域資料庫,採用座標系統 WGS84,提供了所有國家及 其行政區域細分的地圖和空間資料,可免費用於學術和其他非商業用途。未經 GADM事先許可,不得用於商業用途。 下載步驟:點選網址後,在Country欄位選取Taiwan或直接輸入Taiwan,點選 Shapefile 即可下載台灣地圖檔案。 下載的壓縮檔,共有16個檔案,除了 license.txt 外有3種檔案名稱,分別代表3組 不同區域劃分的地圖資料: gadm36_TWN_0 - 無區分縣市,只有台灣形狀, gadm36_TWN_1 - 區分六都區域, gadm36_TWN_2 - 區分22個縣市區域。 從壓縮檔內容,也同時確認繪製地圖必須的3個檔案 .shp、.shx 和 .dbf 檔都有 包含在其中。

(4)

4 2. 鄉鎮市區界線(TWD97經緯度) https://data.gov.tw/dataset/7441 為我國各鄉(鎮、市、區)行政區域界線圖資,採用座標系統 TWD97 由內 政部國土測繪中心提供,在政府資料開放平臺提供免費下載,可直接點網址下 載檔案,或依照以下路徑:政府資料開放平臺->公共資訊->內政部國土測繪中 心->鄉鎮市區界線(TWD97經緯度),找到該檔案的所在位置。除此之外,也有其 它台灣相關的地圖資料可下載。 下載步驟:點選網址頁面的 SHP 下載。

(5)

5

下載的壓縮檔,共有7個檔案,其中也確認有包含繪製地圖必須的3個檔案 .shp、.shx 和 .dbf 檔。

(6)

6 2-2、讀入地圖資料檔案 以下範例使用 GADM 區分22縣市的地圖檔案 gadm36_TWN_1。 [程式碼] - 讀入 shp 檔案 install.packages("rgdal") #安裝 讀取shp檔案的套件 library(rgdal) taiwan_shap2=readOGR("C:/Users/CWTsai/Desktop/plot -繪製地圖/gadm 36_TWN_shp/gadm36_TWN_2.shp") #讀入shp檔案 - output - 顯示此檔案有22筆資料,13個欄位 [程式碼] - 使用 head 函數,查看前6筆資料的內容(位於地圖檔案的data屬性裡) head(taiwan_shap2@data) - output - 欄位名稱NL_NAME_1、VARNAME_2、NL_NAME_2、TYPE_2的資料呈現亂 碼,需進一步處理亂碼問題。

(7)

7

[程式碼] - 使用 iconv 函數,對亂碼欄位執行轉換

taiwan_shap2@data$NL_NAME_1 <- iconv(taiwan_shap2@data$NL_NAME_ 1, from = "UTF-8", to="UTF-8")

taiwan_shap2@data$VARNAME_2 <- iconv(taiwan_shap2@data$VARNAME_ 2, from = "UTF-8", to="UTF-8")

taiwan_shap2@data$NL_NAME_2 <- iconv(taiwan_shap2@data$NL_NAME_ 2, from = "UTF-8", to="UTF-8")

taiwan_shap2@data$TYPE_2 <- iconv(taiwan_shap2@data$TYPE_2, from = "UTF-8", to="UTF-8")

head(taiwan_shap2@data)

- output - 主要的中文縣市名稱欄位 NL_NAME_1 已成功轉換成繁體中文,雖然 VARNAME_2 欄位還是有亂碼情況,但因後續不會用到此欄位,可忽略 亂碼問題。

(8)

8 2-3、繪製地圖

[程式碼] - 使用ggplot() 函數繪製台灣地圖,其中,x 放經度(long),y 放緯度(lat), group 為22個縣市 install.packages("ggplot2") #安裝 讀取shp檔案的套件 library(ggplot2) ggplot(taiwan_shap2, aes(x=long,y=lat,group=group))+geom_path() - output - 區分22個縣市的台灣地圖 2-3、統計資料檔案 使用內政統計指標縣市排名資料https://www.moi.gov.tw/stat/node.aspx?sn=6718, 分別下載100~106年的縣市指標檔案,採用戶政檔案裡的出生登記嬰兒性比例資料 (性比例係指每百個女嬰的男嬰數量),再匯整在一起存成CSV檔。

(9)

9

原始檔案內容:

(10)

10

[程式碼] - 讀入匯整完成的嬰兒性比例統計資料檔案 index.csv

index.data=read.csv("C:/Users/CWTsai/Desktop/plot -繪製地圖/政府開 放資料/縣市指標歷年資料/index.csv", head=T)

View(index.data)

- output - 縣市欄位名稱city,bsex_ratio100 ~ bsex_ratio106 欄位為100~106年嬰兒 性別比例 2-4、地圖資料和統計資料整理合併 地圖資料(taiwan_shap2)的縣市名稱和嬰兒性比例統計資料( index.data )的縣市名 稱要一致,才能將資料合併在一起。 [程式碼] - 地圖資料(taiwan_shap2)的縣市名稱 taiwan_shap2$NL_NAME_2 - output - 地圖資料的縣市名稱欄位 NL_NAME_2 屬性為文字字串,名稱頭尾有雙 引號 [程式碼] - 嬰兒性比例統計資料( index.data )的縣市名稱 index.data$city - output - 嬰兒性比例資料的縣市名稱欄位city屬性為因子變數(factor),Levels顯示 有幾種分類變數

(11)

11 雖然在合併檔案時,factor屬性欄位和字串屬性欄位合併沒有問題,但建議還是將 欄位屬性轉成一致再合併。 [程式碼] - 地圖資料的縣市名稱和嬰兒性比例資料的縣市名稱整理成一致 #另外建立一個 zh_NAME 欄位,存放修正後的地圖資料縣市名稱 taiwan_shap2$zh_NAME=taiwan_shap2$NL_NAME_2 #將地圖和嬰兒性比例資料的5個不一致縣市名稱,重新命名與嬰兒性比例的縣市名稱 一樣 taiwan_shap2$zh_NAME[taiwan_shap2$zh_NAME=="馬祖列島"]="連江縣" taiwan_shap2$zh_NAME[taiwan_shap2$zh_NAME=="台中"]="臺中市" taiwan_shap2$zh_NAME[taiwan_shap2$zh_NAME=="台南"]="臺南市" taiwan_shap2$zh_NAME[taiwan_shap2$zh_NAME=="台北市"]="臺北市" taiwan_shap2$zh_NAME[taiwan_shap2$zh_NAME=="台東縣"]="臺東縣" #將嬰兒性比例資料的縣市名稱欄位city屬性轉成文字字串 index.data$city=as.character(index.data$city) #排序後print修正後的地圖和嬰兒性比例資料,確認縣市名稱是否一致 [程式碼] - print修正後的地圖資料其縣市名稱欄位 zh_NAME,排序sort是為了方便 比對 sort(taiwan_shap2$zh_NAME) - output - [程式碼] - print修正後的嬰兒性比例資料其縣市名稱欄位city sort(index.data$city) - output -

(12)

12

[程式碼] - 合併地圖資料和嬰兒性比例資料

taiwan_shap2@data=merge(taiwan_shap2@data, index.data, by.x="zh_ NAME", by.y="city", sort=F)

head(taiwan_shap2@data) - output - 100~106年的嬰兒性比例資料與地圖資料成功合併在同一個檔案 2-5、繪製合併後的統計地圖資料 2-5-1、方法1 - 使用qtm() 函數 install.packages("tmap") #安裝 qtm 函數案的套件 library(tmap)

qtm(taiwan_shap2,fill="bsex_ratio100", text="zh_NAME", text.size =0.7, fill.title="birth sex ratio", fill.palette="Reds")

(13)

13

- output -

[程式碼] - 100年和106年的嬰兒性比例統計地圖

qtm(taiwan_shap2,fill="bsex_ratio100", fill.title="birth sex rat io(100)")

qtm(taiwan_shap2,fill="bsex_ratio106", fill.title="birth sex rat io(106)")

- output - qtm()函數會根據嬰兒性比例的數值分佈,自動分類比例組別和給定組別 顏色,100年和106年的圖例,雖然性別比例都是以每 5% 分組,但給定 的組別顏色不一樣,同時放一起比較時,會因圖例不一致而不易比較。

(14)

14

2-5-2、方法2 – plot() 函數  調色盤介紹

Palettes {grDevices} 套件裡,有rainbow、heat.colors、terrain.colors、topo.colors、 cm.colors 5種類型的調色盤。

[程式碼] - 呈現5種類型的調色盤顏色 #Display a Color Image

display.col=function(color,title){

image(1:10,1,as.matrix(1:10),col=color,xlab=title, ylab="",xaxt="n",yaxt="n",bty="n") } par(mfcol=c(5,1), mai=c(0.1,0.1,0.5,0)) display.col(rainbow(10),"rainbow"); display.col(heat.colors(10),"heat.colors"); display.col(terrain.colors(10),"terrain.colors"); display.col(topo.colors(10),"topo.colors"); display.col(cm.colors(10),"cm.colors" ) - output - 由上至下分別為rainbow、heat.colors、terrain.colors、topo.colors、cm.colors 5種調色盤的顏色

(15)

15 使用plot() 函數,繪製統計地圖步驟,依序程式碼1~程式碼3, [程式碼 1] - 建立繪圖色盤,100 年~106 年的性別比例介於 85%~135% 之間,每 5 個百分比為一組,共分成 10 組,指定繪圖色盤有 10 種顏色。 my.color=rev(heat.colors(10)) [程式碼 2] - 建立一致性的性別比例資料分類,100 年~106 年的性別比例介於 85%~135% 之間,每 5 個百分比為一組,並從建立的繪圖色盤給每 一種比例組別一個固定的顏色 ratio.color100=my.color[cut(taiwan_shap2@data$bsex_ratio100, bre aks=seq(85,135,5),labels=F)] [程式碼 3] - 繪製性別比地圖 #設定畫布邊界 par(mai=c(0.1,0.1,0.5,0)) #繪製統計地圖 plot(taiwan_shap2, col=ratio.color100) #22個縣市區域標註縣市名稱 text(coordinates(taiwan_shap2)[,1], coordinates(taiwan_shap2)[, 2], taiwan_shap2$zh_NAME, cex=0.6, col="gray5")

#圖例設定

legend("bottomleft",

title="birth sex ratio(100)",

legend=c("85~90", "90~95", "95~100", "100~105", "105~110", "110~115", "115~120", "120~125", "125~130", "130~135"), fill=my.color,

(16)

16

- output -

[程式碼] - 呈現 100~106 年 縣市性別比例地圖 attach(taiwan_shap2@data)

#將繪製統計地圖步驟(程式碼1~程式碼3)寫成function,年度為參數 ratio.plot=function(year){

bsex_ratio=paste("bsex_ratio", year, sep="") #建立繪圖色盤

ratio.color=my.color[cut(get(bsex_ratio), breaks=seq(85,135,5),l abels=F)]

#繪製統計地圖

plot(taiwan_shap2, col=ratio.color, main=paste(year, "年 出生性別比 ", sep="") )

#圖例設定

legend("bottomleft",

title=paste("birth sex ratio", sep="") ,

(17)

17 "110~115", "115~120", "120~125", "125~130", "130~135"), fill=my.color, cex=0.8) } #設定畫布邊界 par(mai=c(0.1,0.1,0.5,0), mfrow=c(2,4)) #將繪製統計地圖步驟的function ratio.plot(year),使用迴圈,畫出100~1 06年的統計地圖 for (year in 100:106){ ratio.plot(year) } - output - 每個年度的性別比例分組和指定的分組顏色都一致,就可直接從顏色觀 察 22 個縣市,從 100 年~106 年的出生嬰兒性別比例的變化。

(18)

18 2-5-3、方法3 – ggplot() 函數 先將 .shp 檔案中的地區屬性資料和經緯度資料分別抓出,再合併在一起 [程式碼] - 抓出 .shp 檔案中的地區屬性資料(包含地區名稱和先前已合併的統計 資料) taiwan_shap2.data=taiwan_shap2@data head(taiwan_shap2.data) output

(19)

-19 [程式碼] - 在taiwan_shap2.data(地區屬性資料)中增加一個 id 欄位(用於合併經緯 度資料) taiwan_shap2.data$id=as.character(0:21) head(taiwan_shap2.data) output -[程式碼] - 使用 fortify() 函數將 .shp 檔案中的經緯度抓出來 taiwan_shap2.latitude_ongitude=fortify(taiwan_shap2) head(taiwan_shap2.latitude_ongitude)

(20)

20 [程式碼] - 合併 taiwan_shap2.data(地區屬性+統計資料) 和 taiwan_shap2.latitude_ongitude (經緯度資料) taiwan_mapdata=merge(taiwan_shap2.data, taiwan_shap2.latitude_on gitude, by="id")

ggplot(taiwan_mapdata, aes(x=long, y=lat, group=group, fill=bsex_ ratio100 )) + geom_polygon(colour="black")

(21)

-21 2-6、使用鄉鎮市區界線(TWD97經緯度)地圖資料,依照前面執行步驟,繪製台灣368 個鄉鎮區域地圖 [程式碼] – 讀入縣市鄉鎮地圖檔案 town_shap=readOGR(dsn="C:/Users/CWTsai/Desktop/plot -繪製地圖/鄉鎮 地圖資料/TOWN_MOI_1071226.shp") - output - 縣市鄉鎮地圖檔案中,有368筆資料(368個鄉鎮),7個欄位 [程式碼] # 中文亂碼轉換 iconv {base} town_shap@data$COUNTYNAME=iconv(town_shap@data$COUNTYNAME , from = "UTF-8", to="UTF-8")

town_shap@data$TOWNNAME=iconv(town_shap@data$TOWNNAME , from = "U TF-8", to="UTF-8")

head(town_shap@data) # 中文正常顯示 output

(22)

-22 [程式碼] #繪製縣市鄉鎮地圖-方法1 qtm() 函數 library(tmap) qtm(town_shap) # 未給定鄉鎮顏色 output

(23)

-23 [程式碼] - #以鄉鎮市區界線(TWD97經緯度)地圖檔案中,現有的22個縣市 COUNTYCODE欄位編 碼,做為fill的指定縣市顏色 qtm(town_shap, fill="COUNTYCODE") output

(24)

-24

[程式碼] -

#以鄉鎮市區界線(TWD97經緯度)地圖檔案中,現有的365個鄉鎮 TOWNCODE 欄位編 碼,做為fill的指定鄉鎮顏色

qtm(town_shap, fill="TOWNCODE") #fill最多只能指定30種顏色,分類超過 30個組別會自動合併分成30

(25)

-25

[程式碼]

#繪製縣市鄉鎮地圖-方法2 plot()函數 plot(town_shap) #未給定鄉鎮區域顏色 output

-plot(town_shap, col=rainbow(365)) #使用rainbow指令給定鄉鎮區域顏色 output

(26)

-26 參考資料 坐標系統相關 • Taiwan datums https://wiki.osgeo.org/wiki/Taiwan_datums#TWD97_vs._WGS84 • 上河文化 – 大地坐標系統漫談 http://www.sunriver.com.tw/grid_tm2.htm • 中華民國內政部國士測繪中心 – 平面控制 https://www.nlsc.gov.tw/Home/MakePage/42?level=42 • 台灣常用坐標系統簡介 https://blog.xuite.net/lwkntu/blog/364536963-%E5%8F%B0%E7%81%A3%E5%B8%B8%E7%94%A8%E5%9D%90%E6%A8%99%E7 %B3%BB%E7%B5%B1%E7%B0%A1%E4%BB%8B • 坐標系統 Coordinate systems http://140.121.160.124/GEO/%E5%BA%A7%E6%A8%99%E7%B3%BB%E7%B5%B1 .pdf shp 檔案 • 維基百科 Shapefile https://zh.wikipedia.org/wiki/Shapefile • Ultimanual 的 ArcGIS https://ultimanualarcgis.wordpress.com/2016/06/03/2-fu-dang-ming-de-yi-yi/ 繪製地圖相關 • qtm function | R Documentation https://www.rdocumentation.org/packages/tmap/versions/2.2/topics/qtm • 深入淺出繪製統計地圖 4 – R https://ariheart2011.wordpress.com/2017/05/10/%E6%B7%B1%E5%85%A5%E6% B7%BA%E5%87%BA%E7%B9%AA%E8%A3%BD%E7%B5%B1%E8%A8%88%E5%9C

(27)

27

%B0%E5%9C%96-4-r/

• 主題式地圖(Thematic map) - 政府開放資料為例

參考文獻

相關文件

利用 Microsoft Access 資料庫管理軟體,在 PC Windows 作業系 統環境下,將給與的紙本或電子檔(如 excel

根據商務活動之舉辦目標及系統需求,應用 Microsoft Office 文書處理 Word、電子試算表 Excel、電腦簡報 PowerPoint、資料庫 Access

利用 Microsoft Access 資料庫管理軟體,在 PC Windows 作業系統環境 下,將給與的紙本或電子檔(如 excel

(軟體應用) 根據商務活動之舉辦目標及系統需求,應用 Microsoft Office 文書處理 Word、電子試算表 Excel、電腦簡報 PowerPoint、資料庫 Access

審查整理呈現資料:蒐集到的資料應先審核 是否完整、正確、合理與一致,然後利用敘

利用 Microsoft Access 資料庫管理軟體,在 PC Windows 作業系 統環境下,將給與的紙本或電子檔(如 excel

D5.1 應用1個具體圖像代表 1個單位,製作象形圖 D5.2

‡ RFID 運作原理是透過一片小型硬體的無線射頻辨識技 術晶片( RFID chips),利用內含的天線來傳送與接