细粒度的中国地理数据分析（四）

本文利用 Delaunay 算法进行地块的邻居识别，并将计算统计的结果进行展示。开源代码可见我的前端笔记本

Neighborhood of land types^[1]

细粒度的中国地理数据分析（四）^[2]
- 转移矩阵的计算方法^[3]
- 转移矩阵的呈现^[4]
- 附录：Delaunay 计算的核心代码和输出^[5]

转移矩阵的计算方法

由于 Delaunay 方法已经在前文有所提及，因此本文不再赘述。本文使用 Python 对全部地块进行计算得到三角关系如附录代码所示。为了简化起见，我们采用其中 simplices 字段的数据进行分析，它代表平面三角形的三个顶点的集合。由 Delaunay 算法的原理可知，这些三角形将各个地块的中心紧密地连接起来，三角形的边可以看作是两个邻居节点之间边，因此这些边就描述了地块之间的邻居关系。

将它们按地块种类进行统计，统计的内容为邻居地块的数量，得到统计结果如下。统计结果用的矩阵所示，每个行、列分别代表种地块之一，其中的数字代表该地块的邻居属于另一种地块的数量有几个。举例来说，512000 的四川省某行政区划内，其全部居住地块（首行）中，其邻居为居住（1）、商业（2）、工业（3）、交通（4）和公共（5）地块的数量分别为由于邻居关系是无顺序的，因此该矩阵为对称矩阵。

该矩阵的意义可以理解为，当你身处某个地块之中，向其他地块“转移”时，新地块的种类的概率分布，因此我将该矩阵称为“转移矩阵”。

{512000: ('四川省',
  array([[292.,  32., 158.,   8.,  30.],
         [ 32.,  16.,  16.,   0.,   0.],
         [158.,  16.,  84.,  10.,  28.],
         [  8.,   0.,  10.,   0.,   2.],
         [ 30.,   0.,  28.,   2.,  20.]])),
 469000: ('海南省',
  array([[1744.,  224.,  868.,   46., 1068.],
         [ 224.,   68.,   94.,    0.,  120.],
         [ 868.,   94., 1260.,   50.,  998.],
         [  46.,    0.,   50.,    8.,   56.],
         [1068.,  120.,  998.,   56., 2216.]])),
 520200: ('贵州省',
  array([[336.,  42., 156.,   0., 126.],
         [ 42.,  12.,   6.,   0.,   6.],
         [156.,   6., 232.,   0.,  90.],
         [  0.,   0.,   0.,   0.,   0.],
         [126.,   6.,  90.,   0.,  88.]])),
 370700: ('山东省',
  array([[10412.,   554.,  2710.,    74.,  2640.],
         [  554.,   176.,   168.,    12.,    66.],
         [ 2710.,   168.,  2884.,    42.,  1548.],
         [   74.,    12.,    42.,     0.,    42.],
         [ 2640.,    66.,  1548.,    42.,  2192.]])),
 430100: ('湖南省',
...
  array([[1800.,  242.,  678.,   10.,  518.],
         [ 242.,  140.,  124.,    2.,   66.],
         [ 678.,  124., 1204.,   16.,  624.],
         [  10.,    2.,   16.,    4.,   20.],
         [ 518.,   66.,  624.,   20.,  764.]]))}

转移矩阵的呈现

针对转移矩阵进行呈现的方式有很多，最直接的保持矩阵的形式，如下图所示。下图中我们比较了北京市（右）和上海市（左）的地块转移概率，其横坐标代表从该地块出发，而纵坐标代表抵达的地块类型。从图中可以看到

两个城市都有大量的居住地块，从居住地块（首列）出发，能够到达的邻居地块大概率（）还是居住地块；
两个城市的不同点在于工业地块和公共地块，上海市的工业地块与居住地块较为相似，仍显示出聚集的特点；而北京市的公共地块则更加容易到达。

转移矩阵，北京市-上海市

我拉天津市进来进行对比，发现天津市和上海市的转移概率矩阵相似性更强。

转移矩阵，天津市-上海市

