求解公共自行車再平衡問題的克隆選擇算法

劉喜梅，潘立軍

湖南工程學院 管理學院，湖南 湘潭411104

1 引言

公共自行車再平衡問題（Bike Sharing Rebalancing Problem，BRP）是一類NP-難問題[1]。其可描述為利用有相同載重的車輛從單車維護補給中心出發完成一定區域內各單車停放點自行車的再分配，使各停放點達到預先設置容量后回到中心，如何使所有的車輛行駛的距離或花費的時間成本最少。20 世紀60 年代，公共自行車系統最早在荷蘭阿姆斯特丹出現，隨著社會經濟的不斷發展，特別是交通擁堵、環境污染等城市問題不斷突顯，公共自行車作為健康、環保、便捷的出行方式逐步在世界各地普及。目前，全世界約有49個國家共計超過500個公共自行車系統，如法國巴黎（20 600 輛）、中國杭州（78 000 輛）、中國武漢（90 000 輛）、中國株洲（20 000輛）[2]。公共自行車系統的普及使各地公共自行車運營部門單車再平衡任務加劇，通常公共自行車運營公司組織自有運輸車輛從維護補給中心出發，按照事先各停靠點分配單車取送量依次訪問各停靠點，再回到維護補給中心，在該過程中，運輸車輛既可以在各停放點間調配自行車，也可在出發前從維護補給中心載入一部分自行車調配到各停放點，或將各停放點多余的自行車運回維護補給中心。

國外有關公共自行車再平衡的研究報道目前主要關注二類BRP，一是靜態BRP[3-9]，即在車輛開始服務各停放點后，各點平衡數量不發生變化，如在夜間（22點至第二天凌晨6點間）用戶幾乎不使用自行車的時段進行再平衡作業；二是動態BRP[10-12]，即在車輛開始平衡作業服務后，各停靠點單車平衡數量依舊發生變化，如在公共自行車使用頻率較高的區域，再平衡期間仍有用戶用車。從BRP的問題特點來看，求解目標有最小化車輛行駛費用[11]、最小化再平衡總費用（該費用包括車輛運輸費用、裝卸單車費用）[10，12]、最小化作業時間[3-8]、最大化客戶滿意度（滿意度表示為站點初始單車保有量與用戶取還車頻率的函數，用戶在站點取還車失敗的可能性越高，滿意度越低，反之則越大）[9-10，12]。從BRP求解方法來看，目前報道較多的主要有三類，一類是精確算法，如分枝定界法[3-5]；二是各類經典啟發式算法，如節約算法[6]、貪婪算法[7]，可變鄰域下降搜索[8]等；三是智能啟發式算法，如禁忌搜索算法[9]。表1 對國外主要的研究報道進行了梳理總結。

我國雖然是公共自行車保有量最大的國家，但目前國內公共自行車再平衡的研究報道比較少，白雪等[13]研究了考慮維修車輛的公共自行車系統再平衡問題，提出了基于動態規劃的精確算法。前期也對求解BRP的遺傳算法[14]進行了研究，運用標準算例測試表明其具有較好的求解效率。此外，葉麗霞[15]、喬曉[16]對我國南京市江寧區、西安市雁塔區的公共自行車再平衡作業進行了案例研究。

綜合分析來看，國外對BRP 的研究起步早，成果較多，對BRP問題特點及精確算法、經典啟發式算法作了較深入的研究，但應用智能啟發式算法求解BRP的研究比較少，尚沒有求解BRP的人工免疫克隆選擇算法的報道，而國內研究多為案例研究，研究成果難以推廣應用，因此均不適合開發具有一定通用性的計算機調度軟件來輔助公共自行車系統再平衡。在實踐中，公共自行車的再平衡是一項周期性的調度工作，用戶的用車行為具有慣性（即一段時間內各停靠點取車、用車數量不會出現劇烈波動），前一天的再平衡調度方案對第二天的調度方案優化具有重要的參考價值。

