传输理论

传输理论（英語：tansport theory、transportation theory），又称为运输理论，是数学、经济学等学科中研究最优运输和资源配置的理论。该问题最早由法国数学家加斯帕尔·蒙日于1781年提出。^[1]

1920年代，A·N·托尔斯泰是最早运用数学方法研究传输问题的学者之一。1930年，他在苏联国家交通部编纂的《运输规划》第一卷中发表了题为《寻找太空货物运输的最小千公里方法》的论文。^[2][3]

第二次世界大战期间，苏联数学家、经济学家列昂尼德·坎托罗维奇在该领域取得了重要进展。^[4]因此，这一问题有时也被称为蒙日-坎托罗维奇运输问题（Monge–Kantorovich transportation problem）。 ^[5]该问题的线性规划形式也被称为希区柯克（英语：Frank Lauren Hitchcock）-库普曼斯运输问题。^[6]

背景

矿山与工厂

假设有${\displaystyle m}$个开采铁矿石的矿山以及${\displaystyle n}$个工厂使用这些矿山生产的铁矿石，且这些矿山和工厂构成欧几里得平面${\displaystyle \mathbb {R} ^{2}}$中两个不相交的子集${\displaystyle M}$和${\displaystyle F}$。同时假设存在一个成本函数${\displaystyle c:\mathbb {R} ^{2}\times \mathbb {R} ^{2}\to [0,\infty )}$， 其中${\displaystyle c(x,y)}$表示将一批矿石从${\displaystyle x}$运送到${\displaystyle y}$的成本。为简化起见，此处忽略运输所需的时间。我们还假定每个矿山只能供应一家工厂（不能拆分运输），并且每家工厂需要恰好一批货物才能运营（工厂不能半负荷或双倍负荷运转）。基于上述条件，一个传输计划可以看作是一个双射${\displaystyle T:M\to F}$ 。换句话说，每个矿井${\displaystyle m\in M}$仅供应一个目标工厂${\displaystyle T(m)\in F}$，而每个工厂也只由一个矿山供货。我们希望找到最优传输计划${\displaystyle T}$，使得总成本

${\displaystyle c(T):=\sum _{m\in M}c(m,T(m))}$

在所有${\displaystyle M}$到${\displaystyle F}$的传输计划中是最小的。该问题是传输问题的一个特例，可看成一个任务分配问题。更具体地说，它等价于在二分图中寻找最小权重匹配。

书籍移动：成本函数的重要性

下面这个简单的例子说明了成本函数在确定最优传输计划中的重要性。假设我们有${\displaystyle n}$本宽度相等的书摆放在书架上（可以看成具像化的实数线），形成连续一排书。我们希望将它们重新排列，在保持其连续性的同时将整体向右移动一本书的宽度。针对这个问题，有两个显而易见的最优传输候选方案：

1. 将所有${\displaystyle n}$本书全都向右移动一本书的宽度（许多小步）；
2. 将最左侧的书向右移动${\displaystyle n}$书本的宽度，其他书则保持不动（一大步）。

如果成本函数与欧几里得距离成正比（即 ${\displaystyle c(x,y)=\alpha \|x-y\|}$，其中${\displaystyle \alpha >0}$ ），那么这两种候选方案都是最优的。而如果我们选择与欧几里得距离的平方成比例的严格凸成本函数（即 ${\displaystyle c(x,y)=\alpha \|x-y\|^{2}}$，其中${\displaystyle \alpha >0}$ ），则“许多小步”的方案则是唯一的最优解。

需要注意的是，上述成本函数仅考虑书籍本身移动的水平距离，而没有考虑拿起每本书并将其移动到位的设备所行进的水平距离。如果考虑后者，那么在两种传输计划中，第二种方案对于欧几里得距离始终是最优的，而第一种方案则对于平方欧几里得距离是最优的（至少有三本书的情况下）。

希区柯克问题

以下传输问题的表述由弗兰克·劳伦·希区柯克提出：

假设有${\displaystyle m}$个供应源${\displaystyle x_{1},\ldots ,x_{m}}$为某一商品供货，且每个供应源${\displaystyle x_{i}}$处有${\displaystyle a(x_{i})}$个单位的供应量。同时有${\displaystyle n}$个需求点${\displaystyle y_{1},\ldots ,y_{n}}$需要该商品，每个需求点${\displaystyle y_{j}}$处有${\displaystyle b(y_{j})}$个单位的需求。若${\displaystyle c(x_{i},\ y_{j})}$表示从${\displaystyle x_{i}}$到${\displaystyle y_{j}}$的单位运输成本，任务是找到一个流量分配方案，在满足供应需求的同时最小化运输成本。这一物流问题由德尔伯特·雷·富尔克森提出^[7]，并在他与小莱斯特·伦道夫·福特合著的《网络流》（Flows in Networks）（1962年）一书中得到了阐述。^[8]

