陳禹伊， 陳 璐

(上海交通大學(xué) 機(jī)械與動(dòng)力工程學(xué)院， 上海 200240)

隨著信息技術(shù)的不斷發(fā)展，電子商務(wù)(簡(jiǎn)稱電商)已成為最受歡迎的消費(fèi)方式之一.而在線交易量的急劇增加為物流網(wǎng)絡(luò)帶來了巨大壓力，因此高效的商品交付是電商公司提高服務(wù)水平的關(guān)鍵因素之一.電商物流中的“最后一公里”交付問題，即商品從配送站到客戶的配送過程是一個(gè)典型的有容量約束的車輛路徑規(guī)劃問題(Capacitated Vehicle Routing Problem，CVRP)，可以描述為由一個(gè)配送中心向n個(gè)客戶派遣K輛車配送貨物，每個(gè)客戶由一輛車服務(wù)且僅被服務(wù)一次.每輛車設(shè)有運(yùn)送容量上限，而每個(gè)客戶有特定的貨物量需求.其目的為得到成本最低的客戶配送服務(wù)路線.

目前，利用CVRP求解算法[1-2]估算現(xiàn)實(shí)生活中各條道路的行駛成本較為困難，得到的優(yōu)化結(jié)果存在較大誤差，實(shí)際應(yīng)用價(jià)值往往較低.Bertsimas等[3]認(rèn)為，各條道路的行駛成本之間存在相互依賴性，因此對(duì)基于最短路徑或最短行駛時(shí)間成本矩陣的CVRP模型進(jìn)行優(yōu)化難以得到符合實(shí)際需求的路徑.此外，經(jīng)驗(yàn)豐富的駕駛員(專家)并不總是遵循最短路徑成本矩陣選擇配送路徑，如優(yōu)先避免易產(chǎn)生擁堵或交通信號(hào)燈多的道路.因此，在路徑規(guī)劃中還需要融入駕駛員的經(jīng)驗(yàn)或偏好.Kovacs等[4]研究考慮一致性的車輛路徑規(guī)劃問題.其中，一致性包括人員一致性，即客戶希望被所熟悉的駕駛員服務(wù)[5]；區(qū)域一致性，即讓駕駛員在更為熟悉的服務(wù)區(qū)域進(jìn)行配送[6-7].然而，駕駛員的偏好并不總需要完全遵循一致性原則.Liu等[8]通過分析駕駛員的連續(xù)數(shù)據(jù)軌跡來揭示駕駛員的操作模式.韋增欣等[9]建立基于多屬性決策方法的路徑優(yōu)化模型，并根據(jù)駕駛員的偏好調(diào)整模型參數(shù).Pahlavani等[10]利用局部線性神經(jīng)模糊模型學(xué)習(xí)駕駛員偏好來預(yù)測(cè)路線.He等[11]根據(jù)市區(qū)出租車司機(jī)的經(jīng)驗(yàn)，利用協(xié)作偏好發(fā)現(xiàn)算法和智能司機(jī)網(wǎng)絡(luò)生成算法計(jì)算在旅行時(shí)間變化的情況下更可靠的行駛路線.黃敏等[12]根據(jù)出租車司機(jī)的經(jīng)驗(yàn), 提出約束深度強(qiáng)化學(xué)習(xí)算法，并在線計(jì)算不同時(shí)間段內(nèi)的最快行駛路線.以上研究多針對(duì)特定應(yīng)用場(chǎng)景進(jìn)行數(shù)據(jù)或?qū)嶒?yàn)分析，因此方法的應(yīng)用價(jià)值有限.

