基于高維點云分布的最優(yōu)傳輸問題的研究

阮一鳴

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

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

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

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

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

1?分片Wasserstein距離理論

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

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

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

1.2?分片Wasserstein距離

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

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

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

2?理論與算法

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

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

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

2.1.1?隨機梯度下降算法

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

2.1.2?標準動量法

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

2.1.3?Nesterov?Momentum動量法

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

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

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

2.2?算法基本理論

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

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

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

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

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

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

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

其中η表示衰減因子，表示對先前方向的依賴程度，通常取值為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距離的計算算法

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

Require：初始點云分布X（0）=X，初始動量v（0）

While?不滿足停止準則?do

1：對X進行隨機梯度法迭代，?初值迭代X（0）=X

2：計算上述投影算子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ù)值實驗

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

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

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

3.2.1?數(shù)值實驗

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

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

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

3.2.2?評價指標

接下來利用標準差、峰值信噪比和平均梯度三個評價指標對結(jié)果圖像進行評價。

（1）標準差

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

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

而圖像標準差反映了圖像像素值與均值的離散程度，標準差越大說明圖像的對比度越高.標準差計算公式為：

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

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

（2）峰值信噪比

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

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

PSNR=10lg10f2maxMSE?（22）

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

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

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

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

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

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

4?結(jié)?論

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

參考文獻

[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].?運籌學(xué)學(xué)報，?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.