人工免疫算法是模仿生物免疫系統的一種智能優化方法。生物免疫系統通過構建具有動態性、自適應性和自組織性的信息防御系統，以此來抵御外部無用、有害和干擾信息的侵入，從而保證系統接受信息的有效性與無害性。代表性的人工免疫算法是巴西人工免疫學專家Castro 等人[17]借鑒生物免疫系統的克隆選擇原理提出的克隆選擇算法，該算法通過模擬生物免疫系統的對抗原初次免疫應答和二次免疫應答機制，能較好地處理優化求解中基于過往方案優化現有方案的場景[18]，作為一種全局性鄰域搜索智能啟發式算法，在處理NP-難問題時有較好效果，Kim[19]與Coello[20]分別運用該方法處理了多目標的數值優化問題，文獻[21]則運用該方法求解帶工作時間與時間窗的開放式車輛路徑問題，均顯示出較好的求解效果。

表1 國外BRP研究代表性文獻

2 BRP數學模型

借鑒文獻[5]，BRP 的數學模型表示如下：設有一完備圖G=(V,E)，V={1,2,…,n}表示頂點集，其中1表示自行車維護與補給中心，V′ {2,3,…,n}表示單車停放點集，E={(i,j)|i,j ∈V,i ≠j}表示弧集或邊集。其符號定義為：M 表示實施再平衡運輸的專用車輛數；Q表示專用運輸車輛的最大載重量，di表示車輛在第i 個單車停放點取送量，di的取值范圍為[-Q,Q]，di小于0，表示該停放點需補充單車，di大于0，表示該停放點要回收單車；cij表示專用車輛訪問弧(i,j)產生的運輸成本或時間；θj表示專用車輛訪問第j 個停放點后的載重量，即車輛線路載重量。決策變量xij=1，表示車輛訪問了弧(i,j)，否則，xij=0。

式（1）為目標函數，為最小化運輸總成本或時間，式（2）、式（3）確保除維護補給中心點外其余單車停放點均被訪問且只訪問一次，式（4）、式（5）確保所有的專用運輸車輛在完成任務后均回到維護與補給中心，式（6）為消除子回路約束，式（7）、式（8）、式（9）為BRP車輛載重約束，確保各專用運輸車輛在實施單車回收與投放過程中，車輛線路載重量不超過額定載重量，其假設在某一停靠點回收的自行車可被投放到任意其他需要的點，且車輛出發或者回到補給中心時，載重可不為零，因為可以事先裝入部分單車，或者帶回部分單車。

3 求解BRP的克隆選擇算法

3.1 算法主要框架

求解BRP 的克隆選擇算法基本步驟如圖1 所示。其中抗體的編碼方式采用多維整數編碼方法[14]，初始解生成借鑒文獻[22]的方法，采用隨機選擇停靠點、插入位置插入。算法的終止條件設定為迭代次數Max_iter或當前最優解連續不改進的次數Max_Not_Improve 達到事先設定的常數。

圖1 克隆選擇算法基本流程圖

算法運行中設定了一定容量的記憶細胞群體Memory，算法每迭代100 次，將抗體群中的最優抗體存入Memory，算法終止后，從當前Memory 輸出最優抗體。當有新的問題求解時，也可直接從Memory 中選取一定數量的記憶細胞，直接計算是否滿足新的問題約束條件，若滿足，可直接計算抗體與抗原的親和力，進入算法后續步驟迭代，其過程如圖1 中虛線部分所示，這部分反映了人工免疫系統中特有的二次免疫應答機制，即當有新的抗原侵入免疫系統時，可直接借助成熟的記憶細胞展開免疫應答，以達到快速形成免疫能力的目標。

在實踐中，由于部分公共自行車用戶需求以通勤出行為主，因此公共自行車各停靠點的再平衡量作業量在相鄰的幾天內變動較小，前一天的再平衡調度方案對形成新的再平衡調度方案有重要的參考價值。同時分析BRP模型來看，約束條件式（9）為兩不等式，表明可行線路中車輛通過每點的裝載量有一定的變動空間。因此，在克隆選擇算法中引入二次免疫應答機制對提升算法求解效率有幫助。

3.2 親和力定義與克隆選擇算子

在本文算法中，求解目標為使目標函數式（1）最小。因此將抗原與抗體的親和力定義為目標函數式（1）的倒數，即目標函數值越小，則該抗體與抗原的親和力越大。

克隆選擇算子設計參考文獻[21]的方法，主要有以下兩步：

步驟1 計算當前抗體群中每個抗體與抗原的親和力，并按親和力降序排列。

步驟2 按照下式對種群中親和力高的抗體進行克隆，得到新的抗體群Nc：