逆向優(yōu)化是指已知一個(gè)優(yōu)化問題的專家解(最優(yōu)解)，從而逆向推斷優(yōu)化模型的參數(shù)[13].Ahuja等[13-14]研究了基于最優(yōu)性條件和對(duì)偶理論的逆線性優(yōu)化問題，并將其擴(kuò)展到逆向網(wǎng)絡(luò)流問題.Zhang等[15-16]利用牛頓法解決了在l∞范數(shù)下的逆向組合問題.Heuberger[17]綜述了不同逆向組合優(yōu)化問題的算法.隨著大數(shù)據(jù)應(yīng)用的發(fā)展，研究者提出在逆向優(yōu)化中結(jié)合數(shù)據(jù)來設(shè)計(jì)開發(fā)算法.You 等[18]對(duì)城市卡車貨運(yùn)量進(jìn)行建模和預(yù)測(cè)，并提出利用觀測(cè)數(shù)據(jù)對(duì)車輛路徑規(guī)劃問題的模型參數(shù)進(jìn)行校準(zhǔn)，使得優(yōu)化結(jié)果更符合實(shí)際需求.牟健慧等[19]將逆向優(yōu)化理論和方法引入車間調(diào)度領(lǐng)域，建立考慮調(diào)度效率和穩(wěn)定性的數(shù)學(xué)模型，以解決多目標(biāo)流水車間的逆調(diào)度問題.B?rmann等[20]結(jié)合逆向優(yōu)化和在線學(xué)習(xí)算法的思想，將模型解決方案與專家決策之間的差異設(shè)為成本矩陣的梯度，并利用迭代法令梯度逐漸降低并收斂.Aswani等[21-22]考慮了在凸優(yōu)化問題最優(yōu)解的測(cè)量結(jié)果被噪聲破壞情況下的逆向優(yōu)化問題.Saez-gallego等[23]利用逆向優(yōu)化法預(yù)測(cè)對(duì)價(jià)格敏感的消費(fèi)者的電力總需求，建立電力總需求的價(jià)格響應(yīng)模型，并利用廣義逆向優(yōu)化法估算模型參數(shù).

綜上分析，提出一種逆向優(yōu)化方法.通過學(xué)習(xí)配送過程中專家的過往路徑?jīng)Q策，逆向優(yōu)化得到能夠表示專家對(duì)道路選擇偏好的成本矩陣，并將該矩陣應(yīng)用于CVRP模型，令配送路徑能夠融入專家經(jīng)驗(yàn).

1 逆向優(yōu)化模型

對(duì)“最后一公里”的配送問題進(jìn)行建模，并根據(jù)該模型的對(duì)偶模型建立逆向優(yōu)化數(shù)學(xué)模型.

1.1 CVRP模型

配送問題描述如下：共有K輛車為n個(gè)客戶進(jìn)行配送，每個(gè)客戶僅由一輛車負(fù)責(zé)服務(wù)，每輛車的容量約束為Q，第i個(gè)客戶的需求量為qi(i=1, 2, …,n), 0≤qi≤Q.定義集合V={0, 1, …,n}，其中0表示配送中心.C={cij,i,j∈V}為成本矩陣，其中cij為從客戶i到客戶j的行駛成本.優(yōu)化目標(biāo)為在滿足車輛容量約束的前提下，得到總行駛成本最小的車輛配送路徑.

模型的決策變量包括xij(i,j∈V)和車輛到達(dá)客戶i時(shí)的負(fù)載ui(i∈V/{0}).如果車輛直接從客戶i行駛到客戶j，則xij=1，否則xij=0.

CVRP的數(shù)學(xué)模型(模型L)表示如下：

min ∑i∈V∑j∈Vcijxij

(1)

s.t. ∑j∈Vxij=1， ?i∈V/{0}

(2)

∑i∈Vxij=1， ?j∈V/{0}

(3)

∑j∈Vx0j=K

(4)

ui-uj+Q(1-xij)≥qi

(5)

?i,j∈V/{0}

ui≤Q， ?i∈V/{0}

(6)

ui≥qi， ?i∈V/{0}

(7)

xij∈{0, 1}， ?i,j∈V

(8)

目標(biāo)函數(shù)式(1)使得總行駛成本最小化；式(2)和(3)確保每一位客戶僅被服務(wù)一次并控制網(wǎng)絡(luò)流守恒；式(4)確保配送車輛數(shù)目為K；式(5)為子回路消除約束和車輛容量約束；式(6)和(7)用于定義ui；式(8)為決策變量的屬性定義.

1.2 逆向優(yōu)化模型

min ∑i∈V∑j∈V(ηij+σij)

(9)

∑i∈V∑j∈Vφij-∑j∈V/{0}∑j∈V/{0}

[αi+βj+(Q-qi)γij+Qλi+qiμi]=0

(10)

δ≤c00+η00

(11)

(12)

βj+δ+φ0j=c0j+η0j-σ0j

(13)

?(0,j)∈J

(14)

αi+φi0=ci0+ηi0-σi0， ?(i, 0)∈J

(15)

αi+βj-Qγij≤cij+ηij

(16)

αi+βj-Qγij+φij=cij+ηij-σij

(17)

?(i,j)∈J,i,j∈V/{0}

∑j∈V/{0}γij-∑j∈V/{0}γji+λi+μi≤0

(18)

?i∈V/{0}

αi,βj,δ∈R， ?i,j∈V/{0}

(19)

λi,φij≤0， ?i,j∈V/{0}

(20)

γij,μi,ηij,σij≥0， ?i,j∈V/{0}

(21)

