DMSP/OLS夜间灯光影像可在较大尺度上分析城市发展的空间变化，是从事大尺度城镇化研究的有效数据手段^[1]。该影像为灰度影像，表达城镇的灯光强度信息和空间分布关系，像元灰度(DN)值代表灯光强度，灰度值范围为0~63，灰度值越大表示城镇建成区的可能性越大^[2]。

目前，利用DMSP/OLS夜间灯光影像进行城镇建成区提取的方法主要包括两大类：基于像元的图像分类法^{[3, 4]}和基于区域阈值的图像分割法。第一类方法，以像元为单位，将大于某一阈值的像元识别为建成区像元，小于该阈值的像元识别为非建成区像元，大时空尺度应用的普适性没有得到有效证实。第二类方法，主要指统计数据比较法^{[5, 6]}，该方法以行政边界为划分依据，借助统计数据逐一确定各行政单元的阈值，当行政单元个数较多时，重复性操作多，效率低下。

由于我国的城市化水平和经济发展程度在单个行政单元内部存在明显的不均衡性^[7]，必须打破行政边界的划分局限，对同一行政单元的不同区域确定不同的最佳阈值^{[8, 9]}。因此，本文提出一套打破单个像元和行政边界的划分局限、兼顾城市化发展水平差异和城镇发展规模的聚类阈值法。该方法采用空间聚类思想，对识别的潜在建成区对象进行空间聚类，划分不同的集合，然后确定各集合的最佳阈值提取城镇建成区。

聚类阈值法最基本的空间单元为所识别的潜在建成区对象而不是单个像元或行政单元。

参考文献[2]的研究，DN值大小与城市化强度成正比，且建成区规模越大，区域的城镇化发展等级越高。即研究对象灯光影像DN值的大小和建成区规模的大小可以作为不同城镇化发展等级的衡量指标。据此，聚类阈值法以各潜在建成区对象的地理重心处DN值大小和各对象所包含像元的数量多少为聚类依据，进行空间聚类，划分不同的城镇化发展等级。

本文对NGDC发布的1992~2011年的F10(1992~1994)、F12(1994~1999)、F14(1997~2003)、F15(2000~2007)、F16(2004~2009)和F18(2010~2011)共32期DMSP/OLS(Version 4)非辐射标定夜间灯光强度影像进行处理，主要技术思路如下：① 数据预处理。利用一元二次回归模型对DMSP/OLS夜间灯光数据进行校正，获取1992~2011年的具有连续性和可比性的灯光影像；② 潜在建成区对象的识别。利用连通区域循环标识算法完成2 728个潜在建成区对象的识别；③ 空间聚类。根据地理重心处DN值大小和对象规模进行空间聚类，对潜在建成区对象划分不同的城镇化发展等级；④ 确定阈值。借助统计数据确定各发展等级的最佳阈值，提取城镇建成区；⑤ 建成区对象识别。再次利用连通区域标识算法对提取结果进行建成区对象识别。⑥ 几何形态优化。对所识别建成区对象进行几何形态优化，消除细小的城镇碎片、内部空洞和边缘锯齿；⑦ 精度评定。对得到的历年中国行政区划(本文中不含南海诸岛)内的城镇建成区信息，从数量尺度和空间格局两方面进行精度评定。

受传感器探测性能、卫星过境时间及成像环境不同等的影响，未经校正的灯光影像不具有连续性和可比性^[10]。本文应用相交分析和叠加分析，去除影像中的非稳定灯光像元，并应用一元二次回归模型对影像DN值进行校正：

式中，DN、DN_C分别为校正前、后像元的灰度值;DN_m=a×DN²+b×DN+c;a、b、c为回归参数。

参考文献[11]的研究：F16传感器2007年的灯光数据具有最高累积灰度值；黑龙江省鸡西市历年灯光影像的扩展情况与城镇建成区的空间发展变化高度一致。故选取F16(2007)灯光影像的鸡西市数据为标定数据集，确定式(1)中各回归参数，进行拟合回归^[12]，然后，根据式(2)对回归校正后的影像DN值进行修正，完成DMSP/OLS影像数据的校正：

式中，DN_(n-1,i)、DN_(n,i)表示DMSP/OLS夜间灯光影像上i像元在第n-1年、n年的DN值。

首先，利用较小阈值去除背景噪声。将灰度值大于DN_T=19的像元识别为前景像元，作为可能的城镇建成区。

然后，基于像元四邻域(上、下、左、右)的连通区域循环标识算法^[13]，将提取的前景影像识别为由多个像元组成的、相互独立的潜在建成区对象并编号，完成2 728个对象(C₁~C₂₇₂₈)的识别。过程为：从影像左上角开始逐行扫描，遇到第一个前景像元时，将之识别为具有一个单独像元的潜在建成区对象编号为C₁，并作为C₁对象的种子像元标记为已识别；然后，对其四邻域进行扫描，若为前景像元，则编入C₁对象，并加入C₁对象的种子像元群，同时标记为已识别；递归此过程，直至所有与此对象相连通的前景像元被识别出来。继续逐行扫描影像像元，判断其是否为前景像元及是否被标记，若为前景像元且没有被标记，则识别为新的潜在建成区对象C₂，并进一步完成C₂对象的标识。以此类推，直至将全部前景像元标识为潜在的建成区对象。

参考国家统计局发布的《中国城市发展研究报告》，本文将全国按城市化水平和经济发展程度，划分为3个等级进行空间聚类。

