黄河流域作为中国重要的生态屏障与资源集聚带, 其生态环境问题备受关注.长期以来, 由于黄河流域生态系统较为脆弱, 生态环境承受了巨大的压力, 植被破坏、水土流失、土地沙化加剧和水源涵养功能减退等生态问题严重^{[1, 2]}.与此同时, 黄河流域以高投入、高消耗为特征的粗放式经济增长模式谋求发展^[3], 致使生态环境与生态系统多样性面临严重威胁, 发展与保护的矛盾日益突出.为缓解区域间生态资源供需矛盾和防范化解黄河流域生态安全风险, 我国政府高度重视黄河流域生态补偿机制建设工作, 相继出台《支持引导黄河全流域建立横向生态补偿机制试点实施方案》《黄河流域生态保护和高质量发展规划纲要》, 并指出要构建全方位推进生态补偿制度及相关领域改革的新发展格局.由此可见, 在黄河流域“共抓大保护, 协同大治理”理念的引领下, 生态补偿机制已成为解决区域间生态资源供需矛盾与推动全流域绿色可持续发展的战略依托^{[4, 5]}.

生态补偿作为一项保护资源环境的市场化政策工具, 能够将损害或保护资源环境行为的外部效应内部化, 从而实现保护生态环境的目的^{[6 ~ 10]}.生态补偿额度的确定是生态补偿研究从理论走向实践的关键环节, 从生态补偿测度方法来看, 学者们主要采用排污量分配法^{[11 ~ 15]}、水足迹法^{[2, 16]}、生态系统服务价值法^{[3, 4, 17, 18]}、机会成本法^[19]和条件价值法^[20]等方法.上述方法为黄河流域生态补偿额度的定量研究提供了借鉴, 但这类方法存在时间与空间范围的局限性, 且未能综合考虑区域间社会经济发展与生态系统资源供需匹配的协同关系.为此, 部分学者采用能够反映生态环境基质变化情况的生态足迹模型以测算区域生态补偿额度^{[21 ~ 24]}.然而, 传统生态足迹模型中的均衡因子、产量因子与全球平均生产力等指标会随时间与地区变化而发生较大波动, 致使测算结果的稳定性较差.随着研究的深入, Odum^[25]将能值分析理论纳入传统生态足迹模型分析中, 提出了能值生态足迹模型, 该模型能够将生态环境系统中不同类型、不同质量、不同量纲与不可比较的能量转换为具有可比性的太阳能值形式.基于能值生态足迹模型, 学者们对生态补偿研究的区域主要在省域^[26]、市域^[27]、草原^[28]和流域^[29]等.刘璐瑶等^[26]选取33项能值生态足迹指标, 采用能值三维生态足迹模型对珠江干流4省区生态补偿进行测度, 发现能值生态足迹在空间分布上表现为：广东（粤） > 广西（桂） > 贵州（黔） > 云南（滇）, 而能值生态承载力反之, 因此上游地区（滇）为生态受偿区, 而下游地区（粤、桂和黔）为生态补偿支付区.郭婷等^[27]以陕南国家重点生态功能区为例, 采用14项生物资源产品与7项能源产品构建能值生态足迹模型, 发现陕南三市均处于生态盈余状态, 且应获得的生态补偿标准均呈增加态势.邓汉琨等^[28]以畜产品作为主要测算对象, 并根据能值生态足迹模型计算的生态外溢价值结果, 对内蒙古自治区牧业盟市的草原生态补奖资金进行测算, 发现鄂尔多斯、呼伦贝尔和乌兰察布等城市长期处于生态补奖资金赤字区, 而巴彦淖尔和阿拉善等城市为生态补奖资金盈余区.王奕淇等^[29]将能值分析思想与水足迹法结合, 发现2013年渭河流域上游地区处于生态盈余状态, 其生态外溢能值（solar emjoules, se, 太阳能焦耳）为1.16×10²² J, 因此上游各市能够获得16.31亿元的生态补偿资金.

