基于GIS和RS的山東日照海岸帶1988—2018年間演化分析

苗苗，謝小平

(曲阜師范大學地理與旅游學院，日照 276826)

0 引言

海岸帶作為海、陸、氣交互作用的地帶，不僅是海洋生態系統向陸地生態系統的過渡區域，也是全球最重要的生態交錯帶，各種過程耦合多變，機制復雜，疊加人類活動，使得海岸帶成為響應全球變化最敏感的區域[1-4]。隨著改革開放的發展，我國的海岸帶開發加速，海岸帶局部侵蝕問題日益嚴重，海岸侵蝕已經與風暴潮和巨浪災害并列為我國三大海洋災害，海岸帶變化的研究也日趨受到研究學者們的重視[5-7]。中國海岸分布主要集中于遼東半島、山東半島和華南海岸3個區域，國內學者對山東半島海岸變化情況研究較多，但是對日照海岸帶時空變化的研究卻很少。山東省日照市作為山東半島藍色經濟區的重要組成部分和我國東部沿海地區的新興港口城市，近年來經濟發展迅速，海岸帶地區面臨日趨嚴峻的生態環境風險。因此，對日照海岸帶進行系統研究，了解其發展趨勢是極其重要的。

遙感(remote sensing,RS)和地理信息系統(geographic information system,GIS)技術作為新的監測手段在海岸帶動態變化監測上得到了廣泛應用[8-9]，國內學者基于RS和GIS技術對我國自北向南多個區域海岸線變化進行了研究，如姚曉靜等[10]利用衛星遙感資料與現場實測資料對海南島30 a海岸線時空變化進行了研究；陳金月[11]結合GIS和RS技術分析了珠江三角洲40 a間海岸線變遷及驅動因素；康波等[12]基于遙感和GIS技術對長島南五島30 a來海岸線時空變遷進行詳細分析。近年來，一些學者對日照海岸帶地貌和地質景觀、海岸線形態、海岸帶景觀格局變化過程進行了研究[13-15]，但只是對日照海岸帶的局部特性做了研究。因此，本文綜合利用RS與GIS技術，對日照海岸帶1988—2018年間近30 a的時空變化特征及其演變的驅動力進行了探究。

1 研究區概況及數據源

1.1 研究區概況

山東省日照市海岸帶位于黃海中部，北起日照市東港區兩城河口，與青島膠南市海岸接壤，南至嵐山區繡針河入?？?，與江蘇贛榆縣相連，西至204國道，東至20 m水深線附近(圖1)。沿岸發育平直的基巖砂礫質海岸，中部的奎山、南部的嵐山兩段為剝蝕丘陵，形成突出的岬角濱海陡崖和水下巖灘；南北部的剝蝕平原和現代海岸線之間為寬廣的沙壩潟湖沉積體系。區域斷裂構造以NE向的日照—青島斷裂和NW向的梭羅樹斷裂為主，海岸呈NE-SW走向。境內地表水水體發育，河流為內陸河入海段，流向多為NW-SE向，自北而南主要有白馬河、兩城河、潮白河、付疃河、巨峰河、龍王河和繡針河，除傅疃河等為永久性河流外，其他多屬季節性河流。由于該地區地理位置優越，形成了眾多的天然良港，海岸線上的主要港灣包括石臼灣、佛手灣與日照港、嵐山港組成日照港群。日照海岸帶是日照經濟、文化發展中心，現已發展成集臨海工業、港口貿易、商貿、金融、交通運輸、旅游度假等多功能、外向型、綜合性、現代化的海岸帶[16-18]。

圖1 研究區位置Fig.1 Location of study area

1.2 數據源及其預處理