式中：αi,βj,δ,γij,λi,μi和φij分別為模型L中式(2)～(8)的對(duì)偶變量.目標(biāo)函數(shù)式(9)使C的變動(dòng)量最小化；式(10)確保x*成為模型L的最優(yōu)解；式(11)～(18)為對(duì)偶約束；式(19)～(21)為對(duì)偶變量的屬性定義.

利用模型Inv-L的迭代算法(Inv-L Algorithm，ILA)，輸入過往專家解，迭代多次求解Inv-L模型，每次求解時(shí)對(duì)C進(jìn)行變動(dòng)，最終得到一個(gè)能夠代表專家經(jīng)驗(yàn)的成本矩陣.然而，該算法在每次迭代時(shí)均需要對(duì)模型Inv-L求最優(yōu)解，因此求解時(shí)間較長(zhǎng).對(duì)此，開發(fā)一種在線學(xué)習(xí)算法以實(shí)現(xiàn)對(duì)過往專家解的學(xué)習(xí).

2 乘性權(quán)重更新算法

乘性權(quán)重更新(Multiplicative Weights Updates，MWU)算法[24]采用迭代更新的方式，在每一次迭代中，通過定義每輪的損失成本，并利用乘法更新規(guī)則對(duì)決策進(jìn)行更新，能夠使得總損失成本的期望最小化.采用改進(jìn)的MWU算法，將原MWU算法中基于列向量的運(yùn)算擴(kuò)展至基于二維矩陣的運(yùn)算，并針對(duì)問題屬性確定算法學(xué)習(xí)率，使得改進(jìn)后的算法具備收斂性.

2.1 算法描述

算法中的相關(guān)參數(shù)如表1所示.

表1 參數(shù)定義Tab.1 Definition of parameters

wt+1=wt-μ(wt⊙mt)

(22)

式中：a⊙b表示維數(shù)相同的兩個(gè)矩陣對(duì)應(yīng)分量相乘.

MWU算法偽代碼

fort=1, 2, …,Tdo

xt=arg min (Ct·xt|xt∈X(pt))

mt=0

else

wt+1=wt-μ(wt⊙mt)

end if

end for

returnC1,C2, …,CT

對(duì)于第t輪，算法首先通過標(biāo)準(zhǔn)化wt得到Ct，并基于Ct和X(pt)正向求得xt.然后計(jì)算mt.最后基于mt和式(22)將wt更新為wt+1.經(jīng)過反復(fù)迭代，得到一個(gè)成本矩陣序列C1,C2, …,CT，其中CT為最終學(xué)習(xí)得到的代表專家經(jīng)驗(yàn)的成本矩陣.由于MWU算法具有收斂性，所以可以保證經(jīng)過足夠多次學(xué)習(xí)迭代后得到的CT無限接近于Ctrue.

在 MWU算法中, 僅存在一個(gè)循環(huán)體，在每個(gè)循環(huán)中均需要調(diào)用數(shù)學(xué)優(yōu)化技術(shù)CPLEX執(zhí)行一次CVRP的求解.由于CPLEX求解線性規(guī)劃的內(nèi)置算法為分支定界算法，所以MWU算法的時(shí)間復(fù)雜度為O(n′2n′)，其中n′為輸入?yún)?shù).

2.2 算法收斂性

(23)

式中：abs(·)為求數(shù)據(jù)絕對(duì)值的函數(shù).

(24)

即MWU算法具有收斂性.

3 算例分析

利用不同維度的算法指標(biāo)評(píng)價(jià)MWU算法的有效性，然后分別利用VRPLIB算例庫中的隨機(jī)算例和電商真實(shí)算例進(jìn)行實(shí)驗(yàn)分析，所需計(jì)算機(jī)配置為Intel?CoreTMi7-8500@3.60 GHz，內(nèi)存為8 GB.程序由MATLAB R2019b進(jìn)行編譯，模型求解調(diào)用CPLEX IBM 12.9.0.

VRPLIB算例庫中的隨機(jī)算例包括成本矩陣、車輛容量和客戶需求等信息.基于給定的成本矩陣，通過更改客戶需求和車輛容量能夠生成多個(gè)不同算例，求得各算例的最優(yōu)解并作為專家解.將同一規(guī)模下的專家解分成兩部分，其中90%的專家解為訓(xùn)練集，10%的專家解為驗(yàn)證集.

3.1 算法指標(biāo)

從算法的有效性、可實(shí)現(xiàn)性和收斂性等方面提出以下評(píng)價(jià)指標(biāo).

平均目標(biāo)函數(shù)值偏差：

其中：M為驗(yàn)證集的數(shù)量.該指標(biāo)用于衡量由逆向優(yōu)化學(xué)習(xí)算法得到的最優(yōu)解與真實(shí)專家解目標(biāo)函數(shù)值之間的平均偏差.