综上可知, 既有研究仍然存在改进之处：一是生态补偿测度方法有待完善, 部分研究虽采用能值生态足迹模型测算生态补偿额度, 但针对黄河流域的研究成果较少, 且能值生态足迹指标有待扩充；二是流域生态补偿具有典型的外部性特征, 不同地区之间的生态环境与资源利用存在着明显的空间关联性, 现有研究依旧停留在对生态补偿额度测算结果的简单概括层面, 对于生态补偿额度的时空格局、区域差异与分布动态演进特征分析的研究较少.为此, 本文基于能值生态足迹模型, 从上中下游、省域和市域这3个维度对2007~2021年黄河流域横向生态补偿额度进行差异化分析, 并结合探索性空间数据分析、核密度估计、Dagum基尼系数与空间Markov链分析方法全面考察黄河流域横向生态补偿额度的时空格局、区域差异及分布动态演进特征, 以期为建立健全差异化、市场化与多元化的生态补偿机制以及推动黄河流域生态环境保护和区域协调发展提供科学依据与理论参考.

1 材料与方法 1.1 研究区概况

考虑到地理单元的连续性与地区经济发展间的关联性, 本文研究的黄河流域包括自然流域所涉及的区域与黄河流经地区以外社会经济发展深受黄河流域影响的区域^[30].鉴于内蒙古东部地区已被纳入东北振兴范围, 本文未将呼伦贝尔、兴安、通辽、赤峰和锡林郭勒纳入研究区域, 同时由于莱芜于2019年划规济南管辖, 故将其数据合并计算至济南^[31].本文最终选取2007~2021年黄河流域90个城市作为数据样本, 并将其划分为上中下游三大区段（图 1）.

1.2 数据来源

原始数据均来自2007~2022年《中国城市统计年鉴》《中国农村统计年鉴》《中国林业和草原年鉴》《中国能源统计年鉴》和各省市统计年鉴及《水资源公报》等, 对于缺失数据采用均值插补法和回归插补法进行填补.

1.3 研究方法 1.3.1 跨区域流域横向生态补偿测度模型

（1）能值生态足迹  杨灿等^[32]认为无论是出口还是留存的生物质资源均占用了本区域样本期内的生态环境资源, 在以生态足迹为基础计算生态补偿额度时, 应以生产性而非消费性生态足迹对生物质产品进行核算.因此, 本文将生物资源与能源资源分别按照生产性生态足迹与消费性生态足迹进行核算.根据黄河流域各类资源的实际生产情况, 选取6类42项能值生态足迹指标, 其中生物资源生产性项目包括农产品、林产品、畜产品与水产品4类33项指标, 能源资源消费性项目包括化石能源与电力2类9项指标, 参考既有研究^{[25, 32 ~ 35]}, 各指标的能值折算系数与能值转换率如表 1所示.

能值生态足迹计算公式如下：

(1)

(2)

式中, EEF_i为i地区能值生态足迹；C_ik为i地区第k种生物资源的能值；ER_i为i地区可更新资源能值；A_i为i地区土地面积.

（2）能值生态承载力  为客观合理地反映能值生态承载力状况, 通常将自然资源划分为可更新资源与不可更新资源两类^[26].由于不可更新资源的再生速度远低于其消耗速度, 且在经济社会发展过程中终将走向枯竭, 因此, 能值生态承载力研究的主要对象是可更新资源.本文在考虑了太阳能、雨水势能、雨水化学能、风能、径流势能、径流化学能与地球旋转能这7种可更新资源的基础上, 根据Odum^[25]的能值分析理论, 将不同类型的能值换算为统一的太阳能值标准进行比较.根据能值分析理论, 相同性质的能量一般只选择其中的最大值.因此, 为避免重复计算, 在太阳能、雨水势能、雨水化学能、风能、径流势能与径流化学能这6种具有相同性质的能量中, 选取最大的某一类能值与地球旋转能值之和作为可更新资源的太阳能值^[32].根据世界环境与发展委员会的报告《我们的未来》建议, 在测度生态承载力时需扣除12%的生物多样性修正系数^[35].能值生态承载力计算公式如下：

