共享单车停车需求识别与地图匹配

在我们日常生活中，大家都骑过共享单车，现在的共享单车app里面都有这样一个功能：在地图上规划出了很多的单车停放栅栏（停放点），用户需要把单车停在停放栅栏内，否则就要扣钱交调度费。但有时候停放栅栏的规划并不合理，城市里有些地区的停放栅栏比较密集，有些地方的停放栅栏却很稀疏，骑车骑到一个较远的地方，却找不到停放点，用户就只能多花钱，同样地，对于共享单车公司来说，停车栅栏规划的不合理，也加大了调度的人力时间。所以根据单车停放需求来就合理规划停车栅栏就很重要。如果我们能识别出来某些停车需求量很大的地区，这些地区的停放栅栏却很稀少，那就可以方便共享单车公司作进一步的业务调整。同时，我们也可以动态地探索停放需求量变化规律，协同调整。

要达到上述目标，需要完成以下任务：

在某个时间点（比如早晨8点），有哪些单车处于停放状态？
这些停放状态的车，哪些停在停放栅栏内部，哪些在停放栅栏外部？
在停放栅栏外部的那些单车，距离它们各自最近的停放栅栏是哪些？单车和最近停放栅栏之间的距离是多少？（寻优匹配问题）
不同时段单车停放点需求量怎么变化？哪些区域的停放需求更高？停放栅栏怎样增加？停车需求量随时间如何变化？

该项目用到两份数据：
1.厦门市某共享单车公司2020年12月25日的共享单车订单数据
2.厦门市共享单车电子栅栏地理信息数据。

数据大致信息大致如下：

共享单车订单数据：

样本量：232806*5

变量	含义	类型
BIKE_ID	单车订单号	object
LATITUDE	单车所在位置纬度	float64
LONGTITUDE	单车所在位置经度	float64
LOCK_STATUS	单车开/锁状态	int64
DATA_TIME	单车改变状态的时间	object

电子栅栏地理信息数据：

样本量：14071*2

变量	含义	类型
FENCE_ID	停车栅栏ID	object
FENCE_LOC	停车栅栏位置	object

用到的python算法包有：
1.pandas
2.numpy
3.matplotlib
4.scipy : kd tree寻优算法
5.geopandas：处理空间地理数据的利器，基本可以替代Arcgis或者QGIS
6.shapely.geometry: 存储点、线、面要素

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import geopandas as gpd
from shapely.geometry import Polygon,Point,LineString
import scipy

地理空间数据可视化：kepler.gl

一、共享单车GPS信息数据处理

共享单车GPS数据是一个dataframe,状况如下：
这份数据的记录逻辑是这样的：共享单车的锁状态改变一次，数据库就会生成一条这辆单车的状态数据。比如有人扫A车开锁了，这时数据库里面就会记录一条A车的信息：单车的ID、所在位置的经纬度信息、锁状态为0(开锁)，状态改变时的时间；用户骑行结束锁了A车，又会生成一条关于A车的数据：单车的ID、所在位置的经纬度信息、锁状态为1(上锁)，状态改变时的时间。

2.1单车运行轨迹处理

我们想看每一辆车的行动记录，用sort函数先对BIKE_ID进行排序，再在每辆单车的记录当中对DATA_TIME进行排序，这样就能看到每辆单车在一天之内行动轨迹、状态变化。

data_bike=data_bike.sort_values(by=['BIKE_ID','DATA_TIME'])

2.2 停放状态单车的提取

由于我们想要的是单车停放信息的数据，所以该dataframe里面处于停放状态的记录是我们需要提取出来的。而且我们希望知道它停放的时段，也就是某一辆车从停放状态转变为开锁状态。即LOCK_STATUS从1变为0的这一段时间。这里需要用.shift()函数对dataframe做一个上移操作，这时处理时序数据经常会用到的技巧，看下面图示：
把原来的dataframe上移一个数据，让每条数据和它自身的下一条数据拼接在一起，就像红框框起来的数据是我们想要的，它刻画了车1从停放（状态1）变为开锁（状态0）的这一段时间。我们要把数据中这样的记录提取出来，用来分析单车停放需求。车1从开锁（状态0）变为停放（状态1）的这一段时间是在骑行的，这一段动态的变化我们不需要，删掉这样的数据。还有前面是两辆不同车辆拼接在一起的记录也要删掉（车1拼接车2），没有什么用处。
代码如下：

