李毅， 程麗娜， 魯瑩瑩， 張博淳， 于森， 賈明明

(1.吉林建筑大學(xué)測(cè)繪與勘查工程學(xué)院，長(zhǎng)春 130118； 2.中國(guó)科學(xué)院東北地理與農(nóng)業(yè)生態(tài)研究所濕地生態(tài)與環(huán)境重點(diǎn)實(shí)驗(yàn)室,長(zhǎng)春 130102； 3.吉林大學(xué)地球科學(xué)學(xué)院,長(zhǎng)春 130061； 4.長(zhǎng)春新區(qū)北湖英才學(xué)校，長(zhǎng)春 130000)

0 引言

全面掌握濱海灘涂變化情況的關(guān)鍵在于灘涂斑塊信息的高精度完整提取。其困難來(lái)源于衛(wèi)星過(guò)境時(shí)潮位的不確定性，一定時(shí)間段內(nèi)的遙感影像并不能提取出準(zhǔn)確完整的灘涂信息。針對(duì)這一問(wèn)題，國(guó)內(nèi)外學(xué)者進(jìn)行了一系列的相關(guān)研究。Murray等[8]利用區(qū)域潮汐模型獲取海水高潮與低潮時(shí)的遙感影像，以確定瞬時(shí)水邊線，并進(jìn)行疊加，完成了東亞濱海的潮間帶濕地地圖的繪制； Wang等[9]基于GEE(Google Earth Engine)使用Landsat時(shí)序影像，通過(guò)計(jì)算像元的年度水體頻率值，確定最高潮位線和低水線的水體頻率閾值范圍來(lái)識(shí)別和劃分潮間帶，然后根據(jù)AOIs(areas of interest)的植物頻率分布閾值進(jìn)行灘涂和植被的分類； Jia等[10]利用最大光譜指數(shù)合成算法(maximum spectral index composite, MSIC)和最大類間方差算法(Otsu algorithm, OTSU)，基于GEE使用Sentinel-2數(shù)據(jù)繪制了中國(guó)10 m空間分辨率灘涂分布圖； 智超等[11]借助GEE構(gòu)建時(shí)序遙感指數(shù)，基于樣本分析遙感指數(shù)特征分類水體與植被，根據(jù)潮間帶不同植物物候及動(dòng)力特征進(jìn)行潮間帶濕地提取； Zhao等[12]通過(guò)計(jì)算長(zhǎng)時(shí)序影像中每個(gè)像元的分位數(shù)，分別選取第5和第95分位數(shù)像素合成低、高潮位基準(zhǔn)面影像，然后使用自適應(yīng)二分類器對(duì)基準(zhǔn)面影像進(jìn)行水體和陸地的二分類，根據(jù)潮灘高潮時(shí)淹沒(méi)、低潮時(shí)露出水面的特性，繪制了中國(guó)南方的潮灘地圖； 李振等[13]建立各類土地利用類型的疊合光譜圖，通過(guò)建立決策樹(shù)法進(jìn)行地物分類，繪制了膠州灣海岸帶土地利用圖； 張康永[14]基于GEE平臺(tái)，通過(guò)分析灘涂光譜信息的周期性變化規(guī)律，選取分類特征進(jìn)行影像合成，綜合隨機(jī)森林分類法和面向?qū)ο蟮膱D像分割技術(shù)，開(kāi)展了杭州灣以北濱海灘涂30 a來(lái)的遙感監(jiān)測(cè)和驅(qū)動(dòng)力分析研究。然而目前，尚缺乏針對(duì)萊州灣地區(qū)灘涂變化的相關(guān)研究。

基于以上分析，考慮到數(shù)據(jù)的獲取性、一致性及空間分辨率，本研究借助GEE平臺(tái)構(gòu)建Landsat高質(zhì)量密集時(shí)序影像堆棧，基于MSIC和OTSU算法實(shí)現(xiàn)1990年、2000年、2010年和2020年萊州灣地區(qū)灘涂信息的高精度提取，結(jié)合灘涂周邊土地覆被變化數(shù)據(jù)集從自然和人為因素方面探討30 a間萊州灣灘涂資源的數(shù)量變化以及時(shí)空分布特征，分析其演變驅(qū)動(dòng)機(jī)制，以期為該地區(qū)科學(xué)開(kāi)發(fā)和合理保護(hù)灘涂資源提供技術(shù)參考和數(shù)據(jù)基礎(chǔ)。