(3)

式中, EEC_i为i地区能值生态承载力；P_i为i地区的区域能值密度；P为全球平均能值密度, 根据Odum^[25]的计算结果, 取值（se）为3.104×10¹⁴ J·hm^-2；ER_i为i地区可更新资源能值；A_i为i地区土地面积.

(4)

式中, EB_ij为i地区生态系统供给的7种能量（由表 2对应参数乘积求得）, ET_j为不同能量对应的能值转换率, j = 1，2，⋯，7, 分别表示太阳能、雨水势能、雨水化学能、风能、径流势能、径流化学能和地球旋转能, 参考现有研究^{[25, 32 ~ 35]}, 各项指标参数见表 2.

量化横向生态补偿额度需要考虑各地区在自然资源禀赋、经济发展水平与生态环境治理投入等方面的差异.为此, 本文在构建横向生态补偿测度模型时, 综合考虑了区域间生态系统服务的外溢性价值、区域间经济发展差异及其资源配置效率造成的机会成本差异等因素.跨区域流域横向生态补偿测度模型的分析框架如图 2所示.可以发现, 横向生态补偿额度是根据生态服务供给与生态服务消费的差额计算得到的, 其中, 若差额为正（即EES > EED）, 则说明该地区生态系统提供的自然资源产品与服务数量远大于其自身对生态环境资源消耗的数量, 该地区能够向其他地区提供剩余的生态系统产品与服务, 属于生态补偿受偿区；若差额为负（即EES < EED）, 则表明该地区生态系统提供的自然资源产品与服务数量无法满足地区社会经济发展所依赖的生态环境资源需求量, 该地区需要消耗其他地区提供的生态系统产品与服务, 属于生态补偿支付区.

生态服务供给金额是指某一地区生态系统因提供特定自然资源产品与服务而应获得的生态补偿资金额度^[36], 公式如下：

(5)

式中, EES_i为i地区生态服务供给金额；M_i为i地区环境污染治理投资额, 囿于数据获取限制, 参考现有研究^[37], 本文以地级市国内生产总值占所在省份国内生产总值的比重为权重, 乘以该省份环境污染治理投资总额, 得到地级市层面的环境污染治理投资总额；α_i为i地区生态服务供给系数；EEC_i为i地区能值生态承载力总量.

生态服务消费金额是指某一地区因消费或使用特定自然资源产品与服务而应支付的生态补偿资金额度^[36], 公式如下：

(6)

式中, EED_i为i地区生态服务消费金额；β_i为i地区生态服务消费系数；rec_i为i地区综合修正系数, 公式如下：

(7)

(8)

(9)

(10)

(11)

式中, U_i为i地区万元GDP能值生态足迹；W_i为i地区支付意愿；U和W分别为i地区所在省份万元GDP能值生态足迹与支付意愿平均水平；y_i为i地区人均GDP；y为i地区所属省份人均GDP；In_i为i地区人均收入；In为i地区所属省份人均收入；I_i¹与I_i²分别为i地区城镇人均可支配收入与农村人口纯收入；m_i与n_i分别为i地区城镇与农村人口数量.

生态补偿额度是指某一地区因提供生态服务价值剩余（或消费生态服务价值剩余）而应受偿（或支付）的生态补偿资金额度, 公式如下：

(12)

式中, EP_i为i地区应受偿（或支付）的生态补偿额度；当EP_i > 0时, 表明i地区生态补偿额度为净流入, 即为生态补偿受偿区；反之, 当EP_i < 0时, 表明i地区生态补偿额度为净流出, 即为生态补偿支付区.

