一種基于概率關聯的局部高斯過程回歸算法

肖文鑫 張文文

高斯過程回歸(Gaussian process regression,GPR)廣泛應用于機器人、控制系統和航空航天飛行器等領域,如機器人的逆動力學模型估計[1-3].它是僅利用簡單的線性代數處理非線性模型的有效工具,提供了一種簡單但有效的方法來表示數據的先驗分布,其特點是使用較少的參數得到準確的學習結果[1].作為應用范圍廣泛的非參數回歸,它在樣本外預測測試中脫穎而出,已成為機器學習的重要手段[4].然而,隨著現實世界中訓練集和測試集的增長,高斯過程回歸模型在面對大數據量時暴露了缺點.當數據量n很大時,標準高斯過程模型的立方復雜度 O(n3)會導致不可接受的計算量,使其不適用于高維系統或在線學習[5].

降低這種復雜度的解決方案主要分為兩種類型,其中一種方法是稀疏高斯過程回歸和高斯回歸混合模型[6-7].另一種是將數據分割并使用多個局部模型代替全局高斯過程模型,這種方法稱為局部高斯過程回歸(Local GPR,LGPR).在稀疏高斯過程中,按照設定的規則減少輸入點的數量,以降低計算復雜度[8].局部高斯過程回歸用多個局部子模型替換整體高斯過程模型,由于每個局部模型具有較少數量的數據點并且被單獨訓練和更新,因此總體計算成本降低[4].對于局部高斯過程回歸,最關鍵的是如何劃分數據并將其分配給不同的模型.

根據最近的文獻,夏嘉欣等[9]提出一種基于噪聲輸入的稀疏高斯算法,并將其應用于解決人體姿態估計問題.喬少杰等[10]利用高斯過程回歸預測移動對象最可能運動軌跡.王傳云等[11]對圖像局部區域建立高斯混合模型解決圖像的背景模型構建問題.Da 等[12]研究了一種基于單任務聚合方案的局部高斯過程回歸模型.Liu 等[13]提出了一種基于流形學習方法的局部高斯過程相似度維持算法.Binois 等[14]使用基于似然的高斯過程回歸方法,將異方差性問題表征為單目標函數,使用完整的閉合形式導數,實現了基于庫的優化.任志剛等[15]利用加權的優秀樣本預估高斯均值改進了傳統高斯分布估計算法.Sarkar 等[16]使用高斯過程回歸構建了一種魯棒算法,能夠處理不確定性和數據中的噪聲,并驗證了相比經典諧波分析方法的優勢.文獻[17]中的工作通過測量距離來選擇合適的局部模型并最終計算加權預測,但該方法在每次預測時仍然考慮所有模型,因此不完全遵循局部支持的概念.

為了進一步降低計算量,本文擬提出一種基于邊界平滑過渡的局部高斯過程回歸模型,即訓練集中的數據不是確定性地插入到單個局部模型中,而是按概率分布分配給若干相鄰的局部模型.這意味著數據點屬于某一模型的概率分布在特定區域內是連續的,而在超出該區域的概率分布則為零.這種分布將每個模型限制在一個定義清晰的區域,通過僅考慮實際相關且可能對所考慮的數據點有貢獻的模型來提高計算效率方面的性能,也就是說對于特定的數據點只需要考慮少量數據點所屬范圍內的局部模型.另外,模型更新過程中使用矩陣塊的求逆方法進行矩陣求逆,大大減少了計算量.為了實現更新與預測階段局部模型的快速檢索,本文還引入并利用了來自計算幾何領域的k-d 樹最近鄰搜索算法[18].

1 局部高斯過程回歸模型

1.1 高斯過程回歸

高斯過程回歸是機器人和控制理論領域中用于監督機器學習的常用技術,它的優點是通過少量的超參數實現高預測精度.但是,當數據點的總數n過大時,高計算復雜度 O(n3)將帶來嚴重問題.對于由n個數據點組成的數據集 {X,y},每個數據點由d維輸入數據x和一維輸出數據y構成.高斯過程回歸假設輸出數據由函數f(x)和附加噪聲?組成