1 研究區(qū)概況及其數(shù)據(jù)源

1.1 研究區(qū)概況

萊州灣位于山東省北部，地理坐標(biāo)為E118°39′～120°35′，N36°55′～38°00′，主要包括東營(yíng)市、濰坊市和煙臺(tái)市部分地區(qū)，內(nèi)有以黃河為主的多條河流以及龍口港等重要港口，坐落著黃河三角洲國(guó)家自然保護(hù)區(qū)，擁有豐富的灘涂資源，區(qū)內(nèi)氣候?yàn)闇貛Ъ撅L(fēng)氣候。本研究依據(jù)海岸帶環(huán)境地質(zhì)調(diào)查規(guī)范(1∶25萬(wàn))，以萊州灣地區(qū)海岸線為中心，向陸海兩側(cè)做10 km緩沖區(qū)作為研究區(qū)邊界(圖1)。研究區(qū)內(nèi)受不規(guī)則半日潮影響，平均潮差為1.01 m，最大潮差為2.97 m。

圖1 研究區(qū)地理位置Fig.1 Geographical location of the study area

1.2 數(shù)據(jù)源及處理

本研究以1990年、2000年、2010年和2020年為研究年份，基于GEE云計(jì)算平臺(tái)，便捷調(diào)用覆蓋研究區(qū)域目標(biāo)年份及其前后1 a(時(shí)間窗口為3 a)所有符合云量要求的Landsat TM和OLI地表反射率影像集合，構(gòu)建初始密集時(shí)間序列影像堆棧。利用對(duì)應(yīng)的去云函數(shù)進(jìn)行云掩模以獲取有效觀測(cè)像素，構(gòu)成研究區(qū)域目標(biāo)年份的高質(zhì)量密集時(shí)間序列影像堆棧。經(jīng)統(tǒng)計(jì)，研究區(qū)域目標(biāo)年份每個(gè)像素的有效觀測(cè)數(shù)量空間分布及整體分布比例如圖2所示。2000年、2010年和2020年濱海灘涂以及其周邊土地覆蓋類型(內(nèi)陸地物)驗(yàn)證樣本均在Google Earth高空間分辨率影像中直接獲取。其中，采集的灘涂驗(yàn)證樣本均選取低潮時(shí)期位于灘涂邊緣的點(diǎn)。1990年驗(yàn)證樣本通過(guò)文獻(xiàn)檢索，查找歷史地圖等方式獲取。最終，1990年驗(yàn)證樣本共1 022個(gè)(其中灘涂232個(gè)，其他土地覆被790個(gè))； 2000年驗(yàn)證樣本共1 136個(gè)(其中灘涂251個(gè)，其他土地覆被885個(gè))； 2010年驗(yàn)證樣本共1 127個(gè)(其中灘涂281個(gè)，其他土地覆被846個(gè))； 2020年驗(yàn)證樣本共1 068個(gè)(其中灘涂231個(gè)，其他土地覆被837個(gè))。

(a) 1990年 (b) 2000年圖2-1 研究區(qū)域目標(biāo)年份有效觀測(cè)數(shù)量空間分布及比例Fig.2-1 Interannual distribution and proportion of each image stack

(c) 2010年 (d) 2020年圖2-2 研究區(qū)域目標(biāo)年份有效觀測(cè)數(shù)量空間分布及比例Fig.2-2 Interannual distribution and proportion of each image stack

2 研究方法

2.1 技術(shù)路線

由于潮汐的不確定性，為了提取完整的灘涂斑塊信息，本研究結(jié)合前人的相關(guān)研究[1,10]。首先，通過(guò)構(gòu)建Landsat高質(zhì)量密集時(shí)序影像堆棧代替單景影像獲取遙感信息的方法； 其次，基于MSIC算法獲取研究區(qū)目標(biāo)年份最高、最低潮影像； 最后，利用OTSU算法進(jìn)行圖像二值化分類，實(shí)現(xiàn)1990年、2000年、2010年和2020年萊州灣地區(qū)灘涂信息的高精度提取。結(jié)合模糊邏輯分析方法(fuzzy-based segmentation parameter，F(xiàn)bSP)和面向?qū)ο蠓治黾夹g(shù)獲取灘涂的周邊土地覆被變化數(shù)據(jù)集，在此基礎(chǔ)上探討近30 a間萊州灣灘涂資源的數(shù)量變化以及空間分布特征，分析該地區(qū)灘涂演變驅(qū)動(dòng)機(jī)制，具體流程如圖3所示。