1.3.2 探索性空间数据分析

探索性空间数据分析方法通过研究不同空间位置的内在关联程度, 进而揭示空间事物的相似集聚和分异特征, 包括全局空间自相关与局部空间自相关^{[1, 31]}.其中, 全局空间自相关用于描述事物在整个区域的空间关联及分异程度, 常用全局Moran's I指数进行衡量^{[38, 39]}.局部空间自相关旨在描述某一区域与邻近单元之间的高-高集聚与低-低集聚的空间关联特征, 即形成的热点区与冷点区, 常用Getis-Ord G_i^*指数^[40].

1.3.3 Dagum基尼系数及其分解方法

Dagum基尼系数能够测度黄河流域生态补偿的相对差异与差异来源^[41].根据Dagum基尼系数及其分解方法, 总体差异来源G可以分解为区域内差异贡献G₁、区域间差异贡献G₂和超变密度贡献G₃, 且三者关系满足G = G₁ + G₂ + G₃^[42].通常, 区域内差异与区域间差异越小, 表明区域内或区域间的生态补偿额度越接近, 且拥有较高的协同性, 反之协同性则较低；超变密度则衡量不同区域交叉重叠部分所具有的差异性.

1.3.4 核密度估计

核密度估计方法能够从密度曲线的分布位置与分布形态深入剖析黄河流域生态补偿额度绝对差异的分布特征^{[30, 42]}.其中, 分布位置反映生态补偿额度的高低；分布形态反映生态补偿额度的区域差异大小与空间极化程度, 区域差异大小由曲线的波峰高度与宽度表示, 空间极化程度由曲线的波峰数量表示.

1.3.5 空间Markov链分析

本文采用空间Markov链分析考察空间效应在黄河流域生态补偿动态演进中的影响^{[43 ~ 45]}, 并将生态补偿额度在相邻类型中发生变动的现象定义为正向（或负向）转移, 如由高水平转移至中高水平为负向转移.同时, 将生态补偿额度跨类型的变动现象定义为正向（或负向）跨越式转移, 如由高水平转移至中水平为负向跨越式转移.

2 结果与讨论 2.1 黄河流域生态补偿的时空演变特征分析 2.1.1 时间演变分析

图 3（a）报告了黄河流域横向生态补偿额度的时间演变趋势.可以发现, 黄河流域生态补偿额度在波动中略有下降, 年均值为0.20亿元, 除2016年外, 其余年份均为正值, 说明黄河流域在经济发展过程中因消耗过多的生态资源产品与服务, 致使生态系统承载压力增加, 因此需要支付的生态补偿额度有所上涨.需要说明的是, 由于2016年太原、大同与朔州等城市的生态足迹大幅增加, 且生态系统资源供给能力与生态系统承载力均较为疲软, 导致中游地区需要支付较多的生态补偿额度, 整体拉低了当年生态补偿均值.为全面梳理黄河流域生态补偿的差异化情况, 本文进一步从上中下游地区、省域和城市3个维度展开分析.

（1）基于不同区段差异维度的分析  如图 3（a）所示, 黄河流域各区段生态补偿额度存在不同的时间演变趋势, 按年均值由大到小依次为：下游（3.164亿元） > 中游（-0.148亿元） > 上游（-2.435亿元）.其中, 黄河上游生态补偿额度整体为负, 自2007年开始波动下降, 但在2014起开始有所上升；中游生态补偿额度自2007年开始波动上升, 但自2016年大幅降低, 且由生态补偿受偿区转变为支付区；下游生态补偿额度始终为正, 区域内生态补偿受偿区较多.各区段生态补偿额度不同是由于经济发展水平差异导致的, 上游地区受自然区位因素的影响, 产业发展相对薄弱, 综合经济发展水平较其他地区低, 生态系统产品与服务的消耗数量较低, 用于生态补偿的资金较少, 而下游地区经济发展水平较高, 在实施生态环境保护与治理过程中, 对环境投入成本较高, 因此下游地区获得的生态补偿资金远高于其他地区.

