本发明属于商业地理分析技术领域，涉及栅格数据处理和商圈划分方法，具体涉及一种基于耗费距离栅格的商圈划分方法。
背景技术：
现代商业经济和信息技术的高速发展，使得商业数据呈几何式增长，商业数据主要表现出海量、多源、异构、不确定性高等特点。据估计商业数据中超过80％的信息与地理位置相关，这使得利用地理信息技术进行空间区位与商业业务信息的综合分析成为可能。
商业地理分析技术的一个重要方面就是进行商圈分析。商圈是以区域内的中心地带为核心开始向外辐射，对消费者产生吸引力而形成的最大辐射范围的空间。商圈分析对于整个城市的商业网点布局与层次结构调整、商业业态创新与结构调整、产业集群与产业园区的规划建设、促进区域经济的发展与城市经济的繁荣都有着重要的意义。
现今，商圈分析的技术主要集中在商圈的界定和划分上，常用的划分商圈的主要方法有顾客点绘法、模型法(reily法则、断裂点模型、huff模型等)、基于gis的商圈划分法。
顾客点绘法(customerspotting)是指将某一商店的顾客的来源地用点表示出来，直接绘制在地图上。其中来源地的信息可以通过问卷调查、顾客记录等方法获得，可以直观的获得顾客的集聚情况及其顾客距离商店的距离。此方法不便于大规模操作。
现有应用中有针对现实的市场功能区域进行空间分析和模拟的模型，如reily的零售重力模型、batty的裂点方程、huff模型、tobler的价格场以及众多的修正式，但是，这些理论模型在应用中存在一个问题，即在分析有多个中心地的市场功能区域时由于缺乏可视性和连续性，其分析结果的成图十分困难。
基于gis的商圈划分方法主要是针对矢量数据处理的voronoi图划分方法。voronoi图具有势力范围的重要性质，能够对地理实体空间影响范围进行划分。voronoi图针对矢量数据的空间剖分，可进行单影响要素的商圈划分，但是，该方法限制了商圈划分过程中的多要素考量，使得商圈划分的准确性不高。
技术实现要素：
为了克服上述现有技术的不足，本发明提供一种基于耗费距离栅格的商圈划分方法，耗费距离栅格的商圈划分是一种针对栅格数据的区域完整剖分方法，选取人口、地价、社交网络签到、回归模型等商业活动要素作为耗费栅格数据，对各种商业活动要素对商圈划分的影响进行计算分析，实现对商圈的划分及其布局优化；能够综合考量多种商业活动要素对商圈进行划分，从而提高商圈划分的准确性。
本发明的技术方案是：
一种基于耗费距离栅格的商圈划分方法，首先选取并构造空间化的商业活动要素，经过对称差处理形成耗费栅格，栅格单元的耗费值代表了穿越该单元的代价，同时栅格以连续场模型来表达地理空间，通过耗费距离分析计算得到的分配栅格是进行商圈划分的基础。包括如下步骤：
a选取并构造空间化的商业活动要素(包括但不限于人口、交通要素、地价、社交网络签到、商业网点)，得到商业活动要素的标准格网栅格图层，包括其像素单元值及对应的坐标点的经纬度值；
在进行商业活动要素研究过程中，考虑到资料获取的难易程度以及影响因素难以量化等原因，本发明主要按照基础性、前瞻性、易获取、可定量、区域全覆盖等原则，筛选对城市商业空间格局形成产生直接影响、而相互之间独立性和可操作性均较强的要素。故而选取了城市常住人口数目、城市交通网络中心性指数、城市建筑物平均房价以及社交网络公众签到(check-in)数量以及以上因素的多元线性回归分析因子，作为本发明方法中的影响商业活动的要素。并将这些商业活动要素作为计算耗费距离栅格的基础。耗费距离栅格的值表示该单元格到最近的源的耗费值。选择计算耗费值的权重即选择对城市商业活动产生重要影响的耗费因子。
a1构造空间化的人口要素：
表征人口要素对商业活动影响的要素包括但不限于常住人口数目、性别比例、民族类别等；所述的常住人口指国家统计局发布的人口普查常住人口数据，本发明具体实施中选取常住人口数目作为人口要素表征指数。首先，将原始人口数据进行预处理，删除部分奇异数据；其次，修正人口数据，把其中建筑物空间数据及其人口数进行合并处理，以求建筑物中常住人口数据与官方公布数据相一致；最后，将上述数据进行面插值处理，运用arcgis软件中fishnet格网生成工具生成500m规则格网，500m为标准的商圈分析空间尺度，按照建筑物单元面积权重将人口分摊到规则格网中，并且将矢量数据转换成栅格数据，由此得到500m常住人口规则格网栅格数据。
a2构造空间化的交通要素：
表征交通要素对商业活动影响的要素包括但不限于道路网络可达性指数、道路网络中心性指数、道路交叉路口数量等；本发明具体实施中选取道路网络中介中心性指数作为反映城市交通要素重要性的指标。
中介中心性指标是用所有节点对之间的最短路径经过给定节点的次数来衡量，这个参数反映出该节点在网络中的中转和衔接功能。
式1中，cb(i)为中介中心性指数；其中，n代表道路网络中节点的个数；njk代表节点j和节点k之间的网络最短路径的总条数；njk(i)代表节点j和节点k之间的最短路径中经过节点i的条数。
根据上述公式计算道路网络中心性指数，形成带有权重的道路网络空间数据，通过空间连接操作，将道路网络的权重信息赋予其最邻近的建筑物空间数据上，使得建筑物空间数据带有道路网络中心性属性。最后，将上述数据进行面插值处理，运用arcgis软件中fishnet格网生成工具生成500m规则格网，按照建筑物单元面积权重将道路网络属性分摊到规则格网中，并且将矢量数据转换成栅格数据，由此得到500m交通规则格网栅格数据。
a3构造空间化的商业网点要素：所述的商业网点要素是指国家进行的经济普查数据，包括每个商业网点的从业人口数、期末总值、经纬度等信息，经过数据预处理，删除部分奇异数据。然后，利用arcgis软件空间连接工具进行500m格网化，按照单元面积权重方法将商业网点属性分摊到规则格网中，并且将矢量数据转换成栅格数据，由此得到500m商业网点规则格网栅格数据。
a4构造空间化的社交网络签到要素：社交网络签到要素的数据是指互联网(包括但不限于新浪微博、大众点评、街旁网、人人网等)上抓取了社交网络签到数据，获取签到地点名、地址、经纬度、城市代码、兴趣点类别、签到次数、照片次数等，经过数据预处理，删除部分奇异数据。然后，经过投影转换和数据纠偏，即得到社交网络签到兴趣点数据。之后，按照建筑物单元面积权重将道路网络属性分摊到500米规则格网中，并且将矢量数据转换成栅格数据，由此得到500m社交网络签到规则格网栅格数据。
a5构造空间化的地价要素：地价要素的数据是指从互联网(包括但不限于链家、房天下等)上抓取了居民小区的房价数据，获取每栋楼的平均房价、小区评价、小区地址、经纬度等，经过数据预处理，删除部分奇异数据。然后，经过投影转换和数据纠偏，即得到房价空间点数据。通过空间连接操作，将房价空间点数据赋值到最邻近的建筑物空间数据上，使得建筑物空间数据带有房价属性。之后，按照建筑物单元面积权重将房价属性分摊到500米规则格网中，并且将矢量数据转换成栅格数据，由此得到500m房价规则格网栅格数据。
a6对以上人口、交通、社交网络、地价、商业网点要素进行多元回归，构造回归模拟要素：选取典型商圈的商业网点期末从业人员数目表征商业网点规模并作为模型因变量；选取人口、交通路网中心性指数、社交网络签到、地价四种商业活动要素作为自变量。在spss统计软件的支撑下，对自变量和因变量进行多元线性回归建模。
y＝-0.044x1-0.009x2+0.059x3+0.285x4(式2)
式2中，y为回归模型因子；x1～x4分别为进行多元回归分析的人口、交通、社交网络签到和地价因子。
经过模型的检验表明，本发明采用的多元线性回归模型满足正态性假设、等方差性假设和独立性假设,具有良好的拟合度和解释能力。由此得到500m回归模拟标准格网栅格数据。
b.计算以上商业活动要素的对称差，得到耗费距离栅格运算所需耗费图层；所述的耗费图层包括但不限于经过对称差处理的人口、交通、社交网络、地价规则格网栅格数据。
由于上述500m规则格网栅格像元值表示的是耗费值，我们选定的栅格图层与耗费的概念是负相关，以人口分布图层为例，也就是说在耗费值大的地方应该是人口少的地方，故在进行计算过程中，需要将选定的栅格图层进行对称差分布处理，即获得耗费图层。
c.选取待划分的商圈的源图层；源图层可以是矢量数据(点、线、面)，也可以是栅格数据,对于栅格数据，要求除标识源以外的单元格为无值。源图层表示的耗费值为零的地方，即是我们选择的将要被划分的商圈中心建筑物的中心点。
d.计算每个耗费距离栅格；计算耗费图层中单元格与源图层中各个源点(即各个商圈中心点)的距离，根据距离数值确定分配给各个源点的像元，生成距离栅格与分配栅格；耗费距离栅格的值表示该单元格到最近的源点的耗费值。最近源是由单元格与源之间的耗费距离决定的。具体计算耗费距离栅格的算法过程：
d1.计算耗费距离栅格；
在图论的表达方法中，各栅格的中心被视为结点，并且各结点通过多条连接线与其相邻结点连接。每条连接线都带有关联的耗费值。耗费值是根据与连接线各端点上的栅格相关联的成本和在栅格中的移动方向确定的，分配给各栅格的成本表示在栅格中移动每单位距离所需的成本；每个栅格的最终值由栅格大小乘以成本值求得。可以用公式3表示：
式中，lcd表示最小累计耗费距离，dij表示地理要素从源点j运动到单元i的空间距离，ci表示单元i对某要素运动的耗费值，∑表示单元i与源j之间穿越所有单元的距离和耗费值的累积，min表示耗费图层的栅格值对于不同的源取累积耗费最小值。
也就是说耗费距离栅格的值表示该单元到最近源的最小耗费值(可以是各种类型的耗费因子，也可以是各感兴趣的耗费因子的加权)。最近源是当前单元到达所有的源中耗费最小的一个源。耗费栅格中为无值的单元在输出的耗费距离栅格中仍为无值。
d2.生成分配栅格；
按照耗费距离栅格中的数值大小将其进行分割，使其分割为不同的区域，按照各个源的值将其填充，即完成分配栅格过程。
e.根据分配栅格计算结果，把整个地理空间划分为以商圈中心点为核的栅格区域，区域个数即为待划分的商圈个数，商圈中心点(源点)内所分配的栅格即为相应商圈的范围，进一步为应用方便，将划分出来的栅格区域矢量化，得到商圈的矢量范围，进而得到商圈的市场域范围。
与现有技术相比，本发明的有益效果是：
本发明的商圈分析主要从城市商圈的角度进行，可应用于整个城市商业格局宏观规划的领域上，聚焦于城市商圈的划分和界定问题上。采用本发明技术方案，一方面可为广大商业投资者及其商业企业在选址开店时，提供所需商圈内的详细资料组分；另一方面可为政府商业管理部门提供商圈发展规划、布局以及商业发展潜力区分析等宏观决策参考。而北京城市商圈在发展过程中经历了进化和升级，取得了长足的进步，但在仍然有一些问题存在。通过本发明提出的方法能够揭示北京城区商圈存在的问题，进而提出未来商圈发展的优化意见，能够为城市商业发展提供技术支持。
附图说明
图1是本发明实施例中基于耗费距离栅格的商圈划分方法的流程框图。
图2为本发明实施例的残差直方图；
其中，因变量为商业网点期末从业人数；标准偏差为0.995。
图3为本发明实施例中人口500m格网数据处理流程图。
图4为本发明实施例中500m常住人口规则格网数据图。
图5为本发明实施例中水平与垂直结点计算示意图。
图6为本发明实施例中累积耗费结点计算示意图。
图7为本发明实施例中斜对角结点计算示意图。
图8为本发明实施例中源栅格与成本栅格；
其中，a为源栅格；b为成本栅格。
图9为本发明实施例中初次迭代的输入栅格与活动累积成本像元列表；
其中，a为输入栅格，b为活动成本像元列表及图例。
图10为本发明实施例中持续迭代过程图；
其中，a为输入栅格，b为活动成本像元列表及图例。
图11为本发明实施例中获得更廉价累积成本路径图；
其中，a为输入栅格，b为活动成本像元列表及图例。
图12为本发明实施例中耗费距离栅格。
图13为本发明实例中的源栅格、耗费距离栅格和分配栅格；
其中，a为源栅格，b为耗费距离栅格，c为分配栅格。
图14为本发明实施例中基于耗费距离栅格的北京区域级商圈划分结果图；
其中，a是道路因子商圈划分效果图；b是社交网络签到因子效果图；c是人口因子商圈划分效果图；d是多元回归模拟因子商圈划分效果图；e是商业网点密度因子商圈划分效果。
图15为本发明实施例中基于耗费距离栅格的北京商圈划分结果对比图；
其中，横坐标为商圈；纵坐标为商圈面积；其结果表明不同商业要素影响下的北京商圈界定结果趋同，证明算法的稳定性与准确性。
具体实施方式
下面结合附图，通过实施例进一步描述本发明，但不以任何方式限制本发明的范围。
本发明属于商业地理分析技术领域，具体涉及一种基于耗费距离栅格的商圈划分方法，并利用该方法表达各种商业活动要素对商圈划分的影响，并应用于北京商圈划分界定中。
该方法的具体操作步骤如图1所示，下面进行详细叙述：
a.源数据层选取
本发明选取待划分的北京市的27个区级商圈作为源图层，如表1所示：
表1北京主要商圈发展概况
这27个区域级商圈遍布于北京城中心城区、规模等级相当。
本发明将上述基于耗费距离栅格的商圈划分方法应用于北京市的商圈划分界定，根据不同的商业活动要素对商圈进行划分界：
a1.道路网络因子的商圈界定表达出更多的城市交通的细节部分，即商圈内部交通路网中心性指数数值越高，耗费距离也就越大，所形成的势力范围也就越大。
a2.社交网络签到因子的商圈界定更多的表现出城市热点程度，即商圈内部社交网络签到的数目越多，所形成的势力范围也就越大。
a3.人口要素是商圈得以建立和发展的基础，基于人口要素的商圈界定在与其他的界定结果相比，虽整体上基本保持一致，但具体商圈内部差异较大。
a4.商业网点规模因子结果所确定的商圈划分与现阶段商业网点的期末从业人数相匹配，表现出其在商圈划分中的优势地位。
a5.回归模型因子是对上述各因子进行多元回归分析后模拟出来的结果，所形成的商圈划分结果与其他因素相比保持了较高的一致性，并且能够综合上述各因子优点。
本发明采用标准回归方程如下：
y＝-0.044x1-0.009x2+0.059x3+0.285x4(式4)
式4中，y为回归模型因子；x1～x4分别为进行多元回归分析的人口、交通、社交网络签到和地价因子。
通过方程显著性检验可知，方差分析表中(如表2所示)回归的f值为8.909，其显著性概率为0.000<0.05，模型的整体拟合效果通过了显著性检验。
表2回归线性模型方差分析表
回归系数表中(如表3所示)进入模型的变量t的显著性检验都通过了显著性检验。
表3回归方程相关系数表
残差直方图中(如附图2所示)，接近于标准的n(0,1)正态分布。
b.计算耗费栅格
由前所述，本发明选取了城市常住人口数目、城市交通网络中心性指数、城市建筑物平均房价以及社交网络公众签到(check-in)数量四种指标，作为计算耗费栅格的权重。下面以人口因子为例，来说明耗费栅格如何生成。
b1.数据选择
采用2008年全国第二次经济普查时制作的北京市建筑物空间单元统计数据作为原始统计区集合，以建筑物为统计单元的人口数据可近似作为人口统计的总体。人口统计指标选取北京城六环内2010年11月1日零时的常住人口数目。
b2.数据网格化
具体过程如附图3所示：
b2.1将原始人口数据进行预处理，删除部分奇异数据；
b2.2修正人口数据，把其中建筑物及其人口数进行合并处理，以求建筑物中常住人口数据与官方公布数据相一致；
b2.3将上述数据进行面插值处理，运用arcgis软件生成500m规则格网，按照建筑物单元面积权重法将人口分摊到规则格网中，由此得到城六区500m常住人口规则格网数据(如附图4所示)。
b3.对称差处理
由于栅格像元值表示的是耗费值，与栅格图层成负相关，所以要对选定的栅格图层进行对称差分布处理。即将上述计算生成的人口因子插值图层进行对称差处理即为本发明所用的耗费栅格图层。
c.计算耗费距离栅格图层
使用图论创建累积成本距离栅格可被视作尝试识别最低成本像元，并将其添加到输出列表。它能通过若干次迭代建立一个累积成本像元列表，创建一个输入栅格，不断扩散累积成本像元列表中最小值周围的像元创建可能的最低成本像元，通过迭代不断更新累积成本像元列表中的最小值并输出到输出栅格中，最终实现每一个单元格到自己最近源的路径值均为最小。这是起始于源像元的迭代过程。每个像元的目标是快速分配到输出成本距离栅格中。下面具体叙述其工作过程：
c1.累积成本计算
c1.1相邻结点成本
从一个像元移动到四个与其直接连接的近邻之一时，跨越连接线移动到相邻结点的成本为用1乘以像元1与像元2的和，然后再除以2，如式5：
a1＝(cost1+cost2)/2(式5)
其中cost1是像元1的成本，cost2是像元2的成本，a1即从像元1到像元2连接线的总成本。如附图5所示：
c1.2累积垂直成本
累积成本由以下公式确定：
accum_cost＝a1+(cost2+cost3)/2(式6)
其中：cost2为像元2的成本，cost3为像元3的成本，a2为从像元2移动到3的成本，accum_cost为从像元1移动到像元3的累积成本。如附图6所示：
c1.3对角结点成本
如果沿对角线移动，则连接线上的行程成本为1.414214(或2的平方根)乘以像元1的成本加上像元2的成本，再除以2：
确定对角线移动的累积成本时，使用以下公式：
计算斜对角结点如附图7所示：
c2.第一次扩展输入栅格
如附图8所示为源栅格图层(图8a)与成本(本步骤中的成本概念等同于耗费)栅格图层(图8b)。图8左所示源栅格图层中数字1与2标识了源点。初次迭代中，在输入栅格中识别出源像元并分配0值，因为它们返回自身不消耗累积成本。接下来，启用全部源像元的近邻，使用c1.1-c1.3计算累积成本公式将成本分配到源像元结点与邻近像元结点之间的连接线。计算结果如附图9中的a输入栅格所示。各邻域像元都可以达到某个源。因此，可以选择它们或将它们分配给输出累积成本栅格。要分配到输出栅格，像元必须具有到达某个源的最低成本路径。同时将扩展的各个像元的累积成本值按由低到高的顺序排列于活动累积成本像元列表中如附图9中的b所示。
c3.扩大邻域
再次选择列表中具有最低成本的活动像元，扩大邻域，计算新的成本，并将新的成本像元添加到活动列表。在上图活动累积成本像元列表中最小值为1.5，则扩展其邻域像元，同时将此像元被输出到输出栅格中。更新累积成本像元列表，此时最小值为2.0，即此时将要扩展2.0邻域。从图9可知，2.0周围像元已全部扩展。
c4.持续迭代
此过程持续进行，如果通过将新像元位置添加到输出栅格创建新的成本较低的路径，则将更新活动列表上的像元。如附图10所示为扩充第3行第4列像元的过程(从第1行第1列算起)。图10a中扩散像元成本4.9的位置即单元格(3，4)时，发现其下面单元格无值，则将无值像元直接输出到输出成本，仍为无值，即无值像元不做处理。如附图10a的输入栅格中所示10.6为新扩充像元成本值。图10b为更新的活动累积成本像元列表。
c5.获得更廉价累积成本路径
当活动列表上出现新的像元路径时会进行此更新，因为更多的像元被分配到输出栅格。当活动累积成本列表上具有最低值的像元被分配到输出栅格时，计算所有累积成本。也会计算新分配的输出像元的相邻像元的成本，即使相邻像元位于其他像元的活动列表上。如果活动列表上的位置的新累积成本大于这些像元当前的累积成本，则忽略该值。如果活动列表上的位置新累积成本小于这些像元当前的累积成本，则使用新值替换掉活动列表上该位置的原有累积成本。此时已具有到达某个源的更廉价和更理想路径的像元在活动选择列表中上移。
比如，将附图10中的a输入栅格中的第3行、第1列的像元位置放在活动列表上时，它达到栅格顶部的源的累积成本为11.0。然而，因为较低的源扩展到此位置，该像元获得了到达其他源的更廉价的累积成本路径。由于存在这一较低的累积成本，因此应早些在活动列表上更新该位置的值，并分配到输出。更新后的输入栅格如附图11中a所示，活动累积成本像元列表如附图11中b所示：
c5.得到成本距离输出栅格
如图12所示为最终迭代生成的成本距离输出栅格，即耗费距离栅格图层。
d.同步分配栅格
按照耗费距离栅格中的值的大小(如图13.b)将源图层(如图13.a)与耗费图层进行划分，即完成分配栅格过程(图13.c)。
划分结果说明
e1.图14中的a为道路网络因子的商圈划分结果：选取道路网络因子为耗费图层的商圈划分结果可以看出，亚奥商圈、崇文门商圈范围远大于其他商业活动要素作为耗费图层的划分结果，体现出该地区道路网路优势使得商圈向竞争域外扩张的态势。其中，亚奥商圈由于奥运会建设工程的影响，商圈势力范围内交通较为便利，使其竞争域范围明显扩大。崇文门商圈内交通主要依靠前门东大街、崇文门东大街，珠市口东大街、广渠门大街，崇文门外大街形成的两横一纵的“王字”交通网络格局，使得商圈辐射能力进一步扩大。
e2.图14中的b为社交网络签到因子的商圈划分结果：从各商业活动要素结果对比可以看出，工体-三里屯商圈和五道口商圈的势力范围明显大于其他图层中该商圈的势力范围。主要原因是社交网络签到主要是在移动终端上进行，并且年轻人居多，同时这两个商圈消费群体也主要是年轻人，使得这两个商圈在以社交网络因子的商圈界定过程中表现出较强的辐射能力。
e3.图14中的c为基于人口要素的商圈划分结果：虽整体上基本保持一致，但具体商圈内部差异较大。商圈内常住人口数量多一定程度上提供了客源，另一方面也使得商圈的核心竞争力不断增强，表现出竞争力高的商圈在压制周边较弱商圈时的绝对优势，当竞争力平衡时表现出一种商圈间的互利共存的稳定的商业生态系统。
e4.图14中的d为商业网点规模因子商圈划分结果：从各商圈所占面积来看，广渠门商圈相比其他结果势力范围明显大于其他图层中该商圈的势力范围。由于广渠门商圈主要还是居住为主的地区，商圈内居民服务业态商业网点及其发达，也就增长了商业网点整体的从业人数，使之在商圈范围上有所外扩。
e5.图14中的e为回归模型因子商圈划分结果：回归模型因子所形成的商圈划分结果与其他因素相比保持了较高的一致性，在各商圈界定上均没有极值的情况出现，回归模型较好地综合各种商业活动要素对商圈空间结构的影响。
从商圈划分结果(图14、图15)可以看出，利用各种商业活动要素的商圈划分总体上趋同(图15曲线形状基本一致)，但局部商圈内部又存在细节差异，从而反映出各种商业活动要素对商圈划分的作用。
需要注意的是，公布实施例的目的在于帮助进一步理解本发明，但是本领域的技术人员可以理解：在不脱离本发明及所附权利要求的精神和范围内，各种替换和修改都是可能的。因此，本发明不应局限于实施例所公开的内容，本发明要求保护的范围以权利要求书界定的范围为准。