论文解读
郭王懿,孙楚天,陈泰劼,张云天
编者按
共乘出行极大地改变了人们的日常出行方式。如何高效运营背后的双边平台是极具挑战性的工作。滴滴出行、Lyft公司是其中的佼佼者。本专题将探讨双边平台运营中的一个关键问题:如何科学度量双边平台的运营效率?
本专题分为上下两篇,递进推出。其中,上篇围绕滴滴出行、上海财经大学发表在统计学“四大”期刊之一《Journal of the American Statistical Association》的一篇文章 [1] 展开;下篇关注解读Lyft公司的一篇优质技术博客 [2]。
引言
近年来,动态多时空供需网络被应用于气候科学、社会科学、神经科学等多个学科,网约车平台的动态需求和供应网络是多时空动态网络的典型应用。网约车司乘匹配问题是在相同的无向(或有向)图 G = ( V , E ) G =(\mathbb{V},\mathbb{E}) G=(V,E)上匹配“供应”和“需求”两个未归一化的测度序列,其中 V \mathbb{V} V 和 E \mathbb{E} E 分别是顶点集和连接顶点组的边集。具体的做法是将一个城市划分为数百个不重叠的网格,作为顶点集 V \mathbb{V} V,边结构 E \mathbb{E} E 由道路网络和位置功能确定。在每个决策时间段内通过观测订单需求和空闲司机的分布,使用订单分发策略将客户请求与可能的周围空闲司机进行匹配。
本文考虑的基本问题是如何量化描述双边市场动态多时空供需网络的空间匹配程度,特别是网约车平台(例如 Uber 和 DiDi)。在基于转移成本的匹配度定义下,计算出最优匹配度就意味着计算出将供应网络分布转移到需求网络分布的最优转移策略,即订单分配策略,因此研究多时空供需网络的匹配度非常关键。为了解决这个问题,首先引入一个加权图结构 ( G , W , C ) (G, W, C) (G,W,C)来描述城市的交通网络和交通成本。具体来说,将每个市场划分为 N 个不相交的区域,并将它们视为顶点,表示为 V = { v 1 , . . . , v N } \mathbb{V} = \{v_1, . . . , v_N\} V={v1,...,vN}。设 E \mathbb{E} E是任意可能的顶点对之间的边集,使得 ( v i , v j ) ∈ E ⊂ V × V (vi, vj) ∈ \mathbb{E} ⊂ \mathbb{V} × \mathbb{V} (vi,vj)∈E⊂V×V是具有非负权重 w i j w_{ij} wij(例如,运输成本)的边。对于所有 w i j w_{ij} wij,设 w i j w_{ij} wij。加权图结构包括一个无向(或有向)图 G = ( V , E ) G =(\mathbb{V},\mathbb{E}) G=(V,E)以及一个权重矩阵 G = ( V , E ) G =(\mathbb{V},\mathbb{E}) G=(V,E),其中 w i j w_{ij} wij是非负权重。从 w i j w_{ij} wij到 v j v_j vj的基于图的运输成本定义为 c i j = m i n K ≥ 0 , ( i k ) k = 0 K : v i → v j { ∑ k w i k , i k + 1 : ∀ k ∈ [ [ 0 , K − 1 ] ] , ( v i k , v i k + 1 ) ∈ E } c_{ij} = min_{K≥0,(i_k)_{k=0}^{K}:vi→vj} \{\sum_k w_{i_k,i_{k+1}} : ∀k ∈[[0, K - 1]],(v_{i_k}, v_{i_k+1}) ∈ \mathbb{E}\} cij=minK≥0,(ik)k=0K:vi→vj{∑kwik,ik+1:∀k∈[[0,K−1]],(vik,vik+1)∈E},其中 ( i k ) k = 0 K : v i → v j (i_k)_{k=0}^{K} : v_i → v_j (ik)k=0K:vi→vj代表以 v i 0 = v i v_{i_0} = v_i vi0=vi开始并以 v i K = v j v_{i_K} = v_j viK=vj结束的 G G G上通过 E \mathbb{E} E的任何路径。因此, c i j c_{ij} cij是从 v i v_i vi到 v j v_j vj的测地距离,即从 v i v_i vi到 v j v_j vj运输一个单位物体的最小成本。因此,可以在 ( G , W ) (G, W) (G,W)上定义一个运输成本矩阵,表示为 C = ( c i j ) ∈ R N × N C = (c_{ij}) ∈ R^{N×N} C=(cij)∈RN×N。
其次,需要引入一个距离来量化每个时间区间上的需求和供应质量间的差异,以及在 ( G , W , C ) (G, W, C) (G,W,C)上跨时间的差异。在给定的时间区间上,定义 ν j = ν ( v j ) ν_j= ν(v_j) νj=ν(vj)和 μ j = μ ( v j ) μ_j= μ(v_j) μj=μ(vj)为顶点 v j v_j vj的点质量两个Lebesgue测度 ν ν ν和 μ μ μ,分别代表了打车平台节点 v j v_j vj内的用户请求和可用司机的数量。每个时间段的供应和需求网络可以建模为 ( G , W , C ) (G, W, C) (G,W,C)上的两个离散的Lebesgue测度 μ μ μ和 ν ν ν。
基于最佳传输理论的 Wasserstein 类度量已被证明是度量复杂空间中对象的有效工具。然而,现有的Wasserstein类度量并不直接适用于在2.1节中详细定义的 ( G , W , C ) (G, W, C) (G,W,C)上的两个不均衡度量的比较。因此,作者引入了一种基于图的均衡度量(GEM),并将其表述为一个不均衡的最优运输问题。
一、基于图结构的均衡度量方法
1. 现有的Wasserstein类距离
有很多经典的度量两个分布之间差异的标准,大致分为两类,第一类是用像素的差异定义距离,第二类是基于转移成本定义的距离。
L
p
L_p
Lp距离、TV距离、KL散度等都是第一类度量准则,Hellinger距离是经典案例:
D
H
2
(
μ
,
v
)
=
0.5
∫
X
(
d
μ
/
d
ν
0
−
d
ν
/
d
ν
0
)
2
d
ν
0
D_H^2(\mu, v)=0.5 \int_X\left(\sqrt{d \mu / d \nu_0}-\sqrt{d \nu/ d \nu_0}\right)^2d\nu_0
DH2(μ,v)=0.5∫X(dμ/dν0−dν/dν0)2dν0,意味衡量两个分布在测度
ν
0
\nu_0
ν0上的均方根差异。但这一类基于像素差异的指标有两大主要缺陷:
- Hellinger类指标没有考虑不同位置的连接关系,即没有考虑 X X X的拓扑分布。
- Hellinger类指标需要进行归一化操作,即要求 μ ( X ) = ∫ X d μ = ν ( X ) = ∫ X d ν \mu(X)=\int_X d \mu=\nu(X)=\int_X d \nu μ(X)=∫Xdμ=ν(X)=∫Xdν。
为了解决如上两个缺陷,部分学者提出了Wasserstein距离等第二类基于转移成本的距离,在组合优化等分配问题领域有很深的研究。Wasserstein距离的定义如下:
D
W
(
μ
,
v
∣
c
)
=
inf
γ
∈
M
+
(
X
×
Y
)
{
∫
X
×
Y
c
d
γ
+
ι
{
=
}
(
P
#
X
γ
∣
μ
)
+
ι
{
=
}
(
P
#
Y
γ
∣
v
)
}
\begin{aligned} & D_W(\mu, v \mid c) =\inf _{\gamma \in M_{+}(X \times Y)}\left\{\int_{X \times Y} c d \gamma+\iota_{\{=\}}\left(P_{\#}^X \gamma \mid \mu\right)+\iota_{\{=\}}\left(P_{\#}^Y \gamma \mid v\right)\right\} \end{aligned}
DW(μ,v∣c)=γ∈M+(X×Y)inf{∫X×Ycdγ+ι{=}(P#Xγ∣μ)+ι{=}(P#Yγ∣v)}
其中 X , Y X,Y X,Y是两个Hausdorff拓扑空间,损失函数c是 X × Y → R ∪ { ∞ } X \times Y \rightarrow R \cup\{\infty\} X×Y→R∪{∞}的映射,转移函数 γ \gamma γ被定义为一个非负测度 γ ∈ M + ( X × Y ) \gamma \in M_{+}(X \times Y) γ∈M+(X×Y),等式约束 ι { = } ( α ∣ β ) \iota_{\{=\}}(\alpha \mid \beta) ι{=}(α∣β)代表当 α = β \alpha=\beta α=β时为0,其余为无穷。该度量的前提条件是 μ ∈ M + ( X ) \mu\in M_{+}(X) μ∈M+(X)和 ν ∈ M + ( Y ) \nu\in M_{+}(Y) ν∈M+(Y)总量一样,即 μ ( X ) = ν ( Y ) \mu(X)=\nu(Y) μ(X)=ν(Y)。 P # X γ 和 P # Y γ P_\#^X\gamma和P_\#^Y\gamma P#Xγ和P#Yγ分别是 γ \gamma γ的两个边际分布。
该方法可以获得最优的分配方案 γ \gamma γ,但当供需总量不一样时,Wasserstein距离无法获得一个可行解,而在实际中通常存在供需分布总量不一致。因此,部分学者提出了基于散度的不均衡运输问题度量方法。首先是 ϕ − \phi- ϕ−散度的定义。首先有熵函数 ϕ \phi ϕ,对于 μ , v ∈ M ( T ) , d μ d v v + μ ⊥ \mu, v \in M(T), \frac{d \mu}{d v} v+\mu^{\perp} μ,v∈M(T),dvdμv+μ⊥是 μ \mu μ对于 ν \nu ν的勒贝格分解。散度 D ϕ D_\phi Dϕ的定义为
D
ϕ
(
μ
∣
ν
)
:
=
∫
T
φ
(
d
μ
d
v
)
d
v
+
φ
∞
′
μ
⊥
(
T
)
D_\phi(\mu|\nu):=\int_T \varphi\left(\frac{d \mu}{d v}\right) d v+\varphi_{\infty}^{\prime} \mu^{\perp}(T)
Dϕ(μ∣ν):=∫Tφ(dvdμ)dv+φ∞′μ⊥(T),因此可以定义广义Wasserstein距离(GWD):
D
φ
1
,
φ
2
(
μ
,
ν
∣
c
)
=
inf
γ
∈
M
+
(
X
×
Y
)
{
∫
X
×
Y
c
d
γ
+
D
φ
1
(
P
#
X
γ
∣
μ
)
+
D
φ
2
(
P
#
Y
γ
∣
v
)
}
.
\begin{aligned} D_{\varphi_1, \varphi_2}(\mu, \nu \mid c) \quad=\inf _{\gamma \in M_{+}(X \times Y)}\left\{\int_{X \times Y} c d \gamma+\mathcal{D}_{\varphi_1}\left(P_{\#}^X \gamma \mid \mu\right)+\mathcal{D}_{\varphi_2}\left(P_{\#}^Y \gamma \mid v\right)\right\} . \end{aligned}
Dφ1,φ2(μ,ν∣c)=γ∈M+(X×Y)inf{∫X×Ycdγ+Dφ1(P#Xγ∣μ)+Dφ2(P#Yγ∣v)}.
GWD用 ϕ − \phi- ϕ−散度来量化输运计划与供应测度和需求测度两个不平衡测度 μ \mu μ和 ν \nu ν之间的偏差。但是由于GWD测度是基于熵函数 ϕ \phi ϕ的度量,其求解结果 γ ∗ \gamma^* γ∗不具有任何物理意义,因此不能提供可靠的调度方案。除此之外,GWD还有以下两个主要问题:
- 很多实际双边平台中,只有一个测度中的质点允许移动, 另一个测度中是固定的,比如打车平台中只有空闲司机的分布可以调度,客户订单是无法移动的。
- 在供需网络中,运输费用通常是一个分布而不是一个常数,并且顶点到自身的运输成本可能不是零。
2.Graph-based Equilibrum Metrics
针对上述research gap,本文提出Graph-based Equilibrum Metrics(GEM),首先在权重图问题
(
G
,
W
,
C
)
(G,W,C)
(G,W,C)中,我们将供应和需求定义为两个
M
+
(
V
)
M_+(\mathbb{V})
M+(V)的测度
μ
,
ν
\mu,\nu
μ,ν,其中
μ
\mu
μ的质点可以移动,
ν
\nu
ν的质点是固定的。
N
i
\mathcal{N}_i
Ni定义节点
v
i
v_i
vi的邻接节点集。其次由于交通单行道等原因,
v
i
∈
N
j
v_i\in\mathcal{N}_j
vi∈Nj并不能保证
v
j
∈
N
i
v_j\in\mathcal{N}_i
vj∈Ni,因此我们定义GEM为:
ρ
λ
(
μ
,
v
∣
G
,
C
)
=
inf
μ
~
∈
M
+
(
V
)
,
γ
∈
M
+
(
V
×
V
)
{
∣
v
−
μ
~
∣
+
λ
∫
V
×
V
c
d
γ
}
\begin{aligned} & \rho_\lambda(\mu, v \mid G, C) \\ & =\inf _{\tilde{\mu} \in M_{+}(\mathbb{V}), \gamma \in M_{+}(\mathbb{V} \times \mathbb{V})}\left\{|v-\tilde{\mu}|+\lambda \int_{\mathbb{V} \times \mathbb{V}} c d \gamma\right\} \end{aligned}
ρλ(μ,v∣G,C)=μ~∈M+(V),γ∈M+(V×V)inf{∣v−μ~∣+λ∫V×Vcdγ}
s . t . ∣ μ ∣ = ∣ μ ~ ∣ , ( P # 1 V γ ) ( v i ) = ∑ v j ∈ N i γ ( v i , v j ) = μ i and ( P # 2 V γ ) ( v i ) = ∑ v i ∈ N j γ ( v j , v i ) = μ ~ i , \begin{aligned} s.t. & |\mu|=|\tilde{\mu}|, \quad\left(P_{\# 1}^V \gamma\right)\left(v_i\right)=\sum_{v_j \in \mathcal{N}_i} \gamma\left(v_i, v_j\right)=\mu_i \text { and } \\ & \left(P_{\# 2}^V \gamma\right)\left(v_i\right)=\sum_{v_i \in \mathcal{N}_j} \gamma\left(v_j, v_i\right)=\tilde{\mu}_i, \end{aligned} s.t.∣μ∣=∣μ~∣,(P#1Vγ)(vi)=vj∈Ni∑γ(vi,vj)=μi and (P#2Vγ)(vi)=vi∈Nj∑γ(vj,vi)=μ~i,
其中 λ \lambda λ是非负的超参数,约束用于限制 γ \gamma γ将 μ \mu μ转移到 μ ~ \tilde{\mu} μ~。GEM本质上是Wasserstein距离和供需差异的 L 1 L_1 L1范数。相比于现有的非均衡最有运输问题解决方法,本文提出的GEM有两大优势:
- 考虑到了实际的拓扑信息,只允许邻接节点之间的运输动作,如图所示。
- 可以通过
λ
\lambda
λ来均衡供需匹配度和运输成本之间的比例。通常通过data-driven的方法来选择模型的
λ
\lambda
λ。
GEM考虑到拓扑结构,而下图的GWD没有考虑拓扑结构
二、计算方法
GEM的计算方法主要分为两步,第一步是通过紧凑形式对原定义重构,并对绝对值进行线性化处理,第二步是对LP问题进行对偶,降低决策变量的维度。为了描述计算求解过程,我们需要定义一些紧凑形式的变量,首先根据GEM的离散形式定义:
ρ
λ
(
μ
,
v
∣
G
,
C
)
=
min
γ
∈
Γ
{
∥
v
−
μ
~
∥
1
+
λ
∑
v
i
∈
V
∑
v
j
∈
V
c
i
j
γ
i
j
}
subject to
∑
v
j
∈
N
i
γ
i
j
=
μ
i
,
∑
v
j
∉
N
i
γ
i
j
=
0
,
∀
i
,
and
μ
~
i
=
∑
v
i
∈
N
j
γ
j
i
for
∀
v
i
∈
V
,
\begin{aligned} & \rho_\lambda(\mu, v \mid G, C)=\min _{\gamma \in \Gamma}\left\{\|\boldsymbol{v}-\tilde{\boldsymbol{\mu}}\|_1+\lambda \sum_{v_i \in \mathbb{V}} \sum_{v_j \in \mathbb{V}} c_{i j} \gamma_{i j}\right\} \\ & \text { subject to } \sum_{v_j \in \mathcal{N}_i} \gamma_{i j}=\mu_i, \quad \sum_{v_j \notin \mathcal{N}_i} \gamma_{i j}=0,\quad\forall i, \quad \text { and } \\ & \tilde{\mu}_i=\sum_{v_i \in \mathcal{N}_j} \gamma_{j i} \text { for } \forall v_i \in \mathbb{V},\end{aligned}
ρλ(μ,v∣G,C)=γ∈Γmin⎩
⎨
⎧∥v−μ~∥1+λvi∈V∑vj∈V∑cijγij⎭
⎬
⎫ subject to vj∈Ni∑γij=μi,vj∈/Ni∑γij=0,∀i, and μ~i=vi∈Nj∑γji for ∀vi∈V,
对于
v
j
∉
N
i
v_j \notin \mathcal{N}_i
vj∈/Ni的
γ
i
j
\gamma_{ij}
γij,其取值一定是0,因此使用
γ
~
=
V
e
c
{
γ
i
j
,
j
∈
N
i
}
∈
R
N
0
×
1
\tilde{\gamma}=Vec\{\gamma_{ij},j\in\mathcal{N}_i\}\in R^{N_0\times1}
γ~=Vec{γij,j∈Ni}∈RN0×1来简化
γ
\gamma
γ,其中
N
0
=
∑
i
=
1
N
n
i
N_0 = \sum_{i=1}^N n_i
N0=∑i=1Nni代表所有节点的相邻节点总数目之和。这种简化方式可以变量的复杂度从
O
(
N
2
)
O(N^2)
O(N2)降低到
O
(
N
0
)
O(N_0)
O(N0)。为了方便后续计算,现在定义
A
1
,
A
2
A_1,A_2
A1,A2两个
N
×
N
0
N\times N_0
N×N0的矩阵,
A
1
A_1
A1的第
∑
j
=
1
i
−
1
n
j
+
1
\sum_{j=1}^{i-1}n_j+1
∑j=1i−1nj+1到
∑
j
=
1
i
n
j
\sum_{j=1}^{i}n_j
∑j=1inj的元素全为1,以外的元素全是0;
A
2
A_2
A2则是
∑
p
=
1
j
−
1
n
p
+
1
\sum_{p=1}^{j-1}n_p+1
∑p=1j−1np+1到
∑
p
=
1
j
n
p
\sum_{p=1}^{j}n_p
∑p=1jnp之间是1,其余是0。
A
1
A_1
A1结构示意图如下所示,
如果所有的路径是双向联通的,即问题是一个无向图即 γ i j = γ j i \gamma_{ij}=\gamma_{ji} γij=γji,则 A 1 A_1 A1与 A 2 A_2 A2完全一样。因为这里 A 1 A_1 A1与 A 2 A_2 A2分别代表 P # 1 V γ P_{\# 1}^V \gamma P#1Vγ和 P # 2 V γ P_{\# 2}^V \gamma P#2Vγ。
我们定义 C ~ ∈ R N 0 × 1 \widetilde{C} \in R^{N_0 \times 1} C ∈RN0×1是成本向量,用于描述每个路径的单位运输成本。因此我们定义:
A = [ A 1 0 0 0 ∥ A 2 I N − I N 0 ∥ A 2 − I N 0 I N ] , b = [ μ T , v T , v T ] T , A=\left[\begin{array}{llllllllll} A_1 & \mathbf{0} & \mathbf{0} & \mathbf{0} \| A_2 & \mathbf{I}_N & -\mathbf{I}_N & \mathbf{0} \| A_2 & -\mathbf{I}_N & \mathbf{0} & \mathbf{I}_N \end{array}\right],\\ \mathbf{b}=\left[\boldsymbol{\mu}^T, \boldsymbol{v}^T, \boldsymbol{v}^T\right]^T, A=[A1000∥A2IN−IN0∥A2−IN0IN],b=[μT,vT,vT]T,
其中, μ = ( μ 1 , … , μ N ) T , A ∈ R 3 N × ( N 0 + 3 N ) , b ∈ R 3 N \boldsymbol{\mu}=\left(\mu_1, \ldots, \mu_N\right)^T, A \in R^{3 N \times\left(N_0+3 N\right)}, \mathbf{b} \in R^{3 N} μ=(μ1,…,μN)T,A∈R3N×(N0+3N),b∈R3N, I N \boldsymbol{\rm{I}}_N IN是单位阵。
将上述紧凑形式的变量引入定义式可以重构为:
min
{
∥
v
−
A
2
γ
~
∥
1
+
λ
C
~
T
γ
~
}
subject to
A
1
γ
~
=
μ
,
γ
~
≥
0.
\min \left\{\left\|\boldsymbol{v}-A_2 \widetilde{\gamma}\right\|_1+\lambda \widetilde{C}^T \widetilde{\gamma}\right\}\\ \text{subject to}\quad A_1 \tilde{\gamma}=\mu, \tilde{\gamma} \geq 0.
min{∥v−A2γ
∥1+λC
Tγ
}subject toA1γ~=μ,γ~≥0.
引入 S ∈ R N × 1 S\in R^{N\times1} S∈RN×1将非线性项线性化以后等价于:
min { 1 T S + λ C ~ T γ ~ } subject to A 1 γ ~ = μ , A 2 γ ~ + S ≥ v , A 2 γ ~ − S ≤ v , γ ~ ≥ 0 , and S ≥ 0 , \begin{aligned} & \min \left\{\mathbf{1}^T S+\lambda \widetilde{C}^T \widetilde{\gamma}\right\} \\ &\text { subject to } A_1 \tilde{\gamma}=\boldsymbol{\mu}, A_2 \widetilde{\gamma}+S \geq \boldsymbol{v}, A_2 \tilde{\gamma}-S \leq \boldsymbol{v}, \tilde{\gamma} \geq 0, \text { and } S \geq 0, \end{aligned} min{1TS+λC Tγ } subject to A1γ~=μ,A2γ +S≥v,A2γ~−S≤v,γ~≥0, and S≥0,
上述线性规划问题可以写为紧凑形式:
min
X
{
B
T
X
}
subject to
A
X
=
b
,
X
≥
0
,
\begin{aligned} &\min _X\left\{B^T X\right\}\\ &\text { subject to } A X=\mathbf{b}, X \geq 0, \end{aligned}
Xmin{BTX} subject to AX=b,X≥0,
其中
B
=
(
λ
C
~
T
,
1
T
,
0
T
,
0
T
)
T
B=\left(\lambda \widetilde{C}^T, \mathbf{1}^T, \mathbf{0}^T, \mathbf{0}^T\right)^T
B=(λC
T,1T,0T,0T)T,$ X=\left(\tilde{\gamma}^T, S^T, \boldsymbol{w}_1^T, \boldsymbol{w}_2T\right)T,$
w
1
T
,
w
2
T
\boldsymbol{w}_1^T, \boldsymbol{w}_2^T
w1T,w2T是不等式约束转为等式约束的松弛变量,将其矩阵拆解可以表示为下图。
该典型LP问题的对偶形式为
max y ∈ R 3 N { b T y } subject to A T y ≤ B \begin{aligned} &\max _{\mathbf{y} \in R^{3 N}}\left\{\mathbf{b}^T \mathbf{y}\right\} \\&\text { subject to } A^T \mathbf{y} \leq B \end{aligned} y∈R3Nmax{bTy} subject to ATy≤B
按照矩阵拆解的示意图如下所示:
对偶之后的LP问题的变量维度可以从 N 0 + 3 N N_0+3N N0+3N缩减到 3 N 3N 3N.
三、GEM在共乘出行平台的三方面应用
首先需要基于每个城市的环境建立动态权重的图结构,首先将城市分为不重合的 N N N个正六边形, N i = ∪ k = 0 2 N i k \mathcal{N}_i=\cup_{k=0}^2 \mathcal{N}_i^k Ni=∪k=02Nik其中 N i k \mathcal{N}_i^k Nik是 ν i \nu_i νi的第 k k k层外部节点,第零层只有它本身。构建 G = ( V , E ) G=(\mathbb{V}, \mathbb{E}) G=(V,E),下一步构建 W t = ( w i j t ) W_t=(w_{ijt}) Wt=(wijt),其中 w i j t w_{ijt} wijt是两个顶点在 t t t时间段的距离,下一步通过 W t W_t Wt计算 C t = ( c i j t ) C_t=(c_{ijt}) Ct=(cijt)。最后得到动态权重图 ( G , W t , C t ) (G,W_t,C_t) (G,Wt,Ct)
1. 应答率预测
指标定义
首先可以测量观测到的动态供应网络与需求网络之间的实时最优距离,定义 O = { ( o i t ) } t \mathbf{O}=\{(o_{it})\}_t O={(oit)}t, D = { ( d i t ) } t \mathbf{D}=\{(d_{it})\}_t D={(dit)}t分别为需求和供应动态系统。定义 μ t = ( d i t ) i , ν t = ( o i t ) i \mu_t=(d_{it})_i,\nu_t=(o_{it})_i μt=(dit)i,νt=(oit)i,其中i是节点 ν i \nu_i νi的索引。然后求解 ρ ( t ) = ρ λ ( μ t , v t ∣ G , C t ) \rho(t)=\rho_\lambda\left(\mu_t, v_t \mid G, C_t\right) ρ(t)=ρλ(μt,vt∣G,Ct)求得t时段的最优运输策略 ( μ ~ t ∗ = ( d ~ i t ∗ ) i , γ t ∗ ) \left(\tilde{\mu}_{t *}=\left(\tilde{d}_{i t *}\right)_i, \gamma_{t *}\right) (μ~t∗=(d~it∗)i,γt∗)。
本文将每个顶点 ν i \nu_i νi的最优供需比定义为需求 o i t o_{it} oit比上“最优供应” d ~ i t ∗ + ι { = } ( d ~ i t ∗ = 0 ) \tilde{d}_{i t *}+\iota_{\{=\}}\left(\tilde{d}_{i t *}=0\right) d~it∗+ι{=}(d~it∗=0),记作 DSr i t \operatorname{DSr}_{it} DSrit,定义中的第二项为了防止分母为0。同时我们可以定义每个顶点在每个时刻的最优的供需差为 DSd i t = o i t − d ~ i t ∗ \operatorname{DSd}_{it}=o_{it}-\tilde{d}_{i t *} DSdit=oit−d~it∗。
在上述两个定义下,我们可以基于GEM测量方法 ( DSr i t , DSd i t ) (\operatorname{DSr}_{it},\operatorname{DSd}_{it}) (DSrit,DSdit)创造时空地图。进一步的,我们将 ( DSr i t , DSd i t ) (\operatorname{DSr}_{it},\operatorname{DSd}_{it}) (DSrit,DSdit)扩展到更宽的时间尺度 T 0 \mathcal{T}_0 T0和一个更大的区域 V 0 ∈ V \mathbb{V}_0\in \mathbb{V} V0∈V,因此本文定义了一个加权平均供需比和一个加权平均供需差:
DSr T 0 ( V 0 ) = ∫ t ∈ T 0 ∑ i ∈ V 0 w i t D S r i t d t ∫ t ∈ T 0 ∑ i ∈ V 0 w i t d t and ADSd T 0 ( V 0 ) = ∫ t ∈ T 0 ∑ i ∈ V 0 w i t ∣ D S d i t ∣ d t ∫ t ∈ T 0 ∑ i ∈ V 0 w i t d t \begin{aligned} \operatorname{DSr}_{\mathcal{T}_0}\left(\mathbb{V}_0\right) & =\frac{\int_{t \in \mathcal{T}_0} \sum_{i \in \mathbb{V}_0} w_{i t} \mathrm{DSr}_{i t} d t}{\int_{t \in \mathcal{T}_0} \sum_{i \in \mathbb{V}_0} w_{i t} d t} \text { and } \\ \operatorname{ADSd}_{\mathcal{T}_0}\left(\mathbb{V}_0\right) & =\frac{\int_{t \in \mathcal{T}_0} \sum_{i \in \mathbb{V}_0} w_{i t}\left|\mathrm{DSd}_{i t}\right| d t}{\int_{t \in \mathcal{T}_0} \sum_{i \in \mathbb{V}_0} w_{i t} d t} \end{aligned} DSrT0(V0)ADSdT0(V0)=∫t∈T0∑i∈V0witdt∫t∈T0∑i∈V0witDSritdt and =∫t∈T0∑i∈V0witdt∫t∈T0∑i∈V0wit∣DSdit∣dt
这种定义方式中的 w i t w_{it} wit可以设置为 o i t o_{it} oit来强调高需求节点的指标。一个良好的网约车市场均衡对应着较小的 ( DSr i t , DSd i t ) (\operatorname{DSr}_{it},\operatorname{DSd}_{it}) (DSrit,DSdit)
算例分析
本文使用DiDi GAIA实验室的2018年四月21日到5月20日的需求数据和空闲的司机数据。将城市分割为了 N = 800 N=800 N=800个六边形子区域,每个区域边长1400m,设置边权重 w i j w_{ij} wij为两个子区域之间的行车距离,并且计算每分钟每个区域内的空闲司机数量和订单数量。
订单应答率通常定义为应答的订单数量除以总的订单数量。本文根据历史数据预测未来10分钟和60分钟的应答率对数值。首先分别计算10分钟间隔的Hellinger距离、L2距离、Wasserstein距离和GEM。其中L2距离和Hellinger距离通过10个连续的一分钟时间间隔的距离计算;Wasserstein距离的计算首先需要将10个连续的一分钟间隔的供需数据进行归一化,来保证可以求得可行的Wasserstein距离。然后我们通过对应的权重,获得每10分钟一个的聚合后的Wasserstein距离度量值。关于GEM相关指标,首先计算每一分钟的 DSr i t \operatorname{DSr}_{it} DSrit,然后每十分钟计算一次权重聚合后的 DSr T ( V ) \operatorname{DSr}_{\mathcal{T}}(V) DSrT(V)。
本文使用百分比绝对值误差(MAPE)和均方根误差(RMSE)来衡量相应指标对预测精度的影响。本文将数据集分为了训练集和测试集,使用线性回归方法,使用先前10个时段的指标和五天前同一时间窗口的历史指标来预测未来6个时刻的订单响应率的对数值。
t
+
10
,
t
+
60
t+10,t+60
t+10,t+60分别代表短时间和长时间尺度的预测,结果如下图所示,可以看出GEM明显优于其他三个指标,传统的指标可能无法充分描述加权图的动态传输和系统平衡。
使用不同指标的应答率预测结果
2. 订单调度策略
通过利用历史数据中的供需两侧信息
{
(
D
S
r
i
t
,
D
S
d
i
t
)
:
(
v
i
,
t
)
∈
V
×
T
0
}
\left\{\left(\mathrm{DSr}_{i t}, \mathrm{DSd}_{i t}\right):\left(v_i, t\right) \in \mathbb{V} \times \mathcal{T}_0\right\}
{(DSrit,DSdit):(vi,t)∈V×T0},我们可以设计大范围的双边平台调度策略。传统的调度策略关注当前时刻需求的满足,比如就近原则以及first-come-first-go策略等,但大多无法利用城市的空间结构,因此无法获得最优的调度。本文利用GEM记录更多的历史供需网络信息,用于描述其对当前订单平均预期收益的影响。
N
o
,
N
d
N_o,N_d
No,Nd分别代表订单数量和空闲司机数量,当吧司机
l
l
l分配给订单$k
时,
时,
时,A(k,l)
作为二分图权重是司机的期望收益。
作为二分图权重是司机的期望收益。
作为二分图权重是司机的期望收益。x_{kl}$为1代表订单的分配。因此全局订单分配策略可以归纳为如下二分图匹配问题:
arg
max
x
k
l
∑
k
=
0
N
d
∑
l
=
0
N
o
A
(
k
,
l
)
x
k
l
,
s.t.
∑
k
=
0
N
d
x
k
l
≤
1
∀
l
;
∑
l
=
0
N
o
x
k
l
≤
1
∀
k
;
x
k
l
=
0
if
c
k
l
>
ϵ
∀
k
,
l
.
\begin{aligned} & \arg \max _{x_{k l}} \sum_{k=0}^{N_d} \sum_{l=0}^{N_o} A(k, l) x_{k l}, \quad \text { s.t. } \sum_{k=0}^{N_d} x_{k l} \leq 1 \quad \forall l ; \\ & \sum_{l=0}^{N_o} x_{k l} \leq 1 \quad \forall k ; x_{k l}=0 \text { if } c_{k l}>\epsilon \quad \forall k, l . \end{aligned}
argxklmaxk=0∑Ndl=0∑NoA(k,l)xkl, s.t. k=0∑Ndxkl≤1∀l;l=0∑Noxkl≤1∀k;xkl=0 if ckl>ϵ∀k,l.
约束用于保证每个订单最多分配给一个司机,同样每个司机只能接到一个订单。在实际距离中距离如果过远则司机无法接单。
订单分配的二分图匹配问题
本文主要对比基于不同
A
(
k
,
l
)
A(k,l)
A(k,l)的三种调度策略。第一种是只考虑当前分配回报
A
(
1
)
(
k
,
l
)
=
α
1
r
k
−
α
2
c
k
l
,
A^{(1)}(k, l)=\alpha_1 r_k-\alpha_2 c_{k l},
A(1)(k,l)=α1rk−α2ckl,这两项分别是收益和距离cost。第二种策略是
A
(
2
)
(
k
,
l
)
=
α
1
r
k
−
α
2
c
k
l
+
α
3
{
η
Δ
t
l
k
V
1
(
s
l
k
′
)
−
V
1
(
s
l
)
}
,
A^{(2)}(k, l)=\alpha_1 r_k-\alpha_2 c_{k l}+\alpha_3\{\eta^{\Delta t_{lk}}V_1(s'_{lk})-V_1(s_{l})\},
A(2)(k,l)=α1rk−α2ckl+α3{ηΔtlkV1(slk′)−V1(sl)},增加的第三项用于考虑当前动作对未来收益的影响,其中
V
1
(
s
)
V_1(s)
V1(s)是从当前时刻开始到当天结束的期望收益。第三种策略是
A
(
3
)
(
k
,
l
)
=
A
(
2
)
(
k
,
l
)
+
α
4
{
η
Δ
t
l
k
V
2
(
s
l
k
′
)
−
V
2
(
s
l
)
}
,
A^{(3)}(k, l)=A^{(2)}(k, l)+\alpha_4\{\eta^{\Delta t_{lk}}V_2(s'_{lk})-V_2(s_{l})\},
A(3)(k,l)=A(2)(k,l)+α4{ηΔtlkV2(slk′)−V2(sl)},新增的一项用于衡量供需匹配程度,其中
V
2
(
s
)
=
ν
t
(
v
)
−
μ
~
t
(
v
)
V_2(s)=\nu_t(v)-\tilde{\mu}_t(v)
V2(s)=νt(v)−μ~t(v)是通过GEM定义中计算的,
V
2
(
.
)
V_2(.)
V2(.)项可以大大提高顾客被临近司机响应的概率。算例测试如下图所示,结果表明,基于GEM的第三种策略可以获得更高的司机收益和应答率。
不同调度策略对比
3. 策略评估
本文提出的GEM指标的一个重要应用是直接比较两个(或多个)乘车平台的调度策略。其关键思想是检测在同一平台环境下,两套竞争性策略的GEM之间是否存在显著差异。考虑到平台中需求和供应的联合分布,GEM越小,许多全局运营指标,如订单应答率、订单完成率和司机的工作时间就越好。与这些全局性的运营指标相比,GEM更直接地衡量了一个乘车平台的运营效率。
本文使用2018年12月3日至12月16日同一城市H的另一个供需数据集进行实验,以比较两种订单调度策略之间的有效性。两种策略在连续的半小时的时间间隔内交替执行。此外,本实验从前半小时的基线政策开始,每半小时改变一次政策,贯穿整个一天,并在另一天颠倒它们的顺序。本实验还包括一个A/A测试,通过使用从11月12日到11月25日获得的历史数据作为直接比较,将基线政策与自身进行比较。
每个30分钟的时间窗口内计算GEM的方法如下。每天总共有 M T = 48 M_T=48 MT=48个时间区间。为了得到每个时间间隔 τ \tau τ的GEM,本文通过使用归一化权重 o i t / ∑ t ∈ τ ∑ v i ∈ V o i t o_{it}/ \sum_{t∈\tau} \sum_{v_i∈V}o_{it} oit/∑t∈τ∑vi∈Voit汇总了30个GEM值,每个GEM值都是在1分钟的时间戳内计算的。
本文用
y
m
(
t
k
)
y_m(t_k)
ym(tk)表示汇总的GEM值,用
x
m
(
t
k
)
x_m(t_k)
xm(tk)表示2×1的预测因子向量,包括第m日的第k个时间间隔内所有司机的需求总数和总供应时间。
a
m
(
t
k
)
=
1
a_m(t_k)=1
am(tk)=1代表使用新策略,否则为
−
1
-1
−1。为了考察策略对GEM的边际影响,本文考虑以下回归模型:
y
m
(
t
k
)
=
β
0
(
t
k
)
+
β
1
(
t
k
)
T
{
x
m
(
t
k
)
−
x
ˉ
(
t
k
)
}
+
β
2
(
t
k
)
a
m
(
t
k
)
+
η
m
(
t
k
)
+
ε
m
(
t
k
)
y_m(t_k) = β_0(t_k) + β_1(t_k)^T \{x_m(t_k) − \bar{x}(t_k)\} + β_2(t_k)a_m(t_k) + η_m(_tk) + ε_m(t_k)
ym(tk)=β0(tk)+β1(tk)T{xm(tk)−xˉ(tk)}+β2(tk)am(tk)+ηm(tk)+εm(tk)
其中,
β
(
t
k
)
=
(
β
0
(
t
k
)
,
β
1
(
t
k
)
T
,
β
2
(
t
k
)
)
T
β(t_k) = (β_0(t_k), β_1(t_k)^T,β_2(t_k))^T
β(tk)=(β0(tk),β1(tk)T,β2(tk))T是
t
k
t_k
tk处的回归系数向量,
x
ˉ
(
t
k
)
\bar{x}(t_k)
xˉ(tk)是
k
=
1
,
.
.
.
,
M
T
k = 1,...,M_T
k=1,...,MT 的所有
x
i
(
t
k
)
x_i(t_k)
xi(tk)的样本平均值。此外,本文假设
η
m
=
(
η
m
(
t
1
)
,
.
.
.
,
η
m
(
t
M
T
)
)
T
η_m = (η_m(t_1),...,η_m(t_{M_T}))^T
ηm=(ηm(t1),...,ηm(tMT))T 和
ε
m
=
(
ε
m
(
t
1
)
,
.
.
.
,
ε
m
(
t
M
T
)
)
T
ε_m = (ε_m(t_1),...,ε_m(t_{M_T}))^T
εm=(εm(t1),...,εm(tMT))T是
M
T
×
1
M_T×1
MT×1的随机误差向量,遵循相互独立的多元高斯分布
N
(
0
,
Σ
η
)
N(0, Σ_η)
N(0,Ση)和
N
(
0
,
σ
ε
2
⋅
I
M
T
)
N(0, σ^2_ε\cdot I_{M_T})
N(0,σε2⋅IMT),其中
Σ
η
Σ_η
Ση是一个
M
T
×
M
T
M_T×M_T
MT×MT矩阵,
σ
ε
2
σ^2_ε
σε2是一个正标度。本文将检验以下的无效假设和备择假设:
H
0
:
∫
0
M
T
β
2
(
t
)
d
t
=
0
v
.
s
.
H
1
:
∫
0
M
T
β
2
(
t
)
d
t
≠
0
H_0:∫^{M_T}_0 β_2(t)dt=0 v.s. H1 : ∫^{M_T}_0 β_2(t)dt\neq0
H0:∫0MTβ2(t)dt=0v.s.H1:∫0MTβ2(t)dt=0
其中, ∫ 0 M T β 2 ( t ) d t ≈ ∑ k = 1 M T β 2 ( t k ) Δ t 0 ∫^{M_T}_0 β_2(t)dt≈ ∑^{M_T}_{k=1} β_2(t_k)Δt_0 ∫0MTβ2(t)dt≈∑k=1MTβ2(tk)Δt0表示每天的平均处理效应,其中 Δ t 0 Δt_0 Δt0是每个时间间隔的长度。本文提出了一个基于广义估计方程(GEE)的联合估计程序,迭代估计所有的未知参数,直到达到一个特定的收敛标准。随后,本文计算与每天的平均处理效应相关的 t t t检验统计量及其相应的单侧(或双侧) p p p值。
此外,本文在上述回归模型中考虑了三个全局运营指标,包括订单应答率、订单完成率和商品总值(GMV)作为 y m ( t k ) y_m(t_k) ym(tk)。本文拟合相应的三个回归模型,以研究新的调度策略是否在运营层面显著改善了骑行平台。
下表总结了A/A和A/B两种实验设计的所有回归分析结果。可以看到,在A/B实验设计中,当用新策略取代旧策略时,平均应答率、完成率和商品总值都存在明显的增加,因为所有与平均处理效应有关的
p
p
p值都小于
1
0
−
3
10^{-3}
10−3。新策略也可以显著降低GEM值(
p
p
p值小于0.05),这与我们的假设一致,即GEM可以充分量化供需关系,并随后影响所考察的平台指标。相反,在A/A实验设计中,所有四个指标在5%的显著性水平上没有显示出明显的处理效果。
A/A和A/B实验的相对改进百分比和平均处理效应的双侧p值
下图(A)显示了2018年12月3日A/B实验设计中总共48个30分钟时间窗口的GEM值。在早上7:00到8:00的时间段内,策略从控制组改为实验组策时,GEM值明显减少。下图(B)分别显示了调度策略在控制组和实验组下同一时间段内的
D
S
d
i
t
DSd_{it}
DSdit的热图。由绿色和紫色圆圈标记的三个选定区域的客户要求被附近区域的司机满足,导致更高的供需一致性,从而使GEM值更小。
(A)H市在30分钟范围内随机选择的一天(2018/12/08)的GEM值 (B)H市的30分钟范围内控制组(7:00 - 7:30)和实验组(7:30 - 8:00)的供需差热力图
四、小结
总的来说,针对动态供需网络均衡度评估以及网约车双边平台的司乘匹配等问题,本文主要贡献总结如下:
- 提出了基于图结构的动态供需网络均衡度量方法(GEM),可以看作是一种受限的广义Wasserstein距离,来量化动态需求和供应网络在加权图结构上的距离。相比于传统方法,GEM不仅考虑了图的拓扑结构,而且考虑了不均衡的供需总量。
- 通过改变每个节点的物理区域大小,可以得到多级GEM及其对应的最优调度策略。在最精细的空间尺度上,GEM可以计算包含许多局部细节的供需分配调度策略。相比之下,在相对粗糙的空间尺度上,它给出了最优传输函数的低频模式粗略表示。
- 在数值计算方法上,GEM的计算可以重构为一个标准的线性规划(LP)问题。作者从理论上研究了GEM的几个理论性质,包括计算GEM的LP算法的收敛性、GEM的期望等。并通过LP对偶大大降低了决策变量的复杂度。
- 在滴滴出行的数据集上应用GEM解决了应答率预测、订单派送策略制定等问题,并达到了很好的效果。
参考文献
[1] Zhou, F., Luo, S., Qie, X., Ye, J., & Zhu, H. (2021). Graph-based equilibrium metrics for dynamic supply–demand systems with applications to ride-sourcing platforms. Journal of the American Statistical Association, 116(536), 1688-1699.
[2] Alex Chin, & Tony Qin. (2023.02.25). Quantifying Efficiency in Ridesharing Marketplaces. Link: https://eng.lyft.com/quantifying-efficiency-in-ridesharing-marketplaces-affd53043db2