（2）基于省域差异维度的分析  由图 3（b）中可知, 青海、甘肃、宁夏、内蒙古、陕西与河南的生态补偿额度在0值附近波动, 说明这6个省份生态补偿额度变化幅度相对较小；山西除2016年大幅降低外, 其余年份均在0值附近波动；四川整体呈波动上升趋势, 且生态补偿额度始终为正；山东自2019年起生态补偿额度由正转负, 从生态补偿资金受偿区转变为支付区.总的来说, 省域间生态补偿额度存在一定差距, 且黄河流域内生态补偿受偿区主要集中在生态本底优越或经济发达的地区.原因在于, 对于生态本底优越的地区而言, 其自然资源禀赋较好, 生态系统承载压力较低, 在剔除自身生态产品服务消费后, 能够向其他地区提供其剩余的生态产品服务价值, 故应获得与生态产品服务外溢性价值相对应的生态补偿资金；对于经济发达地区而言, 由于其在生态保护与环境治理过程中, 需要综合考虑GDP、人口密度和资源耗减等问题以及发展机会成本损失, 对生态环境的资金投入需求更大, 因此获得的生态补偿资金往往也更高.

（3）基于城市差异维度的分析  黄河流域各城市生态补偿额度的测度结果见表 3.可以发现, 上游地区支付与获得生态补偿额度最高的地区分别为鄂尔多斯和阿拉善；中游地区支付与获得生态补偿额度最高的地区分别为太原和商洛；下游地区支付与获得生态补偿额度最高的地区分别为东营和临沂.值得注意的是, 玉树、果洛、海南、黄南、临夏、陇南等城市的能值生态承载力远高于能值生态足迹, 均处于生态盈余状态, 这些城市为维系黄河流域生态安全而放弃了诸多经济发展机会, 但未能获得合理的生态补偿, 反而因生态环境保护陷入“政策性返贫”的尴尬境地, 这不符合生态补偿“谁受益, 谁补偿”的基本原则.因此, 在制定与实施黄河流域生态补偿标准时, 应进一步改进生态补偿资金分配方式, 实施差异化生态补偿标准, 并重点向禁止开发区域、限制开发区域以及黄河上游欠发达地区倾斜.

2.1.2 空间演变分析

根据黄河流域90个城市生态补偿额度的测度结果, 采用自然断点法将其划分为高支付区、低支付区、低受偿区与高受偿区4类, 分类结果如图 4所示.可以发现：①结合图 4（a）~4（d）, 从空间分布特征来看, 黄河流域生态补偿额度呈“自西北向东南”梯度增加的格局分布, 生态补偿支付区集中于黄河北部, 受偿区在黄河南部形成集聚；②结合图 4（e）, 从结构层次分布来看, 横向生态补偿额度类型整体呈中间多、两头少的橄榄型结构, 低支付区与低受偿区数量较多, 这两类区域占比分别为28.33%和34.44%, 其次为高受偿区, 占比为25.83%, 高支付区占比为11.39%；③从发展趋势来看, 黄河流域生态补偿空间格局演进明显, 低支付区与低受偿区数量逐渐减少, 呈现出由外围向中心逐步缩小的集聚布局；高支付区与高受偿区数量不断增加, 并形成了分散分布与组团分布并存的空间演进特征, 其中, 高支付区表现为以西安、兰州、郑州和太原等省会城市以及鄂尔多斯、包头、榆林、大同和淄博等主要工业城市为主的点状分散分布格局, 高受偿区呈现出连片式的集聚布局, 强者愈强的“马太效应”特征越发明显.

2.1.3 空间关联分析

