基于python的transbigdata实现出租车轨迹数据分析与可视化

TransBigData是一个为交通时空大数据处理、分析和可视化而开发的Python包。TransBigData为处理常见的交通时空大数据（如出租车GPS数据、共享单车数据和公交车GPS数据等）提供了快速而简洁的方法。TransBigData为交通时空大数据分析的各个阶段提供了多种处理方法,代码简洁、高效、灵活、易用，可以用简洁的代码实现复杂的数据任务。转自同济大学余庆博士的文章：点击跳转

应专业选修课《数据分析与可视化》的期末考核，本人参考了@小旭学长的文章：点击跳转 

本文根据上述文章，对该项目进行了一些代码上的修改（在一些可能报错的地方进行了修改），以及增加了对代码的解释以及图文介绍，并且增加了上海的出租车轨迹数据。

本项目作为课程设计参考，个人感觉本项目在最后的呈现效果上十分出色。该项目使用jupyter notebook开发环境。

项目地址：https://gitee.com/scaukzh/visualization

话不多说，直接看项目。

数据：

①出租车轨迹数据集

来源：https://people.cs.rutgers.edu/~dz220/data.html

介绍：该数据源自国外大学一名助理教授，包含了深圳市2013.10.22一天中大约600多辆出租车的轨迹数据，数据包含的信息依次为出租车ID，时间，纬度，经度，占用状态，速度;入住状态：1-乘客和0-乘客。数据一共大约54w条。

	VehicleNum	Time	Lng	Lat	OpenStatus	Speed
0	34745	20:27:43	113.806847	22.623249	1	27
1	34745	20:24:07	113.809898	22.627399	0	0
2	34745	20:24:27	113.809898	22.627399	0	0
3	34745	20:22:07	113.811348	22.628067	0	0
4	34745	20:10:06	113.819885	22.647800	0	54
...	...	...	...	...	...	...
544994	28265	21:35:13	114.321503	22.709499	0	18
544995	28265	09:08:02	114.322701	22.681700	0	0
544996	28265	09:14:31	114.336700	22.690100	0	0
544997	28265	21:19:12	114.352600	22.728399	0	0
544998	28265	19:08:06	114.137703	22.621700	0	0

544999 rows × 6 columns

补充：我在项目中补充了一个上海的数据集，该数据集是香港科技大学智慧城市研究小组共享的数据集，包含上海市2007年02月20日单日若干辆出租车的轨迹数据。

②地图数据

地图采用json格式，这个地图数据用于标明深圳区域的范围，在后期进行绘图以及数据处理的时候都有作用。

代码部分：

一、第三方包和数据集的导入

1.第三包的导入

1. import transbigdata as tbd
2. import pandas as pd
3. import geopandas as gpd
4. import matplotlib.pyplot as plt

使用pip命令直接装transbigdata和geopands包可能会失败。不过transbigdata我已经放在项目中去了，大家直接调用即可。geopands比较容易装，这里就不详细讨论了。

2.读入出租车GPS轨迹数据

1. data = pd.read_csv('D:/code_for_school/visualtaxi/gps_data/shenzhen_taxi_gps.csv', header=None)#修改自身路径即可
2. data.columns = ['VehicleNum', 'Time', 'Lng', 'Lat', 'OpenStatus', 'Speed']
3. data.head() #head函数其中一个默认参数为5，故只返回头5条数据

	VehicleNum	Time	Lng	Lat	OpenStatus	Speed
0	34745	20:27:43	113.806847	22.623249	1	27
1	34745	20:24:07	113.809898	22.627399	0	0
2	34745	20:24:27	113.809898	22.627399	0	0
3	34745	20:22:07	113.811348	22.628067	0	0
4	34745	20:10:06	113.819885	22.647800	0	54

3.读入深圳市地图数据