其中,yi是輸出集y中的第i點,xi是輸入集X中的第i點,?i是均值為零、方差為的正態分布噪聲,n是數據點的數量.輸出集滿足下面的正態分布[19]

其中,K(X,X)是通過所有訓練點計算得到的協方差矩陣,I表示單位矩陣.輸出訓練集y和對問詢點x*的預測輸出f(x*)的聯合分布由下式給出

f(x*)的條件分布服從一個新的正態分布

其中,f(x*)的均值和方差的表達式為

因此,對問詢點x*的預測值可以通過下面表達式計算

對于內核函數k(xp,xq)的選擇,本文采用充分靈活的平方指數協方差函數,其定義為

使用文獻[20]中提出的共軛梯度算法最大化對數邊際似然

其中,C=K(X,X)+.由以下對數邊際似然導數的表達式可以推導出對數邊際似然取最大值時的超參數為

1.2 局部高斯過程回歸

由于高斯過程回歸在訓練點數量大時表現不佳,為了在保持高斯過程回歸的優勢的同時降低計算復雜度,提出了局部高斯過程回歸.局部高斯過程回歸可以實現模型的快速更新,因此更適合使用在需要通過實時產生的數據不斷更新模型的領域.局部高斯過程回歸方法背后的基本思想是引入多個回歸模型并減少每個模型的數據點數量,從而降低計算復雜度.該方法將數據點分配到多個局部模型并單獨訓練[1].對應于平方指數協方差函數,點與模型中心之間的距離通過如下距離測量來計算,以將數據點分配給最近的局部模型

其中,ci表示的是第i個模型的中心,wi是模型i和給定x之間的距離.因此,為了正確分配新數據點,這個方法必須為每個數據點計算M次距離,其中M是局部模型的總數.基于此度量,選擇具有最小距離的模型,并且將點添加到該最近的模型.一旦選擇了模型,該方法后續步驟與第1.1 節中的常規高斯過程回歸完全相同.如果沒有足夠接近新數據點的模型,則會創建一個僅包含此新點的新模型.亦即如果新點到所有中心的距離大于某個閾值,則創建一個新模型,其中心是新數據點.

式(12)中的距離度量也用于加權預測.使用式(7)對第i個局部模型的預測值進行計算,然后通過對這些模型的加權預測求和來獲得全局預測值

由于每個局部模型中的數據點的數量比高斯過程回歸中的低,用于更新模型的計算工作量得以減少,并且與高斯過程回歸相比可以實現整體性能改進.在第2.1 節,本文提出的算法對以上方法進行了改進,采用概率數據關聯將數據點在其最近的 2d個局部模型中進行分配.

1.3 高效k 近鄰查詢樹

為了實現快速搜索算法,本文使用k-d 樹來存儲數據.利用k-d 樹進行局部模型檢索,可以有效提高模型更新與對未知數據點預測的速度.k-d 樹是二叉搜索樹(Binary search tree,BST)的多維形式,而二叉搜索樹是用于快速可靠查詢的合適方法[18].文獻[21]中提出了這個k-d 樹的優化版本.

1.3.1 k-d 樹的創建

對于樹的第零層,首先在第1 個維度分割數據,并使用此維度中數據的中值作為分離的邊界.此步驟將數據集分為兩部分:一部分是由小于或等于中值的數據組成的樹的左側;另一部分是剩余數據構成的樹的右側.在每個子樹上重復此過程,直到只剩下包含一個點的節點,將其作為樹的子節點.由于在每個節點中重復劃分并將剩余的邊處理為新的子樹,因此這是一個遞歸定義.分割的維度是根據樹中當前節點的級別來決定的,即在層級i上,選擇維度i+1 進行分割.如果迭代遍歷了所有維度,則在樹的下一層過程將再次從第1 個維度開始分割,直到只剩下具有一個數據點的子集表示樹的子節點.