基于ArcGIS平台^[14]，获取2 728个潜在建成区对象所包含的像元个数N_i、地理重心点的空间位置和DN_i值大小。其中N_i表征建成区规模，DN_i表征城镇化强度，N_i和DN_i值越大，表示区域的城镇化发展等级越高。然后以N_i和DN_i值为约束条件，以对象地理重心点为节点，基于Jenks自然断点法构造Voronoi多边形，实现方差最小化聚类，聚类结果如图 1(国界数据来源：国家基础地理信息系统矢量数据)所示，其中，第一聚类区域具有最低的城镇化发展等级，对应最小的提取阈值；第三聚类区域具有最高的城镇化发展等级，对应最大的提取阈值。

图 1 不同城镇化发展水平空间聚类结果 Fig. 1 Spatial Clustered Results with Different Urban Development Levels

图选项

然后，根据式(3)，利用方差拟合优度R_GVF对聚类结果进行评定，一定范围内，聚类后的方差之和S_DCM越小，R_GVF越大，分类效果越好。

式中，N表示对象的个数；z_i表示原始数据中第i个对象的DN值；z表示原始数据的平均DN值；k表示聚类数，聚类后数据被分为k类；N_j表示第j类对象的个数；z_ij表示第j类第i个对象的DN值，表示第j类对象的平均DN值；R_SDAM表示原始数据的方差；R_SDAM表示聚类后的方差之和。本文中k=3，对聚类结果进行计算得到：R_SDAM=15.88，R_SDAM=165.17，R_GVF=0.903 9，取得了较好的聚类效果。

借鉴基于统计数据提取城镇建成区的方法，求取各聚类区域的最佳阈值。参考文献^{[15, 16]}对DMSP/OLS数据与人口、经济数据的相关分析，设定各聚类区域的初始阈值分别为35、50、57，提取历年全国范围内的城镇建成区。用式(6)计算提取结果与同期统计数据间的差值，并根据差值修正各聚类区域的阈值，直到提取结果与统计数据充分接近为止^[17]，根据式(7)，得到各聚类区域的最佳阈值为42、54、60，并依据此阈值提取历年的建成区栅格影像：

式中，A_i为设定的阈值序列；A₀为初始阈值序列；A₀={35,50,57}，E(A_i)为根据A_i提取的城镇建成区面积；S为同期统计数据建成区面积；Δ(A_i)为提取结果与统计数据间的差值；若Δ(A_i) <0，则增大阈值，若Δ(A_i)>0，则减小阈值。

式中，Δ(A_i-1)、Δ(A_i)、Δ(A_i+1)分别表示在阈值序列A_i-1、A_i、A_i+1下得到的统计数据与提取结果间的差值。

再次利用连通区域循环标识算法对所提取结果进行建成区对象识别并编号，同时，对建成区对象的细小碎片、边缘锯齿和内部空洞等破损现象，进行几何形态优化。

根据对象的面积属性，应用消除运算删除面积小于30 km²的建成区对象；应用填充运算对对象的内部空洞进行填充。由于建成区对象中的某些空洞是实际存在的，需设定填充空洞的像元数阈值，对小于此阈值的空洞进行填充，大于此阈值的空洞则进行保留，本文设定的填充像元数阈值为20。最后，应用PAEK指数多项式平滑算法，设定多边形平滑容差值为0.5 km，对建成区对象计算出一个平滑的边界，消除边缘锯齿、缺口，实现对建成区对象边界的平滑。如图 2所示为北京市2008年的建成区信息提取结果。

图 2 聚类阈值法提取的北京市建成区空间信息 Fig. 2 Extraction of the Urban Construction Information of Beijing Based on Clustered Threshold Method

图选项

利用历年《中国统计年鉴》中的城镇建成区统计数据对DMSP/OLS数据提取的1992~2011年的建成区面积总量进行精度评定。如表 1所示，历年提取结果的相对误差均在10%以内，误差最大的年份为1996年，相对误差为9.30%，误差最小的年份为2003年，相对误差仅为-0.12%。由此可见，在数量尺度上，利用聚类阈值法提取DMSP/OLS灯光影像的城镇建成区信息取得了较好的精度。

利用Google Earth ^[18]提取2011年城镇建成区影像对DMSP/OLS数据的提取结果在空间格局上进行精度评定。具体地，选取北京、郑州和秦皇岛3个城市代表发展程度不同的“大城市”、“中等城市”和“小城市”，进行对比。

表 2和图 3显示，以Google Earth 提取的城镇建成区信息作为参考，3个城市DMSP/OLS提取结果的精度均在92%以上，DMSP/OLS提取的建成区空间格局与Google Earth 提取的建成区格局特征基本吻合。

图 3 Google Earth提取结果与DMSP/OLS影像提取结果对比分析 Fig. 3 Comparative Analysis on the Extraction Results of DMSP/OLS Data and Google Earth Data in Spatial Scale

图选项

本文提出的聚类阈值法，以所识别的建成区对象为基本空间单元，既消除了以单一阈值提取同一行政单元内部不同城镇化发展等级区域的不合理性，又顾及了建成区对象的规模与发展等级划分。以该方法为基础，本文以DMSP/OLS (Version 4)的32期非辐射标定夜间灯光强度影像为实验数据，提取历年中国区域(本文中不含南海诸岛)的城镇建成区信息。通过与统计数据和Google Earth提取结果进行对比分析，验证该方法在数量尺度和空间格局上的有效性及可靠性。该方法可在一定程度上弥补现有统计资料空间信息不足、逐一确定各行政单元阈值效率低下的缺陷，实现快速、有效、低成本地更新城镇建成区数据，为DMSP/OLS灯光影像充分发挥其时效性和经济性起到积极的推动作用。