（1）总体空间关联分析  由表 4可知, 样本期内全局Moran's I指数值均大于0, 且均在1%的水平上显著, 说明黄河流域生态补偿额度具有较为显著的空间正相关性, 总体趋于空间集聚特征.纵观整体趋势, 全局Moran's I指数呈波动上升趋势, 由2007年的0.326上升至2021年的0.377, 说明黄河流域各城市生态补偿额度的空间集聚度有所增强, 即生态补偿高值区域与低值区域在空间分布上均呈现明显集聚特征.

（2）冷热点分析  由于全局空间自相关分析仅能说明区域是否具有空间集聚性, 但无法确定这种区域是高值集聚还是低值集聚分布^[40].为此, 本文通过计算2007、2011、2016和2021年的Getis-Ord G_i^*指数, 并采用自然断点法将其划分为冷点区、次冷点区、次热点区和热点区这4类.由图 5可知, 黄河流域生态补偿整体呈“北冷南热”的空间分布态势, 冷热点空间极化现象较为突出, 且冷热点的组团区域逐步转移. ①从空间分布来看, 生态补偿冷点区与热点区均表现为“多核”分布, 形成了以鄂尔多斯、包头、东营和淄博为中心的冷点区与以中卫、陇南、驻马店和济南为中心的热点区, 且次冷点区与次热点区围绕冷热点的中心区域在空间上展开, 形成了次冷点区与次热点区组团集聚的格局；②从发展趋势来看, 冷点区与热点区城市数量均有所增加, 但在空间演进趋势上各有差异, 冷点区呈现由中心向外围局部扩散的趋势, 而热点区呈现由零散分布向带状组团式分布转变的态势, 说明黄河流域生态补偿存在显著的空间外溢效应.

2.2 黄河流域生态补偿的区域差异及其分解

前文结果表明黄河流域生态补偿具有空间非均衡性特征, 为深入分析黄河整体及上中下游地区的差异来源, 本文采用Dagum基尼系数及其分解方法研究黄河流域生态补偿的空间差异大小及其具体来源, 结果见图 6.可以发现, 整体基尼系数均值为0.101, 介于[0.069, 0.123], 呈“下降-上升”频繁交替的周期性波动态势, 且略微上升, 说明黄河流域生态补偿具有一定的空间非均衡特征, 且总体差距有所增加.

2.2.1 区域内差异

结合图 6（a）可知：①从水平值来看, 上中下游地区的差异处于不同水平, 按均值由高到低依次为：下游（0.128） > 上游（0.098） > 中游（0.067）, 说明黄河流域生态补偿存在较为显著的区域内差异, 且下游地区的非均衡发展问题尤为突出；②从演变趋势来看, 样本期内黄河上游和下游地区的区域内差异呈波动下降趋势, 分别从2007年的0.082和0.144降低到2021年的0.079和0.140, 而中游地区呈波动上升趋势, 说明上游和下游地区生态补偿额度差异逐渐缩小, 而中游地区生态补偿额度差异有所增加.

2.2.2 区域间差异

如图 6（b）所示：①从水平值来看, 区域间基尼系数均值由高到低依次为：上游-下游（0.117） > 中游-下游（0.104） > 上游-中游（0.084）, 说明下游地区整体生态补偿额度较高且生态补偿受偿区城市数量较多, 远高于上、中游地区.结合生态补偿额度测算结果可知, 下游地区生态补偿资金充裕, 样本期内河南与山东两省年均环境污染治理投资额达1 018.70亿元, 占黄河9省（区）总额的44.00%, 生态保护投入远高于黄河流域其他省份, 这与前文时空演变特征基本吻合；②从演变趋势来看, 区域间Dagum基尼系数均呈上下交替波动式变化态势.其中, 上游-下游和上游-中游的基尼系数均有所上升且变动态势较为趋同, 而中游-下游的基尼系数总体有所下降, 说明上游-下游和上游-中游的差距均在逐步扩大, 而中游与下游的差距逐渐缩小.由此可见, 黄河流域区域间生态补偿的非均衡现象普遍存在, 尤其应促进下游与其他地区生态补偿的整体协同性.