佳林·库普曼斯也为运输经济学与资源分配问题的表述作出了贡献。

问题的抽象表述

蒙日和坎托罗维奇形式

由于黎曼几何和测度论的发展，在现代或更加技术性的文献中传输问题的表述有所不同。不过，上述矿山与工厂的简单例子还是可以作为考虑抽象形式时一个有用的参照。此时我们可以考虑并非所有矿山和工厂都营业的情况，并允许一个矿山向多家工厂供货，而一家工厂也可以从多个矿山接受矿石。

假设${\displaystyle X}$和${\displaystyle Y}$为两个可分度量空间，使得${\displaystyle X}$（或者${\displaystyle Y}$ ) 上的每个概率测度都是拉东测度，并假设${\displaystyle c:X\times Y\to [0,\infty )}$是一个博雷尔可测函数。给定${\displaystyle X}$上的概率测度${\displaystyle \mu }$和${\displaystyle Y}$上的概率测度${\displaystyle \nu }$，最优传输问题的蒙日形式是指寻找一个传输映射${\displaystyle T:X\to Y}$，使得下式中的下确界成立：

${\displaystyle \inf \left\{\left.\int _{X}c(x,T(x))\,\mathrm {d} \mu (x)\right|T_{*}(\mu )=\nu \right\},}$

其中${\displaystyle T_{*}(\mu )}$由${\displaystyle T}$推进${\displaystyle \mu }$的向前推进算子。如果一个映射${\displaystyle T}$达到这个下确界，则该映射被称为“最优传输映射”。

最优传输问题的蒙日形式有其局限性，因为有时可能不存在满足${\displaystyle T_{*}(\mu )=\nu }$的映射${\displaystyle T}$。例如，当${\displaystyle \mu }$是狄拉克测度而${\displaystyle \nu }$不是时，就会出现这种情况。

此时可以通过采用最优传输问题的坎托罗维奇形式来克服这一局限，即寻找${\displaystyle X\times Y}$上的一个概率测度${\displaystyle \gamma }$，使得下式中的下确界成立：

${\displaystyle \inf \left\{\left.\int _{X\times Y}c(x,y)\,\mathrm {d} \gamma (x,y)\right|\gamma \in \Gamma (\mu ,\nu )\right\},}$

其中${\displaystyle \Gamma (\mu ,\nu )}$表示${\displaystyle X\times Y}$上所有概率测度的集合并满足边缘分布${\displaystyle \mu }$和${\displaystyle \nu }$。可以证明^[9]，当成本函数${\displaystyle c}$是下半连续并且${\displaystyle \Gamma (\mu ,\nu )}$是紧测度集合（拉东空间${\displaystyle X}$和${\displaystyle Y}$蕴含该条件）时，该问题总是存在最小值（另见沃瑟斯坦度量）。西古德·安格嫩特（英语：Sigurd Angenent）、史蒂文·哈克尔（Steven Haker）与艾伦·坦嫩鲍姆（英语：Allen Tannenbaum）提出了蒙日–坎托罗维奇问题解的梯度下降表述。^[10]

对偶形式

坎托罗维奇问题的最小值等于

${\displaystyle \sup \left(\int _{X}\varphi (x)\,\mathrm {d} \mu (x)+\int _{Y}\psi (y)\,\mathrm {d} \nu (y)\right),}$

其中上确界遍历所有有界且连续、满足

${\displaystyle \varphi (x)+\psi (y)\leq c(x,y)}$

的函数对。

经济学解释

将以上表述中的符号翻转更利于从经济学角度解释这一问题。假设${\displaystyle x\in X}$表示工人特征的向量，${\displaystyle y\in Y}$表示企业特征的向量，${\displaystyle \Phi (x,y)=-c(x,y)}$则表示工人${\displaystyle x}$与企业${\displaystyle y}$配对所创造的经济产出。令${\displaystyle u(x)=-\varphi (x)}$、${\displaystyle v(y)=-\psi (y)}$，蒙日–坎托罗维奇问题可以重新表述为：

${\displaystyle \sup \left\{\int _{X\times Y}\Phi (x,y)d\gamma (x,y),\gamma \in \Gamma (\mu ,\nu )\right\}}$

它的对偶形式为：

