0 引言
河口是河流流向海洋或者湖泊的终点，是海洋和陆地相互作用、物质交换最频繁的地区，是经济发达、人口集居之地，沿海地区集中了世界67%的大中城市和60%的人口［1］，但日益加剧的人类活动对河口海岸地区的生态环境产生了很大影响.而卫星遥感影像是研究岸线变化最为有效的方法，有专家学者就利用遥感影像，进行海岸线提取，开展了大量的研究［2］.从研究内容上包括：基于遥感和GIS技术进行海岸线信息提取［3］；海岸线的时空变化进行了分析［4-6］；海岸线变化的驱动因素进行了分析［7-9］，涉及我国河口区域有辽河口［7］、黄河口［8］、长江口［9］和珠江口［6］等.现有研究主要从各个地区进行海岸线的提取及其变迁分析，对多个河口的研究尺度进行比较分析的较少.且现有的针对海岸线驱动因素分析，主要是以基于海岸线变迁结果上的定性分析为主，对海岸线变迁驱动因素的定量分析较少.本文选取辽河口、黄河口、长江口和珠江口这4 个重要河口地区作为研究区，对比分析了近30年来这4 个河口岸线的变迁特征，并结合灰色关联度模型对引起海岸变迁的驱动因素进行了分析.为实现海岸线生态环境变化监测、科学的海域管理、海岸带资源与环境保护和可持续利用都具有重要的指示作用.
1 研究区与数据源
1.1 研究区概况
中国海岸线狭长，拥有众多河口海岸带，地域跨度大，气候条件呈现显著的区域差异［10］.本文选取辽河口、黄河口、长江口和珠江口4 个重要河口海岸线为研究对象，为了便于统计社会经济因素对岸线变迁的影响，每个河口海岸线以临近的地级市作为统计边界.辽河口选用盘锦市海岸线，该市坐落于辽东湾顶部，濒临渤海湾，地处辽河三角洲核心区域，是辽西辽南的交通要道.黄河口选用东营市海岸线，该市位于山东省东北部，渤海湾南岸和莱州湾西岸，是黄河入海口的三角洲地带，地势沿黄河走向自西南向东北倾斜，低平缓和.长江口选用邻接的上海市、南通市和苏州市海岸线，3 市属长江入海之前呈喇叭形状的冲积平原，形成三级分汊、四口入海的格局［11］.珠江口选取河口邻接的深圳市、东莞市、广州市、中山市和珠海市的海岸线，属广东省的中部沿海地区，毗邻南海，被称为中国的“南大门”［12］.4 个河口地区区位经济优势明显，属于我国3 大经济圈环渤海、长江角、珠三角的核心经济带，交通运输、劳动力资源等社会经济条件优渥.
1.2 数据源及预处理
研究所用数据为1987、1997、2007 和2017年的Landsat 系列遥感影像，数据来源于美国地质调查局数据中心（http：//glovis.usgs.gov）.以同一季节、自然环境条件相似、天气特征相近，云量尽可能<10%为原则，选择覆盖整个研究区的遥感影像共28 景，包括1987、1997 和2007年Landsat⁃5 TM 影 像21 景，2017年Landsat⁃8 OLI 影像7 景.对遥感影像数据的预处理包括大气校正、影像配准、图像镶嵌和裁剪等.
此外，本文用到的数据还包括各河口地区城市的矢量边界数据，用于研究区遥感影像的裁切；基于Google Earth 遥感影像提取的岸线数据集［13］，用于验证提取的海岸线数据；以及统计年鉴和调查数据，用于构建河口海岸线变迁驱动因素指标体系.
2 研究方法
2.1 基于阈值法的半自动提取海岸线
2.1.1 基于阈值法的水边线提取
水体和陆地具有不同反射特性，阈值法就是基于这种差异性来提取水边线，且采用阈值法进行灰度值分割时，具有简单易操作、处理速度较快等优点，对目标与背景对比强烈的图像普适性较强.
常用的阈值法都是通过不同指数的计算来扩大水陆信息差异.本文选择归一化水体指数（nor⁃malized difference water index，NDWI）［14］，阈值分割法作为提取水边线的方法.对经过预处理后的遥感影像计算NDWI 扩大水陆信息差异，选取合适的阈值进行水边线的提取.NDWI 的计算公式为
式中ρGreen 和ρNIR 分别表示绿光和近红外波段的反射率.
进行水边线提取时，首先采用波段运算法计算各期图像NDWI，通过反复试验研究区NDWI 图像灰度直方图，选定阈值0.18，得到河口地区水陆分离的二值化图像，在ArcGIS 软件中通过栅格转矢量工具、面数据转线数据，经平滑处理消除栅格转矢量的锯齿特征，最终得到水边线.基于阈值法提取水边线的流程，如图1 所示.
图1 基于阈值法提取水边线的流程
2.1.2 结合目视解译的岸线修正及精度评价
为了得到更高精度的提取结果，本文结合目视解译，先对2017年海岸线进行修正，再以2017年海岸线为基准修正其他3 期海岸线，最终获得4 个河口地区4 期的海岸线.在海岸线上随机取点，计算其与2015年全国海陆岸线（参考岸线）最短距离的均方根误差（root mean squared error，RMSE）［15］，评价提取海岸线精度.
2.2 海岸线变迁分析方法
2.2.1 海岸线长度变化
海岸线长度变化是由岸线空间位置及弯曲复杂程度等变化造成.本文选取海岸线实际长度和曲折度2 个指标分析海岸线长度变化.其中海岸线曲折度是指海岸线在一定尺度上空间走向的弯曲程度［16］，用海岸线实际长度与海岸线起点至终点的直线距离的比值进行度量，无量纲.
2.2.2 围填海面积变化
围填海面积变化是指由于岸线空间位置发生变化，从而导致的海岸陆地面积的变化.本文对2 个时期的岸线进行叠加分析，得到岸线之间发生变化的多边形，分为面积增加和减少2 种，分别计算并进行统计分析，得到海岸陆地面积的变化情况.
2.2.3 岸线变迁速率
岸线变迁速率是指对各期岸线间的平均变化速率的度量，本文采用应用广泛的端点变迁速率（end point rate，EPR）工具［17］来度量.在向陆一侧做出1 条和4 期海岸线走势一致的基线，沿着岸线基线以300 m 为间隔，自北向南获得断面.在辽河口、黄河口、长江口和珠江口地区，分别生成了561、1 085、2 303 和2 671 个断面并依次编号.运用数字岸线分析系统（digital shoreline analysis system，DSAS）的EPR 模块计算得到4 个河口地区的岸线变迁速率.其中，EPR 的计算公式为
式中XA、XB和YA、YB分别是剖面与岸线A 和B 的交点投影坐标，TA和TB分别是岸线A 和B 的年份.
2.3 基于灰色关联度模型的变迁驱动因素分析方法
由邓聚龙［18］提出的灰色系统理论演变出的灰色关联分析法［19］，适用于在不确定的信息状态下描述各因素间的关联程度.故本文采用灰色关联度模型，研究河口地区海岸线时空变化的驱动机制，将引起河口地区海岸线变迁的自然、社会和经济系统视为一个灰色系统，构建河口地区灰色关联度模型，计算各影响因素与海岸线变迁的关联度并进行排序，从而定量分析海岸线变迁的驱动因素，为海岸线资源开发利用提供客观评价.
根据研究区自然、社会和经济的发展状况，结合资料收集难易程度，并在参考国内外相关研究成果的基础上［20-22］，将各河口海岸线长度作为因变量，自变量可分为4 种因素共7 个指标：自然因素选取代表水文站的年径流量和年输沙量（六间房、利津、大通、博罗水文站分别作为辽河、黄河、长江、珠江的代表水文站）；经济因素选取GDP、渔业总产值和港口吞吐量；社会因素选取年末总人口数及城镇化水平.采用SPSS 软件，以河口地区为统计单元，其中GDP、渔业总产值、港口吞吐量、年末总人口为各个城市的累加和，城镇化水平取平均值.由于指标单位不同，在做关联度分析时，首先进行无量纲化处理，再计算各驱动因子对各河口地区海岸线长度的关联系数，计算公式［18］为
式中：Ri(t)为关联系数，Di(t)=|y(t)−xi(t)|称为第t 点y 与xi的绝对差，y 为因变量，xi为不同自变量；h为驱动因子的个数；Di(t)min 称为第t 点海岸线长度与各因素的最小差绝对值；Di(t)max 则称为第t 点海岸线长度与各因素的最大差绝对值；b 通常取0.5［23］.
由于Ri(t)是在不同时刻的关联程度值，因此，将各个时刻的关联系数的平均值，作为海岸线长度因子和驱动因子间关联程度的度量.计算各驱动因子与海岸线长度的灰色关联度公式为
式中Pi 越大，驱动因子与海岸线长度关联性越好，将h 个Pi 由大至小排列，得到海岸线长度驱动因子关联程度主次排序.
3 海岸线变迁特征结果分析
3.1 海岸线长度变迁分析
本文提取的岸线与参考岸线吻合度较好，RMSE 为27.56 m，小于本文所采用影像的空间分辨率30.00 m，在误差允许范围内.将各河口4 期海岸线叠加，得到各河口海岸线的空间分布，如图2所示.4 个河口地区海岸线总长度及曲折度变化见图3 和4.
图2 4 期河口地区海岸线空间分布
（a）辽河口；（b）黄河口；（c）长江口；（d）珠江口
图3 4 个河口海岸线总长度变化
图4 4 个河口海岸线曲折度变化
总体来说，30年间4 个河口海岸线长度均呈增加趋势，且岸线越来越曲折，向海域方向移动.其中：辽河口地区、黄河口地区、长江口地区和珠江口地区岸线总长度分别增加了34.87、90.91、159.61和83.77 km，年均增加长度分别为1.35、3.03、5.32 和2.79 km/a.长江口地区海岸线长度较其他河口增长速度快，辽河口地区岸线长度相较之变化平缓.海岸线曲折度总体上呈升高趋势.其中，珠江口地区曲折度最大，2017年达到6.19，辽河口地区、黄河口地区次之，长江口地区曲折度最小.
3.2 围填海面积变迁分析
4个河口地区1987—1997、1997—2007和2007—2017年围填海面积变化分布和统计如图5 和6 所示.总体来看，整个岸段围填海活动剧烈，呈现持续增加的趋势，增加趋势快慢略有不同.其中，辽河口地区和珠江口地区围填海面积在1987—1997年期间变化最大，分别达到179.9 和125.4 km2；黄河口地区和长江口地区围填海面积在1997—2007年间变化最大，分别达到271.3和368.2 km2.结合解译结果：辽河口地区陆域面积增加主要是由于淤泥质岸滩上养殖区、盐田区的增加和港口的建设；黄河口地区陆域面积的增加体现在黄河入海口向海延伸突出及盐田业和养殖业等围填海活动；长江口地区的增长主要来自于上海市的港口建设及养殖区和南通市近入海口处的养殖池塘；珠江口地区珠海市的围填海活动最为剧烈，珠海市南部与岛屿相连变化较大，随着农业填海和养殖填海的活动岸线向海突出，陆域面积增加.
图5 4 个河口地区围填海面积变化
（a）辽河口；（b）黄河口；（c）长江口；（d）珠江口
3.3 海岸线变迁速率分析
运用GIS 和DSAS 的EPR 模块定量分析海岸线的时空变迁程度，得到了4 个河口地区的岸线变迁速率如图7 所示.1987—2017年，辽河口、黄河口、长江口和珠江口地区岸线平均变迁速率分别为55.09、56.95、35.35 和25.93 m/a.其中：黄河口地区变迁速率最快，这是由于研究期间大量的养殖区在淤泥质岸线上开发建设；珠江口地区变迁速率最慢，由于岸线在珠江口地区绵长，岸线变迁的部分在整个范围内所占比例相对其他地区较少，故整体变迁速率低.
图7 4 个河口地区岸线变迁速率
（a）辽河口；（b）黄河口；（c）长江口；（d）珠江口
图6 4 个河口地区围填海面积变化统计
辽河口地区岸线经历了复杂的变化，岸线变迁较快区域分布集中，其中1987—1997年变迁最为剧烈，向海推进速率在小台子到二界沟镇范围变化最为明显，最快达到691.28 m/a.由于盘锦新港的建设，2007—2017年此区域的推进速率也较快，最快达到522.19 m/a.黄河口地区岸线经历了剧烈的变化，岸线变迁明显的区域在整个岸线都有分布，其中1997—2007年变迁最为剧烈；1987—1997年间，岸线变化明显的区域为东港村附近，最快变迁速率达到492.41 m/a；1997—2007年，岸线向海曲折推进速率加快，变化明显的区域为老十五村到红光屋子村之间，以围垦养殖用地为主的岸线增长，向海推进最快速率为451.34 m/a；2007—2017年，变化明显的区域为东营港附近，向海推进最大距离为4.05 km.长江口地区岸线经历了不同程度的变化，岸段变迁速率在整个岸段变化曲折，其中1997—2007年变迁最为剧烈，岸线变迁速率在南汇新城镇及洋口镇附近变化最大，城镇建设及围垦使岸线迅速向海推进，岸段变迁最大速率为515.43 m/a；2007—2017年间，研究区岸线变化以围垦为主，从新东镇到丁店乡附近海岸线每间隔一段距离均有养殖围堤增加，该岸段变迁最大速率为437.17 m/a.珠江口地区岸段变迁速率在整个岸段变化较为集中，珠海市变化情况最为复杂，其中1987—1997年变迁最为剧烈，研究区岸线变化以围垦养殖用地及农田为主，向海推进速率在洪奇门附近最快，为286.55 m/a；2007—2017年，研究区岸线变化以围垦为主，岸线变化明显的区域为珠海市汉基村到塔石村之间，由于建设而与高栏岛相连，该地区向海推进最大距离为4.48 km.
4 海岸线变迁驱动因素结果分析
本文通过灰色关联度模型量化了年径流量、年输沙量、GDP、港口吞吐量、渔业总产值、年末总人口和城镇化水平这7 个驱动因子与海岸线长度的关联度，如表1 所示.综合4 个河口数据，按7 个指标与海岸线长度的关联度平均值由大到小排序为年末总人口数、渔业总产值、港口吞吐量、年径流量、GDP、年输沙量、城镇化水平.其中，年末总人口数对海岸线变化的影响最大（平均灰色关联度达0.891），随着沿海经济政策的深入，导致河口地区人口数量逐年增加，用地的需求和城市建设的需求更加强烈，引起大规模的围填海活动热潮，使得岸线长度增加.其次为渔业总产值（灰色关联度为0.848），由于河口地区有着良好的渔业养殖条件，围海养殖业发展迅猛，使得陆域面积增加，岸线向海推进.
表1 河口海岸线长度与影响因素之间的关联度
5 河口海岸线变迁原因探析
（1）辽河口地区.影响辽河口地区海岸线长度的首要因素是渔业总产值，其次是港口吞吐量，关联度分别为0.949 和0.832.由于1989年开始实施盘锦辽河三角洲农业综合开发工作［24］，围海养殖面积的逐年增加，水产品产业成为支柱性产业，使渔业总产值逐年增加.随着辽东湾新区海域资源的围填开发，在大辽河口外已围填盘锦新港，2010年盘锦新港的通航，港口吞吐量2007—2017年从115.51 万t 增加到3 348.00 万t，使岸线向海域推进，长度增加，并且更加曲折.
（2）黄河口地区.在黄河口地区年径流量对海岸线长度影响最大，其次为城镇化水平、渔业总产值，关联度分别为0.997、0.967、0.951.黄河含沙量世界最大，较高的径流量及输沙水平使海岸线迅速淤积扩张.同时，河道的变迁也使岸线发生了变化，1996年黄河入海口清八汊改道［25］，黄河向东北方向注入渤海，使得北汊向海淤进，改变了入海口附近岸滩格局，岸线也随之变化.近年来因城市建设引起人口向城镇聚集，劳动力的增加引起用地的需要，填海需求增加海岸线不断向海洋扩展.在研究期内，黄河口地区拥有宽阔的淤泥质潮滩，在潮间带上发展养殖场，据统计年鉴记录，黄河口地区水产养殖面积呈持续增加趋势，1987—1997年 从137.0 km2增加到624.3 km2，2007年持续增加到884.1 km2，2017年增加到1 412.1 km2，渔业总产值及GDP 在经济上客观反映了这一变化.
（3）长江口地区.在长江口地区年末总人口数对海岸线长度影响最大，其次为年输沙量，关联度分别为0.985 和0.892.在研究期内，统计年鉴显示，长江口地区3 市的年末总人口总体呈增加趋势，1987—2017年从2 551.40 万人增加到4 247.43 万人.人口的逐年增加，伴随着对用地的需求更加强烈，导致大规模的城镇建设、围填海活动，使得岸线长度增加.长江入海的泥沙量大，也使河口部分由于泥沙的淤积使岸线加长.
（4）珠江口地区.在珠江口地区年末总人口数对海岸线长度影响最大，其次为城镇化水平，关联度分别为0.974 和0.821.在研究期间，珠江口地区在经济特区建设、沿海城市开放政策等对外开放经济政策的推动下，成为中国经济发展的强力引擎和龙头，聚集了大量的劳动力和资本［26］.劳动力的聚集使人口的增长迅速，推动城市建设用地快速扩张增加土地资源，围填海活动急剧增加，海岸线迅速向海扩张.
研究所选取的7 个因子与岸线长度的灰色关联度区间为0.678～0.985，以关联度>0.900 的指标来看，辽河口岸线变迁以渔业总产值这一经济因素为主导，黄河口岸线变迁以自然、经济和社会多种因素共同作用，长江口、珠江口岸线变迁则以社会因素驱动为主.
6 结论
本文对近30年来我国辽河口、黄河口、长江口和珠江口海岸线变迁及其驱动因素进行了深入研究，主要结论如下：
（1）1987—2017年间，4 个河口地区海岸线长度均保持增长趋势，辽河口、黄河口、长江口和珠江口地区岸线总长度分别增加了34.87、90.91、159.61 和83.77 km；海岸线曲折度指数也逐年增加.1987—2017年间，辽河口、黄河口、长江口和珠江口地区陆域面积分别增加了318.7、696.1、738.7 和83.8 km2，岸线变迁速率分别为55.09、56.95、35.35 和25.93 m/a.
（2）构建了海岸线长度变迁的灰色关联度模型，总体而言，海岸线长度与驱动因子的关联度从大到小依次为年末总人口数、渔业总产值、港口吞吐量、年径流量、GDP、年输沙量、城镇化水平.结合定量和定性分析，4 个重要河口地区近30年岸线长度增长、陆域面积增加的主要原因包括：入海泥沙量引起的自然的淤泥质岸线堆积，经济政策推动的劳动力集聚引起城镇建设，以及养殖业、盐田业的迅猛发展和港口的规模化建设等围填海活动.
综上所述，近30年来，4 个河口的海岸线变化显著，驱动因子以人为活动因素为主导.频繁的围填海活动或加剧了海岸带的环境压力，未来河口地区岸线的开发和景观格局的建设应严格遵循国家生态文明建设战略要求，平衡环境资源与开发利用之间的关系，同时保障自然岸线的保有率、海域的后备空间，对海岸带实施阶段性开发限制，以期促进海岸带资源的保护与可持续利用.我国的河流入海口除了以上4 个地区还有海河、钱塘江和闽江等地区，今后可增加这些河口地区海岸线变迁对比分析.
［1］ 陈吉余，陈沈良.中国河口海岸面临的挑战［J］.海洋地质动态，2002（1）：1-5.
［2］ 梁立，刘庆生，刘高焕，等.基于遥感影像的海岸线提取方法综述［J］.地球信息科学学报，2018，20（12）：1745-1755.
［3］ 侯西勇，毋亭，侯婉，等.20 世纪40年代初以来中国大陆海岸线变化特征［J］.中国科学（地球科学），2016，46（8）：1065-1075.
［4］ 李飞，金茹，王在峰.基于多源遥感的中国沿海地区填海造地特征研究［J］.海洋通报，2019，38（4）：429-437.
［5］ 王冰洁，梁璐，惠凤鸣，等.基于Landsat 数据的1975—2015年中国海岸线时空变化分析（英文）［J］.北京师范大学学报（自然科学版），2019，55（1）：83-100.
［6］ 张怡，李晓敏，马毅，等.基于遥感的珠江口海岸线变迁分析［J］.海洋测绘，2014，34（3）：52-55.
［7］ 王建步，张杰，陈景云，等.近30 余年辽河口海岸线遥感变迁分析［J］.海洋环境科学，2015，34（1）：86-92.
［8］ 王集宁，蒙永辉，张丽霞.近42年黄河口海岸线遥感监测与变迁分析［J］.国土资源遥感，2016，28（3）：188-193.
［9］ 李亚方.基于3S 技术的长江口岸线演变规律及成因分析［D］.北京：中国地质大学（北京），2016.
［10］ 全国海岸带办公室.中国海岸带和海涂资源综合调查专业报告集［M］.北京：气象出版社，1991：1-203.
［11］ 杨桂山，翁立达，李立锋.长江保护与发展报告［M］.武汉：长江出版社，2007.
［12］ 中国海湾志编纂委员会.中国海湾志·第十四分册·重要河口［M］.北京：海洋出版社，1998.
［13］ 刘闯，石瑞香，张应华，等.基于Google Earth 遥感影像全球多尺度海陆（岛）岸线数据（2015）［DB/OL］.［2020⁃11⁃01］.http：//www.geodoi.ac.cn/WebCn/doi.aspx？ID=1265&paperID=66c34073⁃0b8f⁃4302⁃ac31⁃e57df4e329b7.
［14］ MCFEETERS S K.The use of the normalized difference water index（NDWI）in the delineation of open water features［J］.International Journal of Remote Sensing，1996，17（7）：1425-1432.
［15］ ZHANG T，YANG X M，HU S S.Extraction of coastline in aquaculture coast from multispectral remote sensing images：object⁃based region growing integrating edge detection［J］.Remote Sensing，2013，5（9）：4470-4487.
［16］ 寇征.海岸开发利用空间格局评价方法研究［D］.大连：大连海事大学，2013.
［17］ MUJABAR P S，CHANDRASEKAR N.Shoreline change analysis along the coast between Kanyakumari and Tuticorin of India using remote sensing and GIS［J］.Arabian Journal of Geosciences，2013，6（3）：647-664.
［18］ 邓聚龙.社会经济灰色系统的理论与方法［J］.中国社会科学，1984（6）：47-60.
［19］ 刘丽丽.农村产业结构灰色系统分析与预测：以北京市怀柔县汤河口镇为例［J］.首都师范大学学报（自然科学版），1994，25（2）：89-96.
［20］ 李明昱.基于RS 和GIS 的辽宁省海岸线时空变化及驱动因素分析［D］.大连：辽宁师范大学，2011.
［21］ 陈路遥，蒋卫国，陈云浩，等.基于灰色模型GM（1，1）的海岸线变化趋势预测：以辽河口三角洲地区为例［J］.应用基础与工程科学学报，2008（4）：511-517.
［22］ TIAN B，WU W T，YANG ZH Q，et al.Drivers，trends，and potential impacts of long⁃term coastal reclamation in China from 1985 to 2010［J］.Estuarine，Coastal and Shelf Science，2016，170：83-90.
［23］ 隋燕.基于遥感的海南岛近30年海岸线时空变迁监测与分析［D］.青岛：山东科技大学，2018.
［24］ 方汉隆，孟祥哲.辽河三角洲农业综合开发成效显著［J］.垦殖与稻作，1998（1）：44-45.
［25］ 王安东.黄河入海口1996年清八汊改道以来河口段岸线变迁遥感监测［J］.山东林业科技，2018，48（1）：37-39.
［26］ 刘旭拢，邓孺孺，许剑辉，等.近40年来珠江河口区海岸线时空变化特征及驱动力分析［J］.地球信息科学学报，2017，19（10）：1336-1345.