2.2.3 区域差异来源分解

由图 6（c）可知, 区域内差异、区域间差异和超变密度的均值分别为32.74%、6.25%和61.01%, 且3条线的相对位置始终未发生变化.由此可见, 超变密度是导致黄河流域生态补偿额度总体差异的主要来源, 说明区域差异产生的核心根源在于不同区域间交叉重叠问题, 即尽管下游地区受偿区数量较多且获得生态补偿额度较高, 但并非所有城市均属于受偿区并能获得较高的生态补偿资金, 且这类跨群城市偏离地区共同属性的趋势基本保持不变, 致使超变密度对总体差异的贡献值居高不下.因此, 为解决黄河流域横向生态补偿的区域差异问题, 不仅需要重视区域间的整体协同性, 而且还应重视区域内部发展的协调性.

2.3 黄河流域生态补偿的分布动态演进 2.3.1 黄河流域生态补偿分布的时间动态演进

由于Dagum基尼系数描述的是黄河流域横向生态补偿的相对差异规律, 未能体现绝对差异的分布特征及其演进规律^[42], 因此本文借助核密度估计方法进行分析.黄河流域整体生态补偿额度的三维核密度和密度等高线见图 7（a）.可以发现, 黄河流域整体生态补偿水平不断上升, 存在轻微的多级分化现象, 且地区间差异有所扩大.其中, ①从分布位置来看, 核密度曲线主峰整体呈先左移后右移的态势, 说明黄河流域生态补偿额度整体有所增加, 但部分城市在快速的城市化与工业化进程中形成了对生态产品与能源产品的过量需求, 导致生态补偿额度有所下降, 使得曲线在生态补偿低值区域形成聚集, 进而导致主峰在部分年份出现左移趋势；②从分布形态来看, 波峰高度不断下降, 波峰宽度逐渐变宽, 说明横向生态补偿的绝对差异渐趋扩大, 存在“优者更优”效应；③从分布延展性来看, 核密度曲线存在明显的双拖尾现象且呈收敛态势, 说明生态补偿额度较低的区域有所提高, 趋向于接近区域内部生态补偿额度的平均水平；④从极化趋势看, 曲线由一个主峰和多个侧峰构成, 侧峰峰值相对较低, 表明生态补偿具有多级分化趋势.

黄河不同区段生态补偿额度的三维核密度图如图 7（b）~7（d）所示, 可以发现, 三大区段核密度曲线的分布位置和波峰形态演化过程基本相似, 但分布延展性与极化趋势存在一定差异.具体来看：①上游地区对应的核密度曲线位置右移, 波峰高度降低, 宽度略微扩大, 左右拖尾均长且薄, 整体呈多峰形态且侧峰峰值较低, 说明上游地区生态补偿额度有所提高, 绝对差异逐步增加, 且存在一定的多级分化趋势；②中游地区对应的核密度曲线主峰位置不断右移, 波峰形态表现为高度下降且宽度不断拓宽, 左拖尾特征随时间推移逐步扩大, 整体呈一个主峰与多个侧峰构成的多级分化趋势, 说明中游地区生态补偿额度呈上升趋势, 多级分化特征较为明显, 且绝对差异有所增加；③下游地区对应的核密度曲线主峰位置整体向右移动, 波峰高度有所降低, 波峰宽度逐步拓宽, 左拖尾现象明显且随时间推移有所扩大, 具有明显的多级分化趋势, 说明下游地区生态补偿额度逐步提升, 绝对差异逐渐增加, 极化程度有所增加.

2.3.2 黄河流域生态补偿分布的空间动态演进