for i in ['BIKE_ID','LATITUDE','LONGTITUDE','LOCK_STATUS','DATA_TIME']:#每一列上移一位data_bike[i+'1']=data_bike[i].shift(-1)
#提取出前后车辆ID一致的单车记录
data_bike=data_bike[data_bike['BIKE_ID']==data_bike['BIKE_ID1']]
#提取出锁状态从停放（1）变为开锁（0）的单车记录
data_bike=data_bike[data_bike['LOCK_STATUS']==1]
#取出有用的列：单车ID、单车停放位置的经纬度、单车开始停车的时间、单车结束停车（开锁）的时间
data_bike=data_bike[['BIKE_ID','LATITUDE','LONGTITUDE','DATA_TIME','DATA_TIME1']]
#DATA_TIME是开始停车的时间，DATA_TIME1是停车结束的时间

得到的datdaframe如下：
提取出来的单车记录就只剩25000左右的数量了。

2.3 某一个时间点处于停放状态单车的提取

比如要提出在早晨8：00处于停放状态的单车记录，这样操作：

timestamp='2020/12/25 8:00:00'
parking_points=data_bike[(data_bike['DATA_TIME']<=timestamp)&(data_bike['DATA_TIME1']>=timestamp)]

先生成一个时间戳，早晨八点处于停放状态的单车应该是：它的开始停止时间早于8：00，结束停止时间晚于8：00。提取出来的数据记为parking_points，等会就用它来进行停车需求量的分析，（当然可以用所有时间点的数据来分析，这样更全面，效果更好），parking_points状况如下：

2.4 dataframe数据转为Geodataframe数据形式

地理信息数据要用Geodataframe来存储，说白了，Geodataframe就是在DataFrame的基础上增加了一列geometry，这一列里面包含了shapely的图形要素：多边形（Polygon）、线形（LineString）、点（Point）对象。
之后会继续更新关于Geopandas的教程，先把关于该项目的重要关键步骤运用熟练即可。

把DataFrame转换为GeoDataFrame的要点：

为DataFrame加一列geometry，其中每一行都要添加shapely的图形要素。
把DataFrame的类型转换为GeoDataFrame。

在这里我们要把parking_points转换成GeoDataFrame，就是给它加一列deometry,里面的点要素（Points）是它的经纬度信息。

parking_points=gpd.GeoDataFrame(parking_points)
parking_points['geometry']=gpd.points_from_xy(parking_points['LONGTITUDE'],parking_points['LATITUDE'])

这时候我们对parking_points做个可视化来看一下分布状况，直接.plot()：

parking_points.plot()

这就是整个厦门岛在8：00时刻的共享单车停车分布状况

二、电子停车栅栏信息的处理

我们希望在地图上同时显示停车点和停车栅栏，显示每一量单车ID是在停车点外还是停车点内。如下图所示：绿色部分是电子停车栅栏区域，蓝色的是单车停放位置点。接下来处理电子栅栏的地理信息数据，很多步骤和上面一部分其实也有相似部分。
停车栅栏数据如下：FENCE_ID是停车栅栏的编号，FENCE_LOC里面的点是停车栅栏的顶点。现在的这份停车栅栏数据，还是dataframe格式，把它转换成Geodataframe的时候，就是让FENCE_LOC作为geometry这一列，这一列里面存储电子栅栏的Polygon（面要素），它就是由FENCE_LOC里面的这些顶点组成的。
看一下FENCE_LOC里面的数据（一块区域的四个顶点），存储的是一列字符串信息，我们要把字符串里面的顶点提取出来。

data_fence['FENCE_LOC'].iloc[0]

‘[118.103198,24.527338],[118.103224,24.527373],[118.103236,24.527366],[118.103209,24.527331],[118.103198,24.527338]’

怎么提取顶点呢？用python里面的exec()命令,就可以把字符串里面的点提取出来（以第一行数据为例）：

exec('points = ['  + data_fence['FENCE_LOC'].iloc[0]+']')
points