${\displaystyle \inf \left\{\int _{X}u(x)\,d\mu (x)+\int _{Y}v(y)\,d\nu (y):u(x)+v(y)\geq \Phi (x,y)\right\}}$

其中下确界遍历所有有界且连续的函数${\displaystyle u:X\to \mathbb {R} }$和${\displaystyle v:Y\to \mathbb {R} }$。如果对偶问题有解，则有：

${\displaystyle v(y)=\sup _{x}\left\{\Phi (x,y)-u(x)\right\}}$

可以将${\displaystyle u(x)}$解释为${\displaystyle x}$类型工人的均衡工资 ，并将${\displaystyle v(y)}$解释为${\displaystyle y}$类型企业的均衡利润。^[11]

问题求解

一维连续情形

最优传输矩阵

连续最优传输

对于${\displaystyle 1\leq p<\infty }$， 假设${\displaystyle {\mathcal {P}}_{p}(\mathbb {R} )}$表示${\displaystyle \mathbb {R} }$上所有${\displaystyle p}$阶矩有限的概率测度的集合。设${\displaystyle \mu ,\nu \in {\mathcal {P}}_{p}(\mathbb {R} )}$且${\displaystyle c(x,y)=h(x-y)}$，其中${\displaystyle h:\mathbb {R} \to [0,\infty )}$是一个凸函数。

1. 如果${\displaystyle \mu }$没有原子（英语：Atom (measure theory)），即${\displaystyle \mu }$的累积分布函数${\displaystyle F_{\mu }:\mathbb {R} \to [0,1]}$是一个连续函数，则${\displaystyle F_{\nu }^{-1}\circ F_{\mu }:\mathbb {R} \to \mathbb {R} }$是一个最优传输映射。如果${\displaystyle h}$是严格凸的，则它是唯一的最优映射。
2. 可以得到

${\displaystyle \min _{\gamma \in \Gamma (\mu ,\nu )}\int _{\mathbb {R} ^{2}}c(x,y)\,\mathrm {d} \gamma (x,y)=\int _{0}^{1}c\left(F_{\mu }^{-1}(s),F_{\nu }^{-1}(s)\right)\,\mathrm {d} s.}$

拉切夫（Rachev）与吕申多夫（Rüschendorf）于1998年给出了对此的证明。^[12]

离散情形与线性规划

在边缘分布${\displaystyle \mu }$和${\displaystyle \nu }$是离散的情形下，令${\displaystyle \mu _{x}}$和${\displaystyle \nu _{y}}$分别是分配给${\displaystyle x\in \mathbf {X} }$和${\displaystyle y\in \mathbf {Y} }$的概率质量，而${\displaystyle \gamma _{xy}}$则是${\displaystyle xy}$的分配概率。原始坎托罗维奇问题中的目标函数为

${\displaystyle \sum _{x\in \mathbf {X} ,y\in \mathbf {Y} }\gamma _{xy}c_{xy}}$

并满足约束条件

${\displaystyle \sum _{y\in \mathbf {Y} }\gamma _{xy}=\mu _{x},\forall x\in \mathbf {X} }$

${\displaystyle \sum _{x\in \mathbf {X} }\gamma _{xy}=\nu _{y},\forall y\in \mathbf {Y} .}$

为了将这一问题作为线性规划问题处理，我们需要将矩阵${\displaystyle \gamma _{xy}}$向量化，可以通过堆叠其行或列来完成，我们用${\displaystyle \operatorname {vec} }$来表示这一操作。在列主序（英语：Row- and column-major order）的情况下，上述约束条件可改写为

${\displaystyle \left(1_{1\times |\mathbf {Y} |}\otimes I_{|\mathbf {X} |}\right)\operatorname {vec} (\gamma )=\mu }$

${\displaystyle \left(I_{|\mathbf {Y} \|}\otimes 1_{1\times |\mathbf {X} |}\right)\operatorname {vec} (\gamma )=\nu }$

其中${\displaystyle \otimes }$是克罗内克积，${\displaystyle 1_{n\times m}}$是一个大小为${\displaystyle n\times m}$、所有元素为1的矩阵，而${\displaystyle I_{n}}$是大小为${\displaystyle n}$的单位矩阵。设${\displaystyle z=\operatorname {vec} (\gamma )}$，该问题的线性规划形式为

${\displaystyle {\begin{aligned}&{\text{Minimize }}&&\operatorname {vec} (c)^{\top }z\\[4pt]&{\text{subject to:}}&&z\geq 0,\\[4pt]&&&{\begin{pmatrix}1_{1\times |\mathbf {Y} |}\otimes I_{|\mathbf {X} |}\\I_{|\mathbf {Y} |}\otimes 1_{1\times |\mathbf {X} |}\end{pmatrix}}z={\binom {\mu }{\nu }}\end{aligned}}}$