1.3.2 k-d 樹最近鄰搜索

k-d 樹最近鄰搜索是一種快速且可靠的搜索算法,其基于預定義的相似性或距離度量來判斷距查詢點最近的k個局部模型.對于如何找到最佳匹配和k個最近鄰項存在多種方法[22-23],這項工作將搜索算法制定為遞歸策略,其中算法的輸入是每次遞歸對應的節點和其子節點.首先,包含所有數據的根節點作為首次的輸入傳遞到算法.當前最近鄰項的候選項被存儲在優先隊列中,其中查詢點和近鄰項之間的相似性是判斷的依據.接下來會將每個節點與優先隊列中的條目進行比較,如果該節點比隊列中的一個條目更接近查詢點,則相應地更新優先隊列使該節點成為新的候選項.如果節點不是終端,即不是葉子節點,則查詢點所在的子樹再次遞歸地傳遞給算法.如果節點是葉子節點,則此子樹上的遞歸結束.k-d 樹最近鄰搜索的計算復雜度是O(k×n×logn),其中k是數據的維度,n是數據點的個數.

2 基于概率關聯進行模型選擇的局部高斯過程回歸

本研究提出一種無硬邊界局部模型概率關聯方法,這種方法確保了在僅有局部支持的條件下實現模型之間的平滑過渡.模型建立階段,通過定義模型相鄰部分的概率分布,數據點被根據其坐標分配到局部模型中.對未知數據預測階段,由距離查詢點最近的若干模型計算加權平均值,從而求和得到全局預測值.因此,該算法包括3 個階段:1)初始模型的建立;2)局部模型的更新;3)測試集的加權預測.第2.1 節將討論局部模型的概率型數據關聯,第2.2 節和第2.3 節描述了模型的建立和模型的更新,第2.4 節提出了對未知數據的預測方法.

2.1 邊界約束的局部高斯過程回歸模型

緊致支持的邊界約束局部模型是本在線學習模型的基本組成部分,本節介紹如何構建本地模型.

根據文獻[1,17]局部高斯過程回歸方法通常使用模型內數據點的平均值作為每個模型的中心.因此,對于每個新輸入的數據點,這些中心都需要更新,即局部模型需要持續地改變.為了減少每次更新造成的計算量,本文提出一種固定的模型網絡,同時為了實現局部模型與數據之間的概率關聯,我們提出了一種激活函數.激活函數φ(ci|x)表示數據點x屬于局部模型ci的概率,滿足概率分布的條件.距數據點最近的 2d個局部模型對x的激活函數可以組成一個離散概率密度函數f(c),其中,d是輸入數據點x的維度,c=c1,···,ci,···,,ci是局部模型i的中心,本文以模型中心代指具體的模型.最后,根據f(c)中的概率分布產生一個偽隨機數,以決定數據點 (x,y)分配到哪個局部模型中.

所有局部模型組成一個覆蓋數據集空間的模型網格,在維度為一維時,該網格由每個中心之間的距離w定義.當維度k增加時,在所有維度中應用這種設置來延伸網格,從而生成由各局部模型組成k維立方體.

相鄰的局部模型在共享邊界處做出不同的預測,導致邊界處的預測是不連續的[24].實際上,對于幾個模型的交匯處的點,屬于每個模型的概率可能都大于 0,這意味著它不是絕對屬于最近的模型.然而,當局部模型被設計為具有硬邊界時,預測值僅由最近的模型確定,這可能導致大的誤差.為了實現相鄰局部模型之間的平滑過渡,在創建局部高斯過程回歸模型時,本文中使用緊湊支持的邊界約束來構建局部模型.這里采用激活函數定義模型的邊界約束,以表示某一坐標的數據點屬于這個模型的概率.在模型中心位置,定義只有這個模型是“激活”的,這意味著在這個位置模型ci的激活函數φ(ci|x)等于 1,其余模型的激活函數等于 0,位于這個位置的數據點只可能屬于這個模型.通過這樣的設置,局部模型與其相鄰的各模型共享邊界區域.由于本算法面向的是工程中在線更新的場景,主要處理的是多維數據,因此激活函數應當盡可能考慮對稱的形式,以避免多維空間中模型構造過于復雜.基于通過簡潔的結構實現更好的性能的理念,將共享區域的激活函數定義為范圍從 0 到 1 的斜坡函數.這樣,所有局部模型使用相同形式的激活函數.如圖1 所示,激活函數可表示為