[[118.103198, 24.527338],
[118.103224, 24.527373],
[118.103236, 24.527366],
[118.103209, 24.527331],
[118.103198, 24.527338]]

对每一行遍历进行上面的操作，对于得到的这些点，用Polygon()转化成面要素，把每一个电子栅栏的Polygon存储到一个空列表里面，操作如下：

geometry=[]
for i in range(len(data_fence)):exec('points = ['  + data_fence['FENCE_LOC'].iloc[i]+']')geometry.append(Polygon(points))geometry

现在的geometry这一列里面都是Polygon要素。

接着，dataframe转换成Geodataframe：加一列geometry,再gpd.GeoDataFrame（)转换格式。

data_fence['geometry']=geometry
data_fence=gpd.GeoDataFrame(data_fence)
data_fence=data_fence.drop('FENCE_LOC',axis=1)

现在的电子栅栏信息数据如下：

三、共享单车与电子栅栏匹配

3.1 空间连接——识别车是在栅栏内，还是在栅栏外

接下来，对GPS停车数据和停车栅栏数据做空间连接（确定每辆车是否停在停车栅栏里面），要用到一个sjoin函数，和pd.merge函数功能很像，只不过前者是对地理要素进行merge。

parking_points=gpd.sjoin(parking_points,data_fence,how='left')

将单车停放数据和电子栅栏数据连接在一起，以单车停放数据为基准，如果有些单车不在任何栅栏之中，和它连接起来的FENCE_ID会是空值。
可以看到表里面index_right、FENCE_ID这两列有好多的空值，说明8：00这时期停放的共享单车，许多都没有停放在规定停车点内部。

接下来分别提取出8：00时，停在栅栏里面的单车，和停在栅栏之外的单车。分别存储到out_fence和in_fence里面。

#停在栅栏外部的车
out_fence=parking_points[parking_points['FENCE_ID'].isnull()]
out_fence

#停在栅栏内部的车
in_fence=parking_points[-parking_points['FENCE_ID'].isnull()]
in_fence

3.2 最邻近电子栅栏寻优搜索

对于栅栏外部的单车，我们希望找到离它们距离最近的电子栅栏，并计算距离大小，如果某个区域停放在栅栏外部的单车，它们离最近电子栅栏的距离非常远，就说明这个区域需要合理增设一些电子栅栏，来缓解停放需求量。

3.2.1 寻找离栅栏外单车最近的电子栅栏

对于每一辆在栅栏外的单车，寻找离它最近的电子栅栏，基于计算单车与电子栅栏之间的距离，然后比较距离大小。但是当单车数据和电子栅栏数据很多时，两两计算距离这个计算量实在是太大了。假设有一万量单车，一万个停放点，每次匹配要计算一亿次！怎样快速寻优匹配呢？

这里用到的是Kd-tree 算法。学习过机器学习中KNN 算法的小伙伴一定不会陌生的，KNN算法其实就是基于Kd-tree的一种算法。Kd-tree是一种高维索引树形数据结构，常用于在大规模的高维数据空间进行最近邻查找。Kd-Tree是一棵二叉树，树中存储的是一些K维数据。在一个K维数据集合上构建一棵Kd-Tree代表了对该K维数据集合构成的K维空间进行的一个划分，即树中的每个结点就对应了一个K维的超平面。
暴力搜索（线性搜索）寻优的时间复杂度是O（N），但是Kd-tree的时间复杂度是O（logN），在算法性能上是一个很大的提升。

3.2.2.1 Kd-tree原理

其实Kd-tree相当于BST（二叉查找树 Binary Search Tree）在高维空间上的推广，而这大体步骤都差不多，都是先构造查找树，再进行数据查询。先看二叉查找树的原理，它适合对一维数据进行搜索寻优。

构造查找树

根据BST的性质来创建二叉树，将数据点一个一个插入到BST树中，插入后的树仍然是BST树，即根结点的左子树中所有结点的值均小于根结点的值，而根结点的右子树中所有结点的值均大于根结点的值。每一个节点，其实都是对二维平面的一个划分，如下图所示：

数据查询
构建好了查找树之后，只需要将查询数据与结点值进行比较然后选择对应的子树继续往下查找。

