R與SAS軟件在繪制統計地圖方面的應用研究

天津醫科大學公共衛生學院流行病與衛生統計學教研室(300070)

鄧曉偉 陳云峰 劉紅偉 段同慶 馬 駿△

統計地圖是用不同的顏色和花紋表示統計數值在地理分布上的變化，適宜描述研究指標的地理分布[1]。尤其是在流行病學研究中，統計地圖用來描述疾病在地區、時間的分布再合適不過。

地圖的繪制、處理經常出現在地理信息系統(geographical information systems)中。隨著學科交融和統計學的發展，地圖的繪制，特別是統計地圖的繪制已經不僅限于使用固定的地理信息系統專用軟件，如ArcGIS等。SAS作為統計分析商用軟件有Graph模塊可以實現統計制圖；R作為開源統計編程語言，Base中有plot函數和maptools、sp、sf等許多包可以實現相應的統計制圖，其中最廣泛應用的當屬ggplot2和sf包了。本研究旨在應用R與SAS繪制統計地圖，并比較其優缺點。

本研究應用地圖數據分別是SAS軟件SAS/GRAPH中自帶的中國地圖，和GADM(global admin-istrativeareas)開源數據庫[2]的數據，包含中國(除港澳臺地區)地圖、中國香港特別行政區地圖、中國澳門特別行政區地圖和中國臺灣地區地圖。應用這兩種軟件，可以實現除港澳臺地區中國省級行政區劃衛生機構數量[3]統計地圖的繪制。所有制圖均保持基本參數調節一致，不進行過度修飾。

程序介紹

1.SAS數據結構和程序

SAS邏輯庫下MAP路徑下有China1和China2兩個文件，其中China1中包含繪圖所需要的經緯度和省份的ID，China2中包含各省份ID對應的省份的信息(包含現用名、別名、曾用名等)。我們只需要在自己的數據中將省份數據對應到唯一識別的ID即可(圖 1)。需要注意的是，沒有數據的省份也需要對應ID，否則繪制的地圖將會缺少無ID的省份。代碼介紹如下：

圖1 SAS數據結構(左為China2，右為個人數據)

/*此處省略導入數據的過程，我們默認將數據導入到WORK下數據集名叫China*/

/*定義圖例標題*/

LEGEND1 LABEL=(“衛生機構數”);

/*定義圖例分段和其文字內容*/

PROC FORMAT;

VALUE MK LOW-10000=“0～1萬”

10000-20000=“1～2萬”

20000-30000=“2～3萬”

30000-40000=“3～4萬”

40000-HIGH=“4萬以上”;

RUN;

QUIT;

/*繪制地圖*/

PROC GMAP MAP=MAPS.CHINA DATA=CHINA; /*MAP為地圖數據，DATA為統計數據*/

FORMAT V MK.; /*將統計指標V按照上面定義的格式進行格式化*/

ID ID;/*標識出ID變量*/

CHORO V/CDEFAULT=RED CEMPTY=GREEN LEGEND=LEGEND1 LEVELS=ALL; /*CHORO為可選的幾種圖形之一，其他還包括BLOCK，PRISM等；CDEFAULT用來標識無數據區域的外邊顏色；CEMPTY用來標識無數據區域內部填充色；LEGEND連接自定義圖例標題；LEVELS=ALL顯示全部圖例分級*/

RUN;

QUIT;

如果希望將圖例按照不同顏色而非漸變色表示，讀者可通過將需要繪制的變量變成字符型變量，然后進行上述操作。

2.R數據結構和程序

我們先調用包，如果沒有安裝請先安裝。

library(sf)

library(ggplot2)

library(dplyr)

然后導入從GADM下載的地圖(CHN為中國除港澳臺地區；TWN為中國臺灣地區；HKG為中國香港地區；MAC為中國澳門地區)和數據。地圖數據中包含唯一識別的ID，變量名為GID，GID后面的數字編碼為區劃編碼，這里用省級GID_1；同時地圖數據中還包含名字、區劃等級等信息(圖 2)。