圖1 一維局部模型激活函數示意圖Fig.1 Schematic diagram of one-dimensional local model activation function

其中,x是數據點在當前維度的值,ci是局部模型的中心在當前維度的坐標,w是相鄰局部模型的間距.由圖中可以看出,每個局部模型在當前維度的寬度為 2w.對于處于相鄰模型中心ci和ci+1之間的x,滿足

由此該激活函數實現了數據點與鄰近局部模型的概率關聯.

圖1 為一維情況下局部模型的初始化示意圖.基于這種一維方法,本文開發了具有更高維度的數據的解決方案.對于模型布局和激活功能,基本原理保持不變.但是,每個數據點可能被分配進的模型最多為 2d個.將一維網格布局自然延伸至d維,得到多維模型網絡.

圖2 展通過在二維空間內構建c1,c2,c3,c4為中心的4 個局部模型,展示了訓練集數據為二維時的局部模型布局和可能激活的模型的情況.示意圖中的c1,c2,c3,c4分別為4 個局部模型的中心.圖中的區域A由于不與其他局部模型相交,此區域內的點只可能屬于模型c1,區域B與一維時的情況相同,最多只有兩個模型處于激活狀態,而這些激活中心c1和c2共享相同的第2 維分量.對于區域A中的點,可能的激活局部模型僅是c1模型;對于區域B,可能的激活模型是c1和c2模型.C區域則有4 個可能激活的模型.另外,該圖還顯示了各局部模型在二維立方體的等距網格中的位置.

圖2 二維局部模型分布示意圖Fig.2 Schematic diagram of two-dimensional local model distribution

對于多維數據,使用輸入數據點在每個維度上的最大值和最小值來計算各局部模型的中心.局部模型的數量和中心可表示為

其中,numi是局部模型在維度i上的個數,和分別是輸入數據在維度i上的最大值和最小值,是相鄰模型在維度i上的間距,ci,n表示的是第i維第n個模型中心的位置.維度i上模型中心的最小值設置為,以保證所有數據點在各維度屬于各局部模型的概率和始終為1.因此,圖1 中的最外側的兩個局部模型均只利用了靠內的一半.擴展到圖2 的情況,數據點范圍是集中在C區域中的,其他區域由于處在整個模型的邊緣部分,未加利用.當模型網格擴大時,始終未使用整個模型的邊緣區域.計算出局部模型網絡中所有模型中心的坐標后,用這些坐標值建立k-d 樹,每個模型都用一個葉節點表示.通過這個數據結構,可以使用文獻[22]中的快速搜索算法來執行有效的最近鄰居查詢,以便后續更新模型和預測未知輸出時搜索模型時使用.

初步確定模型的位置后,根據式(14)的激活函數將訓練集中的所有數據點分配到各個局部模型中.在每個維度上模型中心之間的距離都是w,使得局部模型的范圍覆蓋邊長為w的d維立方體.這些立方體的中心即為局部模型的中心.在整個模型的非邊緣區域,每個位置均由 2d個立方體覆蓋,也就是說可以有 2d個激活模型.為了得到d維情況下的概率分布,首先使用式(14)在每個維度上計算局部模型對應的激活函數φt(ci|x),它們僅取決于該特定維度中的位置,使得它們符合獨立分布的概率.然后通過這些一維值的乘積來計算整體激活函數φ(ci|x),局部模型ci的全局激活函數由其所有子概率的乘積得到,即

2.2 局部模型的初始化

