隨機梯度估計及其在系統管理中的應用

張公伯 ,李海東 ,彭一杰

(北京大學a.光華管理學院;b.工學院,北京 100871)

隨著信息技術的快速發展,系統管理的技術方法正在發生深刻變革。系統仿真研究是一個經典的、傳統的領域,其優點在于可以對一般的復雜系統做高置信度建模,受到學界與業界的廣泛關注。隨著我國產業升級不斷推進,國際競爭程度日益加深,企業對系統管理精益化程度的要求不斷提高,仿真技術在系統管理中應用的重要性日益顯現。大數據、工業物聯網、云計算、人工智能和數字孿生技術等作為實踐系統管理的新興技術在仿真分析中起著重要作用。復雜系統通常具有動態性和隨機性,準確估計系統的平均表現需要做大量仿真,且仿真耗費的計算時間長。復雜隨機系統仿真中通常包含影響系統結構的控制參數與不連續樣本表現,導致系統表現的優化分析困難。因此,需要從理論與方法論層面研究如何將仿真與優化有效結合,應用于復雜隨機系統的管理。有關仿真優化方法可參考文獻[1-7]。

對于具有連續可控參數的優化問題,梯度蘊含了目標函數隨參數變化的信息,在迭代搜索優化方法中有著重要的應用,使得結果逐步逼近局部或全局最優。復雜隨機系統優化問題的目標函數的解析表達通常未知,只能通過蒙特卡洛仿真得到的系統樣本表現進行估計。由于蒙特卡洛仿真存在噪聲,因而觀測到的系統樣本表現具有隨機性。通過蒙特卡洛仿真估計的隨機系統中目標函數關于參數的梯度稱為隨機梯度估計。隨機梯度估計在仿真分析和優化問題中有著廣泛的應用,如敏感性分析、仿真模型的不確定性量化和隨機優化及機器學習中的梯度下降算法。

隨機梯度估計可應用于排隊系統、庫存管理系統、統計質量控制[8-14]、機器維護系統[15-18]和隨機活動網絡[19-21]等系統的優化和調控決策。現實中的系統通常結構復雜并具有隨機性,例如,生產過程中產品訂單的到達、機器加工時間及是否發生故障都是隨機的。這些復雜系統通常不具有解析形式的輸出分布,但具有部分結構信息可以分析,被稱為灰箱系統。恰當地利用這部分結構信息可以在仿真估計系統表現的同時得到其梯度信息,從而實現求解優化問題時高效搜索的目標。例如,質量控制圖中希望在真實生產環境下運行系統的同時實現控制邊界的在線優化設計,但質量控制圖的控制邊界屬于結構參數,它的表現關于控制邊界不連續,使得質量控制圖的結構分析非常困難;生產系統中通過排隊系統對產品的加工流程進行建模,系統中通常存在重加工、批量加工和優先處理等結構,因而使得排隊網絡系統的結構分析非常困難。

采用隨機梯度估計方法對分位數敏感性進行計算和分析,也是仿真和復雜系統管理領域關注的重要問題。相較于系統的期望表現,分位數能夠提供更多隨機系統表現尾部分布的特征,從而更好地刻畫系統在極端環境下的表現,提高極端事件下系統的穩健性。例如,生產訂單在排隊系統中的等待時間過長就可能失效,因此,等待時間的分位數過大將有損訂單的數量,降低企業的效益。

隨機梯度估計方法主要分為間接和直接估計。間接梯度估計通常具有兩個特征:①它只估計真實梯度的一個近似值,如在標量情況下通過割線近似;②它僅使用來自原始(未修改)系統的函數評估(績效指標輸出樣本)[22-23]。直接梯度估計利用仿真模型中隨機變量的分布信息或樣本表現函數的結構信息(如輸入的概率分布或系統表現與樣本路徑的函數關系)給出真實梯度的無偏梯度估計量。

間接梯度估計方法主要有有限差分法(Finite Differences,FD)和同時擾動法(Simultaneous Perturbation,SP)[22-23]兩種。這兩種方法均不需要了解仿真模型隨機性的產生或樣本表現函數,即仿真模型可以被視為黑箱。有限差分法通過仿真兩個臨近輸入點的函數值并用微積分基本原理來近似導數,由此得到的估計量既有導數近似帶來的偏差,也有函數評估噪聲帶來的方差。對于高維問題,有限差分法需要獨立地估計每一個維度,因此,估計一次梯度所需要的仿真樣本的數量與維度相關。同時擾動法的優點是一次仿真同時擾動自變量的多個維度,因此,估計一次梯度所需要的仿真樣本的數量與自變量的維度無關。

