阮一鳴

摘?要：最優(yōu)傳輸問(wèn)題是尋求相對(duì)于給定的代價(jià)函數(shù)，把一種分布轉(zhuǎn)化為另一種分布的最有效的方式，其具有較為深遠(yuǎn)的價(jià)值，其中基于點(diǎn)云的最優(yōu)傳輸問(wèn)題更是得到廣泛的關(guān)注。針對(duì)高維的點(diǎn)云分布的最優(yōu)傳輸問(wèn)題，利用了分片最優(yōu)傳輸?shù)睦碚摷捌鋵?duì)應(yīng)的優(yōu)化模型，提出一個(gè)改進(jìn)梯度迭代的算法，以二維、三維點(diǎn)云分布為例進(jìn)行數(shù)值實(shí)驗(yàn)，并將其與經(jīng)典的方法進(jìn)行比較，驗(yàn)證模型與算法的有效性。實(shí)驗(yàn)結(jié)果表明，本文模型及其算法在標(biāo)準(zhǔn)差、峰值信噪比、平均梯度等指標(biāo)均有更高的性能。

關(guān)鍵詞：最優(yōu)傳輸;Wasserstein距離;點(diǎn)云分布;隨機(jī)梯度迭代

中圖分類號(hào)：O212.4??????文獻(xiàn)標(biāo)識(shí)碼：A

Research?on?Optimal?Transport?Problem?Based?

on?High?Dimensional?Point?Cloud?Distribution

RUAN?Yiming

（School?of?Science，?Hohai?University，?Nanjing，Jiangsu?211100，China）

Abstract：The?problem?of?optimal?transport?is?to?find?the?most?effective?way?to?transform?one?distribution?into?another?relative?to?a?given?cost?function.?It?has?farreaching?value，?among?which?the?problem?of?optimal?transport?based?on?point?cloud?has?gained?wide attention.?Aiming?at?the?problem?of?optimal?transport?of?highdimensional?point?cloud?distribution，?using?the?theory?of?piecewise?optimal?transport?and?its?corresponding?optimization?model，?this?paper?proposes?an?improved?gradient?iterative?algorithm.?Taking?twodimensional?and?threedimensional?point?cloud?distribution?as?examples，?numerical?experiments?are?carried?out?and?compared?with?classical?methods?to?verify?the?effectiveness?of?the?model?and?algorithm.?Experimental?results?show?that?the?model?and?its?algorithm?have?higher?performance?in?Standard?Deviation，?Peak?Signal?to?Noise?Ratio，?Average?Gradient?and?other?indicators.

Key?words：optimal?transport;?Wasserstein?distance;?point?cloud?distribution;?random?gradient?iteration

最優(yōu)傳輸問(wèn)題是尋求相對(duì)于給定的代價(jià)函數(shù)，把一種分布轉(zhuǎn)化為另一種分布的最有效的方式，其具有較為深遠(yuǎn)的價(jià)值，其中基于點(diǎn)云的最優(yōu)傳輸問(wèn)題更是得到廣泛的關(guān)注，例如運(yùn)輸問(wèn)題、顏色遷移等，其中最為常見(jiàn)的便是二維、三維點(diǎn)云分布的應(yīng)用問(wèn)題，在計(jì)算機(jī)圖形學(xué)、計(jì)算機(jī)視覺(jué)、醫(yī)學(xué)圖像處理以及深度學(xué)習(xí)等方面都有顯著的效果[1-4]。