圖3 技術(shù)流程Fig.3 Technical flowchart

GEE是谷歌公司推出的一個(gè)基于云計(jì)算的全球尺度遙感數(shù)據(jù)處理平臺(tái)，其API庫(kù)中提供了MSIC方法。該方法是對(duì)影像堆棧中同一位置的所有像元，逐像元計(jì)算指定光譜指數(shù)值，并利用內(nèi)置的排序函數(shù)進(jìn)行排序，以最大值像元作為合成影像中對(duì)應(yīng)位置的像元。對(duì)所有位置的像元重復(fù)該過(guò)程即可得到指定光譜指數(shù)最大值合成影像。

OTSU算法又名大津算法，是Otsu[15]提出的一種自動(dòng)選取最優(yōu)閾值進(jìn)行圖像分割的方法。該方法根據(jù)圖像灰度值，自動(dòng)選取使類間方差最大、類內(nèi)方差最小的最優(yōu)分割閾值，對(duì)圖像進(jìn)行二值化分割，實(shí)現(xiàn)圖像自動(dòng)化的二分類。

面向?qū)ο蟮姆诸愂峭ㄟ^(guò)對(duì)影像分割，基于分割得到的對(duì)象為基本元素進(jìn)行分類的方法。在分割得到的單個(gè)對(duì)象中，所有的像元均可認(rèn)為是同質(zhì)像元[16]。該方法不僅考慮像元的光譜特征，還考慮到影像的紋理和質(zhì)地等因素。因此，面向?qū)ο蟮姆诸惤Y(jié)果可以很好地避免噪聲的影響，具有良好的整體性[17]。在本研究中，為了提高分割精度和效率，采用FbSP自動(dòng)確定多尺度分割的最佳分割參數(shù)。

依據(jù)CCUS技術(shù)的研發(fā)進(jìn)展，預(yù)計(jì)2035年，第一代捕集技術(shù)的成本及能耗與目前相比降低15%～20%；2045年前后，第二代捕集技術(shù)實(shí)現(xiàn)商業(yè)應(yīng)用，成本及能耗將比第一代技術(shù)降低10%～15%；2050年，CCUS技術(shù)實(shí)現(xiàn)廣泛部署，建成多個(gè)CCUS產(chǎn)業(yè)集群。

2.2 灘涂及周邊土地覆蓋類型提取

2.2.1 灘涂提取

本文將灘涂定義為漲潮時(shí)被海水浸沒(méi)，退潮時(shí)露出海面的砂質(zhì)、淤泥或軟泥質(zhì)的無(wú)植被海岸淺灘。由于水體像素的改進(jìn)型歸一化差異水體指數(shù)(modified normalized difference water index，mNDWI)值明顯高于非水體像素，因此選取mNDWI指數(shù)進(jìn)行MSIC合成，一方面能夠突出水體與非水體像素的特征差異，同時(shí)能夠更好地獲取目標(biāo)年份水體信息的最大覆蓋區(qū)域，即最高潮影像； 另外，水體像素的歸一化差分植被指數(shù)(normalized difference vegetation index，NDVI)值低于非水體像素(植被、灘涂)，使用NDVI指數(shù)獲取最低潮影像可顯著增強(qiáng)植被、灘涂與水體間的特征差異，還可有效避免因潮汐的不確定性對(duì)灘涂提取產(chǎn)生的干擾[1]。因此，本研究采用MSIC算法，分別選用mNDWI指數(shù)、NDVI指數(shù)從Landsat高質(zhì)量密集時(shí)序影像堆棧中合成目標(biāo)年份的最高、最低潮影像。利用OTSU算法對(duì)最高潮影像進(jìn)行二值分割，實(shí)現(xiàn)水體和非水體的分離，通過(guò)保留最大面積斑塊的方法獲取目標(biāo)年份的最大海水面，對(duì)最低潮影像進(jìn)行掩模提取得到最大潮水淹沒(méi)區(qū)影像。為了去除淹沒(méi)區(qū)植被的影響，再次采用OTSU算法進(jìn)行二值分類，實(shí)現(xiàn)灘涂信息的高精度提取。