盡管間接的梯度估計具有通用性,但直接的梯度估計具有如下優點:①它得到的梯度的估計量是無偏的,以無偏的梯度估計量作為迭代方向可以加快仿真優化算法(如隨機逼近算法[24])的收斂速度;②它避免了間接的梯度估計方法中確定合適的擾動大小的需要,而擾動大小的選擇直接影響梯度估計量的準確性;③所得的估計量在計算上更高效,在高維問題中直接的梯度估計方法達到給定的估計精度所需要的仿真樣本的數量通常要比間接的梯度估計方法更少[22-23]。常用的無偏隨機梯度估計方法有無窮小擾動分析法(Infinitesimal Perturbation Analysis,IPA)、似然比法(Likelihood Ratio,LR)(也稱為分數函數法(Score Function,SF))、弱導數法(Weak Derivative,WD)(也稱為測度求導法(Measured-Value Differentiation,MVD))和廣義似然比法(Generalized Likelihood Ratio,GLR)。無窮小擾動分析方法可以用于解決帶有結構參數(即出現在樣本路徑與系統表現度量的函數關系中的參數)和分布參數(即出現在輸入隨機變量的概率分布中的參數)的問題,但要求樣本表現具有連續性;似然比方法雖然能夠處理不連續的樣本表現,但僅適用于系統表現度量關于分布參數的梯度估計而不適用于結構參數;弱導數方法能降低似然比方法的方差,但需要仿真額外的樣本軌道。Gong 等[25]提出了平滑擾動分析法(Smoothed Perturbation Analysis,SPA),對無窮小擾動分析法進行了拓展,通過取條件期望來解決樣本表現不連續的問題,但所選取的條件通常和問題本身相關,而且條件期望估計量的計算可能需要函數取逆和額外的仿真。Rubinstein等[19]提出了推出式似然比法(Pushout LR)來解決同時帶有結構參數和不連續的樣本表現的問題,但是一般需要解析的函數轉換將結構參數從樣本表現中移出。Peng等[11]提出的廣義似然比方法既適用于帶有結構參數和不連續的樣本表現的問題,同時具有無偏性和單樣本軌道性,也不需要取條件期望和函數轉換。

系統管理領域信息化和數字化的轉型,對基于梯度的仿真優化方法的設計提出了挑戰。復雜系統的優化設計需充分考慮問題情境特點及其中的模型假設,以設計簡潔有效的基于梯度的仿真優化方法。本文將首先介紹常用的隨機梯度估計方法,進而考慮4個隨機系統的梯度估計問題,給出其在系統管理中的應用場景。

1 隨機梯度估計方法

其中,當Y是連續隨機變量時,α=P（Y（θ）≤qα（θ））。隨機梯度估計問題研究的是通過對隨機模型進行仿真,得到系統表現度量關于參數的導數的估計。隨機模型中的隨機梯度估計方法可以用于仿真優化的輸入和對參數進行敏感性分析。仿真優化中常用的隨機逼近方法的參數迭代更新過程就采用梯度估計作為其迭代方向,第(n+1)步迭代的計算式為

1.1 間接的梯度估計方法

1.1.1有限差分法 有限差分法不需要獲取復雜系統的信息,僅通過分析系統仿真得到的輸出樣本,就可得到有偏的梯度估計量。有限差分法將復雜系統視為“黑箱”,分別對參數向量θ中每一個元素進行攝動,可得到系統期望表現梯度估計量[22-23]。單邊前向有限差分估計量是有限差分方法中常用的估計量,其第i個元素為:[Y（θ+ci ei）-Y（θ）]/ci。其中:ci為第i個方向的攝動值;ei為第i個方向的單位向量;有限差分法中另一種常用估計量是雙邊對稱有限差分估計量,其第i個元素為

當參數向量的維度很高時,有限差分法需要的計算量很大。此外,有限差分法中需要選取攝動值的大小,如果攝動值的選取過小,得到的系統表現梯度估計量的偏差小、方差大。因此,有限差分法攝動值的選取需要平衡表現梯度估計量的偏差與方差。在實際應用中,雙邊對稱有限差分估計量比單邊前向有限差分估計量的精度高,但在每一次梯度估計時,雙邊對稱有限差分估計量需要消耗更多的仿真資源(前者需要2d個仿真資源,而后者僅需要(d+1)個仿真資源)。對于二階連續可微函數,雙邊對稱有限差分估計量的均方誤差的最優收斂階為n-2/3,所對應的攝動值的階為n-1/6,其中n為仿真資源的數量。得到最優攝動值(即包括n-1/6前的常數)對于小樣本情況下雙邊對稱有限差分估計量的表現十分重要,但往往非常困難。

由于不需要獲取模型的信息,故有限差分法普遍適用于復雜系統表現函數的梯度估計。但采用有限差分法估計梯度需要多次重復仿真,運算成本大,因此不適用于高維參數問題。在基于梯度的優化方法中,有限差分法比直接的梯度估計法的表現差。