最優(yōu)傳輸問(wèn)題的核心是Wasserstein距離，近年來(lái)，關(guān)于Wasserstein距離的計(jì)算得到了蓬勃的發(fā)展，但關(guān)于Wasserstein距離的計(jì)算尤其是高維點(diǎn)云分布仍然具有一定的局限性，目前可以進(jìn)行顯式計(jì)算的只有一維分布或者高斯分布的情形[4]，因此如何準(zhǔn)確快速求解高維的點(diǎn)云分布的最優(yōu)傳輸問(wèn)題是最優(yōu)傳輸領(lǐng)域的一個(gè)難點(diǎn)。關(guān)于一維情形點(diǎn)云分布的最優(yōu)傳輸問(wèn)題，已經(jīng)涌現(xiàn)了許多算法及理論，例如匈牙利算法[5]、拍賣算法[6]、經(jīng)典的線性規(guī)劃算法[7]、熵正則化算法[8]等，但對(duì)于較為常見(jiàn)的高維的點(diǎn)云分布的最優(yōu)傳輸問(wèn)題的研究仍然較少，Rabin[9]等首次提出了用分片Wasserstein距離來(lái)求解高維的概率分布及Wasserstein重心，結(jié)果也表明了算法的有效性，但其收斂性能有時(shí)仍然較不穩(wěn)定;Paulin[10]等利用分片最優(yōu)傳輸重構(gòu)相應(yīng)的點(diǎn)云樣本任意維樣本分布，以提高蒙特卡羅積分的精度;Bonneel[11]等將分片Wasserstein距離與拉東變換應(yīng)用于求解任意測(cè)度的Wasserstein重心;?Xiao[12]等探討了在三維顏色空間中的點(diǎn)云分布的顏色傳輸，其主要考慮對(duì)三維像素點(diǎn)云進(jìn)行顏色遷移的計(jì)算;Reinhard[13]等探討了三維點(diǎn)云分布在不同色彩空間的傳輸過(guò)程，其所考慮在經(jīng)典的顏色空間中進(jìn)行顏色的遷移探究，雖然輸入輸出是RGB像素值，但是中間進(jìn)行算法處理時(shí)很少直接更改RGB值，因此對(duì)結(jié)果仍然具有一定的影響。文獻(xiàn)[14]、[15]則介紹主流的梯度算法及其改進(jìn)算法的基本理論，給出相應(yīng)的介紹及其理論的分析證明。

本文則針對(duì)高維情形下點(diǎn)云分布的最優(yōu)傳輸問(wèn)題進(jìn)行探究，具體以二維、三維點(diǎn)云分布為例，受文獻(xiàn)[9]的啟發(fā)，同時(shí)基于文獻(xiàn)[14]、[15]的基本理論，利用分片Wasserstein距離的相應(yīng)理論，將其投影成一維分布，其積分形式選用相應(yīng)的離散和進(jìn)行近似，考慮利用隨機(jī)梯度下降迭代更新相應(yīng)的分布，同時(shí)，鑒于傳統(tǒng)的隨機(jī)梯度下降法容易收斂到局部最優(yōu)[16]，并且容易被困在鞍點(diǎn)，收斂性能較差，故本文對(duì)傳統(tǒng)的隨機(jī)梯度下降法進(jìn)行了改進(jìn)，在其中加入Nesterov?Momentum項(xiàng)[14-15]以保證更新的時(shí)候在一定程度上保留之前更新的方向，同時(shí)利用當(dāng)前部分的梯度微調(diào)最終的更新方向以在一定程度上增加穩(wěn)定性，并且還有一定擺脫局部最優(yōu)的能力，最后利用數(shù)值實(shí)驗(yàn)驗(yàn)證本文算法及其模型的有效性。

1?分片Wasserstein距離理論

1.1?一維最優(yōu)傳輸

考慮如下的兩個(gè)點(diǎn)云分布μ，v，其所對(duì)應(yīng)的數(shù)據(jù)X=（x1，x2，…，xn），Y=（y1，y2，…，yn），假定其滿足下式：

即每個(gè)點(diǎn)云的概率權(quán)重相等，不失一般性。假定每個(gè)點(diǎn)云的具體位置距離已預(yù)先做好排序即x1≤x2≤…≤xn，y1≤…≤yn，即可得到如下簡(jiǎn)單的形式，點(diǎn)云分布所對(duì)應(yīng)的s維Wasserstein距離為：

1.2?分片Wasserstein距離