式中，Scale 代表抗體群的規模，Pos 表示抗體群在降序排列后該抗體在群中的序位，Pos 為整數且滿足：0 ＜Pos ＜；Nc代表克隆后新的抗體群，Nc ＞Scale；β 為克隆系數，其取值區間為：0.2 ＜β ＜0.3。

3.3 變異算子

抗體經過克隆后則經歷高頻變異，高頻變異是克隆選擇算法進行鄰域搜索的主要手段。本文采用混合鄰域結構，即隨機選擇點對操作對當前解的鄰域進行搜索，設計了以下五類變異算子。為了便于說明，設S1=（1356 147928），編號1 代表維修保養中心，選取進行變換的點對為i=5,j=9，具體變換方法如下：

（1）頂點前向重新指派

將所挑選的頂點i 從線路上當前的位置中取出，并將其插入到所挑選頂點j 的位置之前。例如：S1=（1356 147928）→S2=（136 1475928）。

（2）頂點后向重新指派

將所挑選的頂點i 從線路上當前的位置中取出，并將其插入到所挑選頂點j 的位置之后。例如：S1=（1356 147928）→S2=（136 1479528）。

（3）頂點交換

將所挑選的頂點i、j 交換位置。例如：S1=（1356 147928）→S2=（1396 147528）。

（4）尾巴交換

若所挑選的兩個頂點i、j 位于不同的線路上，則將所挑選的兩個頂點后面的“尾巴”（從被選頂點至線路末尾）互換。例如：S1=（1356 147928）→S2=（13928 14756）。

（5）2-opt

若所挑選的兩點i、j 位于同一條線路上，則將i、j兩點間所有頂點的排列順序逆轉。若所挑選的兩點i、j 位于不同線路上，則執行以下四種變換：變換1，S1=（1356 147928）→S2=（13829 14765）；變換2，S1=（1356 147928）→S2=（139741 6528）；變換3，S1=（1356 147928）→S2=（97416 53128）；變換4，S1=（1356 147928）→S2=（8296 147531）。

以上幾種變異方法隨機采用，且為兼顧算法搜索的質量與速度，算法前期（當前迭代次數＜Max_iter/2）時采取隨機變異，即變異后的抗體可行，且與原抗體目標值不同即接受變異，后期（當前迭代次數≥Max_iter/2）時采取尋優變異，即每種變異取不同的點對執行5 次，若變異后的抗體可行，且親和力優于原抗體，則接收變異。

3.4 抗體相似性定義與抗體抑制

為進行抗體抑制，定義抗體相似性Lxy的度量方法如下：

設x 與y 為同一抗體種群中兩個不同抗體，p 為抗體邊矩陣（n 維方陣，n 為完備圖G 中頂點數量），若抗體x 中有邊E(i,j)，則=1，否則為0(i,j ∈V,i ≠j)。

式（11）中，函數sum(p)的功能為計算抗體邊矩陣中各元素之和，若抗體x 與y 完全相同，則其邊矩陣px=py，則Lxy=0；若抗體x 與y 完全不相同，則Lxy=1。

抗體群經歷克隆與高頻變異后，產生了一些親和力低或與其他抗體結構相近的抗體，通常抗體抑制過程即將這些抗體從種群中刪除并補充新的抗體到原抗體群中，該過程有利于增加算法的搜索區域，避免算法陷入局部收斂。本文抗體抑制的具體步驟如下：

步驟1 定義抗體相似度系數λ ，其取值范圍為[0.3，0.6]。

步驟2 計算抗體與抗原的親和力，按親和力降序排列，得抗體群Pop(a1,a2,…,aNc)，Nc 為抗體群規模。

步驟3 將a1加入新的抗體群Popnew。

步驟4 從a2開始，依次將Pop 中的抗體與Popnew中當前抗體比較，若相似度小于等于λ，則將其加入到Popnew中，直到Popnew的規模達到Scale，若Pop 中的抗體數量不足，則從記憶細胞庫中選擇記憶細胞加入，若還不足，則生成新抗體加入到Popnew中。

4 算法測試

算法采用matlab2014 編程實現，運行于CPU（core i3，3.1 GHz）、ROM（4 GB）的PC機上。

4.1 初次應答求解效率測試

