杜 晓,周玉秀,周 琦,郑 义,张俊辉,林尚纬,陈利军
(1.国家基础地理信息中心,北京 100830; 2.青海省基础测绘院,青海 西宁 810001)
热力图(heat map)一词最初是由软件设计师Cormac Kinney于1991年提出并创造的,即用一个二维图片描述显示实时金融市场信息。最开始的热力图是矩形色块加上颜色编码,经过多年的演化,习语上的热力图被大多数人理解为经过平滑模糊处理的热力图谱。在基于地理信息的大数据应用环境中,热力图一般是指通过密度函数进行可视化,用于表示地图中点的密度的热图,反映较大空间范围内观测量的不同。
随着我国自主可控高分辨率卫星在轨数量的极大增加,以及地理数据获取、处理和建库能力的不断提高,开展全球地理信息资源建设的条件已经具备,全球地理信息数据正在逐步完善和丰富。在已有成果的基础上,实现基于地理信息数据的统计分析,开展基于现有成果数据的深度挖掘,对政府决策和公共服务显得尤为重要[1]。现有的全球地理信息资源成果主要包含数字正射影像、数字表面模型、核心矢量要素、地名、地表覆盖等数据。在信息提取和数据挖掘方面,国内的研究多集中在影像目标识别、区域统计、聚类分析等传统的遥感影像与测绘分析手段,而对如何将地理信息蕴含的本质特征和空间分布与社会公众的直观感受结合起来的探索略显不足[2]。对于核心矢量要素数据,其空间特征、要素属性、表达方式等包含了大量的与人类活动息息相关的信息,可以作为反映经济发展水平、评价社会发展状况的有效数据源。本文基于核心矢量要素数据中的道路网,尝试建立道路热力分布计算规则,以直观新型的热力方式反映世界部分国家交通情况的发展水平。
传统的热力计算以带有空间坐标的点位为基础,通过其属性值、权重和影响范围,计算图面内每个栅格点位的密度重叠等信息,可以准确直观地识别关注信息的空间分布。基于离散点的热力图计算已经在各级地理信息数据中得以全面应用。但在基于道路网的道路发达水平热力计算中,以下3个问题亟需解决:①道路为线状要素不能直接参与散点空间计算;②道路等级的高低对热力计算的影响应远大于其本身的空间分布密度;③搜索半径与输出格网大小应根据图面范围实时调整。部分学者也对线状地物的缓冲算法进行了深入的研究和探索[3],同时Esri公司也提供了相应的核密度分析工具[4],都可以作为良好的借鉴。本文拟建立基于全球道路网的道路发达程度热力计算技术流程,构建相关规则体系,实现道路发达水平信息的自动化处理与直观展示。
道路发达程度热力计算技术流程主要包括3部分:道路赋权、道路线曲面化、整体优化协调。具体如图1所示。
图1 道路发达程度热力计算技术流程
图1中虚线图框所示为在不同流程阶段对应的处理规则,主要包括路线类型(铁路/公路/地铁等)的权重赋值、道路技术等级权重赋值、分段道路自动化融合、线到面的反距离权重计算、分布密度/叠加密度协调、全图综合调整等数十条规则项。
热力计算的重点是建立单位面积内道路发达程度与格网热力值之间的对应关系,而道路的类型、等级、密度、碎片化、缓冲半径、输出格网大小等因素与道路发达程度的计算结果均有相关性,其相关性为综合影响且非简单的线性叠加对应关系。热力计算规则的合理制定可有效避免单一因素对道路发达程度判断的过度影响,全面考虑各因素对热力值的影响,提高整体的计算效率和准确性。下面以道路赋权、曲面化和优化协调的相关规则类进行重点阐述。
当前全球矢量要素数据中的道路网数据基本以“数据层-要素-属性项”为层次组织存储,为综合考虑其影像因子,需建立数据层、要素、属性三位一体的分级对应赋权规则。通过分别设定权重,经样本训练,得到合理的权重设定规则。道路赋权规则主要包括:
(1)道路层赋权,针对不同类型道路如铁路、公路、地铁、轻轨,设定相应的发达程度影响权重。初始权重种子值设定为铁路∶公路∶地铁∶轻轨=5∶3∶4∶4,各图层的权重为WLayer。
(2)技术等级赋权,公路按照技术等级分为高速公路、主要公路、次要公路、三级公路、其他道路、未分级等多种类型;铁路基本分为单线路和双线路两种类型。因此只对公路和铁路进行按照权重赋值,公路技术等级的初始分值种子值设定为:高速∶主要∶次要∶三级∶其他∶未分级=9∶6∶3∶2∶1∶1,铁路技术等级的初始分值种子值设定为:双线∶单线=8∶5。各技术等级道路的分值为VGrade。
根据道路发达程度与道路线位置及分值的关系,分值在线所在位置处最大,随着与线的距离的增大而逐渐减小,在与线的距离等于指定的搜索半径的位置处分值为零。以Silverman的四次核函数为基础,定义线的核函数[*],道路曲面化的方式如图2所示。
图2 道路线曲面化计算方式
图2为道路线与覆盖在其上方的核曲面。核曲面(道路分值面)与下方的平面所围成的空间体积等于道路长度与道路分值的乘积。
综合考虑上述赋权要求,假设某输出栅格点在其搜索半径R内存在n个类型道路,每个类型道路又有m个道路等级,栅格点与各道路的距离为dij,则该栅格点的综合分值计算方程式为
(1)
式中,Vpoint为输出栅格点的综合分值;WLayeri为路网中第i类型道路的权重;VGradeij为第i类型道路中第j种技术等级的分值。
基于全球道路数据为多段线路的情况,为减少道路碎片化对综合分值的影响,在曲面化之前需做好分段道路的拼接实体化处理,确保同等级同类型道路形成完整的线。
2.3.1 热力值均衡化
上述计算得到的Vpoint为道路分值,需要转换为可以查看的热力值Hpoint。热力值是反映道路发达程度的直观度量,其值应在道路分值的基础上作适当均衡调整。调整应限制在当前可视区域内,根据道路分值的直方图分布作热力值的调整。由于热力值为无量纲的道路发达程度度量值,对其高值区间和低值区间内的变化情况一般不作关注,调整方式可采用直方图规定化后再做归一化赋值,一般采用如图3所示的直方图拉伸方式。
图3 热力值直方图拉伸曲线
2.3.2 叠加密度补偿
道路发达程度体现的是道路空间密度和道路等级的综合,虽然已对高等级路设置了较高权重,但为避免大量密集的低等级道路对热力值的影响,需要进行基于道路等级的空间密度值调整。计算不考虑道路分值的线密度分布图Mdensity,与上述计算的热力值分布图Mheat进行叠加。新的热力分布图计算公式为
Mheat_new=Mheat+K×Mdensity
(2)
式中,K为调整系数,一般设置为-0.2。
2.3.3 参数迭代调整
搜索半径、输出格网大小对热力值的表现具有较大影响,其值的选择应确保既能真实反映实际道路发达程度,又不使结果过于模糊。根据当前视图的范围,设定初始搜索半径为20 km,格网大小为500 m。经过反复迭代,根据不同国家的实际情况,较合理的搜索半径可设定为15~70 km;根据表现结果的精细程度,格网大小可设定为100~900 m。
笔者以上文提出的道路网热力计算规则为核心,完成相应函数和软件模块开发。结合当前已有的全球道路网数据的实际情况,通过参数设定—程序迭代—目视判读的方式进行了数百次调整。初步应用于全球地理信息资源成果应用系统中,并以亚洲的不丹、伊朗和非洲的摩洛哥为例。如图4—图6所示。
图4 不丹道路及发达程度热力图对照
图5 伊朗道路及发达程度热力图对照
图6 摩洛哥道路及发达程度热力图对照
经过反复迭代,不丹和摩洛哥的搜索半径设定为18 km,伊朗由于国土面积较大,搜索半径设定为60 km。经过多次迭代,较为合理的道路类型权重设定为铁路∶公路∶地铁∶轻轨=4∶3∶6∶6,道路技术等级权重设定为高速∶主要∶次要∶三级∶其他∶未分级=11.3∶5.5∶2.8∶1.7∶0.8∶0.8,铁路技术等级权重设定为双线∶单线=7.3∶6.1。从热力图结果来看,不同区域内的热力值基本反映了重点城市发达程度和主要道路的分布状况,对宏观了解全球大范围区域的道路建设发展程度具有较好的指导意义。
本文利用热力图简洁直观的优势,修正了传统热力计算中过分依赖空间密度的缺点,依据道路类型和道路等级,提出了分级赋权重确定分值的思路,并以试验进行了充分验证,效果较好。但由于道路发达程度还与多类其他因素相关,如何充分利用车站、码头等道路附属设施的位置和属性,进一步提高道路热力图的客观真实性,相关的关键技术、规则制定和系统实现有待后续研究完善。