在實(shí)踐應(yīng)用中，對(duì)Wasserstein距離的計(jì)算仍然有較高的要求，且Ws在優(yōu)化包含Wasserstein距離的函數(shù)點(diǎn)云的問(wèn)題中較難處理，鑒于這些原因，我們現(xiàn)在考慮分布之間的替代度量即基于一維投影之間的運(yùn)輸成本，同時(shí)為了將高維的點(diǎn)云分布轉(zhuǎn)化為一維點(diǎn)云分布從而更好地計(jì)算Wasserstein距離，這里我們引入分片Wasserstein距離的概念，其基本的思想是針對(duì)兩個(gè)高維分布隨機(jī)找不同的方向進(jìn)行投影，然后計(jì)算兩個(gè)一維的投影點(diǎn)云的Wasserstein距離，最后將從不同投影方向所得到的一維分布的Wasserstein距離進(jìn)行求和，即得到了相應(yīng)的分片Wasserstein距離。首先給出分片Wasserstein距離的定義：

其中θ表示一個(gè)模長(zhǎng)為1的線性投影，?它將原本高維（d維）的X或者Y投影到了一維空間。根據(jù)上述定義敘述，可以將分片Wasserstein距離用一維的最優(yōu)分配來(lái)表示：

其中SW2指代上述的度量為分片Wasserstein距離。?

2?理論與算法

2.1?基于梯度優(yōu)化的算法

在許多機(jī)器學(xué)習(xí)任務(wù)中，為了求解方便可將機(jī)器學(xué)習(xí)算法模型簡(jiǎn)化為經(jīng)驗(yàn)風(fēng)險(xiǎn)，即求解最優(yōu)化模型：

是第l個(gè)樣本關(guān)于參數(shù)x的損失函數(shù)。當(dāng)樣本容量足夠大時(shí)，經(jīng)驗(yàn)風(fēng)險(xiǎn)最小化能保證有良好的學(xué)習(xí)效果，因此在實(shí)際應(yīng)用廣泛中被采用[14-15]。

2.1.1?隨機(jī)梯度下降算法

隨機(jī)梯度下降法（SGD）算法及其眾多變種算法都是機(jī)器學(xué)習(xí)中應(yīng)用最廣泛的優(yōu)化算法，是求解大規(guī)模學(xué)習(xí)問(wèn)題的一種有效方法。SGD算法以基本的梯度下降形式更新迭代，在每輪更新參數(shù)時(shí)，僅選擇一個(gè)或幾個(gè)樣本計(jì)算其梯度，并以此梯度作為全局梯度的估計(jì)值，極大降低了算法復(fù)雜度。其參數(shù)更新公式為：

2.1.2?標(biāo)準(zhǔn)動(dòng)量法

隨機(jī)梯度下降法可廣泛應(yīng)用于解諸多優(yōu)化問(wèn)題，但其依然存在很多缺陷，例如該算法有時(shí)達(dá)到最優(yōu)解的速度很慢。對(duì)上述問(wèn)題進(jìn)行優(yōu)化，標(biāo)準(zhǔn)動(dòng)量方法增加一個(gè)動(dòng)量項(xiàng)的方式使得目標(biāo)函數(shù)加速收斂，通過(guò)累計(jì)參數(shù)之前梯度指數(shù)衰減的平均值，防止梯度計(jì)算時(shí)過(guò)度震蕩，使目標(biāo)函數(shù)更快達(dá)到局部最優(yōu)。其算法更新迭代規(guī)則如下：

2.1.3?Nesterov?Momentum動(dòng)量法

帶有Nesterov?Momentum項(xiàng)的隨機(jī)梯度下降算法是一階優(yōu)化方法，用于提高常規(guī)下降的穩(wěn)定性和收斂性。Nesterov動(dòng)量是Momentum的變種，即在計(jì)算參數(shù)梯度之前，前瞻一步，超前一個(gè)動(dòng)量單位處xl+ηvl，可理解為往Momentum算法中加入一個(gè)校正因子。迭代過(guò)程中更新規(guī)則如下所示：

其中，xl表示模型參數(shù)，νl表示動(dòng)量，η∈[0，1]表示動(dòng)量系數(shù)，τl表示當(dāng)前迭代的學(xué)習(xí)率，f（x）是目標(biāo)函數(shù)。

與標(biāo)準(zhǔn)動(dòng)量項(xiàng)不同的是，Nesterov動(dòng)量是對(duì)某個(gè)位置的參數(shù)向量xl更新，取決于當(dāng)前參數(shù)位置ηνl的最后一次動(dòng)量更新。Nesterov動(dòng)量相比標(biāo)準(zhǔn)動(dòng)量項(xiàng)可以更有效地修正每次迭代較大且不適當(dāng)?shù)乃俣?，這使得它比標(biāo)準(zhǔn)動(dòng)量法更快，并且可以進(jìn)一步降低損失函數(shù)的誤差率。

