距离上次更新已经过去了一个半月之久，通过与各位读者朋友交流，发现有不少地理和gis的朋友关注我的公众号，可能是之前写的文章多与gis有关
这次回归本行，写一篇关于交通的文章，欢迎大家后台私信我与我讨论，尤其是针对技术及idea的讨论，十分欢迎！同时也希望大家在直接开口要数据前有一些自己的思考，毕竟与最终的数据相比，分析的过程与思路才是最重要的。下面开始正文。
这个数据是2015SODA大赛公开的上海公交公司的一卡通数据集，具体的介绍和获取方法网上应该有很多（因此原始数据我不提供，源代码都在文章里，复制粘贴即可），简单的看一下，包括卡的id，线路站点，费用，优惠，刷卡时间几个字段（hour是我后面自己加的）。根据常识，我们进出地铁站要打两次卡，进站不要钱（cost==0），出站时收费，因此我们可以根据这个规则把一个人的进出站的刷卡记录对应起来，找到出行的od站点。用下面这样一行代码，我们对用户和时间进行排序，看看基本情况：
代码语言：javascript复制df['timestamp'] = pd.to_datetime(df['timestamp'])
df.sort_values(['id','timestamp'])
可以看到，id为1的用户，一天出行了一次（从大连路到书院，花了9块），id为4200000172的用户，这天出行了两次（张江高科-人民广场、人民广场-张江高科）。理想的情况应该是，一个人的打卡记录是偶数次，并且一次cost为0（进站），一次cost不为0（出站）。然而，理想很丰满显示很骨感，通过下面这行代码：
代码语言：javascript复制df['id'].value_counts()[df['id'].value_counts()%2 ==1]
我们发现有很多人的打卡次数是奇数次，这可能包括了在前一天开始在今天结束的行程、在今天开始明天结束的行程、和一些可能的没有进站或出站的记录，比如：
代码语言：javascript复制df[df['id']==2102265408]
这个老哥4.2号第一次打卡就是出站（第一列），以及在下午4点多来了一次霸王单（非优惠并且cost为0），对于红线这种数据，都是我们需要清洗的（为了方便清洗规则，这里把霸王单也清洗了）。
还比如这个：
我估计是地铁的员工进出站点，都不要钱。
所以我们要的就是那种上车刷卡cost==0，下车刷卡cost！=0的、并且同一个id，且上下车的刷卡时间挨着的数据，比如这种：
然后把上车和下车合并成一行，就是一个人一次地铁出行的信息。
具体怎么操作的话，最开始我写了一个傻瓜版循环：
代码语言：javascript复制#def get_trip(df):
# df['index'] = list(df.index)
# df['trip'] = -1
# cardholder = df['id'].unique()
# print('共有{}名用户'.format(len(cardholder)))
# trip = 0
# for x,i in enumerate(cardholder):
# if (x000) == 0:
# print('正在处理第{}个用户的数据'.format(x))
#
# df_sub = df[df['id']==i]
# df_sub = df_sub.sort_values('timestamp')
# df_sub.reset_index(inplace=True,drop=True)
# for i in range(len(df_sub)-1):
# if (df_sub.loc[i,'cost'] == 0) & (df_sub.loc[i 1,'cost'] > 0):
# df.loc[df_sub.loc[i,'index'],'trip'] = trip
# df.loc[df_sub.loc[i 1,'index'],'trip'] = trip
# trip = trip 1
# else:
# continue
# df.drop(columns='index',inplace=True)
# return df总的来说，就是循环提取每个人这一天的出行信息，然后进行筛选，为正常数据赋予trip编号（每次出行上下车的trip编号相同），并把脏数据trip的字段为-1。然而，由于用了双循环（python里for循环的速度你懂的
），程序跑起来十分地慢，900万条数据跑完大概需要5个多小时。。。。。。这样肯定是不行的，于是改写了一下代码，增加了几个列用来做关键的判断（前一行后一行的id和cost），利用pandas的apply函数，具体如下：
代码语言：javascript复制### 增加用来判断的列
def get_shift(df):
df['timestamp'] = pd.to_datetime(df['timestamp'])
df = df.sort_values(['id','timestamp'])
df['id_shift_after'] = df['id'].shift(-1)
df['cost_shift_before'] = df['cost'].shift(1)
df['cost_shift_after'] = df['cost'].shift(-1)
df['id_shift_before'] = df['id'].shift(1)
return df### 用来判断数据是否是脏数据
def get_trip_apply(df):
if (df['cost'] == 0) & (df['cost_shift_after'] > 0) & (df['id'] == df['id_shift_after']):
trip = 1
elif (df['cost'] > 0) & (df['cost_shift_before'] == 0) & (df['id'] == df['id_shift_after']) & (df['id'] == df['id_shift_before']):
trip = 1
elif (df['cost'] > 0) & (df['cost_shift_before'] == 0) & (df['id'] != df['id_shift_after']) & (df['id'] == df['id_shift_before']):
trip = 1
else:
trip = -1
return trip### 主函数，完成数据的清洗与整理，并计算行程时间
def get_trip(df):
df = get_shift(df)
df['trip'] = df.apply(get_trip_apply,axis=1)
df = df[df['trip']==1]
df = df.drop(columns=['id_shift_after','id_shift_before','cost_shift_after','cost_shift_before'])
df['trip'] = np.arange(len(df))
df['trip'] = df['trip']//2
df = df.set_index(['trip',df.groupby('trip').cumcount() 1]).unstack().sort_index(level=1,axis=1)
df.columns = ['ori_cost','ori_discount','ori_hour','ori_id','ori_route','ori_timestamp','des_cost','des_discount','des_hour','des_id','des_route','des_timestamp']
df['ori_station'] = df.apply(lambda x: x['ori_route'].split('线')[1],axis=1)
df['ori_route'] = df.apply(lambda x: x['ori_route'].split('线')[0] '线',axis=1)
df['des_station'] = df.apply(lambda x: x['des_route'].split('线')[1],axis=1)
df['des_route'] = df.apply(lambda x: x['des_route'].split('线')[0] '线',axis=1)
df['travel_time(minute)'] = df.apply(lambda x: round((x['des_timestamp']-x['ori_timestamp']).seconds/60,3),axis=1)
order = ['ori_cost','ori_discount','ori_hour','ori_id','ori_route','ori_station','ori_timestamp','des_cost','des_discount','des_hour','des_id','des_route','des_station','des_timestamp','travel_time(minute)']
df = df[order]
return df%%time
df_clean = get_trip(df)
测试一下，可以得到清洗后的数据(440多万条出行记录，包括od的线路与站点、进出站时间以及费用，还有从进站到出站花费的时间），并且相比双循环速度快了很多。（大家如果有更优的方案可以提出来后台私信我，我这个其实也挺慢的
）。
有了这个数据，其实可以分析的东西就很多了。这里推荐大家一篇赵鹏军老师用这个数据写的文章：
故事引人入胜，读完绝对回味无穷。
https://doi.org/10.1016/j.tranpol.2020.03.006
安利完毕，回到主题，今天的主题是可视化客流的特征（主要是od之间的客流特征），这几年有种图特别火，就是一个圈流来流去（学名叫和弦图(chord diagram)），类似这种：
这个是在R语言里画的，python的话也凑合能画（没有R画的好看），实现的具体方法是用holoviews这个库。
先导入holoviews
代码语言：javascript复制import holoviews as hv
from holoviews import opts, dim
import matplotlib.pyplot as plt
plt.rcParams['font.sans-serif']=['Arial Unicode MS']
plt.rcParams['axes.unicode_minus']=False
hv.extension('bokeh')然后画一个各个线路之间客流的和弦图
代码语言：javascript复制station = df_clean.iloc[:,4]
station = station.drop_duplicates()
station = station.reset_index(drop=True).reset_index()
station.columns = ['index','route']
od = df_clean.groupby(['ori_route','des_route'])['ori_cost'].count().to_frame().reset_index()
od['ori_route'] = od['ori_route'].map(station.set_index('route').to_dict()['index'])
od['des_route'] = od['des_route'].map(station.set_index('route').to_dict()['index'])
nodes = hv.Dataset(station, 'index', 'route')
chord = hv.Chord((od, nodes), ['ori_route', 'des_route'], ['ori_cost'])
chord.opts(
opts.Chord(cmap='glasbey', edge_color=dim('ori_route').str(),
labels='route',node_color=dim('index').str(),width=1000,height=1000,node_size=8,edge_alpha=0.4,label_text_font_size='12pt'))
还行吧，可以看到，1号线和2号线还是大哥，无论是进站客流还是出站客流都非常的大，除此之外我们还可以进行站点等级的客流od分析，这里选取进站客流最大的前30个站点之间的流量进行可视化：
代码语言：javascript复制station = df_clean.iloc[:,5]
station = station.drop_duplicates()
station = station.reset_index(drop=True).reset_index()
station.columns = ['index','station']
od = df_clean.groupby(['ori_station','des_station'])['ori_cost'].count().to_frame().reset_index()
od['ori_station'] = od['ori_station'].map(station.set_index('station').to_dict()['index'])
od['des_station'] = od['des_station'].map(station.set_index('station').to_dict()['index'])
nodes = hv.Dataset(station, 'index', 'station')
chord = hv.Chord((od, nodes), ['ori_station', 'des_station'], ['ori_cost'])
top30 = list(od.groupby('ori_station')['ori_cost'].sum().to_frame().sort_values('ori_cost').iloc[-30:].index.values)
top30station = chord.select(ori_station=top30, selection_mode='nodes')
top30station.opts(
opts.Chord(cmap='glasbey_light', edge_color=dim('ori_station').str(),
labels='station',node_color=dim('index').str(),width=1000,height=1000,node_size=8,edge_alpha=0.4,label_text_font_size='8pt'))
封面图就做出来了，看着挺酷，但是乱乱的不好解释
。