本文提出的在線學習模型由初始化和更新階段兩個階段建立,本節解釋了如何使用訓練集構建初始局部模型.

為了初始化模型,本研究使用了一種對訓練集數據點進行聚類從而創建局部模型的方法.首先,通過訓練集中的數據計算各個局部模型中心的坐標,構建第2.1 節中所述的模型網絡.然后將所有數據點分配到相應的局部模型中.為了確定最終分配數據點的模型ci,必須在具有非零概率的激活模型中確定一個模型.如果直接選擇最有可能的模型,會在模型相鄰處形成硬邊界,導致預測在這個區域的不連續性.為了解決這個問題,對于訓練集中的數據點 (x,y),通過第1.3.2 節介紹的的k-d 樹實現k-最近鄰搜索,找到距離數據點最近的 2d個局部模型的編號.通過式(18)計算數據點對應這些局部模型的激活函數φ(ci|x)的值,可得到至多 2d個非零概率,構造由這些概率值組成的離散概率密度函數f(c),其中的概率值分別對應一個局部模型.根據各模型對應的概率從f(c)中抽取隨機數i,確定激活的局部模型ci,將數據點分配到這個模型.每個數據點僅分配給一個模型,但與其他局部高斯過程回歸方法不同,這種關聯是用概率方法完成的.完成訓練集中所有數據點的分配后,通過第2.1 節中的方法計算每個局部模型的預測矢量α,并創建相應的LGPR 模型.

2.3 局部模型的更新

本文提出的緊致支持的邊界約束局部模型最突出的特點是支持在線學習,這意味著局部模型能夠通過實時數據不斷進行更新.模型的在線更新可以應對一些現實世界的工程問題,例如機械特性和機器人交互隨時間的變化[25].在更新階段,依據第2.2節中分配數據點的方法,新輸入點 (x,y)分配到一個局部模型,需要對這個模型的參數進行更新,供預測階段調用.對未知點的預測是通過式(7)計算得到,其中協方差函數K(x*,X)只能在查詢點x*已知的情況下計算,而預測矢量α=(K(X,X)+)-1y可以在模型更新階段利用已知數據更新,供預測階段直接調用.

為了更新預測矢量,首先,將局部模型的數據集更新為

其次,預測矢量α中更新后的協方差矩陣Kn(Xn,Xn)+項可以表示為

其中,K(X,X)+是更新前的協方差矩陣,K(X,x)是由X和x的協方差組成的矢量,K(x,x)+由x的自協方差和模型的超參數相加得到.使用式(20)對協方差矩陣求逆,得到

其中

根據以上方法得到 (Kn(Xn,Xn)+))-1和yn后,可以通過下式計算更新的預測矢量,即

由式(21)可知,每次在計算新的協方差矩陣的逆時,都利用了它在更新前的值以減少計算量.因此,在模型的更新階段需要更新并儲存數據集、協方差矩陣的逆和預測矢量.

如果在使用新數據點更新模型階段檢測到數據點不在任何已有局部模型范圍內,即在儲存局部模型的k-d 樹中沒有搜索到足夠接近數據點的模型,則會創建一個僅包含此新點的新模型,其中心是新數據點,并補齊模型網格的其余部分.最后使用第1.1 節中的方法計算這個模型的超參數θ和預測矢量α,并創建LGPR 模型,計算新建模型的中心坐標和間距參數.

2.4 局部模型對未知數據的預測

基于輸入數據預測輸出數據是本局部高斯過程回歸模型的一個重要功能.在多個模型處于激活狀態的區域,可以通過激活函數計算查詢點屬于各局部模型的概率,再通過各模型的預測結果計算加權預測.首先,通過k-d 樹搜索查找到距查詢點x*最近的 2d個局部模型,使用激活函數φ(ci|x)算出各模型對應的概率,這些概率值組成離散概率密度函數f(c).基于這些概率值和這 2d個局部模型各自的預測值,可以計算得到全局加權預測值.如第1.1 節所述,局部模型ci的預測遵循高斯分布,其平均值為

