2017年第七届MathorCup高校数学建模挑战赛

B题 共享单车

原题再现：

  共享单车是指企业在校园、地铁站点、公交站点、居民区、商业区、公共服务区等提供自行车单车共享服务，是一种分时租赁模式。共享单车是一种新型共享经济。共享单车已经越来越多地引起人们的注意，由于其符合低碳出行理念，政府对这一新鲜事物也处于善意的观察期。
  很多共享单车公司的单车都有 GPS 定位，能够实现动态化地监测车辆数据、骑行分布数据，进而对单车做出全天候供需预测，为车辆投放、调度和运维提供指引。
  请根据下面附件给出数据及结合根据需要自己收集的数据，完成以下问题：
  （1）根据附件 1 中共享单车的骑行数据，估计共享单车的时空分布情况。如从某地点 A 出发，到达不同地点的分布情况。可分时间段讨论。
  （2）假如根据调查，得到人们的骑行需求估计数据，见附件 2。
  根据问题 1 的估计结果，建立数学模型解决如何优化共享单车的调度问题。
  （3）根据附件 1 的骑行数据和附件 2 的需求数据，判断各区域所需共享单车的满足程度，给出你的度量指标。若增加 100 辆单车，如何进行投放更优。
  （4）附件 3 是某地区投入不同数量共享单车后打车人次的数据。据此分析研究共享单车的投入对该地区打车市场的影响。同时请你收集实际数据进行量化研究。

整体求解过程概述(摘要)

  随着共享单车的出现和普及，共享单车系统在城市公共交通中的地位越来越重要，由于其快捷、方便、环保的特点，共享单车已成为居民解决出行“最后一公里”问题的重要选择。
   首先，我们利用 VBA 编程对数据进行处理，整理成为包含以下五项指标的数据表:出发时间、到达时间、出发区域、到达区域、骑行时长，对骑行时长进行方差分析，从而假设从 A 地出发到达 B 地的骑行时长以平均值代替，并假设骑行速度都统一为正常人的平均速度 15km/h，从而解出各个区域相对距离矩阵。
   针对问题一，我们分为时间分布和空间分布讨论。从时间分布上看，各个区域单车在 420—900 分钟流动频繁、分布密集，其他时段逐渐下降，1380 分钟以后分布趋于静止。从空间分布来说，我们通过 VBA 编程统计分析得到从出发区域至到达区域车次的OD 矩阵，可以明显看出，共享单车在 1 到 2、2 到 3、3 到 4、……、9 到 10 之间的分布较为密集、其他地区基本呈现均匀分布。
  针对问题二，建立了单调度中心调度模型和动态调度优化模型。根据附件 1、附件2 的数据计算出现有单车情况下的各个地区的单车需求量，确定不同地区单车的需求时间和可接受时间。根据问题 1 结果估计不同地区间的最短距离，建立单调度中心调度模型并利用 MATLAB 软件和遗传算法求解得出初始调度方案，最优调度路线为2-5-6-3-7-10-9-2。以单调度中心软时间窗调度模型为基础建立动态需求调度优化模型，向初始静态优化解中不断插入新的调度需求，应用“初始静态优化+实时动态优化”的方法对多个连续静态调度问题进行求解不断优化调度路线，初次优化结果为5-2-3-4-8-10-9，其中 5、2、3 已完成调度。
  针对问题三，要判断各区域所需共享单车的满足程度。首先我们定义了从出发区域i 到到达区域 j 的需求量、满足量以及需求满足比作为衡量各区域所需共享单车满足程度的指标。根据问题一得到的 OD 矩阵以及附件 2 所给的需求矩阵，即可计算出从出发区域 i 到到达区域 j 的需求满足比，进而判断各区域所需共享单车的满足程度。由结果可知，各区域的满足程度从 78%到 99%不等，总体来看，满足程度较好。按照问题二的需求比例即可求得 100 辆单车的投放方案。
  针对问题四，要分析研究共享单车的投入对该地区打车市场的影响。根据附件 3 的数据，我们将单车投放量作为自变量，打车人次作为因变量，作出散点图。观察得到在单车投放量小于 4000 辆时呈现线性关系，大于 4000 辆时呈现非线性关系。因此我们以x=4000 为分段点建立一元线性回归模型、一元非线性回归模型，利用 MATLAB 求解，得出方程。通过观察方程，我们得出共享单车投放量与打车人次呈负相关，即共享单车投放量越大，打车人次越少，从而对打车市场产生负面影响。并且搜集到上海市 4 月份单车使用次数及打车人次的数据，从时间维度分析我们可以得出，固定日期分时段单车使用次数与打车人次数的图像呈现近似的趋势，而固定时段分日期单车使用次数与打车人次数呈现线性负相关的关系。