假设这里有一个一维数组A=[3,0,4,6,11,8,7]。x=2，要找在A中最接近x的元素。

step1：把数组A进行升序排列，得到[0, 3, 4, 6, 7, 8, 11]
step2：数组中间的元素是6，2小于6，所以2只可能存在于6的左边，我们只需要在数组[0, 3, 4]中继续查找
step3：左边的数组中间的元素是3，2小于3，所以2只可能存在于3的左边，即数组[0]

构建Kd-tree
Kd-tree可以看成是二叉查找树推广到高维空间上，但是也有一些不同，比如：

二叉查找树：数据存放在树中的每个结点（根结点、中间结点、叶子结点）中。
Kd-Tree：数据只存放在叶子结点，而根结点和中间结点存放一些空间划分信息（例如划分维度、划分值）；而且Kd-tree是选择方差最大的特征作为分割特征。

构建Kd-tree的步骤：

1. 建立根节点；2. 选取方差最大的特征作为分割特征；3. 选择该特征的中位数作为分割点；4. 将数据集中该特征小于中位数的传递给根节点的左儿子，大于中位数的传递给根节点的右儿子；
5. 递归执行步骤2-4，直到所有数据都被建立到KD Tree的节点上为止。

3.2.2.2 Kd-tree应用在单车与栅栏匹配

利用电子栅栏信息来构建Kd-tree,然后用构建好的Kd-tree对单车停放点GPS数据进行查询。这个过程可以用scipy包当中的cKDTree函数，直接可以构建好树，并对停车点GPS数据进行查询，而且返回离查询数据最近的点和它们之间的距离。

    #为表A的点建立KDTreebtree = cKDTree(A)#在表A的KDTree中查询表B的点,dist为距离,idx为A中离B最近的坐标点dist, idx = btree.query(B, k=1)

3.2.2 最近的电子栅栏离单车距离多远？

编写了两个函数来对geodataframe中数据进行寻优计算。
1.ckdnearest_point(gdA, gdB)：输入两个geodataframe，gdA、gdB均为点要素，用gdB来进行Kd-tree的构建，然后对gdA中的每一个元素，在Kd-tree中找到离它最近的元素，并返回它们之间的距离。该函数会为gdfA表连接上gdfB中最近的点，并添加距离字段dist，返回含有最近邻点的ID字段与距离字段的gdf。

2.ckdnearest_line(gdfA, gdfB)：输入两个geodataframe，其中gdfA为点，gdfB为线，用gdB来进行Kd-tree的构建，然后对gdA中的每一个元素，在Kd-tree中找到离它最近的元素，并返回它们之间的距离。该函数会为gdfA表连接上gdfB中最近的线，并添加距离字段dist，返回含有最近邻点的ID字段与距离字段的gdf。

直接附上代码，可以把这两段代码保存成.py文件，用的时候直接调用。
注释已经写在代码里面了。

在该项目中，我们要用第二个函数ckdnearest_line(gdfA, gdfB)，因为计算单车停车点到电子栅栏的最近距离，看得是停车点到栅栏边上的距离。

首先导入要用到的这些库。

import numpy as np
from scipy.spatial import cKDTree
import itertools
from operator import itemgetter

def ckdnearest_point(gdA, gdB):'''输入两个geodataframe，gdfA、gdfB均为点，该方法会为gdfA表连接上gdfB中最近的点，并添加距离字段dsit'''#提取gdA中的所有点要素nA = np.array(list(gdA.geometry.apply(lambda x: (x.x, x.y))))#提取gdB中的所有点要素nB = np.array(list(gdB.geometry.apply(lambda x: (x.x, x.y))))#为gdB表的点建立KDTreebtree = cKDTree(nB)#在gdB的KDTree中查询gdA的点,dist为距离,idx为gdB中离gdA最近的坐标点dist, idx = btree.query(nA, k=1)#构建匹配的结果gdf = pd.concat([gdA.reset_index(drop=True), gdB.loc[idx, gdB.columns != 'geometry'].reset_index(drop=True),pd.Series(dist, name='dist')], axis=1)return gdf