2.2.2 灘涂周邊土地覆蓋類型提取

針對(duì)萊州灣地區(qū)潮間帶的景觀特點(diǎn)，從研究的目的和遙感制圖的角度出發(fā)，將研究區(qū)中除灘涂以外的地物類型分為： 內(nèi)陸植被、濱海植被、養(yǎng)殖池/鹽田、人工表面、農(nóng)田、裸地、內(nèi)陸水體、海水。

本研究基于eCognition9.2平臺(tái)，通過(guò)FbSP方法自動(dòng)確定目標(biāo)年份單景最優(yōu)影像多尺度分割的最佳分割參數(shù)，在獲取最佳分割尺度下基于面向?qū)ο蠹夹g(shù)目視解譯獲取灘涂周邊土地覆被圖。并以2020年的灘涂周邊土地覆被圖作為底圖參考制作1990年、2000年、2010年的灘涂周邊土地覆被圖。

2.3 精度驗(yàn)證

本研究基于混淆矩陣，計(jì)算總體分類精度和Kappa系數(shù)，用以評(píng)價(jià)灘涂和周邊地物的分類精度。混淆矩陣是根據(jù)像元的地表真實(shí)情況與該像元的分類情況進(jìn)行比較計(jì)算得出的。總體分類精度是被正確分類的像元數(shù)與總像元數(shù)的比值，是分類結(jié)果和檢驗(yàn)數(shù)據(jù)一致性的評(píng)價(jià)指標(biāo)。Kappa系數(shù)采用一種多元技術(shù)，綜合考慮混淆矩陣所有因素的一種精度評(píng)價(jià)指標(biāo)。通常認(rèn)為介于0.61～0.80則具有高度一致性； 高于0.80則幾乎完全一致[18]。

3 結(jié)果與分析

3.1 灘涂及周邊地物分類精度

利用選取的驗(yàn)證點(diǎn)對(duì)4期分類結(jié)果分別進(jìn)行精度驗(yàn)證，計(jì)算得到了各研究年份的總體分類精度和Kappa系數(shù)。表1是根據(jù)選取的驗(yàn)證樣本生成的研究區(qū)2020年灘涂及周邊地物遙感解譯結(jié)果的混淆矩陣，結(jié)果表明2020年灘涂及周邊土地覆被分類結(jié)果總體精度為94%，Kappa系數(shù)為0.92。1990年、2000年和2010年分類結(jié)果的總體分類精度分別為88%，91%和90%，Kappa系數(shù)分別為0.85，0.89和0.88，表明分類結(jié)果與地表真實(shí)情況具有高度的一致性，各年份的分類精度均能滿足進(jìn)一步分析的要求。

表1 2020年灘涂及周邊地物土地覆被分類結(jié)果的混淆矩陣Tab.1 Confusion matrix of land cover classification results of tidal flats and surrounding land covers in 2020

3.2 萊州灣濱海灘涂時(shí)空變化分布特征