1.1.2同時擾動法 同時擾動法由Spall[26]提出。同時擾動估計量消耗的仿真資源與參數向量的維度無關,適用于求解高維參數問題,它的第i個元素為:［Y（θ+CΔ）-Y（θ-CΔ）］/（2ciΔi）。其中:Δ=(Δ1,Δ2,…,Δd)是d-維擾動向量,通常假設各元素Δi獨立;C是對角線為攝動值ci的對角矩陣。在每一次梯度估計時,同時擾動法需要生成d個隨機數來構造擾動向量Δ。與雙邊對稱有限差分估計式(1)相比,同時擾動估計量的每個元素的分子都相同,只需要兩個樣本表現函數擾動估計,與參數向量的維度d無關。由于生成樣本表現函數擾動估計的運算成本大,因而同時擾動法比有限差分法的計算效率高。為保證同時擾動法幾乎處處收斂,要求擾動向量中每一項獨立,均值為0且逆二階矩有限。均勻分布和正態分布不滿足上述條件,常用的Δi的分布是成功概率為1/2的伯努利分布。類似于擬蒙特卡洛方法的核心思想,文獻[27-28]中提出了采用確定的擾動向量代替隨機擾動向量Δ的方法。

文獻[29]中介紹了采用同時擾動法作為迭代方向的隨機逼近方法。同時,擾動法求解高維參數問題的優越性使其在排隊系統[30]、交通系統[31-35]和統計質量控制[36]等系統的優化和調控決策中有廣泛的應用。Fu等[30]采用同時擾動法估計排隊系統的系統停留時間關于服務臺服務速率的梯度;朱承元等[34]優化設計珠三角地區多機場系統的航班時刻,采用同時擾動法估計系統中航班總延誤時間的梯度,降低了航班總延誤成本;齊心[35]優化設計城市交叉口排隊系統的信號燈時長,采用同時擾動法估計車輛在交叉口的通行時間的梯度,提升了城市交叉口的通行能力。

1.2 直接的梯度估計方法

1.2.1擾動分析法 擾動分析法(Perturbation Analysis,PA)是一種通過分析單個樣本路徑對隨機模型的期望表現梯度進行估計的方法。文獻[37-38]中最早對該方法進行了研究。這種方法通過樣本路徑對表現度量的影響以及系統參數對樣本路徑產生的擾動來衡量系統表現梯度。擾動分析法中最常用的是無窮小擾動分析法,該方法通過樣本表現關于參數的導數來估計系統期望表現梯度。系統期望表現為

如果期望與導數是可以交換的,則系統期望表現的梯度為

無窮小擾動分析法可以解決樣本表現帶有結構參數的隨機梯度估計問題,但不能解決樣本表現不連續的問題。為克服這一理論缺陷,文獻[25,42-43]中提出利用取條件期望處理樣本表現不連續的平滑擾動分析法,但該方法需要解析地計算一個條件期望。在實際應用中,該方法需要視實際問題情況而選取適當的條件。Hong[44]首次提出利用無窮小擾動分析法解決分位數敏感性估計問題,但得到的梯度估計量是漸近無偏,且需要對數據進行批處理才可以保證估計量的弱相合性。Jiang等[45]在文獻[44]的基礎上提出一種不需要進行批處理的無窮小擾動分析估計量,但該估計量僅在特定的條件下具有強相合性。有關無窮小擾動分析法和平滑擾動分析法的詳細介紹,可以參考文獻[46]。在復雜系統的優化設計中,基于擾動分析法的隨機逼近方法可以實現生產參數的估計。排隊系統是對生產工藝流程建模的一種常用方法,可以刻畫生產訂單到達和裝配生產線加工等生產流程。優化工藝流程的資源配置(如訂單到達速率、生產線機器數量和緩沖區數量等)可以提高生產效率、降低生產成本。Cheng[47]考慮多階段存貨式生產系統,采用無窮小擾動分析法估計樣本路徑的梯度,實現減少訂單等待時間和提高生產加工效率的權衡。Fu等[48]采用平滑擾動分析法估計多服務臺單隊列的排隊系統的穩態系統時間的梯度。有關擾動分析法在排隊系統優化設計中的更多應用,可以參考文獻[49-53]。企業需要管理生產原料和成品的庫存,根據客戶的需求確定補充庫存的時間和訂貨數量。優化設計企業的庫存管理策略,可以避免企業超儲或缺貨,實現降低庫存成本的目標。Fu[54]考慮帶有隨機客戶需求和延遲交貨的(s,S)策略庫存系統,應用擾動分析法估計表現函數的樣本路徑梯度。Bashyam 等[55]采用擾動分析法估計每周期的期望費用關于參數s和S的梯度。擾動分析法在其他復雜系統管理中也有廣泛的應用。Heidergott[16]采用平滑擾動分析法估計維護系統修復時間關于系統閾值的梯度,優化設計維護系統的修復時間。Donohue等[56]應用擾動分析法得到收益最大化的生產線產能配置方案。Yan等[57]基于擾動分析法提出了一種求解近似最優閾值的算法,用于控制帶有兩臺串聯機器的不可靠柔性制造系統。Liberopoulos等[58]采用無窮小擾動分析法估計生產流的一階和二階梯度,用于生產流的設計和控制。Brooks等[59]采用無窮小擾動分析法估計異步傳輸網絡最小平均延遲的漸近無偏梯度。何寧等[60]采用擾動分析法估計衛星網絡信息流的延時關于帶寬的梯度,實現最優的信息流帶寬分配。