def ckdnearest_line(gdfA, gdfB):'''输入两个geodataframe，其中gdfA为点，gdfB为线，该方法会为gdfA表连接上gdfB中最近的线，并添加距离字段dsit'''#提取gdA中的所有点要素A = np.concatenate([np.array(geom.coords) for geom in gdfA.geometry.to_list()])#把gdfB的几何坐标提取到B，此时B为一个大list中包含多个小list，每个小list代表一个几何图形，小list中为坐标#B=[[[要素1坐标1],[要素1坐标2],...],[[要素2坐标1],[要素2坐标2],...]]B = [np.array(geom.coords) for geom in gdfB.geometry.to_list()]#B_ix代表B中的每个坐标点分别属于B中的哪个几何图形B_ix = tuple(itertools.chain.from_iterable([itertools.repeat(i, x) for i, x in enumerate(list(map(len, B)))]))#把B表展开，B=[[要素1坐标1],[要素1坐标2],...,[要素2坐标2],[要素2坐标2],...]B = np.concatenate(B)#为B表建立KDTreeckd_tree = cKDTree(B)#在B的KDTree中查询A的点,dist为距离,idx为B中离A最近的坐标点dist, idx = ckd_tree.query(A, k=1)#由坐标点对应到几何要素idx = itemgetter(*idx)(B_ix)#构建匹配的结果gdf = pd.concat([gdfA.reset_index(drop=True), gdfB.loc[idx, gdfB.columns != 'geometry'].reset_index(drop=True),pd.Series(dist, name='dist').reset_index(drop=True)], axis=1)return gdf

接下来我们调用ckdnearest_line(gdfA, gdfB)来为那些停放在停车栅栏之外的单车，寻找距离它们最近的栅栏边框，以及二者之间的距离。但是在调用函数之前，还要对停车栅栏数据中的Polygon要素提取边框。毕竟我们的函数是计算点线之间的距离，所以要提取Polygon的边框。
先来可视化看一下停车栅栏的Polygon以及边框：
提取前十个停车栅栏的Polygon看一下：

data_fence[:10].plot()

一条一条的带状停车栅栏：
提取栅栏边框（data_fence.boundary）可视化一下：

data_fence[:10].boundary.plot()

开始计算：

data_fence_boundary=data_fence.copy()
#copy栅栏信息做一个备份，把备份数据中的geometry设置为栅栏的边框数据
data_fence_boundary['geometry']=data_fence.boundary
#停车栅栏外面的车，每一辆都匹配最近的停车栅栏
out_fence=ckdnearest_line(out_fence,data_fence_boundary)

看寻优结果：停在栅栏外部的单车，算法都给他们匹配上了离得最近的停车栅栏，计算出了二者之间的距离dist。
发现这份数据里面有两个FENCE_ID，要重新设置一下Columns，只保留我们需要的字段，然后把在栅栏外数据与栅栏内数据拼接在一起（用concat函数）：

out_fence.columns=['BIKE_ID', 'LATITUDE', 'LONGTITUDE', 'DATA_TIME', 'DATA_TIME1','geometry', 'index_right', 'FENCE_ID1', 'FENCE_ID', 'dist']
in_fence=in_fence[['BIKE_ID', 'LATITUDE', 'LONGTITUDE', 'DATA_TIME', 'DATA_TIME1','geometry','FENCE_ID' ]]
out_fence=out_fence[['BIKE_ID', 'LATITUDE', 'LONGTITUDE', 'DATA_TIME', 'DATA_TIME1','geometry','FENCE_ID', 'dist' ]]
parking_points=pd.concat([in_fence,out_fence])

最终的数据就整理好啦！结果如下：
dist为0的记录是那些停放在栅栏里面的单车，dist不为0的记录是那些停放在栅栏外面的单车。

3.3 matplotlib可视化栅栏与单车停放分布

以停车栅栏’蔡岭路 _R_3’为例，想看看这个停车栅栏附近的单车停车分布如何，对其进行可视化。

fig=plt.figure(1,(10,10),dpi=200)
ax=plt.subplot(111)
fence_id='蔡岭路 _R_3'
parking_points[parking_points['FENCE_ID']==fence_id].plot(ax=ax)
data_fence[data_fence['FENCE_ID']==fence_id].plot(ax=ax,color='g',alpha=0.5)
plt.show()

结果如下：

四、用Kepler进行空间可视化