模型假设：

  ❖ 假设 0 分钟为凌晨 0 点；
  ❖ 假设各单车的行驶速度均匀固定；
  ❖ 假设出发区域 i 到到达区域 j 的相对距离=骑行时长×骑行速度。
  ❖ 假设区域 1 到 10 只有一个调度中心且位于区域 1，配备足够多的同型号调度车辆，能为地区 1 到 10 内所有地点提供调度服务；
  ❖ 假设调度车辆必须从调度中心出发，完成调度任务后需返回调度中心；
  ❖ 假设各地区的调度需求不能大于调度车辆容量；
  ❖ 假设各地区间的最短路网距离已知，调度车辆在任意两地区间行驶时必须选择两地区间的路网最短路径；
  ❖ 假设调度车辆在路网中的平均行驶速度已知且固定，不受其他因素影响；
  ❖ 假设前一△t内完成调度的最后一个地点是下一△t内新调度路径的起点；
  ❖ 假设所收集到的上海市单车使用人次数和出租车使用人次数真实可靠。

问题分析：

  针对问题一，在分析共享单车的时空分布情况时，考虑到时空既可以指不同时间下的空间分布，又可以指不同区域单车的时间分布。我们将其综合考虑，对于时间分布，直接将各个点做出散点图并不直观，所以我们将其划分成每 60 分钟一个时间段，统计时间段内的单车出发、到达车次，绘制图表观察。对于空间分布，我们通过建立 OD 矩阵，描述区域间交通中一定时间范围内所有出行起点和终点之间的出行交通量，间接反映空间分布情况。
  针对问题二，题目要求根据问题 1 的共享单车的时空分布情况的估计结果，建立数学模型解决如何优化共享单车的调度问题。首先，假设骑行人的骑行距离一定，由平均骑行时间可计算出各地区间的相对距离。其次，根据客户需求是否实时变动，调度问题可分为静态和动态调度问题。对动态调度问题求解时可将动态调度问题转化为多个连续静态调度问题进行求解，即“初始静态优化+实时动态优化”。先根据系统内各地点在某一时刻的调度需求，形成初始的调度方案，然后不断更新各地点的调度需求，同时更新调度方案，实现调度需求与调度路径的不断反馈，最终达到动态优化的目的。
  针对问题三，要判断各区域所需共享单车的满足程度。首先我们应该衡量各区域所需共享单车满足程度的指标，建立指标评价体系，但是由于题目所给数据不全面，无法得到指标体系的具体数据，所以我们只用需求满足比来作为衡量各区域所需共享单车满足程度的指标。
  针对问题四，要分析研究共享单车的投入对该地区打车市场的影响。根据附件 3 的数据，我们可以将单车投放量作为自变量，打车人次作为因变量，作出散点图，观察单车投放量与打车人次的关系，探究是线性还是非线性，是正相关还是负相关。再用MATLAB 软件得到具体的函数关系。进而得到单车投放量对打车市场的影响关系。

模型的建立与求解整体论文缩略图

全部论文请见下方“ 只会建模 QQ名片” 点击QQ名片即可

程序代码：(代码和文档not free)

The actual procedure is shown in the screenshot

import pandas as pd
from math import radians, cos, sin, asin, sqrt,ceil
import numpy as np
import geohash