2.2?算法基本理論

對(duì)于高維的情形，要得到相應(yīng)的正交投影，對(duì)一些常用的算法例如線性規(guī)劃算法、熵正則化算法來(lái)說(shuō)則較難處理。受文獻(xiàn)[9]、[14]、[15]的啟發(fā)，本文則考慮利用隨機(jī)梯度下降法來(lái)計(jì)算分布之間的分片Wasserstein投影，即隨機(jī)尋找相應(yīng)的投影方向?qū)⒏呔S概率分布投影下來(lái)，將其轉(zhuǎn)化為一維點(diǎn)云分布，然后計(jì)算兩個(gè)一維情況下概率分布的Wasserstein距離。

首先計(jì)算關(guān)于μY局部最小的分片Wasserstein距離：

E（t）=W2（μt，μY）（10）

利用計(jì)算結(jié)果逼近，其中上式與分布X足夠逼近且E是Cl連續(xù)函數(shù)。為得到逼近X的點(diǎn)云分布，數(shù)據(jù)考慮對(duì)E使用梯度下降法，從X開(kāi)始進(jìn)行迭代：

X（0）=X，X（l+1）=X（l）-τl

EΘ（X（l）+η×v（l））

X（l+1）=X（l）+v（l+1）（17）

其中η表示衰減因子，表示對(duì)先前方向的依賴程度，通常取值為0.9。且根據(jù)式（15），則式（17）迭代更新規(guī)則如下所示：

v（l+1）=η×v（l）-τl×（X（l）+η×v（l）-（l））

X（l+1）=X（l）+v（l+1）（18）?

綜上，可得本文算法。

算法1?基于NMSGD的分片Wasserstein距離的計(jì)算算法

Require：點(diǎn)云分布數(shù)據(jù)X，Y，學(xué)習(xí)率τl，衰減因子η

Require：初始點(diǎn)云分布X（0）=X，初始動(dòng)量v（0）

While?不滿足停止準(zhǔn)則?do

1：對(duì)X進(jìn)行隨機(jī)梯度法迭代，?初值迭代X（0）=X

2：計(jì)算上述投影算子P

3：?更新v（l+1）∶v（l+1）←η×v（l）-τl×EΘ（X（l）+η×v（l））

4：更新X（l+1）∶X（l+1）←X（l）+v（l+1）

end?while

Output：Xl+1

3?數(shù)值實(shí)驗(yàn)

3.1?二維點(diǎn)云分布的最優(yōu)傳輸問(wèn)題

本文考慮對(duì)二維與三維點(diǎn)云分布數(shù)據(jù)進(jìn)行數(shù)值實(shí)驗(yàn)。首先隨機(jī)給出兩組二維點(diǎn)云數(shù)據(jù)即d=2，點(diǎn)云數(shù)據(jù)個(gè)數(shù)為200，點(diǎn)云分布見(jiàn)圖1。根據(jù)上述算法，對(duì)上述兩組點(diǎn)云進(jìn)行迭代計(jì)算P直至滿足收斂。其中niter=1，3，10，100傳輸過(guò)程與最終的結(jié)果見(jiàn)圖2與圖3。

3.2?三維點(diǎn)云分布的最優(yōu)傳輸問(wèn)題

3.2.1?數(shù)值實(shí)驗(yàn)

對(duì)于三維情形的離散分布的最優(yōu)傳輸問(wèn)題即d=3，本文考慮對(duì)應(yīng)的三維彩色圖像，則其像素對(duì)應(yīng)的離散分布即是三維的點(diǎn)云分布。給出兩組兩幅128*128的三維彩色圖像見(jiàn)圖4。

故對(duì)于上述三維點(diǎn)云分布數(shù)據(jù)，基于上述算法，利用隨機(jī)梯度下降法最小化分片Wasserstein距離計(jì)算迭代次數(shù)niter=1，3，10，100的目標(biāo)傳輸圖像，具體迭代傳輸過(guò)程與結(jié)果見(jiàn)圖5。