本文選取的數據包括：研究區1988年、1998年、2008年和2018年4期遙感影像數據，來自地理空間數據云網站(http://www.gscloud.cn/)，其中2018年影像數據為Landsat8 OLI，2008年、1998年和1988年影像數據為Landsat5 TM；研究區各年潮汐數據來源于中國港口網(http://www.chinaports.com)；參考數據包括日照市行政區劃、1∶50 000基礎地理信息數據、數字高程模型(digital elevation model,DEM)數據，以10 a為間隔來研究日照海岸帶時空變化特征。采用ENVI5.2軟件數據處理平臺，對研究區4期遙感影像依次進行輻射校正、幾何糾正、鑲嵌、裁剪等數據預處理操作，以消除空間上和輻射物理量上的畸變或噪聲，得到研究區各時相遙感影像圖[19-20]。

2 研究方法

2.1 海岸線類型的劃分及提取

海岸線是指平均大潮高潮時水陸分界的痕跡線，本文根據海岸形態、成因、物質組成和發展階段等特征將日照市海岸線劃分為基巖海岸線、砂質海岸線、淤泥質海岸線、人工海岸線4類，使用紅光、近紅外和短波紅外1波段合成遙感影像，根據不同海岸線的紋理、形狀、色調等特征建立各海岸線的解譯標志(表1)[21]。

表1 日照市海岸線類型Tab.1 Coastline type in Rizhao City

在提取海岸線前，先利用濾波器突出各類型海岸線的水陸分界線[22-23]。采用修正歸一化差異水體指數(modified normalized difference water index，MNDWI)提取Landsat5遙感影像的水陸分界線，表達式為：

MNDWI=(G-MIR)/(G+MIR)，

(1)

式中：MIR為中紅外波段的反射率；G為綠光波段的反射率。采用自動提取水體指數(automatic extraction water index，AWEI)提取Landsat8遙感影像的水陸分界線，表達式為：

AWEI=B+2.5G-1.5(NIR+S1)-0.25S2，

(2)

式中：B為藍光波段的反射率；NIR為近紅外波段的反射率；S1和S2分別為2個短波紅外波段的反射率。運用式(1)—(2)獲得4個時相的水體指數影像，分別設定 MNDWI和AWEI的取值范圍來提取水陸分界線[24-26]。由于衛星過境時刻獲取的海岸線影像不能保證其為平均大潮高潮線，而海岸線的位置又易受潮汐等因素的影響，故需要對瞬時水邊線進行潮位校正[27-29]。人工海岸和基巖海岸水陸分明、坡度大而陡，潮汐影響小，不需要做潮位校正，僅對砂質海岸和淤泥質海岸做潮位校正。本文采用外牽引法[30]進行潮位校正，根據衛星過境時潮位和日最高潮位計算出潮位差，再由潮位差和海岸平均坡度計算外牽引距離(L)，計算公式為：

L=(H2-H1)/tanα，

(3)

式中：H1為瞬時潮位高度；H2為最高潮位高度；α為海岸平均坡度(°)。海岸線位置計算的原理如圖2所示，將需要校正的水邊線向前移動距離L就可以得到真正意義的海岸線位置。

圖2 海岸線位置計算原理Fig.2 Theory of counting coastline position

2.2 海岸帶土地利用時空變化研究方法

海岸帶是指海陸之間相互作用的地帶，是海岸線向陸海兩側擴展一定寬度的帶狀區域，包括陸域與近岸水域，對于其范圍，至今尚無統一的界定。本文規定的海岸帶調查范圍為：由海岸線向陸地方向延伸10 km，向海方向延伸2 km。利用各時相海岸線向陸區域做距離為10 km的緩沖區，向海區域做距離為2 km的緩沖區，以緩沖區為掩模裁剪各時相遙感影像。根據《土地利用現狀分類標準》的中的一級分類體系和研究區特征，將研究區土地利用類型分為水體、林地、耕地/草地和建設用地(主要為鹽田和養殖池)。在ENVI5.2軟件支持下，對各時相影像采用最大似然法與傳統目視解譯相結合的方法進行監督分類，以Google地圖與2018年9—10月實地考察數據為驗證數據，對分類結果進行精度驗證[31-33]。

2.3 研究方法

在ArcMap10.2軟件支持下，基于研究區DEM數據，運用水文分析工具，依次進行洼地填充、流向分析、匯流累積量計算、提取河流網絡柵格，進行流域分析，生成河網。把集水閾值設為10 000的河網圖層與研究區各時相海岸線圖層疊加，分析河口區海岸線變化特征。再結合各時相入海河流遙感影像解譯標志，分析1988—2018年30 a間研究區入海河流河口區演化特征。

3 結果與分析

3.1 海岸線演變分析

3.1.1 海岸線潮位校正及精度驗證

海水在天體引潮力作用下產生的潮汐現象具有規律性，日照海岸帶潮汐屬于半日潮型，即一日有2次高潮和低潮[34]。本文獲取了研究區2017—2019年6月的大潮日的潮高數據，并用2017—2019年各個時刻潮位數據擬合潮位曲線(圖3)。

圖3 2017—2019年潮位擬合曲線Fig.3 The fit curve of tidal level in 2017—2019

由圖3可知，不同年份潮位曲線具有相似的走勢，對高低潮時刻和潮高做相關性分析，相關系數均達到0.99以上。計算圖4中任意2條曲線同一時刻的潮位高度差，可以得到2017—2019年時刻的平均潮位高差為30.23 cm，根據三角函數關系計算此誤差在水平方向上為7.5 m，而衛星分辨率為30 m，誤差在合理范圍內，故可利用潮汐規律可以對缺失數據的年份做潮位修正。利用研究區DEM圖像生成坡度圖(圖4)，計算可得從水邊線向陸地推80 m的范圍內平均坡度為2.3°，通過計算可得潮位校正外牽引距離(表2)。

圖4 研究區坡度圖Fig.4 Slope map of the study area

表2 潮位數據以及計算結果Tab.2 Tidal level data and calculation results

由表2可知，為得到研究區各時相正確的海岸線，需將獲得的1988年、1998年、2008年、2018年的未利用的淤泥質海岸和砂質海岸的瞬時水邊線分別向前移動88 m，83 m，8 m和70 m，獲得各時相完整的海岸線。為了驗證提取效果，首先在各時相原圖像不同類型海岸線上分別隨機取40個像素點，然后在提取出的海岸線圖像上找出相對應的像素點，比較二者之間的位置變化。如果2個像素點的位置重疊則表示該點海岸線提取效果好；如果2個像素點相鄰則表示該點海岸線提取效果在誤差允許范圍之內；否則，則表示該點海岸線提取失敗。經計算可獲得不同時相海岸線像素點位置的檢驗結果(表3)[35-36]。

表3 像素點位置檢驗結果Tab.3 The veracity of pixel position

從檢驗結果可以看出基巖海岸和人工海岸的提取效果最好，砂質海岸提取效果較差，這是因為砂質海岸海岸線平直，易產生點的位移；而基巖海岸、砂質海岸和未開發淤泥質海岸的海岸線相對曲折，海岸線位置變化不大。

3.1.2 海岸線動態演變分析

把各時相海岸線疊加顯示(圖5)，通過ArcMap的數據統計功能提取出各類型海岸線的長度(表4)，綜合分析不同類型海岸線時空變化情況。結合圖5和表4分析可知，1988—2019年30 a間日照海岸線變化劇烈，海岸線演變整體趨勢為向海推進。海岸線長度及增速均呈現逐年增加的趨勢，海岸線長度總體增加達49.17 km，其中，1988—1998年增加12.1 km，年變化率為1.21 km/a；1998—2008年增加22.4 km，年變化率為2.24 km/a；2018年比2008年增加14.8 km，年變化率為1.48 km/a，3個時期均保持較快增長，1998—2008年增速最快。1988—1998年間，除日照港附近的海岸線向海推進外，其他地區海岸線均后退；1998—2008年，海岸線變遷顯著，隨著日照港到傅疃河口沿岸建設幅度增大，海岸線向海方向增加明顯；2008—2018年，韓家營子到繡針河口段為嵐山港建設開發區，此段海岸線變化劇烈。在不同類型海岸線中，人工海岸線所占比重呈逐年上升趨勢，30 a內凈增加85.72 km，到2018年比重已達70.04%，海岸工程、養殖池是人工海岸線的主要形式，如中部的奎山和南部的嵐山兩段目前改造為碼頭、養殖池，南部的沙壩外普遍修筑了人工堤壩，沙壩間的潟湖區被改造為養殖池；砂質海岸線分布較多，受海水侵蝕作用小，使沿海地區向陸面積變化小，但是1988年以后，砂質海岸線長度呈逐年遞減趨勢，30 a間減少12.84 km；基巖海岸線相間分布所占比例很少，主要分布在任家臺咀、石臼咀、嵐山咀和奎山咀，近幾年基巖海岸線變化很??；淤泥質海岸線保持緩慢下降趨勢，主要分布在河流入?？谔?。

(a)1988—1998年 (b)1998—2008年

表4 各時期各類型海岸線長度統計Tab.4 Length and proportion statistics of coastline in each period

3.2 海岸帶土地利用變化過程

3.2.1 海岸帶土地利用分類精度驗證

在ENVI5.2軟件支持下，以Google地圖與2018年9—10月實地考察數據為驗證數據，對進行監督分類后的分類結果進行精度驗證，以總體分類精度和Kappa系數來表示分類精度(表5)。由于本文分類系統簡單，主要分為了4類，且各地類在光譜上有明顯區別，故而分類精度較高，從表5中可以看出其Kappa系數均在0.9以上，表明分類結果可以滿足研究目的。利用Post Classification工具消除噪聲，獲得研究區各時相海岸帶土地利用圖(圖6)[37-38]。

表5 分類精度Tab.5 Classification accuracy

3.2.2 海岸帶土地利用變化結果分析

得到研究區各時相土地利用分類結果后，利用分類后處理工具統計并計算各地類所占面積及其所占百分比(表6)，分析其結構組成的時空變化。由表6可知，1988—2018年30 a間建設用地面積呈現逐年增加的趨勢，由1988年的213.77 km2增加到2018年的413.93 km2，其所占總地類面積的比重由21.73%增加到39.43%；30 a間耕地和草地面積則逐年減少，由1988年的445.50 km2減少到2018年的287.03 km2，其所占總地類面積的比重從45.29%減少到27.34%，嵐山區、東港區附近的大部分耕地轉換為建設用地；水體和林地分布時空變化起伏不大。

(a)1988年 (b)1998年

表6 1988—2018年日照海岸帶各土地利用類型面積及比例Tab.6 Land cover type area and proportion from 1988 to 2018 in Rizhao coastal zone

用2018年海岸帶緩沖區作為感興趣區裁剪1988年和2018年遙感影像，進行監督分類，制作1988—2018年轉移矩陣表(表7)，分析日照海岸帶30 a間各地類轉化情況。由表7可知，1988—2018年各地類相互轉化情況如下：有60.82 km2建設用地轉化為耕地/草地，占其轉化總面積的27.69%，相對應的有185.51 km2耕地/草地轉化為建設用地，占其轉化總面積的比重最大，為39.60%，兩者間凈轉移量為124.69 km2；其次，水體轉化為建設用地的面積占其轉化為其他地類的總面積的比重最大，轉為建設用地的面積84.63 km2，占其轉化總面積的12.97%；林地轉為建設用地的面積也最大，為33.8 km2，占其轉化總面積的31.97%?？傮w來看，建設用地在各地類使用中的重要程度逐漸加大。

表7 1988—2018年日照海岸帶土地利用面積轉移矩陣Tab.7 Transition matrix of land over types in Rizhao coastal zone from 1988 to 2018

3.3 入海河流對海岸帶的影響

研究區沿岸的主要入海河流有兩城河、付疃河和繡針河，根據入海河流河口區特征(表8)和河口區海岸線位置變化(圖7)對日照海岸帶河口區演化進行分析。結合表8和圖7分析可知，1988—1998年間，兩城河口段海岸線處于蝕退狀態，1988年后兩城河河口南側的沙壩潟湖水域出現了大面積圍海和向陸開挖養殖池的活動，導致潟湖消失；1998—2008年間，在兩城河下游和河道內建設蝦池、育苗廠，大面積占用河道；2008—2018年間，對占用河道的設施進行清除，并建立地下水庫攔蓄和調節地下水流，向海輸出的泥沙持續淤積，在河口區重新堆積砂體[39]。

表8 河口區特征Tab.8 The characteristics of the estuary

(a)1988—1998年兩城河河口區 (b)1988—1998年付疃河河口區 (c)1988—1998年繡針河河口區

(g)2008—2018年兩城河河口區 (h)2008—2018年付疃河河口區 (i)2008—2018年繡針河河口區

1988—1998年間，付疃河有一半泥沙被日照水庫攔淤，河流供沙大幅度減少，淤進速度減緩，河口附近海岸線還被用于建設養殖池，轉變為人工海岸；1998—2008年間，在付疃河流域境內分布著多處采沙點，嚴重影響了付疃河的泄洪能力；2008年以來，對付疃河進行了綜合治理，使其向海輸出的泥沙持續淤積，向海伸出形成突灘[40-41]。

1988—1998年間，嵐山港相繼建成防波堤和多個噸級泊位，攔阻了大量泥沙，使河口區蝕退；1998—2008年間，又在河口灣潟湖中修建了人工砌石護岸，岸灘泥沙被阻擋，人工開采海灘沙和興修水庫也使入海泥沙大為減少，導致河口東岸的長條狀沙嘴消失；2008—2018年間，由于政策保護河口濕地，河口三角洲處轉蝕為淤，河口處形成向陸微凹的小型三角洲沖積平原[42-43]。

4 討論

研究區海岸線變化與沿海土地利用之間是有直接聯系的，海岸線劇烈變化的區域也是海岸帶土地利用變化劇烈的區域。30 a間，研究區海岸線變化劇烈，北岸海岸線較為完整，變化不大；中岸和南岸海岸線變化較大，海岸線破碎。海岸線長度持續增加，自然海岸線持續減少，人工海岸線持續增加，其中變化較大的區域主要集中在日照港、嵐山港及河流入?？诟浇?，如日照港和嵐山港相繼建成防波堤和多個噸級港作泊位和生產泊位，作為建設開發區；在河口區出現了大面積圍海和向陸開挖養殖池的活動、分布著多處采沙點。30 a間，研究區主要土地利用類型為建設用地、耕地/草地，建設用地面積呈現逐年增加的趨勢，養殖池和鹽田增加面積所占比重較大，區域內嵐山區、東港區附近的大部分耕地轉換為建設用地，且建設用地擴張呈現中心擴張趨勢；隨著入?？诤恿鲾y帶的泥沙逐漸減少及大面積圍海和向陸開挖養殖池、采沙活動的增加，沿岸水體減少，相應建設用地增多；在海岸線變化較小的任家臺咀、石臼咀、嵐山咀、奎山咀地區，沿岸土地利用變化也相應較小。由此可知，在海岸線變化劇烈的時間段內沿岸土地利用速率也相應變化。

5 結論

本文系統總結了日照海岸帶影像解譯標志，并結合前人研究成果，利用RS和GIS技術，分析了日照海岸帶時空變化特征及其驅動力，得出如下結論：

1)1988—2018年的30 a間日照海岸線變化劇烈，海岸線演變整體趨勢為明顯的向海推進。海岸線長度及增速均呈現逐年增加的趨勢，3個時期均保持較快增長。在不同類型海岸線中，人工海岸線所占比重呈逐年上升趨勢；砂質海岸線逐年減少；基巖海岸線相間分布于砂質海岸之間，所占比例很少，幾乎不發生變化；淤泥質海岸線始終變化不大，保持緩慢下降趨勢。

2)30 a間日照海岸帶各土地利用類型中，建設用地面積及其所占比重逐年遞增，且養殖池和鹽田所占比重較大；耕地和草地面積及其所占比重逐年減少；林地和水體面積時空變化不明顯。轉移矩陣表顯示各地類轉化以耕地/草地轉化為建設用地為主要趨勢。

3)30 a間河口區海岸線動態變化明顯，1988—1998年間河口區進入蝕退階段，沙嘴前端變淺，甚至部分解體，河口附近岸段基本上轉為人工海岸線；1998—2008年間各入海河流河口地貌保持相對穩定，說明河口區水動力對沉積物的沖蝕和淤積基本保持動態平衡；2008—2018年間河口區恢復以堆積狀態為主，在河流的兩側分布有規模較小的帶狀沖積平原。