1.2.2似然比法 似然比法,也稱為分數函數法。該方法將系統期望表現對參數的依賴全部轉移至輸入分布的密度函數中,即參數θ是輸入分布的密度函數中的參數(分布參數)。似然比法最早在文獻[61]中進行討論,隨后,文獻[62-65]也獨立地發現并對該方法進行了研究。系統期望表現為

如果期望與導數可以交換,則系統期望表現的梯度為

稱h(x)為似然比估計量。似然比法可以解決樣本表現不連續的問題,但不能解決樣本表現帶有結構參數的問題。文獻[46,66]中聲稱在同一模型表示下,無窮小擾動分析法比似然比法方差小,該結論在仿真分析中被視為經驗準則。Cui等[67]給出了該經驗準則成立的充分條件,也給出了該經驗準則不成立的反例。

為克服似然比方法無法解決樣本表現帶有結構參數問題的方法局限,文獻[19]中提出了推出似然比法(Push-out LR),該方法利用變量替換將不連續的樣本表現中的結構參數推出至輸入分布,再采用似然比法;但采用變量替換推出的似然比法依賴于解析形式的變量替換,在實際應用中這種方法需要視實際問題情況而定,這給該方法的理論分析和實際應用帶來了困難。與無窮小擾動分析法相比,似然比法的無偏性條件通常更容易滿足。L’Ecuyer[68]在樣本表現帶有結構參數且連續的情況下,采用導數中的乘積法則,將無窮小擾動分析法與似然比法相結合,提出了IPA-LR 估計量。該估計量具有單樣本軌道的性質,但在處理樣本表現不連續的問題時是有偏的。

在復雜系統的優化設計中,似然比法不適用于估計系統表現關于系統結構參數的梯度。例如,在生產訂單隨機到達的(s,S)策略庫存系統中,似然比法不適用于估計庫存量關于參數s和S的梯度,但可以估計庫存量關于訂單隨機到達時間分布中的參數的梯度。Fu[1]采用似然比法估計了GI/G/1排隊系統停留時間關于隨機服務時間分布中的參數的梯度。Fu等[8]采用似然比法估計質量控制圖的平均運行長度關于控制邊界的梯度,用于控制圖控制邊界的敏感性分析,更加準確地判定系統是否處于受控狀態。

1.2.3弱導數法 弱導數法,也稱為測度求導法。該方法最早由Pflug[69-70]進行討論。與似然比法相似,系統的期望表現梯度可以表示為

稱c（θ）（h（X(+)（θ））-h（X(-)（θ）））為弱導數估計量,其中,c（θ）是常數,隨機變量X(+)（θ）的概率密度函數為,隨機變量X(-)（θ）的概率密度函數為由于概率密度(或概率質量)函數的分解不是唯一的,因而弱導數估計量不是唯一的。

與無窮小擾動分析法和似然比法相比,弱導數法的無偏性條件通常更容易滿足,但要仿真額外的樣本軌道。當參數向量的維度很高時,弱導數法需要的計算成本大。Heidergott等[71]采用弱導數法研究了分位數敏感性估計問題,解決了文獻[44]中估計量無法處理樣本表現不連續的問題,并將其應用于排隊系統中,優化了設計考慮系統停留時間分位數的系統服務速率。

在復雜系統的優化設計中,弱導數法不適用于估計系統表現關于系統結構參數的梯度,且不適用于高維生產參數的估計問題,因而有關弱導數法在復雜系統中的應用的文獻較少。Heidergott[17]考慮帶有老化更換策略的多組件維護系統的優化設計問題,應用弱導數法估計維護成本關于閾值參數的梯度。

1.2.4廣義似然比估計法 復雜隨機系統中經常存在不連續樣本表現與結構參數,例如,質量管理中控制圖的控制邊界屬于結構參數,它的表現函數關于控制邊界不連續。存在不連續樣本表現與結構參數時,經典梯度估計方法的無偏性無法得到保障,如何系統地解決該問題是仿真領域的重要難題之一,很多學者對這類問題進行了研究。Wang等[72]提出一種混合無窮小擾動分析與似然比的方法,但該方法處理不連續樣本表現的手段非常依賴于問題的特定結構,適用范圍有限(主要在金融領域)。Liu等[73]提出了一種核函數方法來處理不連續樣本表現,但該方法得到的隨機梯度估計是有偏的。文獻[71,74-75]也為解決該問題做了大量研究。這些方法在應用時需要視問題的特定情況分析,并且通常不具有單樣本軌道。