萊州灣濱海地區(qū)1990—2020年土地覆蓋變化情況如圖4所示，疊加分析可知，2020年萊州灣地區(qū)的灘涂總面積為822.38 km2，較之1990年的1 365.84 km2減少了543.46 km2，縮減幅度約為40%。近30 a期間，該地區(qū)的灘涂面積呈現(xiàn)持續(xù)減少的趨勢(shì)。1990—2000年灘涂面積減少了150.94 km2； 2000—2010年間的減少幅度最大，共減少304.78 km2； 2010—2020年灘涂的縮減速率明顯減緩，10 a間僅減少了87.74 km2。由于萊州灣地區(qū)涵蓋范圍較廣，包含濰坊和煙臺(tái)等多個(gè)市，各市對(duì)于灘涂濕地的保護(hù)力度和開(kāi)發(fā)程度均不同，從而導(dǎo)致了萊州灣地區(qū)各岸段的灘涂變化狀況之間也有所差異。根據(jù)圍墾程度和岸線侵占情況，本研究將研究區(qū)分為東營(yíng)市岸段、濰坊市岸段和煙臺(tái)市岸段。從整體上來(lái)看，東營(yíng)市岸段灘涂面積占比最大，濰坊市岸段次之，煙臺(tái)市岸段最??； 縮減幅度最大的是濰坊市岸段，為323.56 km2，最小的是煙臺(tái)市岸段，僅82.01 km2。東營(yíng)市岸段的灘涂在3個(gè)岸段中是最廣闊的，占萊州灣灘涂總面積的一半以上，該岸段在1990—2000年的減少幅度較小，為51.82 km2，在2000—2010年的減幅最大，達(dá)到了131.40 km2，在2010—2020年有小幅度的增加，增幅為45.32 km2。該岸段除黃河口地區(qū)附近灘涂海一側(cè)邊界整體上不斷向海推進(jìn)，且有較大的推進(jìn)幅度外，其余地區(qū)灘涂海一側(cè)邊界均不斷后移，陸側(cè)邊界持續(xù)向海推進(jìn)。濰坊市岸段的灘涂整體上寬幅變化不大，呈現(xiàn)較規(guī)則的條狀分布，在1990—2000年、2000—2010年、2010—2020年灘涂的面積減幅先增大后減小，同樣在2000—2010年達(dá)到峰值，為144.31 km2，灘涂面積呈持續(xù)減少的趨勢(shì)，灘涂海一側(cè)邊界緩慢的向陸推進(jìn)的同時(shí)，陸一側(cè)邊界不斷向海推移。煙臺(tái)市岸段的灘涂寬幅變化較大，呈現(xiàn)西寬北窄的形式，是灘涂面積最小的一個(gè)岸段，在1990—2000年、2000—2010年和2010—2020年灘涂面積的減少幅度逐漸升高，分別為10.65 km2，29.09 km2和42.27 km2，灘涂面積同樣呈持續(xù)減少的趨勢(shì)，整體邊界變化與濰坊市岸段相同。

(a) 1990年 (b) 2000年

(c) 2010年 (d) 2020年圖4 1990—2020年萊州灣地區(qū)潮間帶土地覆被Fig.4 Land cover of intertidal zone in Laizho Bay from 1990 to 2020

3.3 灘涂與周邊土地覆被的相互轉(zhuǎn)化

1990—2020年灘涂與周邊地物相互轉(zhuǎn)化的情況如圖5所示，灘涂面積呈凈流失狀態(tài)，灘涂持續(xù)被侵占，最主要的侵占類型為養(yǎng)殖池/鹽田，萊州灣沿岸除黃河口地區(qū)以外，侵占灘涂的情況均十分嚴(yán)重。其次是濱海植被，雖總體上不斷減少，但仍有相當(dāng)一部分灘涂轉(zhuǎn)化為了濱海植被。灘涂與海水、內(nèi)陸水體(河流)的相互轉(zhuǎn)化量雖然較大，但面積總體變化上幾乎持平，一方面是海平面上升使灘涂邊界縮退，另一方面是河流泥沙輸入使灘涂邊界向海擴(kuò)張，二者對(duì)灘涂的影響呈均勢(shì)，整體上對(duì)灘涂面積變化影響不大。人工表面對(duì)灘涂的侵占雖然對(duì)比養(yǎng)殖池/鹽田類型少，在1990—2000年間幾乎沒(méi)有侵占情況，然而，在2000—2010年和2010—2020年2個(gè)時(shí)段中侵占面積劇烈增加，有持續(xù)增長(zhǎng)的趨勢(shì)，空間上主要發(fā)生在濰坊市岸段，主要表現(xiàn)為圍海圈地。其他地物如裸地和內(nèi)陸植被等，雖然與灘涂也存在相互轉(zhuǎn)化，但轉(zhuǎn)化量較小。

圖5 1990—2020年灘涂與周邊地物相互轉(zhuǎn)化Fig.5 Mutual transformation tidal flats and surrounding features from 1990 to 2020