其中,x*是輸入的查詢點,Xi是分配到模型ci中的輸入數據點集,αi是第i個局部模型的預測矢量.

通過每個局部模型對應的概率對其預測值進行加權求和,計算對查詢點x*的預測輸出y*.如前所述,激活模型ci的概率P(ci|x)通過激活函數φ(ci|x)表示,可以導出y*的期望值的表達式為

因此,最終的全局預測值表示為

3 實驗結果及分析

本節通過實驗,在更新時間、預測時間和誤差三個方面評估邊界約束的局部高斯過程回歸算法的性能.第3.1 節通過人工生成的低維數據,將本文算法與全局高斯過程回歸方法和具有硬邊界的局部高斯回歸過程方法進行比較,驗證了低維數據時本算法的優勢.第3.2 節使用來自F16 飛機的真實多維工程數據進行實驗并將本算法應用于機械手臂上,分析了本算法在多維數據和大數據量的場景下,實現在線學習的可行性.

3.1 低維數據實驗

為了驗證該基于邊界約束的局部高斯過程回歸算法的有效性,本文通過模擬帶有噪聲的高斯過程對該算法進行測試,并與全局高斯過程回歸方法和具有硬邊界的局部高斯回歸過程方法進行比較[17,19].對于訓練和測試,數據集由正弦函數的組合產生,超參數通過第1.1 節中的方法計算.本算法的性能取決于局部模型間距參數的選擇,本文首先評估它們對預測結果的影響.通過下面正弦函數組合生成二維數據集

并為輸出y附加σn=0.1 的噪聲.首先,將各維度的輸入值都限定在 0 到 40 之間,取 6 400 個數據點作為總的數據集.從中隨機選取 80% 作為訓練集,剩余部分用于測試集.使用這組數據重復模擬邊界約束的高斯局部模型 50 次,得到了不同參數情況下的平均預測誤差與運行時間,如圖3 所示.

圖3 局部模型參數對邊界約束模型性能的影響Fig.3 Influence of local model parameters on performance of boundary constraint model

圖中展示了局部模型間距以8 為間隔由4 增加到44 時,用單個數據點更新模型和預測未知輸出所用時間及預測誤差值的變化情況.小的局部模型寬度意味著總的模型數量多,即每個局部模型分配的數據點較少,對應著較大的誤差和較短的模型更新和預測時間.當局部模型寬度增加時,每個局部模型分配的數據點也隨之增加,使得預測結果更為準確.但同時由于求協方差逆矩陣的計算工作量隨著更多數據點而上升,消耗的時間也隨模型寬度增加而延長.因此,模型設計過程中需要對預測精度和計算時間做出權衡.由圖中可見,當模型寬度增加到20 以后,模型預測的準確度變化不大,但更新和預測時間仍然隨模型寬度增加而大幅增加,因此20 是這個數據集范圍內較好的局部模型寬度參數.

為了直觀地比較3 種方法的預測結果,本文使用以下等式生成一維輸入與輸出數據集

圖4 展示了3 種方法的預測情況,其中本文提出的算法的預測誤差為0.03292,全局回歸方法的誤差為0.03239,硬邊界的局部回歸方法的誤差為0.50499.在這個數據集中,本算法實現了與全局高斯過程回歸同一量級的精度,比硬邊界局部高斯過程回歸方法的精度高一個數量級.

圖4 對一維測試集的預測結果Fig.4 Prediction results for one-dimensional test sets

為了進一步比較3 種算法的性能,本文使用式(24)生成的二維數據,將3 種方法的模型各重復50 次,分別得到它們的平均預測誤差、單個數據點的更新和預測時間,如表1 所示.其中,邊界約束的局部高斯過程模型的模型寬度參數設置為20.

表1 3 種方法的性能對比Table 1 Performance comparison of three methods