接下来，我们显然不能满足于矩阵式的呈现方式，因为它看不到全局的情况分布。因此，我提供了交互式的呈现方案。用户可以选择感兴趣的对比维度，将全国所有城市的对应维度信息都展示出来。如下图所示，图中横坐标代表该种类地块的面积占比，纵坐标代表转移矩阵中该维度的数值。从图中可以看到，北京市与全国其他城市之间的差异均较大。

接下来，再将全部个维度的信息展开，形成如下图宫格的形式。从图中可以看到不同种类的地块具有不同的转移概率分布模式

首先，居住、工业和公共地块呈现出地块面积越大，其内部转移概率也越大的分布规律，这代表它们的分布相对集中，（见图中对角线的子图）；
其次，第 2、4 列的子图表明，商业和交通地块的面积虽然不大，但它们转移到其他地块的概率却较大，这代表了商业和交通地块被其他地块所包围；
另外，第 5 列的公共地块除了具有以上特性之外，其向居住和工业地块转移的概率呈现“面积越大，转移概率越小的趋势”，这暗示了这些地块呈现出所谓“高内聚、低耦合的特点”；
高内聚、低耦合的特点在居住地块向工业地块的转移中也有所体现。

附录：Delaunay 计算的核心代码和输出

'''
The output of Delaunay module in python
'''

from scipy.spatial import Delaunay

position = table[['Lat', 'Lon']].to_numpy()
delaunay = Delaunay(position)
delaunay.__dict__

# The output of delaunay.__dict__ is following

{'_qhull': None,
 'paraboloid_scale': 0.00903880731831826,
 'paraboloid_shift': -65.40899681217851,
 'simplices': array([[331453, 398246, 155552],
        [155570, 155566, 155552],
        [129221,   1667,   1594],
        ...,
        [410183, 410194, 410186],
        [410198, 410194, 410192],
        [410200, 410194, 410198]], dtype=int32),
 'neighbors': array([[   297,      4,     18],
        [    -1,    132,    289],
        [   301,   2203,      3],
        ...,
        [764724, 757687, 767210],
        [767210, 762046, 767213],
        [767212, 762054, 764725]], dtype=int32),
 'equations': array([[  0.11548202,   0.93673225,  -0.33044908, -92.44259623],
        [  0.29560449,   0.88340894,  -0.3636023 , -89.78548361],
        [  0.78081188,   0.54314449,  -0.30875048, -77.37581242],
        ...,
        [  0.21400833,   0.8794097 ,  -0.42525171, -81.09379807],
        [  0.21401068,   0.87940924,  -0.42525148, -81.09382386],
        [  0.21401068,   0.87940924,  -0.42525148, -81.09382386]]),
 'coplanar': array([[398247,      0, 398246],
        [155553,      0, 155552],
        [ 12616,     14,  12614],
        ...,
        [388812, 767120, 388814],
        [388813, 767120, 388814],
        [388809, 767120, 388814]], dtype=int32),
 'good': array([1, 1, 1, ..., 1, 1, 1], dtype=int32),
 'nsimplex': 767214,
 '_transform': None,
 '_vertex_to_simplex': None,
 '_vertex_neighbor_vertices': None,
 '_vertices': array([[331453, 398246, 155552],
        [155570, 155566, 155552],
        [129221,   1667,   1594],
        ...,
        [410183, 410194, 410186],
        [410198, 410194, 410192],
        [410200, 410194, 410198]], dtype=int32),
 '_points': array([[ 33.24462, 104.24732],
        [ 33.24613, 104.24452],
        [ 33.24644, 104.24623],
        ...,
        [ 27.23314, 107.8998 ],
        [ 27.23314, 107.8998 ],
        [ 27.24416, 107.86895]]),
 'ndim': 2,
 'npoints': 440798,
 'min_bound': array([18.20895, 75.24152]),
 'max_bound': array([ 50.80089, 133.00325]),
 'furthest_site': False}

参考资料

[1]

Neighborhood of land types: https://observablehq.com/@listenzcc/neighborhood-of-land-types

[2]

细粒度的中国地理数据分析（四）: #细粒度的中国地理数据分析四

[3]

转移矩阵的计算方法: #转移矩阵的计算方法

[4]

转移矩阵的呈现: #转移矩阵的呈现

[5]

附录：Delaunay 计算的核心代码和输出: #附录delaunay-计算的核心代码和输出