这一问题可以很容易地通过大规模线性规划求解器计算。^[13]

半离散情形

在半离散情况下，令${\displaystyle X=Y=\mathbb {R} ^{d}}$，且${\displaystyle \mu }$是${\displaystyle \mathbb {R} ^{d}}$上的连续分布，${\displaystyle \nu =\sum _{j=1}^{J}\nu _{j}\delta _{y_{i}}}$则是分配概率质量${\displaystyle \nu _{j}}$到${\displaystyle y_{j}\in \mathbb {R} ^{d}}$的离散分布。此时，坎托罗维奇问题的原始形式为：^[14]

${\displaystyle \inf \left\{\int _{X}\sum _{j=1}^{J}c(x,y_{j})\,d\gamma _{j}(x),\gamma \in \Gamma (\mu ,\nu )\right\}}$

其中${\displaystyle \gamma \in \Gamma (\mu ,\nu )}$满足${\displaystyle \int _{X}d\gamma _{j}(x)=\nu _{j}}$和${\displaystyle \sum _{j}d\gamma _{j}(x)=d\mu (x)}$ 。

而其对偶形式则为

${\displaystyle \sup \left\{\int _{X}\varphi (x)d\mu (x)+\sum _{j=1}^{J}\psi _{j}\nu _{j}:\psi _{j}+\varphi (x)\leq c(x,y_{j})\right\}}$

还可以写为：

${\displaystyle \sup _{\psi \in \mathbb {R} ^{J}}\left\{\int _{X}\inf _{j}\left\{c(x,y_{j})-\psi _{j}\right\}d\mu (x)+\sum _{j=1}^{J}\psi _{j}\nu _{j}\right\}}$

这是一个有限维凸优化问题，可以通过梯度下降等方法求解。

当${\displaystyle c(x,y)=|x-y|^{2}/2}$时，可以证明分配给特定${\displaystyle j}$的${\displaystyle x\in \mathbf {X} }$集合是一个凸多面体，而得到的配置称为幂图（英语：Power diagram）。^[15]

二次正态情形

假设一个特殊情形${\displaystyle \mu ={\mathcal {N}}(0,\Sigma _{X})}$、${\displaystyle \nu ={\mathcal {N}}(0,\Sigma _{Y})}$ 且${\displaystyle c(x,y)=|y-Ax|^{2}/2}$，其中${\displaystyle A}$是可逆矩阵。此时有

${\displaystyle \varphi (x)=-x^{\top }\Sigma _{X}^{-1/2}\left(\Sigma _{X}^{1/2}A^{\top }\Sigma _{Y}A\Sigma _{X}^{1/2}\right)^{1/2}\Sigma _{X}^{-1/2}x/2}$

${\displaystyle \psi (y)=-y^{\top }A\Sigma _{X}^{1/2}\left(\Sigma _{X}^{1/2}A^{\top }\Sigma _{Y}A\Sigma _{X}^{1/2}\right)^{-1/2}\Sigma _{X}^{1/2}Ay/2}$

${\displaystyle T(x)=(A^{\top })^{-1}\Sigma _{X}^{-1/2}\left(\Sigma _{X}^{1/2}A^{\top }\Sigma _{Y}A\Sigma _{X}^{1/2}\right)^{1/2}\Sigma _{X}^{-1/2}x}$

加利雄（Galichon）于2016年证明了该情况下的解。^[16]

可分希尔伯特空间

令${\displaystyle X}$是一个可分希尔伯特空间，定义${\displaystyle {\mathcal {P}}_{p}(X)}$为${\displaystyle X}$上所有${\displaystyle p}$阶矩有限的概率测度的集合，${\displaystyle {\mathcal {P}}_{p}^{r}(X)}$则表示其中高斯正则的测度集合，即如果${\displaystyle g}$是${\displaystyle X}$上任何严格正的高斯测度（英语：Gaussian measure）且满足${\displaystyle g(N)=0}$ ，则${\displaystyle \mu (N)=0}$也成立。

假设${\displaystyle \mu \in {\mathcal {P}}_{p}^{r}(X)}$、${\displaystyle \nu \in {\mathcal {P}}_{p}(X)}$，并且${\displaystyle c(x,y)=|x-y|^{p}/p}$，其中${\displaystyle p\in (1,\infty ),p^{-1}+q^{-1}=1}$ 。则坎托罗维奇问题存在一个唯一解${\displaystyle \kappa }$ ，并且该解对应一个最优传输映射：即存在一个博雷尔映射${\displaystyle r\in L^{p}(X,\mu ;X)}$，使得