Peng等[11]首次提出了廣義似然比估計法,該方法是目前仿真優化領域解決復雜隨機系統包含結構參數和不連續樣本表現的單樣本軌道無偏隨機梯度估計問題中適應性最強的方法之一。擾動分析法與似然比法均可以放到廣義似然比法的理論框架中分析。“廣義”的含義是當樣本表現沒有結構參數時,廣義似然比估計量可以退化為似然比估計量。系統的期望表現梯度可以表示為

稱φ（h（X;θ））w（X;θ）為廣義似然比估計量,其中,φ是幾乎處處連續的可測函數,

d（x;θ）是與函數h、函數h的雅可比矩陣和樣本表現函數有關的量。廣義似然比估計法不但適用于絕大多數帶有結構參數的問題,而且適用于現有方法通常很難處理的不連續且非線性的樣本表現問題。該方法同時解決了無窮小擾動分析法不能處理不連續樣本表現,與似然比法不能處理結構參數的問題。該方法為無偏隨機梯度估計提出了新的估計方法或途徑,提供了擾動分析法與似然比法由于樣本表現對于結構參數不連續導致的偏差的解析表示,并提出了可積性條件用來補償光滑性條件從而得到無偏梯度估計的方法創新。

廣義似然比法在復雜系統的優化設計中具有廣泛的應用前景。廣義似然比法可以解決排隊系統分析、統計過程控制、供應鏈管理、制造工廠和服務系統等一系列梯度估計問題。Peng等[11-12,76]討論了廣義似然比法在質量控制圖、隨機活動網絡、維護系統和具有積壓的單一產品、離散時間及多期庫存模型中的應用。Zhang等[13]和張公伯等[14]提出了一種質量控制圖的經濟優化設計方法,采用廣義似然比法估計平均運行成本關于控制圖控制邊界的梯度。Zhou等[77]考慮在給定當前時刻信息的條件下,隨機系統在未來的預期性能,采用廣義似然比法估計了樣本表現條件期望的梯度。Peng等[78]采用廣義似然比法研究了復雜隨機模型的分布敏感性,導出了仿真模型輸出分布函數關于參數與變量各階導數的估計,可以解決生產服務與管理等復雜系統數據驅動的隨機建模,減小模型結構假設不準確導致的誤導性輸出結果。Glynn 等[79]采用廣義似然比法研究了基于分布敏感性的分位數敏感性估計。但在采用廣義似然比法進行分位數敏感性估計的問題中,仿真樣本有限時估計量存在偏差,因而在梯度下降算法中,計算估計量需要的仿真樣本數量隨著迭代次數增加,從而消耗大量的仿真資源。如何在分位數梯度優化中逐步消除分位數敏感性估計量的偏差值得進一步研究。

表1比較了幾種梯度估計方法及其在復雜系統管理中的應用場景。

表1 隨機梯度估計方法的比較

2 系統管理中的應用

隨機梯度估計在復雜系統管理中有著廣泛的應用。復雜系統的輸入和輸出的關系未知,輸入通常受到測量誤差、信息缺失等不確定性因素的影響。敏感性分析可以通過識別模型輸入來減少不確定性,提高系統的穩健性。隨機梯度估計方法可以估計復雜系統表現函數關于生產參數的敏感性,實現優化復雜系統結構、控制風險的目標。由于間接的梯度估計方法在實際問題總是可用,本節主要介紹直接的隨機梯度估計方法在隨機活動網絡、質量控制圖、排隊系統和(s,S)策略庫存系統4個系統中的應用。

2.1 隨機活動網絡

網絡計劃技術是制定項目進度計劃的一種常用技術。在項目實施過程中,受天氣、資金和材料供應等隨機因素的影響,生產活動持續時間具有隨機性。項目管理人員采用隨機活動網絡表示項目中計劃完成的各項生產活動的先后順序和邏輯關系,進而可以預測項目的工期和費用,優化資源配置,實現總工期最短的目標。應用隨機梯度估計方法,可以估計項目的工期和費用關于某項生產活動的參數的梯度。如下考慮一種簡單的隨機活動網絡模型[23],并介紹隨機梯度估計方法在該模型中的應用。

隨機活動網絡由具有M個節點、N條邊的有向非循環圖表示,其中,每條有向邊表示生產活動,生產活動時間為隨機變量Xi,其分布為Fi,i=1,2,…,N,X1,X2,…,XN相互獨立。假定節點1為任務起始點,節點M為任務終止點,P為構成起始點至終止點路徑的有向邊的集合,為所有路徑的集合,P*∈為最優路徑(活動時間最短或最長),最優路徑的時間可以表示為:Y(X)=問題的目標是估計d E[Y]/dθ,其中,θ為Xi的分布參數。應用隨機梯度估計方法,可得無窮小擾動分析估計量為