為驗(yàn)證本文算法，接下來(lái)將對(duì)其與圖像遷移領(lǐng)域中幾種有代表性的算法進(jìn)行比較。文獻(xiàn)[9]中，Rabin等人提出隨機(jī)梯度算法，并將其應(yīng)用于Wasserstein?重心，本節(jié)則將其應(yīng)用于三維點(diǎn)云分布的最優(yōu)傳輸問(wèn)題;文獻(xiàn)[12]中Xiao等提出了一種在三維顏色空間中對(duì)像素進(jìn)行傳輸?shù)乃惴?文獻(xiàn)[13]中?Reinhard等人根據(jù)lαβ顏色空間中各通道互相不關(guān)聯(lián)的特點(diǎn)，提出了一組適用于各顏色分量的色彩遷移算法。文獻(xiàn)[12]、[13]所考慮在經(jīng)典的顏色空間中進(jìn)行顏色的遷移探究，本文則將其考慮在Wasserstein空間進(jìn)行應(yīng)用，同時(shí)基于文獻(xiàn)[9]進(jìn)行改進(jìn)。以上方法及其所得的最終結(jié)果見(jiàn)圖6。

3.2.2?評(píng)價(jià)指標(biāo)

接下來(lái)利用標(biāo)準(zhǔn)差、峰值信噪比和平均梯度三個(gè)評(píng)價(jià)指標(biāo)對(duì)結(jié)果圖像進(jìn)行評(píng)價(jià)。

（1）標(biāo)準(zhǔn)差

圖像的均值表示圖像的亮度，均值越大圖像越亮，其中均值計(jì)算公式為：

u=1W×H∑Wi=1∑Hj=1f（i，j）（19）

而圖像標(biāo)準(zhǔn)差反映了圖像像素值與均值的離散程度，標(biāo)準(zhǔn)差越大說(shuō)明圖像的對(duì)比度越高.標(biāo)準(zhǔn)差計(jì)算公式為：

δ=1W×H∑Wi=1∑Hj=1[f（i，j）-u]2（20）

其中，W×H表示圖像中像素點(diǎn)的個(gè)數(shù)，f（i，j）為（i，j）點(diǎn)處的灰度值。

（2）峰值信噪比

峰值信噪比（PSNR）是用于評(píng)價(jià)圖像失真度的客觀指標(biāo)，峰值信噪比數(shù)值越大表示失真度越小，圖像質(zhì)量越好，其主要通過(guò)最大可能像素值fmax和均方誤差（MSE）來(lái)定義，其中均方誤差和峰值信噪比的計(jì)算公式為：

MSE=∑Wi=1∑Hj=1-f0（i，j）]2W×H（21）

PSNR=10lg10f2maxMSE?（22）

其中，W×H表示圖像中像素點(diǎn)的個(gè)數(shù)，f（x，y）代表目標(biāo)圖像的像素點(diǎn)，f0（x，y）代表結(jié)果圖中的像素點(diǎn)。

（3）平均梯度（Average?gradient）

平均梯度反映圖像的清晰程度，平均梯度越大，圖像越清晰，圖像層次越豐富，圖像的平均梯度（AG）計(jì)算公式為：

AG=1W×H∑Wi=1∑Hj=1fx2+fy22（23）

其中，W×H表示圖像中像素點(diǎn)的個(gè)數(shù)，fx表示水平方向的梯度，fy表示垂直方向的梯度。不同方法在兩組實(shí)驗(yàn)中的不同指標(biāo)結(jié)果及其比較見(jiàn)表1與表2。

從兩組實(shí)驗(yàn)可以看出本文算法評(píng)價(jià)指標(biāo)結(jié)果。MSE略低于文獻(xiàn)[9]方法、PSNR略低于文獻(xiàn)[13]方法外，其余指標(biāo)結(jié)果均優(yōu)于其他算法。綜上，本文算法可較好地應(yīng)用于求解點(diǎn)云分布的最優(yōu)傳輸問(wèn)題。

4?結(jié)?論