首先運用BRP標準測試算例，測試人工免疫克隆選擇算法初次應答求解效率，主要參數設置為：變異概率取0.9，種群規模Scale=50,β=0.3,λ=0.4 ，最大循環次數Search_Max 為4 000，Not_Improve 為800。每個算例運行10 次，取最好解與平均運行時間。將本文算法與分支定界法（B&C）[5]、遺傳算法[14]進行比較（三種算法的測試硬件環境基本相同）。

表2、表3 中，Si代表各算法的求解結果，CPUi表示各算法求解CPU消耗，N/C 分別代表算例的停靠點數量與運送卡車額定載重量，g_B&C=(S2-S1)/S2、g_GA=(S3-S1)/S3，表示本文算法相較于B&C、GA 在求解質量上的差距，g_t_B&C=(CPU2-CPU1)/CPU2、g_t_GA=(CPU3-CPU1)/CPU3，為本文算法相較于B&C、GA在求解時間上的差距。（1）在規模小于50個點算例上，本文算法能找到所有算例的當前最好解，平均CPU 消耗為9.6 s，比B&C 的快96.80%，但比GA 慢，如表2所示；（2）在規模為50至100個點的算例，本文算法找到6 個算例的最優解（表3 中加粗表示），平均求解質量比B&C 低7.43%，但比GA 高1.02%，平均CPU 消耗相較于B&C快96.8%。

4.2 算法二次應答的求解效率測試

為了測試人工免疫克隆選擇算法二次應答的求解效率，本文對BRP 標準測試問題Roid(N=55/C=20)的停靠點取送車量（d）進行修改，生成3組新的測試問題，用以模擬不同場景下公共自行車運營系統連續兩天的再平衡需求。修改方法為：分別取p=1、2、3，隨機選取Roid 中80%的停靠點，對原停靠點取送單車數量按d±p 進行修改（p 的正負號隨機選取，但需確保新的停靠點取送車數量不等于0），每組算例生成4個新的問題。算法測試的參數設置除Not_Improve=400 外，其他參數設置與第一部分測試相同。二次應答測試方法為算法先運行Roid問題，并將求解過程中的記憶細胞保存，再分別運行改進后的測試問題，記錄運行結果及CPU消耗。

測試算例Roid中各點再平衡車輛數的取值區間為[-7，7]，表4～6分別模擬測試連續兩天內，再平衡量變化幅度為14.3%(p=1)、28.6%(p=2)、42.9%(p=3)的情況下，算法二次應答效率。測試結果表明：二次應答的求解質量分別比一次應答的高0.65% (p=1) 、0.24%(p=2) 、0.16% (p=3) ，CPU 消耗比一次應答分別快39.78%(p=1)、46.18%(p=2)、40.66%(p=3)。表4～6中：g_32=(S1-S2)/S1,t_32=(CPU1-CPU3)/CPU3。

表2 小于50個點的算例測試結果

表3 大于50個點小于100個點的算例測試結果

5 結束語

公共自行車系統在我國普及度高，系統的再平衡調度方法研究有著非常廣泛的應用價值。本文在已有研究基礎上，設計了求解BRP的人工免疫克隆選擇算法，算法采用多維整數編碼方法，結合問題特點設計了新的抗體相似性定義方法及抗體抑制策略，并結合BRP周期性調度特點引入二次應答求解機制。運用標準算例測試表明：（1）一次應答測試中，本文算法在求解規模小于50個點的問題中，均能找到最優解，求解質量與分枝定界法相當，但速度更快；（2）在求解規模為50個點到100個點間的算例上，算法求解質量相較于GA提高1.02%；（3）二次應答測試中，二次應答的求解質量比一次應答的略高，但CPU消耗比一次應答快39%以上。

表4 二次應答測試(p=1)

表5 二次應答測試(p=2)

表6 二次應答測試(p=3)

結合本文算法一次應答、二次應答的求解測試結果，本文認為，求解BRP的人工免疫克隆選擇算法的求解速度較精確算法更有優勢，求解質量與遺傳算法相差不大，其二次應答機制能加快周期性問題的求解速度，在實現公共自行車再平衡實時調度方面能體現更好的實用性。公共自行車再平衡是一項系統工程，即與停靠點的規劃布置有關，也與用戶取用車習慣有關，還與運營平臺再平衡調度管理策略有關，將這些因素融入BRP模型，設計新的問題類型與對應求解方法將是后續BRP問題的重要研究方向之一。