其中,1{·}是示性函數。似然比估計量和弱導數估計量分別為:

數值算例。如圖1所示,考慮包含5個節點、6條邊的隨機活動網絡。假定節點5為任務終止點,所有路徑的集合為假設生產活動時間Xi服從均值為i的指數分布,即令θ=θ1=1為X1中的分布參數,P*為使生產活動時間最短的最優路徑。應用隨機梯度估計方法,可得無窮小擾動分析估計量、似然比估計量和弱導數估計量分別為:

圖1 隨機活動網絡

實驗分別在指數分布和正態分布的設定下,比較3種無偏隨機梯度估計量和單邊前向有限差分估計量FD(c)。基于106次獨立的仿真實驗,得到各估計量的數值(平均值±標準誤差),其中,標準誤差由估計的標準偏差除以樣本數的平方根計算得到。

由表2可以觀察到,當攝動值為0.1時,有限差分估計量的偏差較大;當攝動值為0.01時,有限差分估計量的方差較大。無窮小擾動分析估計量和弱導數估計量的方差比似然比估計量的方差小。基于108次獨立的仿真實驗,在指數分布的設定下,FD(0.001)的梯度估計量的數值為0.585±0.103 8;在正態分布的設定下,FD(0.001)的梯度估計量的數值為0.843±0.126 2(平均值±標準誤差)。這表明,無窮小擾動分析估計量,似然比估計量和弱導數估計量具有無偏性。

表2 隨機活動網絡最優路徑時間的梯度估計量(平均值±標準誤差)

2.2 質量控制圖

在生產過程中,受生產人員、機械設備、材料、工藝方法和環境等因素的影響,生產過程可能處于失控狀態,導致產品的質量特性值偏離規定的質量特性值。質量控制圖是統計過程控制的一種重要方法,它利用控制圖提供的信息判定生產過程是否處于統計控制狀態,有助于將生產過程維持在受控狀態,減少生產階段因產品異常帶來的損失。生產過程輸出的質量特性值有隨機性,應用隨機梯度估計方法,可以估計生產過程平均運行長度關于質量控制圖的控制邊界的梯度。如下考慮指數加權移動平均控制圖 (Exponentially Weighted Moving Average,EWMA)[13-14],并介紹隨機梯度估計方法在該模型中的應用。此外,通過一個考慮最小化平均運行成本的控制邊界優化設計的數值算例,以期為隨機梯度估計方法在復雜系統管理中的具體應用、改善的指標和得到的結論提供范例。

在統計過程控制中,生產過程的平均運行長度被定義為

其中:N為系統失控的時間,即N=min{i:Yi?[θ1,i,θ2,i]};Yi為在得到第i個樣本后觀測到的檢驗統計量,對于指數加權移動平均控制圖:Yi=αXi+(1-α)Yi-1,i＞1,Y1=X1。Xi為第i個樣本輸出,θ是一類參數,可以是θ1,i或θ2,i,分別為第i個檢驗統計量的控制下界和控制上界。當系統處于受控狀態和非受控狀態時,系統將輸出具有不同分布的樣本。取條件于Z=z,其中,Z=R/Δ,R為系統受控的隨機持續時間(不可觀測),服從密度為q0(·)的分布,Δ為采樣間隔,則Xi的條件密度為

其中,f1(·)和f2(·)分別為系統受控和不受控時輸出Xi的密度函數。假定系統失控時間服從參數為λ的指數分布,受控樣本輸出服從均值為μ0、方差為σ2的正態分布,失控樣本輸出服從均值為μ1、方差為σ2的正態分布。令θ1,i=1,θ2,i=2,問題的目標是估計d ARL（θ）/d2。由于結構參數和不連續樣本表現的存在,無窮小擾動分析估計量,似然比估計量和弱導數估計量均不存在,由此可得廣義似然比估計量為

2.2.1數值算例 考慮如下參數設定的指數加權移動平均控制圖:α=0.5;系統失控時間服從均值為20的指數分布,即λ=0.05;采樣間隔Δ=1;控制下界1=-1,控制上界2=1;受控或失控樣本方差σ=1;受控樣本均值μ0=0,失控樣本均值μ1=1,μ1=3。通過實驗比較廣義似然比估計量和單邊前向有限差分估計量FD(c)。分別基于104和106次獨立的仿真實驗,得到各估計量的數值(平均值±標準誤差),其中,標準誤差由估計的標準偏差除以樣本數的平方根計算得到。

由表3可觀察到,當攝動值為0.1時,有限差分估計量的偏差較大;當攝動值為0.01時,有限差分估計量的方差較大。基于108次獨立的仿真實驗,在失控樣本均值μ1=1的設定下,FD(0.001)的梯度估計量的數值為10.2±1.12;在失控樣本均值μ1=3的設定下,FD(0.001)的梯度估計量的數值為7.9±1.02(平均值±標準誤差)。這表明,廣義似然比估計量具有無偏性。