平均收斂誤差：

該指標(biāo)用于衡量算法的收斂誤差.

平均編輯距離：

其中：Edit為求解兩個(gè)變量編輯距離的函數(shù).編輯距離為將一個(gè)字符串轉(zhuǎn)換為另一個(gè)字符串所需要的最小操作數(shù).該指標(biāo)用于衡量驗(yàn)證集中的專家解和由逆向優(yōu)化學(xué)習(xí)法得到的最優(yōu)解之間的平均編輯距離.

平均相似度：

其中：|H|為算例規(guī)模.判斷兩個(gè)解的相似度需要同時(shí)考慮E和|H|.此外，評(píng)價(jià)指標(biāo)還包括算法的學(xué)習(xí)時(shí)間(CPU).

3.2 基于隨機(jī)算例的實(shí)驗(yàn)分析

3.2.1學(xué)習(xí)效果分析 分別利用ILA和MWU算法求解具有10個(gè)客戶點(diǎn)的算例.當(dāng)訓(xùn)練集數(shù)(N1)和驗(yàn)證集數(shù)(N2)逐漸增大時(shí)，兩種算法的性能分析如圖1所示.可知，在指標(biāo)ε1，ε2，S和CPU方面，MWU算法均優(yōu)于ILA.隨著訓(xùn)練集數(shù)量的增多，兩種算法的學(xué)習(xí)效果均有所提升，提升趨勢(shì)均為先激增后趨于平緩，且MWU算法的提升幅度更顯著.因此，在求解大規(guī)模算例時(shí)，可以考慮權(quán)衡訓(xùn)練時(shí)間和訓(xùn)練效果.雖然MWU算法需要花費(fèi)較長(zhǎng)時(shí)間學(xué)習(xí)CT，但在后續(xù)過程中，可以直接利用CT求解新問題.而其在求解新問題時(shí)的效率等同于求解一個(gè)CVRP問題時(shí)的效率.因此，前期訓(xùn)練完成后，MWU算法求解新問題的效率較高.

圖1 兩種算法的效果對(duì)比Fig.1 Comparison of two learning algorithms

根據(jù)上述結(jié)論，選擇MWU算法求解后續(xù)算例，結(jié)果如表2所示.可知，與驗(yàn)證集中已知專家解相比，MWU算法可以提供S=88.05%和ε1=0.16%的配送路徑，表明了算法在學(xué)習(xí)專家經(jīng)驗(yàn)方面的有效性.此外，MWU算法的ε2=1.40，表明算法的收斂性較好.

表2 MWU算法求解不同規(guī)模算例的結(jié)果

3.2.2網(wǎng)絡(luò)圖結(jié)構(gòu)的敏感性分析 上述隨機(jī)算例均為非對(duì)稱網(wǎng)絡(luò)，而在較多配送場(chǎng)景中，網(wǎng)絡(luò)結(jié)構(gòu)呈對(duì)稱性.因此，需要分析MWU算法應(yīng)用于兩種不同網(wǎng)絡(luò)結(jié)構(gòu)時(shí)的性能，結(jié)果如圖2所示.保留非對(duì)稱網(wǎng)絡(luò)成本矩陣的上半部分不變，可以將其轉(zhuǎn)化為對(duì)稱成本矩陣.

圖2中，與求解非對(duì)稱網(wǎng)絡(luò)的算例相比，MWU算法在求解具有對(duì)稱網(wǎng)絡(luò)結(jié)構(gòu)算例時(shí)的ε1值略低，而ε2值相對(duì)一致.表明算法的收斂性不受網(wǎng)絡(luò)對(duì)稱性影響.由于對(duì)稱網(wǎng)絡(luò)中存在多個(gè)最優(yōu)解，會(huì)對(duì)算法學(xué)習(xí)造成干擾，所以對(duì)其進(jìn)行學(xué)習(xí)后所得求解結(jié)果的相似度較差.此外，隨著算例規(guī)模增大，對(duì)稱網(wǎng)絡(luò)的訓(xùn)練時(shí)間呈逐漸增長(zhǎng)趨勢(shì).

圖2 網(wǎng)絡(luò)結(jié)構(gòu)對(duì)稱性敏感性分析Fig.2 Sensitivity analysis of two networks

圖3 研究區(qū)域的客戶聚類分布Fig.3 Clustered customer network of studied district

3.3 基于電商真實(shí)數(shù)據(jù)的算例分析