從表中可以看出,全局高斯過程回歸方法的預測精度最佳,而本算法的預測誤差與它處于同一數量級,硬邊界局部模型方法的預測誤差約為全局高斯過程回歸方法的80 倍.考慮到全局高斯過程回歸使用全部數據構建單個模型,本算法的精度稍低于它是可以接受的.本算法通過與硬邊界局部模型方法接近的局部模型數量(各約100 個),實現了更好的預測精度.比較更新模型的時間時,全局高斯過程回歸方法的弊端顯現出來,其使用單個數據點更新模型的時間為132.753 ms,是局部模型方法的100 倍以上.當數據量繼續增大時,計算量將迅速增大直至無法被接受,而本算法和硬邊界局部模型方法的模型更新時間處于同一數量級.另外,3 種方法對單個數據點的預測時間為同一數量級,其中本文算法的預測時間最短,為1.342 ms.綜上可見,本文提出的算法在保證高精度預測的同時,能大量減少局部模型的更新時間.

3.2 多維工程數據實驗

在機器人控制和航空航天等工程領域,輸入數據通常是多維的,離線學習的數據集很難覆蓋所有情景.即時能夠使用大量數據集提前進行離線學習,計算量也會是巨大的.針對這些工程場景,本文提出的在線學習模型能夠使用實時數據更新模型,實現預測功能.

為了分析該模型處理真實多維工程數據時的性能,本文采用葡萄牙波爾圖大學公開的來自F16 型號飛機升降電梯的數據集進行實驗[26].F16 型號飛機升降電梯的數據集共包括9 517 個數據點.輸入數據有6 個屬性:爬升率、高度、轉輪角速度、轉輪角度、爬升加速度和爬升加速度的變化率,這6 個屬性的數值組成了6 維輸入數據.輸出是控制信號Se,它是一個處于-0.014 到0.013 間的值.為了分析在線學習過程中,邊界約束回歸模型的性能是如何隨著數據量的增長變化的,將9 000 個數據點作為訓練集,以數據流的形式不斷輸入模型進行更新.每增加1 000 個數據點,記錄當前數據量對應的更新時間、預測時間和誤差.這3 項是判斷模型是否滿足在線學習要求的關鍵指標,下面對這3 項指標進行分析.

首先,用一個新數據點更新模型所需的時間Tu是衡量在線學習模型性能的重要指標,只有模型能夠隨數據流快速更新,才能實現在線學習.圖5 展示了Tu隨數據量增長的變化趨勢.

圖5 更新時間隨數據量增長的變化趨勢Fig.5 The trend of update time as data increases

圖5 中N是用于訓練模型的數據點的數量,m是每個維度的局部模型數量,通過m和每個維度輸入數據的范圍可以確定模型寬度.可以看到,當m=3 時,Tu增長的趨勢過快,而m=5 時Tu由N=1 000 到N=9 000 只增加了0.75 倍,m=7 時Tu僅小幅度波動.因此m=5 或 7 對于此數據集是更為合適的參數.

圖6 展示了對一個未知數據點進行預測所需的時間Tp隨數據量增長的變化趨勢.同樣顯示了m=5 或 7 時預測時間增長平緩,分別只增加了0.75 倍和0.28 倍,相比m=3 是更合適的參數.

圖6 預測時間隨數據量增長的變化趨勢Fig.6 The trend of prediction time as the data increases

圖7 顯示了平均預測誤差隨數據量增長的變化趨勢.該圖的趨勢表明,隨著訓練集中數據點的增加,模型的預測誤差減少.也就是說隨著數據流的輸入,在線學習模型能夠越來越精確.

圖7 預測誤差隨數據量增長的變化趨勢Fig.7 The trend of prediction error as data increases

綜合圖5～7,對于這一特定數據集m=5 是最合適的參數.在m=5 時,數據點數量N=9 000雖然是N=1 000 時的9 倍,但它的更新時間0.541 ms只是N=1 000 時的0.309 ms 的1.75 倍,能夠實現快速更新模型,滿足在線學習的要求.預測時間由N=1 000 時的2 ms 增加到N=9 000 時的3.5 ms,只增加了0.75 倍,同樣滿足在線學習的實時性要求.本模型的預測誤差由N=1 000 時的0.0000026989降低到N=9 000 時的0.0000020318,隨著數據量的增加實現了更高的精度.