表3 指數加權移動平均控制圖的梯度估計量(平均值±標準誤差)

2.2.2控制上界的優化設計 生產過程的單位時間內的平均運行成本被定義為

考慮一個指數加權移動平均控制圖,假設失控樣本均值μ1=3,方差σ=2,其余參數設定與例2相同。令損耗成本c=2,系統維修成本C=10,隨機逼近方法中控制圖的初始控制上界為θ(0)=1,迭代步長an=1/n,第n步迭代使用(1 000+n)次獨立的仿真實驗估計梯度ΔE[N]和ΔE[(N -R/Δ)+]。采用隨機逼近方法迭代10 000次得到的控制上界的樣本軌道如圖2(a)所示。基于108次獨立的仿真實驗估計系統的平均運行成本,得到的平均運行成本的樣本軌道如圖2(b)所示。

求解得到的控制圖最優控制上界的值為θ*=1.139。由圖2(b)觀察得到結論:通過控制上界的最優設計,生產過程的單位時間內的平均運行成本從1.52降低為1.32。采用基于廣義似然比梯度估計的隨機逼近方法,解決了控制圖控制上界的優化設計問題,降低了生產過程的單位時間內的平均運行成本,有助于提升企業經濟效益。

圖2 樣本軌道:(a)控制上界(b)平均運行成本

2.3 單服務臺排隊系統

排隊論又稱隨機服務系統理論,是通過分析隨機服務系統中排隊等待現象的穩態概率特征進行系統優化設計和控制的理論。排隊論廣泛應用于制造業產品生產、運輸、庫存等資源分配管理的隨機服務系統。例如,在產品加工生產線流程設計問題中,可以構建加工設備服務待加工制品的排隊系統;在企業設備維修問題中,可以構建維修員工服務故障設備的排隊系統;在產品運輸的集裝箱港口資源配置問題中,可以構建堆場服務集裝箱的排隊系統。排隊系統中存在著大量隨機性,如產品到達時間間隔、工作站服務時間、系統等待時間和系統停留時間等,應用隨機梯度估計方法,可以估計排隊模型表現函數關于模型參數的梯度。如下考慮一種簡單的單服務臺排隊系統[23],并介紹隨機梯度估計方法在該模型中的應用。

考慮一個先到先服務的G/G/1 的排隊網絡。記Ai為生產訂單的間隔到達時間,Xi為生產加工時間,Ti為生產訂單在排隊系統的停留時間。前N個訂單在排隊系統中的平均停留時間可表示為:令θ為Xi的分布Fi中的參數,問題的目標是估計假定系統初始狀態是空的,生產訂單到達過程與服務過程獨立,且到達過程互相是獨立的。根據Lindley方程,一個生產訂單在排隊系統的停留時間可表示為:Ti+1=Xi+1+(Ti-Ai+1)+。應用隨機梯度估計方法,可得無窮小擾動分析估計量為

其中:M為觀測到的忙期數;nm為在第m個忙期中的最后一個訂單;似然比估計量和弱導數估計量分別為:

2.4 (s,S)策略庫存系統

庫存管理是企業物流和供應鏈管理中重要的研究問題。企業的庫存是為了保持持續生產、滿足客戶需求而存儲的原材料、在制品和制成品。庫存管理是通過對客戶需求的預測,在保證企業生產活動的正常進行和產品供應的前提下,加速物資流動速度,減少庫存量,實現降低庫存成本、提升企業經濟效益和提高生產柔性的目標。供應鏈的庫存管理通過企業間的協調配合和信息共享,減少牛鞭效應帶來的需求誤差,實現按需生產和降低庫存總成本的目標。庫存問題的隨機性來自供應者、生產者和客戶的不確定性。應用隨機梯度估計方法,可以估計庫存量和庫存成本關于庫存策略中的參數的梯度。如下考慮一種簡單的(s,S)策略庫存系統[22],并介紹隨機梯度估計方法在該模型中的應用。

考慮一種生產原材料的定期檢查(s,S)策略庫存問題,在每個周期檢查庫存量,確定是否需要補貨。(s,S)訂購策略為如果當前原材料庫存量低于s,則補貨到庫存量S。假設所有的超額需求被積壓,會在下一周期被滿足,原材料的訂單沒有延遲,在每一周期結束時決定是否補貨。令Dn為第n周期的需求量,,Vn為第n周期結束時原材料的庫存量,可以表示為

其中:E[D]為平均需求量;ξn=Yn-s;Yn為滿足需求前的原材料庫存量。

3 結論