采用某電商在2018年9月至11月的日常家用電器配送路線進(jìn)行算例分析.服務(wù)區(qū)域?yàn)樯虾Ｊ行靺R區(qū)，收集到包括422輛廂型貨車在內(nèi)的20萬條訂單的配送路線.數(shù)據(jù)集包括配送的車次、時(shí)間和地址等.根據(jù)企業(yè)反饋，選擇工作效率較高的駕駛員，收集其配送路徑，并剔除異常解和存在客戶時(shí)間窗約束等不屬于研究范圍的解，但此時(shí)所得專家解并不等同于基于最短路徑矩陣的CVRP模型最優(yōu)解.將駕駛員的配送路徑作為專家解進(jìn)行學(xué)習(xí)，以得到能夠代表駕駛員經(jīng)驗(yàn)的成本矩陣，并用于求解新問題.

3.3.1數(shù)據(jù)預(yù)處理 數(shù)據(jù)預(yù)處理包括缺省值、異常值和地址的經(jīng)緯度轉(zhuǎn)換處理.為了增加客戶點(diǎn)的重復(fù)度，減少模型訓(xùn)練難度，根據(jù)客戶的地理位置采用K-means算法進(jìn)行聚類.從算法的可實(shí)現(xiàn)性、聚類合理性及其聚類指標(biāo)等方面綜合分析，選擇聚類數(shù)為55，客戶聚類情況如圖3所示.數(shù)據(jù)預(yù)處理后共得到185個(gè)真實(shí)算例.

3.3.2結(jié)果分析 各聚類質(zhì)心的歐氏距離成本矩陣為C1，算法學(xué)習(xí)過程為從185個(gè)算例中隨機(jī)選擇184個(gè)并利用MWU算法進(jìn)行訓(xùn)練，將剩余的一個(gè)算例作為驗(yàn)證算例.將學(xué)習(xí)得到的CT代入路徑規(guī)劃模型L，求解驗(yàn)證算例，并對(duì)比所得解與已知專家解.按上述過程訓(xùn)練185次，求每一次訓(xùn)練效果的平均值.由于真實(shí)場(chǎng)景下無法得到真實(shí)的成本矩陣，所以本節(jié)不考慮指標(biāo)ε1和ε2.此外，利用C1求解算例，并對(duì)比采用算法前后的優(yōu)化效果，如表3所示.

表3 真實(shí)算例結(jié)果對(duì)比Tab.3 Comparison of real instances

可知，當(dāng)基于CT求解算例時(shí)，求解結(jié)果明顯優(yōu)于基于C1的求解結(jié)果，在指標(biāo)E和S上分別優(yōu)化了37.38%和24.67%.同時(shí)，所得配送路徑與專家解之間的平均相似度為71.84%，表明當(dāng)真實(shí)算例的數(shù)量增加時(shí)，算法的學(xué)習(xí)效應(yīng)也將隨之加強(qiáng).

4 結(jié)語與展望

對(duì)電商物流中的“最后一公里”交付問題進(jìn)行研究.在現(xiàn)實(shí)生活中，經(jīng)驗(yàn)豐富的駕駛員(專家)并不總是遵循最短路徑成本矩陣選擇配送路徑.因此，采用逆向優(yōu)化法學(xué)習(xí)專家的過往路徑?jīng)Q策方案，訓(xùn)練得到代表專家對(duì)路徑選擇偏好的成本矩陣，將其應(yīng)用到路徑規(guī)劃模型中求解.其中，利用MWU算法實(shí)現(xiàn)對(duì)專家經(jīng)驗(yàn)的學(xué)習(xí)，并設(shè)計(jì)隨機(jī)算例和電商真實(shí)物流算例進(jìn)行驗(yàn)證.結(jié)果表明：① MWU算法在收斂性、有效性和可實(shí)現(xiàn)性方面表現(xiàn)優(yōu)異，利用學(xué)習(xí)所得成本矩陣可以得到近似專家的配送路徑；② 算法收斂性不受網(wǎng)絡(luò)結(jié)構(gòu)非對(duì)稱性和對(duì)稱性的影響；③ 在求解真實(shí)算例時(shí)，MWU算法將所得解與專家解的平均相似度提升至71.84%，表明該算法能夠有效利用過往經(jīng)驗(yàn)，提供與專家解相似度極高的配送路徑.

在后續(xù)研究中，可以將逆向優(yōu)化法應(yīng)用于具有不確定性車輛路徑規(guī)劃等更復(fù)雜的優(yōu)化問題.同時(shí)，現(xiàn)實(shí)生活中的專家解集往往來源于多個(gè)專家，如何對(duì)專家解進(jìn)行分類并提供全局逆向優(yōu)化方案，同樣值得探討.