#数据读取
data = pd.read_csv("./mobike_shanghai_sample_updated.csv")
print(data.head(10))
print(data.info())
data['start_time'] = pd.to_datetime(data['start_time'])
data['end_time'] = pd.to_datetime(data['end_time'])
print(data.info())
data["lag"] = (data.end_time - data.start_time).dt.seconds/60
def geodistance(item):
    lng1_r, lat1_r, lng2_r, lat2_r = map(radians, [item["start_location_x"], item["start_location_y"], item["end_location_x"], item["end_location_y，"]]) # 经纬度转换成弧度
    dlon = lng1_r - lng2_r
    dlat = lat1_r - lat2_r
    dis = sin(dlat/2)**2 + cos(lat1_r) * cos(lat2_r) * sin(dlon/2)**2
    distance = 2 * asin(sqrt(dis)) * 6371 * 1000 # 地球平均半径为6371km
    distance = round(distance/1000,3)
    return distance
#data按行应用geodistance()得到distance列的数值
data["distance"] = data.apply(geodistance,axis=1)
#通过摩拜单车的踪迹获取每次交易骑行的路径
def geoaadderLength(item):
    track_list = item["track"].split("#")
    adderLength_item = {}
    adderLength = 0
    for i in range(len(track_list)-1):
        start_loc = track_list[i].split(",")
        end_loc = track_list[i+1].split(",")
        adderLength_item["start_location_x"],adderLength_item["start_location_y"] = float(start_loc[0]),float(start_loc[1])
        adderLength_item["end_location_x"],adderLength_item["end_location_y"] = float(end_loc[0]),float(end_loc[1])
        adderLength_each = geodistance(adderLength_item)
        adderLength = adderLength_each + adderLength
    return adderLength
data["adderLength"] = data.apply(geoaadderLength,axis=1)
data['weekday'] = data.start_time.apply(lambda x: x.isoweekday())
data['hour'] = data.start_time.apply(lambda x: x.utctimetuple().tm_hour)
data['cost'] = data.lag.apply(lambda x: ceil(x/30))
#因数据集仅包含八月份发起的订单数据，故以9月1日为R值计算基准
data['r_value_single'] = data.start_time.apply(lambda x: 32 - x.timetuple().tm_mday) 
# 按每个用户id所有订单日期距9/1相差天数的最小值作为r值
r_value = data.groupby(['userid']).r_value_single.min()    
f_value = data.groupby(['userid']).size()    # 按每个用户id八月累积订单数量作为f值
m_value = data.groupby(['userid']).cost.sum()    # 按每个用户id八月累积消费金额作为m值
#把r值、f值、m值组合成DataFrame
rfm_df = pd.DataFrame({'r_value':r_value,'f_value':f_value,"m_value":m_value})
rfm_df["r_score"] = pd.cut(rfm_df["r_value"],5,labels=[5,4,3,2,1]).astype(float)
rfm_df["f_score"] = pd.cut(rfm_df["f_value"],5,labels=[1,2,3,4,5]).astype(float)
rfm_df["m_score"] = pd.cut(rfm_df["m_value"],5,labels=[1,2,3,4,5]).astype(float)
#后面*1是为了把布尔值false和true转成0和1
rfm_df["r是否大于均值"] = (rfm_df["r_score"] > rfm_df["r_score"].mean())*1
rfm_df["f是否大于均值"] = (rfm_df["f_score"] > rfm_df["f_score"].mean())*1
rfm_df["m是否大于均值"] = (rfm_df["m_score"] > rfm_df["m_score"].mean())*1
#把每个用户的rfm三个指标统合起来
rfm_df["class_index"] = (rfm_df["R是否大于均值"]*100) + (rfm_df["f是否大于均值"]*10) + (rfm_df["m是否大于均值"]*1) 

def transform_user_class(x):
    if x == 111:
        label = "重要价值用户"
    elif x == 110:
        label = "消费潜力用户"
    elif x == 101:
        label = "频次深耕用户"
    elif x == 100:
        label = "新用户"
    elif x == 11:
        label = "重要价值流失预警用户"
    elif x == 10:
        label = "一般用户"
    elif x == 1:
        label = "高消费唤回用户"
    elif x == 0:
        label = "流失用户"
    return label
rfm_df["user_class"] = rfm_df["class_index"].apply(transform_user_class)
data = data.merge(rfm_df["user_class"], on = 'userid', how = 'inner')

全部论文请见下方“ 只会建模 QQ名片” 点击QQ名片即可