隨著新一代信息技術與先進仿真技術的集群式創新、融合發展與突破,大規模復雜系統的仿真數據日益可測可獲,不僅有效提升企業生產質量、效率和效益,也孕育著企業管理決策理論與方法的重大變革,推動管理決策研究向仿真優化領域轉變。仿真技術為企業提供了能夠覆蓋產品論證、研發、設計、采購生產、銷售和服務等供應鏈各階段的全生命周期解決方案,是企業數字化轉型升級的重要技術手段。復雜隨機系統的梯度估計包含比系統表現估計更多有用的信息,梯度估計的準確性常常影響著系統表現的優化分析及后續管理決策的進行。為提升仿真及優化效率,需結合問題情境特點選取高效準確的隨機梯度估計方法。

在實際應用中,根據復雜系統的特定情況選擇合適的隨機梯度估計方法對后續優化分析及管理決策的制定至關重要。在進行方法選擇時,可從以下幾方面考慮。首先,如果沒有獲取到復雜系統的信息,優先考慮間接的梯度估計方法。進一步,如果獲取到復雜系統的一些信息,如系統輸入的分布或樣本表現的特征等信息,則可以考慮應用直接的梯度估計方法。對于不連續型樣本表現,如示性函數或非平滑函數,無窮小擾動分析法不適用。如果樣本表現帶有結構參數,則似然比法和弱導數法不適用。平滑擾動分析法需要視實際問題情況選取適當的條件計算條件期望,會消耗額外的仿真量。推出似然比法依賴于解析形式的變量替換來推出結構參數。對于樣本表現不連續且帶有結構參數的問題,優先考慮廣義似然比法。其次,選擇不同隨機梯度估計方法的另一影響因素是參數向量的維度。有限差分法和弱導數法在處理高維參數向量的梯度時會消耗很多額外的仿真量。在基于梯度的參數估計方法中,選擇合適的隨機梯度估計方法需要同時考慮估計量的偏差和方差。有限差分法得到的梯度估計量是有偏差的。一般情況下,無窮小擾動分析法比似然比法方差小。

基于梯度的參數估計方法的效果提升及其在復雜系統管理實踐中的應用是仿真優化領域關注的方向。

(1) 隨機梯度估計方法與強化學習和人工智能算法的結合可以改善大型裝配線在復雜環境中的質量控制任務,提高系統智能決策的抗干擾能力,實現人機共融的協同決策。深度學習與強化學習中訓練人工神經網絡所采用的反向傳遞算法屬于仿真優化中研究的隨機梯度估計算法。開發AlphaGo的谷歌公司Deep Mind 研究團隊在文獻[80]中對隨機梯度估計方法在機器學習中的應用做了長篇綜述。Peters等[81]討論了有限差分法和似然比法在強化學習策略梯度估計中的應用。Peng等[82]采用推出似然比方法訓練人工神經網絡,該方法通過在神經元傳遞信號中加入噪聲,不但可以處理離散型激活函數和損失函數的問題,也提高了人工神經網絡在對抗攻擊時的魯棒性。

(2) 隨機梯度估計通常可以用于參數的敏感性分析和基于梯度的仿真優化算法的輸入,實現優化復雜系統結構的目標,進而度量和控制復雜系統的風險。隨機梯度估計方法在項目計劃制定、產品生產、運輸、庫存、設備維修和生產過程控制等系統管理領域有廣泛的應用。經典的隨機梯度估計的目標是隨機模型樣本表現的期望,而期望刻畫的是系統在隨機環境下的平均表現,未反映系統在極端環境下的表現。考慮隨機模型樣本表現的分位數的隨機梯度估計可以得到在極端環境下更穩健的系統管理策略。在供應鏈系統管理中,考慮報酬收益分位數敏感性的報童模型可以得到在極端市場需求出現時表現更穩健的訂購策略。Hu等[83]研究了以分位數為目標的仿真優化問題,基于廣義似然比法提出了單樣本軌道多尺度的隨機逼近方法,解決了分位數敏感性估計中仿真樣本軌道的不連續性以及由比值引起偏差的困難,建立了算法的收斂性并分析了其漸進收斂到平衡點的穩定性。Jiang等[84]研究了累計回報的分位數敏感性的強化學習問題,基于似然比法提出了兩尺度的隨機逼近方法,得到了在極端事件下表現更穩健的策略。近兩年疫情導致的極端事件對全球供應鏈造成了巨大影響,極端事件下的系統管理將會是后疫情時代的重點研究方向之一。

綜上所述,本文在系列相關研究的基礎上對系統管理領域中隨機梯度估計的應用進行了梳理,同時還介紹了常用的隨機梯度估計方法和系統管理模型。針對常用的隨機梯度估計方法,分別總結了其研究狀態和應用場景,并通過模型介紹了這些方法在企業項目進度計劃制定、生產過程控制、生產調度流程優化和供應鏈庫存管理領域的應用。后續研究進一步探索隨機梯度估計方法的求解路徑和優化設計,如有限差分法最優攝動值的選取,并進一步思考相關方法在復雜系統管理實踐中的具體應用。