1990—2020年間灘涂與周邊地物的轉(zhuǎn)化量如表2所示?？梢?jiàn)，養(yǎng)殖池/鹽田的面積急劇擴(kuò)張，在30 a間增長(zhǎng)了多達(dá)672.27 km2，增長(zhǎng)率到達(dá)了45%以上。其中侵占灘涂面積為414.20 km2，灘涂面積的減少量對(duì)于養(yǎng)殖池/鹽田面積增長(zhǎng)的貢獻(xiàn)率達(dá)到了62%，是最主要的貢獻(xiàn)者，在1990—2000年、2000—2010年和2010—2020年間養(yǎng)殖池/鹽田分別侵占灘涂135.42 km2，163.16 km2和95.04 km2，在2010年侵占量達(dá)到峰值之后減少。人工表面直接侵占灘涂的面積較少，30 a間共侵占24.68 km2。萊州灣地區(qū)互花米草等濱海濕地植被擴(kuò)張侵占灘涂的面積較大[19]，近30 a間共有100.30 km2的灘涂轉(zhuǎn)化為濱海植被。灘涂分別向海水和內(nèi)陸水體轉(zhuǎn)化了129.08 km2和46.06 km2，同時(shí)亦有155.58 km2海水和10.17 km2內(nèi)陸水轉(zhuǎn)化為灘涂，整體上對(duì)灘涂面積影響較小。

表2 1990—2020年土地利用轉(zhuǎn)移矩陣Tab.2 Land use transition matrix from 1990 to 2020 (km2)

3.4 驅(qū)動(dòng)力分析

灘涂的時(shí)空變化受到自然因素和人為因素的雙重影響。自然因素包括海平面上升、入海河流泥沙輸入以及潮流變化等； 人為因素則為圍海造田和養(yǎng)殖池開(kāi)發(fā)等一系列開(kāi)發(fā)活動(dòng)。

海平面上升是灘涂面積變化重要的自然原因之一。根據(jù)2020年中國(guó)海平面公報(bào)，山東沿海海平面較常年高65 mm，比2019年高11 mm，預(yù)計(jì)未來(lái)30 a山東沿海海平面還將上升55～165 mm[20]。自然狀態(tài)下，海平面上升勢(shì)必會(huì)直接導(dǎo)致灘涂海側(cè)邊界被海水侵蝕，在無(wú)大量泥沙輸入的情況下，引起灘涂面積減少，這也是萊州灣除黃河口地區(qū)以外灘涂海側(cè)邊界縮退，灘涂面積減少的重要誘因。此外，泥沙輸入也是灘涂生長(zhǎng)、變化的重要原因。萊州灣地區(qū)河流資源豐富，擁有以黃河為主的多條河流，這些河流帶來(lái)了巨量的泥沙輸入，其中黃河的泥沙輸入量占總量的99%以上[21]。黃河口地區(qū)附近落潮歷時(shí)大于漲潮歷時(shí)的潮流特性[22]，也使得入海的泥沙在黃河口附近淤積成陸，灘涂海側(cè)邊界向海延伸。正是在上述因素的共同影響下，黃河口地區(qū)灘涂面積總體呈現(xiàn)出與萊州灣整體上相反的增長(zhǎng)趨勢(shì)，僅在2000—2010年期間有短暫下降，后恢復(fù)增加。2020年相較于1990年灘涂面積凈增長(zhǎng)了8.58 km2，在一定程度上減緩了整個(gè)萊州灣地區(qū)灘涂被侵占的形勢(shì)。而萊州灣其他地區(qū)缺少足夠的泥沙輸入，并且落潮歷時(shí)小于漲潮歷時(shí)的潮流特性不利于泥沙淤積成陸[23]，在海平面上升的驅(qū)動(dòng)下，導(dǎo)致灘涂的海側(cè)邊界不斷被侵蝕，呈現(xiàn)出逐漸向陸移動(dòng)的特征。