kepler.gl是由Uber的可视化团队创建的一个基于web的高性能工具，用于大规模地理空间数据集的可视化探索。在我们拿到地理空间数据的Geojson文件时，就可以放进Kepler进行探索分析，相当好用酷炫的一个工具。
下面是对上面处理的单车停放需求量与停车栅栏数据的一个可视化成果展示。
反转一下颜色，使得离停车栅栏越远的单车GPS点颜色越深。点击颜色条，打开Reversed.
也可以看立体图，高度代表该栅格里面的车辆数目，颜色深浅依旧是单车离停车栅栏的距离。

另外一种可视化方式：使用热力图，在地图上离停车栅栏越远的单车点，颜色越深。对应点图里面的Color based,热力图中颜色深浅映射是设置weight这个参数，把weight设置为dist,同样地点击颜色条，进行颜色反转。颜色越深的地带中的单车离停车栅栏越远，这些地带有很强的停车需求，说明这里需要多增设一些停车栅栏。
好，以上是停车需求的图层，接下来我们把停车栅栏的Geojson文件数据也加入进来，新添加一个停车栅栏的图层。

点击Add Data，添加电子停车栅栏的json文件，加载进来之后就成了这样，密密麻麻的黄色网状就是电子停车栅栏：
放大看一下：
可以添加筛选条件，比如我们以DATE_TIME作为筛选条件，可以动态地看到不同时间段停车需求量的变化：

可以使用增量移动窗口，如下面动态图，看随着时间增加，总的停车需求量变化。很明显，从早晨6：00到8：00，随着停车需求量不断地增加，有些地区的停车需求量越发增大。
也可以使用滑动时间窗口，来看不同的时段内，停车需求量的变化。

地理空间数据项目实战——共享单车停车需求识别与地图匹配相关推荐

Python实战案例：我们对共享单车的需求有多大？
现如今,共享单车在生活中可谓处处可见,那么它的租赁需求是多少呢?今天我们就基于美国华盛顿共享单车的租赁数据,利用Python和可视化对租赁需求进行预测. 01 数据来源及背景该数据集是美国华盛顿共享 ...
python共享单车案例分析_python分析数据分析项目：共享单车租用情况影响因素探索分析...
python分析数据分析项目:共享单车租用情况影响因素探索分析
随机森林分析共享单车的需求
随机森林相关API: import sklearn.ensemble as se # 随机森林回归模型 (属于集合算法的一种) # max_depth:决策树最大深度10 # n_estimators ...
自动驾驶系统进阶与项目实战（四）自动驾驶高精度地图构建中的三维场景识别和闭环检测
自动驾驶系统进阶与项目实战(四)自动驾驶高精度地图构建中的三维场景识别和闭环检测闭环检测(loop-closure detection)是SLAM中非常关键的一部分,当然也是自动驾驶高精度地图(HD ...
项目实战——文档扫描OCR识别
扫描全能王的实现,maybe 目录一.文档扫描 1.引入所需要的库 2.图像的读取与预处理读取图像图像reszie, 图像灰度化.滤波.边缘检测. 3.轮廓检测 4.透视与二值变换二.文字识别 ...
基于openCV的项目实战1：信用卡数字识别
目的:识别信用卡卡号: 方法:基于模板匹配: 模板匹配思想:拿4和左边模板进行一一匹配,算一下平方项的差异,恰好4和左边模板中的4差异最小,所以就知道当前的数字是4: 具体步骤: 第一步:找到与你当前 ...
计算机视觉项目实战-图像特征检测harris、sift、特征匹配
【实战】OpenCV+Python项目实战--文档扫描OCR识别
文章目录 1 准备工作(python) 1.1 np.diff用法 1.2 tesseract和pytesseract安装 2代码实现 2.1 文档提取与摆正 2.2 OCR扫描 1 准备工作(pyt ...
摩拜共享单车数据分析项目报告
文章目录项目背景数据探索数据挖掘数据分析时间维度空间维度用户维度项目背景随着智能手机的普及和手机用户的激增,共享单车作为城市交通系统的一个重要组成部分,以绿色环保.便捷高效.经济环保 ...

地理空间数据项目实战——共享单车停车需求识别与地图匹配