為了進一步驗證本模型在機器人控制場景下的可用性,本文將模型應用到具有兩個自由度的機器人的機械手臂上進行實驗,測試其在機器人控制方面應用的可行性.該機器人手臂由文獻[27]中提出的基于PD 反饋的計算扭矩控制方法進行控制.我們為機器人手臂設置期望的二維平面活動軌跡,通過模型的在線學習,預測控制機械臂所需的關節扭矩,并以此控制機器人手臂使得它能按照期望的活動軌跡運動.本文通過下式為機器人手臂設定期望的活動軌跡

其中,θ1和θ2分別是兩段機器人手臂與水平線的夾角,t是從0 開始計時的時間.將θ1和θ2對t求導,得到兩段機器人手臂的角速度ω1和ω2.由θ1,θ2,ω1和ω2組成四維輸入矢量x=[θ1θ2ω1ω2].我們的任務是使用這些輸入數據通過本文模型預測機器人手臂的關節扭矩,并使用該關節扭矩控制機械臂,使得機械臂能夠跟隨預設軌跡,同時將四維輸入數據和對應的真實關節扭矩輸入模型進行在線學習.為了更清晰地展現控制效果,本文截取了機械臂10 s內的運動軌跡參數,共計2 000 個預測點,即每0.005 s 預測一次,對模型有實時性要求.

圖8 顯示了預設和實際運動軌跡分別對應的參數隨時間的變化情況,其中,4 條由星號繪制的粗曲線分別對應的是實際運動軌跡的角度和角速度,4條虛線分別對應的是預設運動軌跡的角度和角速度.從圖中可以看出,在開始階段,實際運動軌跡和預設運動軌跡的曲線相差較大,但隨著模型的運行,兩者的曲線趨于重合.

圖8 預設和實際運動軌跡的參數Fig.8 Preset and actual motion trajectory parameters

為了更直觀地展示圖8 中預設和實際運動軌跡間的誤差變化情況,使用軌跡誤差矢量的范數||xr－x||表示實際運動軌跡和預設運動軌跡之間的誤差

其中,θr,1,θr,2,ωr,1,ωr,2是實際運動軌跡的參數,θ1,θ2,ω1,ω2是預設運動軌跡的參數.由此得到圖9,它顯示了運動軌跡誤差隨著時間振蕩下降的趨勢.圖8 和圖9 驗證了本文模型在高實時性要求的場景下的可用性.

圖9 運動軌跡誤差Fig.9 Motion track error

以上各實驗從不同角度驗證了本文模型能夠適用于工程中的大數據量場景,并且滿足在線學習的實時性需求.

4 結束語

為了在保持全局高斯過程回歸預測精度高等優點的情況下,減小其立方復雜度帶來的過大的計算量,本文提出了一種基于邊界約束的概率相關局部高斯過程回歸模型.本文在該模型中使用一種基于貝葉斯原理的數據關聯方法,定義了用于數據點的分配和預測的激活函數.同時,本文將局部模型放置在空間網格中,實現空間重疊的局部模型設置,提高了模型的精度.在算法方面,利用對計算過程的優化和k-d 樹最近鄰搜索,以及更新模型預測矢量時采用已儲存數據求矩陣的逆,降低了模型更新階段的計算量.評估結果表明,與其他局部高斯過程回歸方法相比,本文方法可以顯著縮短計算時間,同時實現與全局高斯過程回歸方法相似的預測精度,能夠更好地處理大量且持續更新的數據.在真實工程場景下,能夠很好地完成在線學習功能.該模型具有能夠快速更新和預測的特點,滿足在線學習的需求,對于機器人控制、航空航天等實時產生大量多維數據的領域,具有一定的工程意義.