圖2 CHN_1地圖數據結構

x <-st_read(“e：/gadm36_CHN_shp/gadm36_CHN_1.shp”)

y1 <-st_read(“e：/gadm36_TWN_shp/gadm36_TWN_0.shp”)

y2 <-st_read(“e：/gadm36_HKG_shp/gadm36_HKG_0.shp”)

y3 <-st_read(“e：/gadm36_MAC_shp/gadm36_MAC_0.shp”)

y <-read.csv(“e：china.csv”)

同時把統計數據和中國(除港澳臺地區)地圖按照GID_1做左連接。

total <-left_join(x，y，by=“GID_1”)

連接完成后，按照分組新建一個group變量。

total[total$v<=10000，“group”]<-“0～1萬”

total[total$v<=20000 &total$v> 10000，“group”] <-“1～2萬”

total[total$v<=30000 &total$v> 20000，“group”] <-“2～3萬”

total[total$v<=40000 &total$v> 30000，“group”] <-“3～4萬”

total[total$v>40000，“group”] <-“4萬以上”

最后畫圖(圖例為漸變色)，其中geom_sf中data是數據，fill后面是統計數據變量的名字；之后繪制中國港澳臺地區的部分，其中colour代表輪廓線顏色，fill代表填充顏色。

ggplot()+ geom_sf(data=total，aes(fill=v))+ geom_sf(data=y1，colour=“green”，fill=“red”)+ geom_sf(data=y2，colour=“green”，fill=“red”)+ geom_sf(data=y3，colour=“green”，fill=“red”)+ labs(fill=“衛生機構數”)

可改變將fill=group，使統計地圖以不同顏色圖例替換漸變色顯示。

結果與討論

應用以上程序將中國除港澳臺地區的31個省級行政區劃的衛生機構總數繪制成統計地圖，詳見圖3。

圖3中，①③由R繪制，②④由SAS繪制。港澳臺地區由于沒有數據，全部用紅色填充綠色描邊。①②的統計量為數值型(FORMAT并不會影響到變量本身的類型)，R和SAS都默認用漸變色填充，其中②是使用FORMAT將漸變色水平定為5，否則SAS默認按比例劃分圖例區間。我們還觀察到ggplot2默認數值由小到大顏色由深到淺，SAS PROC GMAP與此相反。③④的統計量為字符型變量，SAS和R都采用了不同顏色填充的形式，在默認的情況下，ggplot2顏色更鮮艷。SAS的圖例默認位于底部，ggplot2的圖例默認位于圖外右側。ggplot2默認繪制橫縱參考線，底色灰色且保留橫縱軸及坐標，SASPROC GMAP默認模板底色白色，并無參考線橫縱軸及坐標。

圖3 SAS與R繪制的統計地圖

從布局來看，SAS默認參數下繪制的統計地圖更符合一般要求，不需要額外去除參考線和坐標軸，同時參數設置簡單，比較適合對作圖精細程度要求不高的情況。但是SAS/GRAPH模塊雖然功能強大，可語法并不緊湊，PROCGMAP中的參數只是可調整參數的冰山一角，細致調節需要翻閱PROCGMAP的幫助文件，對于SAS/GRAPH模塊的新手十分不友好。R的ggplot2包功能強大，語法連續性和可讀性非常強。但是，對于完全沒有基礎的人，去框線作修飾的工作很難入門，SAS PROC GMAP默認參數更符合基本需求。如果有一定基礎，再學習的成本較SAS/GRAPH低得多。

小 結

SAS和R在繪制統計地圖上的表現難分伯仲。對于熟悉ggplot2制圖的讀者，R的學習成本低；對于有SAS繪圖基礎的讀者，PROC GMAP是不二之選。針對全無基礎的讀者，PROC GMAP方便快捷，再結合 Adobe Photoshop、Adobe Illustrator或者Corel DRAW，很容易做出令人滿意的統計地圖。