人類活動(dòng)也是灘涂面積變化的主要因素，萊州灣灘涂整體上在海一側(cè)邊界向陸推進(jìn)的同時(shí)，人類活動(dòng)不斷使陸一側(cè)邊界向海推進(jìn)。在這種前追后趕的形式逼迫下，灘涂面積必然持續(xù)性的減少。根據(jù)第七次全國(guó)人口普查公報(bào)，到2020年山東省人口為1億[24]，在人口增加的壓力下，向海要地成為了沿海居民的主要途徑[25]。自1991年山東省政府提出建設(shè)海上山東以來(lái)，實(shí)施耕海牧魚(yú)戰(zhàn)略，重視開(kāi)展灘涂養(yǎng)殖，沿海建立大量鹽田，灘涂開(kāi)發(fā)力度逐年增加[26]。受其影響，僅1990—2000年10 a間，萊州灣地區(qū)灘涂面積就減少了150.94 km2。2000—2010年更是達(dá)到了304.78 km2的峰值，到2010年山東省的鹽田面積在沿海各省市中位居榜首，沿海養(yǎng)殖面積也在劇烈增加[27]。在灘涂保護(hù)刻不容緩的形勢(shì)下，2015年山東省海洋局發(fā)布《中共中央國(guó)務(wù)院關(guān)于加快推進(jìn)生態(tài)文明建設(shè)的意見(jiàn)》，明確濕地保護(hù)目標(biāo)，嚴(yán)控圍填?？偭?，控制海洋開(kāi)發(fā)強(qiáng)度； 2018年國(guó)務(wù)院發(fā)布《關(guān)于加強(qiáng)濱海濕地保護(hù)，嚴(yán)格管控圍填海的通知》，原則上不再接受省級(jí)政府申請(qǐng)的萊州灣地區(qū)的圍填海項(xiàng)目； 2019年山東省政府發(fā)布《關(guān)于加強(qiáng)濱海濕地保護(hù)嚴(yán)格管控圍填海的實(shí)施方案》，嚴(yán)格管理圍填海項(xiàng)目，控制灘涂資源開(kāi)發(fā)?；谶@一系列保護(hù)措施，使得2010—2020年的灘涂侵占面積較之上一個(gè)10 a減少了約71%，侵占面積銳減為87.74 km2，侵占量降低至低谷值。

4 結(jié)論

本研究基于GEE平臺(tái)，采用MSIC-OTSU和面向?qū)ο蠓椒ㄟM(jìn)行萊州灣灘涂高精度提取，研究萊州灣地區(qū)灘涂資源數(shù)量和時(shí)空分布特征，構(gòu)建灘涂及其周邊土地覆被數(shù)據(jù)集，分析1990—2020年萊州灣濱海灘涂時(shí)空演變與驅(qū)動(dòng)機(jī)制，得出了如下結(jié)論：

1)萊州灣地區(qū)1990—2020年30 a間灘涂面積整體呈持續(xù)減少的趨勢(shì)，2020年灘涂面積總量為822.38 km2，比1990年減少了543.46 km2，縮減幅度約40%。在2000—2010年間變化最為劇烈，灘涂面積急劇減少了304.78 km2，在2010—2020年減小幅度變緩，僅減少了87.74 km2。濰坊市岸段灘涂面積縮減量最大，30 a間減少了323.56 km2，煙臺(tái)市岸段灘涂面積縮減量最小，僅為82.01 km2。

2)萊州灣地區(qū)灘涂變化最主要的因素是養(yǎng)殖池/鹽田的擴(kuò)張。2020年養(yǎng)殖池/鹽田的總面積達(dá)到了1 493.18 km2，30 a間擴(kuò)展了約45%，占灘涂減少量的62%。其次是由于植被在灘涂中的擴(kuò)張，使灘涂轉(zhuǎn)化為濱海植被，轉(zhuǎn)化量為100.30 km2。灘涂與海水、內(nèi)陸水體的相互轉(zhuǎn)化較少。此外，雖然當(dāng)前人工表面對(duì)于濱海灘涂的侵占面積較小，但未來(lái)有繼續(xù)增加的趨勢(shì)。

3)萊州灣地區(qū)灘涂變化同時(shí)受到自然和人為驅(qū)動(dòng)力影響。萊州灣黃河口以外其他地區(qū)由于不利于泥沙淤積成陸使灘涂面積增加，而黃河口地區(qū)因局部落潮歷時(shí)大于漲潮歷時(shí)的潮流特性使輸入的大量泥沙淤積成陸，灘涂海側(cè)邊界向海推進(jìn)，一定程度上減緩了萊州灣整體上灘涂面積縮減的形勢(shì)。