研究了高維點(diǎn)云分布數(shù)據(jù)的最優(yōu)傳輸問(wèn)題及其相關(guān)應(yīng)用。基于一維點(diǎn)云分布的最優(yōu)傳輸問(wèn)題的基本理論，給出高維的點(diǎn)云分布的傳輸問(wèn)題的計(jì)算理論及算法，以經(jīng)典的二維、三維情形為例進(jìn)行相應(yīng)的數(shù)值實(shí)驗(yàn)，對(duì)于三維的情形也將其與經(jīng)典的方法進(jìn)行比較，驗(yàn)證模型與算法的有效性。

參考文獻(xiàn)

[1]?PEYRЁ?G，?CUTURI?M.?Computational?optimal?transport：?with?applications?to?data?science?[J].?Foundations?and?Trends?in?Machine?Learning，?2019，?11（5-6）：?355-607.?

[2]?PANARETOS?V?M，?ZEMEL?Y.?Statistical?aspects?of??Wasserstein?distances[J].?Annual?review?of?statistics?and?its?application，?2019，?6（1）：?405-431.

[3]?馬麗濤，?邊?偉.?最優(yōu)傳輸理論及其在圖像處理中的應(yīng)用[J].?運(yùn)籌學(xué)學(xué)報(bào)，?2019，?23（3）：?109-125.

[4]?VILLANI?C.?Optimal?transport：?old?and?new[M].?Berlin：?Springer，?2009.?

[5]?KUHN?H?W.?The?Hungarian?method?for?the?assignment??problem[J]?.Naval?Research?Logistics?Quarterly，?1995，?2（1-2）：83-97.

[6]?BERTSEKAS?D?P.?The?auction?algorithm：?a?distributed?relaxation?method?for?assignment?problem[J].?Annals?of?Operations?Research，1998，?14（1）：?105-123.?

[7]?AURICCHIO?G，?BASSETTI?F，?GUALANDI?S，?et?al.?Computing?Wasserstein?barycenters?via?linear?programming[C].?International?Conference?on?Integration?of?Constraint?Programming，?Artificial?Intelligence，?and?Operations?Research.?Springer，?Cham，?2019：?355-363.

[8]?CUTURI?M.Sinkhorn?distances：?lightspeed?computation?of?optimal?transport[J].?Advances?in?Neural?Information?Processing?Systems，?2013，?26：?2292-2300.

[9]?RABIN?J，?PEYRЁ?G，?DELON?J，?et?al.Wasserstein?barycenter?and?its?application?to?texture?mixing[C]//In?International?Conference?on?Scale?Space?and?Variational?Methods?in?Computer?Vision.?Springer，?Berlin，?Heidelberg，?2011，?435-446.?

[10]PAULIN?L，?BONNEL?N，?COEURJOLLY?D，?et?al.?Sliced?optimal?transport?sampling[J].?ACM Trans.?Graph，?2020，?39（4）：?99-115.

[11]BONNEEL?N，?RABIN?J，?PEYRЁ?G，?et?al.?Sliced?and?radon?wasserstein?barycenters?of?measures[J].?Journal?of?Mathematical?Imaging?and?Vision，?2015，?51（1）：?22-45.

[12]XIAO?X，?MA?L.?Color?transfer?in?correlated?color?space[C]//Proceedings?of?the?2006?ACM?International?Conference?on?Virtual?Reality?Continuum?and?Its?Applications.?2006：?305-309.

[13]REINHARD?E，?ADHIKHMIN?M，?GOOCH?B，?et?al.?Color?transfer?between?images[J].?IEEE?Computer?Graphics?and?Applications，?2001，?21（5）：?34-41.

[14]SUTSKEVER?I，?MARTENS?J，?DAHL?G，?et?al.?On?the?importance?of?initialization?and?momentum?in?deep?learning[C]//International?Conference?on?Machine?Learning.?PMLR，?2013：?1139-1147.

[15]LOIZOU?N，?RICHTRIK?P.?Momentum?and?stochastic?momentum?for?stochastic?gradient，?newton，?proximal?point?and?subspace?descent?methods[J].?Computational?Optimization?and?Applications，?2020，?77（3）：?653-710.

[16]BOTTOU?L.?Online?learning?and?stochastic?approximations[J].?Online?Learning?in?Neural?Networks，?1998，?17（9）：?142-176.