1. # Read the GeoDataFrame of the study area
2. sz = gpd.read_file(r'D:/code_for_school/visualtaxi/geo_data/sz.json')#修改自身路径即可
3. sz.crs = None
4. sz.head()

	centroid_x	centroid_y	qh	geometry
0	114.143157	22.577605	罗湖	POLYGON ((114.10006 22.53431, 114.10083 22.534...
1	114.041535	22.546180	福田	POLYGON ((113.98578 22.51348, 114.00553 22.513...
2	114.270206	22.596432	盐田	POLYGON ((114.19799 22.55673, 114.19817 22.556...
3	113.851387	22.679120	宝安	MULTIPOLYGON (((113.81831 22.54676, 113.81948 ...
4	113.926290	22.766157	光明	POLYGON ((113.99768 22.76643, 113.99704 22.766...

sz.plot()

二、数据处理

1.剔除掉超出深圳地图范围的数据

1. # Data Preprocessing
2. # Delete the data outside of the study area
3. data = tbd.clean_outofshape(data, sz, col=['Lng', 'Lat'], accuracy=500) #剔除超出研究区域的数据
4. data.head(20)

	VehicleNum	Time	Lng	Lat	OpenStatus	Speed
0	34745	20:27:43	113.806847	22.623249	1	27
1	27368	09:08:53	113.805893	22.624996	0	49
2	22998	10:51:10	113.806931	22.624166	1	54
3	22998	10:11:50	113.805946	22.625433	0	43
4	22998	10:12:05	113.806381	22.623833	0	60
5	22396	12:14:28	113.806602	22.623949	1	45
6	32436	08:15:34	113.806465	22.624166	0	50
7	32436	08:15:04	113.807053	22.626350	0	0
8	32436	08:15:19	113.806252	22.625933	0	22
9	32436	08:14:49	113.807053	22.626333	0	0
10	32436	18:10:57	113.806885	22.622967	0	50
11	32436	18:10:27	113.806450	22.626083	0	36
12	32436	18:35:15	113.806686	22.623632	1	29
13	32436	18:10:42	113.805984	22.624683	0	41
14	27373	18:58:43	113.805969	22.624983	0	42
15	27373	08:42:34	113.807114	22.623211	0	61
16	31910	09:49:56	113.806801	22.626200	0	25
17	35161	00:43:41	113.807030	22.622833	1	51
18	35161	00:43:25	113.807098	22.624184	1	54
19	23411	20:20:06	113.806549	22.623550	0	56

如此一来，那些不在深圳市版图中的数据就被我们剔除掉了。

2.对数据按照车辆编号和时间顺序排好

1. # Delete the data with instantaneous changes in passenger status
2. data = tbd.clean_taxi_status(data, col=['VehicleNum', 'Time', 'OpenStatus'])
3. #删除出租车数据中载客状态瞬间变化的记录，这些记录的存在会影响出行订单判断。
4. #判断条件为:如果对同一辆车，上一条记录与下一条记录的载客状态都与本条记录不同，则本条记录应该删去
5. #实际上就是按照编号和时间顺序将数据排好
6. data.head(30)

	VehicleNum	Time	Lng	Lat	OpenStatus	Speed	LONCOL	LATCOL
452072	22396	00:00:29	113.996719	22.693333	1	20	81	65
444077	22396	00:01:01	113.995514	22.695032	1	34	81	66
444078	22396	00:01:09	113.995430	22.695766	1	41	81	66
444075	22396	00:01:41	113.995369	22.696484	1	0	81	66
444079	22396	00:02:21	113.995430	22.696650	1	17	81	66
444073	22396	00:03:01	113.994934	22.697884	0	23	81	66
444074	22396	00:03:41	113.992683	22.697350	0	18	81	66
444076	22396	00:04:21	113.992348	22.696733	0	36	81	66
452073	22396	00:05:01	113.992996	22.693632	0	2	81	65
446704	22396	00:05:41	113.989250	22.693083	0	56	80	65
422019	22396	00:06:21	113.986366	22.691000	0	48	79	65
422218	22396	00:07:41	113.988747	22.685450	0	16	80	63
426967	22396	00:08:21	113.989586	22.681749	0	43	80	63
443922	22396	00:09:01	113.986015	22.681583	0	0	79	63
443923	22396	00:09:41	113.986481	22.680332	0	47	79	62
422292	22396	00:10:21	113.988297	22.676416	0	31	80	61
422295	22396	00:11:01	113.991600	22.677517	0	41	80	62
422291	22396	00:11:09	113.992149	22.677317	0	33	80	62
443991	22396	00:11:41	113.993401	22.674368	0	31	81	61
443987	22396	00:11:49	113.993446	22.674217	0	0	81	61

可以看到编号相同的数据是相邻并且是以时间推进排序的。

3.将数据按照经纬度进行栅格化

栅格化的意思就是把每个数据按照经纬度分到不同位置的格子上

1. # Data gridding
2. # Define the bounds and generate gridding parameters
3. bounds = [113.6, 22.4, 114.8, 22.9]
4. params = tbd.area_to_params(bounds, accuracy=500)
5. params
6. #(113.6, 22.4, 0.004872390756896538, 0.004496605206422906)
7. #栅格参数(lonStart,latStart,deltaLon,deltaLat)，分别为栅格左下角坐标与单个栅格的经纬度长宽
8. # Mapping GPS data to grids
9. data['LONCOL'], data['LATCOL'] = tbd.GPS_to_grid(data['Lng'], data['Lat'], params)
10. data.head()

	VehicleNum	Time	Lng	Lat	OpenStatus	Speed	LONCOL	LATCOL
452072	22396	00:00:29	113.996719	22.693333	1	20	81	65
444077	22396	00:01:01	113.995514	22.695032	1	34	81	66
444078	22396	00:01:09	113.995430	22.695766	1	41	81	66
444075	22396	00:01:41	113.995369	22.696484	1	0	81	66
444079	22396	00:02:21	113.995430	22.696650	1	17	81	66

这一部分完成，我们就已经将数据划分在不同的栈格里了。

4.计算每个栅格中的数据个数

1. # Aggregate data into grids
2. #python中groupby函数主要的作用是进行数据的分组以及分组后地组内运算！
3. datatest = data.groupby(['LONCOL', 'LATCOL'])['VehicleNum'].count().reset_index()
4. # Generate the geometry for grids
5. datatest['geometry'] = tbd.grid_to_polygon([datatest['LONCOL'], datatest['LATCOL']], params)
7. # Change it into GeoDataFrame
8. # import geopandas as gpd
9. datatest = gpd.GeoDataFrame(datatest)
10. datatest.head()
11. #计算出每个栅格里的出租车数量（不按时间，是一天数据中所有的点）

	LONCOL	LATCOL	VehicleNum	geometry
0	36	63	3	POLYGON ((113.77297 22.68104, 113.77784 22.681...
1	36	64	2	POLYGON ((113.77297 22.68553, 113.77784 22.685...
2	36	65	1	POLYGON ((113.77297 22.69003, 113.77784 22.690...
3	36	66	1	POLYGON ((113.77297 22.69453, 113.77784 22.694...
4	36	67	8	POLYGON ((113.77297 22.69902, 113.77784 22.699...

LONCOL代表经度上第几个栈格，LATCOL代表维度上第几个栈格，VehicleNum表示了这个栈格中数据的数量，可用于后面进行GPS轨迹热力图的绘图。

5.OD订单的处理

①将各个数据转化成OD数据（每一个数据代表状态）

载客状态无非就两种状态，即（载客-空载-载客）和（空载-载客-空载）。

1. # Extract taxi OD from GPS data
2. #OD形成就是每一个订单都算一个OD行程，每一行代表一个OD行程
3. #stime代表出发时间，slon和slat是出发的经纬度
4. #etime代表结束时间，elon和elat是结束的经纬度
5. # ID是订单顺序号
6. oddata = tbd.taxigps_to_od(data,col = ['VehicleNum', 'Time', 'Lng', 'Lat', 'OpenStatus'])
7. oddata

	VehicleNum	stime	slon	slat	etime	elon	elat	ID
427075	22396	00:19:41	114.013016	22.664818	00:23:01	114.021400	22.663918	0
131301	22396	00:41:51	114.021767	22.640200	00:43:44	114.026070	22.640266	1
417417	22396	00:45:44	114.028099	22.645082	00:47:44	114.030380	22.650017	2
376160	22396	01:08:26	114.034897	22.616301	01:16:34	114.035614	22.646717	3
21768	22396	01:26:06	114.046021	22.641251	01:34:48	114.066048	22.636183	4
...	...	...	...	...	...	...	...	...
175519	36805	22:49:12	114.114365	22.550632	22:50:40	114.115501	22.557983	5083
212092	36805	22:52:07	114.115402	22.558083	23:03:12	114.118484	22.547867	5084
119041	36805	23:03:45	114.118484	22.547867	23:20:09	114.133286	22.617750	5085
224103	36805	23:36:19	114.112968	22.549601	23:43:12	114.089485	22.538918	5086
170962	36805	23:46:14	114.091217	22.540768	23:53:36	114.120354	22.544300	5087

5088 rows × 8 columns

②将od行程分为载客状态下和空载状态下的

1. data_deliver, data_idle = tbd.taxigps_traj_point(data,oddata,col=['VehicleNum',
2. 'Time',
3. 'Lng',
4. 'Lat',
5. 'OpenStatus'])

data_deliver #载客状态的od路径

	VehicleNum	Time	Lng	Lat	OpenStatus	Speed	LONCOL	LATCOL	ID	flag
427075	22396	00:19:41	114.013016	22.664818	1	63.0	85.0	59.0	0.0	1.0
427085	22396	00:19:49	114.014030	22.665483	1	55.0	85.0	59.0	0.0	1.0
416622	22396	00:21:01	114.018898	22.662500	1	1.0	86.0	58.0	0.0	1.0
427480	22396	00:21:41	114.019348	22.662300	1	7.0	86.0	58.0	0.0	1.0
416623	22396	00:22:21	114.020615	22.663366	1	0.0	86.0	59.0	0.0	1.0
...	...	...	...	...	...	...	...	...	...	...
64411	36805	23:53:09	114.120354	22.544300	1	2.0	107.0	32.0	5087.0	1.0
64405	36805	23:53:15	114.120354	22.544300	1	1.0	107.0	32.0	5087.0	1.0
64390	36805	23:53:21	114.120354	22.544300	1	0.0	107.0	32.0	5087.0	1.0
64406	36805	23:53:27	114.120354	22.544300	1	0.0	107.0	32.0	5087.0	1.0
64393	36805	23:53:33	114.120354	22.544300	1	0.0	107.0	32.0	5087.0	1.0

209250 rows × 10 columns

data_idle #空载状态的od路径

	VehicleNum	Time	Lng	Lat	OpenStatus	Speed	LONCOL	LATCOL	ID	flag
416628	22396	00:23:01	114.021400	22.663918	0	25.0	86.0	59.0	0.0	0.0
396136	22396	00:23:41	114.023148	22.665100	0	21.0	87.0	59.0	0.0	0.0
396129	22396	00:24:21	114.024414	22.665367	0	35.0	87.0	59.0	0.0	0.0
401744	22396	00:25:01	114.027115	22.662100	0	25.0	88.0	58.0	0.0	0.0
394630	22396	00:25:41	114.024551	22.659834	0	21.0	87.0	58.0	0.0	0.0
...	...	...	...	...	...	...	...	...	...	...
170971	36805	23:45:42	114.091217	22.540768	0	0.0	101.0	31.0	5086.0	0.0
170954	36805	23:45:57	114.091217	22.540768	0	0.0	101.0	31.0	5086.0	0.0
170955	36805	23:46:12	114.091217	22.540768	0	0.0	101.0	31.0	5086.0	0.0
64391	36805	23:53:36	114.120354	22.544300	0	0.0	107.0	32.0	5087.0	0.0
64396	36805	23:53:51	114.120354	22.544300	0	0.0	107.0	32.0	5087.0	0.0

332941 rows × 10 columns

可以看出载客状态下数据有20w条，而空载状态下数据有33w条，说明在一天的过程中，总体上处于空载状态的时间是大于载客状态的时间的。
 

三、绘图以及可视化展示

1.GPS轨迹的热力图绘制

1. # Plot the grids
2. fig = plt.figure(1, (16, 6), dpi=300)
3. ax1 = plt.subplot(111)
4. #ax用于调整大小
5. # tbd.plot_map(plt, bounds, zoom=10, style=4)
6. datatest.plot(ax=ax1, column='VehicleNum', legend=True)
7. plt.xticks([], fontsize=10)
8. plt.yticks([], fontsize=10)
9. plt.title('Taxi GPS track point count', fontsize=12);

1. # Plot the grids
2. fig = plt.figure(1, (16, 6), dpi=300) # 确定图形高为6，宽为8；图形清晰度
3. ax1 = plt.subplot(111)
4. datatest.plot(ax=ax1, column='VehicleNum', legend=True, scheme='quantiles')
5. # plt.legend(fontsize=10)
6. plt.xticks([], fontsize=10)
7. plt.yticks([], fontsize=10)
8. plt.title('Taxi GPS track point count', fontsize=12);

1. # Plot the grids
2. fig = plt.figure(1, (16, 6), dpi=150) # 确定图形高为6，宽为8；图形清晰度
3. ax1 = plt.subplot(111)
4. datatest.plot(ax=ax1, column='VehicleNum', legend=True, cmap='OrRd', scheme='quantiles')
5. # plt.legend(fontsize=10)
6. plt.xticks([], fontsize=10)
7. plt.yticks([], fontsize=10)
8. plt.title('Taxi GPS track point count', fontsize=12);

 2.载客状态和空载状态下的od订单图

1. traj_deliver = tbd.points_to_traj(data_deliver)
2. traj_deliver.plot()

1. traj_idle = tbd.points_to_traj(data_idle)
2. traj_idle.plot()

 3.依赖keplergl第三方包绘制可交互式的可视化动画

keplergl是由Uber开源的一款地理数据可视化工具，我们可以在Jupyter notebook中使用该工具让我们的数据可视化做到美观。使用前需要先安装，安装教程可参考该篇文章：点击跳转

keplergl官网：点击跳转

①可交互式热力图

1. # 可视化数据点分布
2. tbd.visualization_data(data,col = ['Lng','Lat'],accuracy=40,height = 500) #accuracy是栅格大小，数值越大表示一个栅格涵盖地图面积越大

②可交互式OD订单图

1. # 可视化数据点分布
2. tbd.visualization_od(oddata,accuracy=2000,height = 500)

 注：同时使用Ctrl和鼠标左键可移动观看视角

③可交互式的视频动图

tbd.visualization_trip(data_deliver) #一天当中的所有载客轨迹

 总结

大家按照上述步骤进行项目，大问题是没有的，之前我遇到的主要问题就是transbigdata的包找了半天，其他的一些如各个包版本之间不适配的问题，可以直接百度解决。本项目的主要难点就是要去理解数据处理的一些函数，这些函数的大概功能我已经在注释中和内容中阐述了。如果有小伙伴想更加深入的了解函数的原理，可以去看transbigdata的源代码，这样可以更加直观的了解这些函数。

补充

如有其他问题，请私信我，谢谢。