空间Markov链分析方法能够探究黄河流域各地区在生态补偿类型分布中的相对位置变动与动态转移特征, 是对核密度估计的较好补充^[43].基于邻接空间权重矩阵, 本文进一步采用空间Markov链分析方法研究空间因素在黄河流域横向生态补偿长期发展趋势中的影响.Q统计量能够判定空间Markov转移概率的显著性水平.不同时长下, Q值分别为62.37（t₁）、62.82（t₂）和74.76（t₃）, 且均通过1%显著性水平检验, 说明黄河流域生态补偿额度的分布演进具有空间效应, 即某一城市生态补偿水平的动态转移会受到周边城市的影响.本文遵循每种类型所涵盖的城市数量相近原则^[45], 按照四分位点将黄河流域生态补偿额度划分为4种类型：低值域、中低值域、中高值域和高值域, 并用N = 1，2，3，4分别表示这4种类型, 进而得到黄河流域生态补偿额度在时间跨度t₁、t₂和t₃下的空间Markov转移概率矩阵（表 5）.

由表 5可知：①不同时长下, 对角线的概率值始终大于非对角线的概率值, 说明在空间因素影响下, 不同值域的研究单元维持原有状态的概率较大；②不同时长下, 低值域与高值域两种类型在对角线上的概率值仍然明显高于中低值域与中高值域, 说明在空间因素影响下, 低值域与高值域存在较高的稳定性；③不同时长下, 4种类型在对角线上的概率值均呈“t₁ > t₂ > t₃”的特征, 说明随时间的推移, 黄河流域各城市生态补偿维持原有水平的概率有所下降, 流动性渐趋增强；④不同邻域对生态补偿额度转移的影响存在较大区别, 这种空间溢出效应容易在一定的地理空间范围内形成“空间俱乐部收敛”现象, 即若某城市与生态补偿高值域空间滞后类型的城市相邻, 则会增加该城市向上转移的概率, 说明邻域城市生态补偿水平较高, 有利于提高周边城市的生态补偿水平.

3 建议与展望

根据研究结果, 建议如下：加大黄河流域生态补偿投入力度, 改进生态补偿资金管理方式, 重点向禁止开发区域、限制开发区域以及黄河上游欠发达地区倾斜；因地制宜制定并实施差异化的生态补偿机制, 推动生态补偿模式由政府主导型向排污权交易、水权交易、绿色债券和产业引导基金等市场化模式转变, 建立健全黄河全流域横向生态补偿机制与生态效益受益者直接补偿机制；充分发挥黄河流域内中心城市辐射带动作用与空间溢出效应, 以黄河流域整体性治理为基础, 建立上下游、左右岸和地区间生态保护与环境治理协同机制, 推动跨区域生态补偿高质量协同发展.

本文基于能值生态足迹模型测度黄河流域90个城市的横向生态补偿额度, 并对其时空格局、区域差异与分布动态演进特征进行分析, 但并未考虑生态系统服务价值、污染物排放和工业产品生产情况等因素, 且本文未对其他流域的适用性进行研究, 因此未来需增加并优化能值生态足迹指标, 完善本文构建的跨区域流域横向生态补偿模型, 并对更多流域进行研究.此外, 本文未开展生态补偿影响机制方面的研究, 这也将是下一步的研究重点.

4 结论

（1）黄河流域生态补偿额度从时间维度上整体略有下降, 在空间分布上整体自西北向东南梯度增加, 生态补偿高值区域集中在生态本底优越或经济发达的地区；

（2）黄河流域生态补偿存在显著的正向空间集聚现象, 整体呈“北冷南热”的空间分布态势；

（3）黄河流域生态补偿总体差异有所扩大, 超变密度是造成总体差异的主要来源, 缓解区域内与区域间差异的着力点位于下游地区；

（4）黄河流域生态补偿存在轻微的多级分化现象, 且地区间差异有所扩大；生态补偿类型具有较强的“空间俱乐部趋同”特征, 但随时长扩大, 各城市维持原有水平的概率有所降低, 流动性渐趋增加.