${\displaystyle \kappa =(\mathrm {id} _{X}\times r)_{*}(\mu )\in \Gamma (\mu ,\nu ).}$

此外，如果${\displaystyle \nu }$具有有界支撑，那么对于${\displaystyle \mu }$-几乎所有的${\displaystyle x\in X}$，存在局部利普希茨、 ${\displaystyle c}$-凹和最大坎托罗维奇势${\displaystyle \varphi }$，使得

${\displaystyle r(x)=x-|\nabla \varphi (x)|^{q-2}\,\nabla \varphi (x)}$

其中${\displaystyle \nabla \varphi }$表示${\displaystyle \varphi }$的加托导数。

熵正则化

考虑上述离散问题的一个变体：在原始问题的目标函数中添加一个熵正则化项

${\displaystyle {\begin{aligned}&{\text{Minimize }}\sum _{x\in \mathbf {X} ,y\in \mathbf {Y} }\gamma _{xy}c_{xy}+\varepsilon \gamma _{xy}\ln \gamma _{xy}\\[4pt]&{\text{subject to: }}\\[4pt]&\gamma \geq 0\\[4pt]&\sum _{y\in \mathbf {Y} }\gamma _{xy}=\mu _{x},\forall x\in \mathbf {X} \\[4pt]&\sum _{x\in \mathbf {X} }\gamma _{xy}=\nu _{y},\forall y\in \mathbf {Y} \end{aligned}}}$

相应的对偶问题为

${\displaystyle \max _{\varphi ,\psi }\sum _{x\in \mathbf {X} }\varphi _{x}\mu _{x}+\sum _{y\in \mathbf {Y} }\psi _{y}v_{y}-\varepsilon \sum _{x\in \mathbf {X} ,y\in \mathbf {Y} }\exp \left({\frac {\varphi _{x}+\psi _{y}-c_{xy}}{\varepsilon }}\right)}$

相较于不含正则化项的问题，原先对偶问题中的硬约束（ ${\displaystyle \varphi _{x}+\psi _{y}-c_{xy}\geq 0}$ ）被替换为了软约束，即惩罚项${\displaystyle \varepsilon \exp \left((\varphi _{x}+\psi _{y}-c_{xy})/\varepsilon \right)}$。对偶问题的最优条件可以表示为

式5.1: ${\displaystyle \mu _{x}=\sum _{y\in \mathbf {Y} }\exp \left({\frac {\varphi _{x}+\psi _{y}-c_{xy}}{\varepsilon }}\right)~\forall x\in \mathbf {X} }$

式5.2: ${\displaystyle \nu _{y}=\sum _{x\in \mathbf {X} }\exp \left({\frac {\varphi _{x}+\psi _{y}-c_{xy}}{\varepsilon }}\right)~\forall y\in \mathbf {Y} }$

令${\displaystyle A}$为${\displaystyle |\mathbf {X} |\times |\mathbf {Y} |}$的矩阵，其中元素${\displaystyle A_{xy}=\exp \left(-c_{xy}/\varepsilon \right)}$。此时对偶问题的求解等价于寻找两个对角正矩阵${\displaystyle D_{1}}$和${\displaystyle D_{2}}$，它们的大小分别为${\displaystyle |\mathbf {X} |}$和${\displaystyle |\mathbf {Y} |}$ ，使得${\displaystyle D_{1}AD_{2}1_{|\mathbf {Y} |}=\mu }$、${\displaystyle (D_{1}AD_{2})^{\top }1_{|\mathbf {X} |}=\nu }$。矩阵${\displaystyle D_{1}}$和${\displaystyle D_{2}}$的存在性是辛克宏定理的推广，可以使用辛克宏-诺普算法进行求解。^[17]该算法通过迭代求解式5.1中的${\displaystyle \varphi _{x}}$和式5.2中的${\displaystyle \psi _{y}}$实现。因此，辛克宏-诺普算法相当于对偶正则问题的坐标下降法。

应用

蒙日-坎托罗维奇运输问题已广泛运用于许多领域，例如：

• 图像配准与变形^[18]
• 反射面设计^[19]
• 从阴影照相与质子射线成像中提取信息^[20]
• 地震层析成像和反射地震学^[21]
• 涉及总替代性质的经济模型

参见

• 沃瑟斯坦度量
• 匈牙利算法
• 交通规划
• 推土机距离
• 齐次蒙日-安培方程