高騰飛 劉勇琰 湯云波 張 壘 陳 丹

(武漢大學(xué)計(jì)算機(jī)學(xué)院 武漢 430072)

從地球、經(jīng)濟(jì)到人腦，在瞬態(tài)頻率和復(fù)雜諧波的演化過程中形成了復(fù)雜系統(tǒng).對(duì)于復(fù)雜系統(tǒng)的工作機(jī)理的探索一直是多學(xué)科研究的核心[1].隨著數(shù)據(jù)采集技術(shù)和采集設(shè)備的飛速發(fā)展，復(fù)雜系統(tǒng)研究中對(duì)于數(shù)據(jù)處理的要求越來越高，所需處理數(shù)據(jù)的規(guī)模、維度及復(fù)雜性也在迅速增加.海量規(guī)模的時(shí)間序列大數(shù)據(jù)記錄了復(fù)雜系統(tǒng)在時(shí)間和空間尺度上的演化.挖掘時(shí)間序列數(shù)據(jù)的潛在深度信息對(duì)于復(fù)雜系統(tǒng)過往數(shù)據(jù)分析和未來行為預(yù)測(cè)具有非常重要的意義[2].時(shí)間序列大數(shù)據(jù)包含了復(fù)雜系統(tǒng)中各個(gè)“個(gè)體單位”之間以及多個(gè)復(fù)雜系統(tǒng)間的相互聯(lián)系.其外部表征具有數(shù)據(jù)相關(guān)性和多重?cái)?shù)據(jù)屬性，這一特性也決定了時(shí)間序列大數(shù)據(jù)具有較高的維度.

提取時(shí)間序列大數(shù)據(jù)潛在因子對(duì)研究復(fù)雜系統(tǒng)的整體機(jī)制起著至關(guān)重要的作用[3].高維科學(xué)數(shù)據(jù)中的一般潛在低維變量揭示了被觀測(cè)系統(tǒng)動(dòng)力學(xué)不斷變化的整體機(jī)制.然而數(shù)據(jù)元素之間復(fù)雜的相互依賴關(guān)系、數(shù)據(jù)的高維度和海量的數(shù)據(jù)規(guī)模阻礙了其對(duì)整體機(jī)制的揭示.在小數(shù)據(jù)時(shí)代，時(shí)間序列大數(shù)據(jù)中低維特征的提取被認(rèn)為是張量分解的數(shù)學(xué)問題.然而，直接采用傳統(tǒng)的分析方法——獨(dú)立成分分析[4](independent component analysis, ICA)及主成分分析[5-6](principal component analysis, PCA)——無法完整準(zhǔn)確地提取出高維張量數(shù)據(jù)的潛在信息.這些方法使得維度之間的信息丟失和相關(guān)性破壞[7]，所提取的線性特征信息物理解釋性較差.平行因子分析(parallel factor analysis， PARAFAC[8])和Tucker模型分解[9]適用于高維張量分解，其中因子化更新通常使用交替最小二乘法(alternating least squares， ALS)[10]或其他相關(guān)的替代算法[11-12]來解決.時(shí)間序列大數(shù)據(jù)具有動(dòng)態(tài)增長的數(shù)據(jù)特征.數(shù)據(jù)規(guī)模不斷增長，小數(shù)據(jù)時(shí)代的張量因子分解方法難以高效地應(yīng)用于大規(guī)模、高維度的時(shí)間序列大數(shù)據(jù).為保證原始關(guān)鍵信息提取的高效和完整，迫切需要開發(fā)一種面向大規(guī)模高維度的時(shí)間序列大數(shù)據(jù)分解分析方法，以實(shí)現(xiàn)快速有效的因子信息提取.

時(shí)間序列大數(shù)據(jù)的數(shù)據(jù)規(guī)模是巨大的，并且在空間中快速增長，無法保證對(duì)問題域有足夠穩(wěn)定的先驗(yàn)知識(shí)，這就導(dǎo)致無法確定張量分解的初始條件.良好的因子化分解初始條件和快速迭代的因子更新過程在大規(guī)模高維時(shí)間序列大數(shù)據(jù)分解過程中極其重要的.如何在缺失先驗(yàn)知識(shí)的前提下，適應(yīng)大規(guī)模和超高維的數(shù)據(jù)特征，快速獲取因子成為了時(shí)間序列大數(shù)據(jù)分析的重要研究挑戰(zhàn).為了解決上述挑戰(zhàn)：1)基于貝葉斯因子化分解方法，本文提出一種在沒有先驗(yàn)信息域和無需精確初始化過程即可快速獲取收斂因子的海量并行貝葉斯因子化分析方法(the massively parallel Bayesian factorization approach， G -BF)；2)G -BF通過NVIDIA通用并行計(jì)算架構(gòu)(compute unified device architecture， CUDA)[13]進(jìn)行核函數(shù)映射，并行加速因子化迭代更新過程；3)G -BF使用高維數(shù)據(jù)映射方法完成對(duì)任意維度的時(shí)間序列大數(shù)據(jù)的分解要求；4)此外，G -BF結(jié)合粗粒度與細(xì)粒度2種并行方式應(yīng)用于子因子融合框架(hierarchical-parallel factor analysis, H-PARAAC)[14]完成對(duì)超大規(guī)模張量的因子化分解.

為了評(píng)估G -BF的性能，本研究進(jìn)行了一系列高維張量數(shù)據(jù)分析實(shí)驗(yàn)并與經(jīng)過GPU加速的G -HALS算法對(duì)比探究G -BF的性能效率.實(shí)驗(yàn)結(jié)果表明：1)G -BF在Maxwell架構(gòu)的GPU運(yùn)行性能優(yōu)于G -HALS；2)G -BF在處理高維數(shù)據(jù)規(guī)模、秩和維度等方面具有出色的可擴(kuò)展性，可處理任意維度數(shù)據(jù)；3)此外，G -BF應(yīng)用于現(xiàn)有的H-PARAFAC架構(gòu)在處理超大規(guī)模高維張量時(shí)表現(xiàn)了優(yōu)異的處理速度和性能，打破了可分析數(shù)據(jù)的規(guī)模限制.

本研究方法主要有3個(gè)方面優(yōu)勢(shì)：

1) G -BF能夠在無需問題域的先驗(yàn)知識(shí)情況下快速準(zhǔn)確地分解大規(guī)模時(shí)間序列大數(shù)據(jù)，獲得低維特征因子信息；

2) G -BF針對(duì)高維時(shí)間序列大數(shù)據(jù)具有普適性，可分解數(shù)據(jù)維度可變，不限于3維或其他限定維度；

3) G -BF嵌入到H-PARAFAC架構(gòu)中作為子張量因子分解過程的解決方案，在分解大規(guī)模張量時(shí)具有優(yōu)越的性能，打破了處理數(shù)據(jù)規(guī)模大小的限制.

1 相關(guān)工作

針對(duì)大規(guī)模張量因子化分解問題，許多成功的方法和框架已經(jīng)得到發(fā)展.大型矩陣及張量分解的性能已經(jīng)成為研究的一大挑戰(zhàn).針對(duì)此問題的相關(guān)工作大致從優(yōu)化分解過程和提高分解性能2個(gè)方面入手.

Bae等人[15]提出了一種高效多維并行可視化處理大規(guī)?？茖W(xué)數(shù)據(jù)算法(SMACOF).該算法利用集群消息傳遞接口并行矩陣操作和矩陣劃分，實(shí)現(xiàn)了高維矩陣到低維歐氏空間的轉(zhuǎn)換.

在線性分析問題中，特征值問題、奇異值問題以及QR分解(QR decomposition)問題極其重要.Agullo等人[16]提出了一種將矩陣分解為正交矩陣和上三角矩陣乘積的高效分解方法.該方法利用加速設(shè)備節(jié)點(diǎn)，采用CPUGPU混合架構(gòu)來加速Q(mào)R分解過程.

Tan等人[17]提出了一個(gè)矩陣分解庫(cuMF)來優(yōu)化交替最小二乘法(ALS)方法以求解超大規(guī)模矩陣分解.cuMF在單個(gè)GPU節(jié)點(diǎn)上使用多種策略優(yōu)化ALS中內(nèi)存訪問，同時(shí)將數(shù)據(jù)與模型并行化用于最小化GPU通信開銷.

Zou等人[18]提出了一種并行張量分解算法(GPUTENSOR)加速大數(shù)據(jù)集的張量分解.該方法將大張量分成小塊，利用圖形處理單元(GPU)固有并行性及高存儲(chǔ)帶寬并行相關(guān)操作，同時(shí)使用零散方式優(yōu)化計(jì)算策略，以避免中間數(shù)據(jù)爆炸而犧牲性能.

最近，Shin等人[19]提出了用坐標(biāo)下降法更新參數(shù)的大規(guī)模高階張量分解算法CDTF.該方法具有良好的數(shù)據(jù)可擴(kuò)展性，運(yùn)行在40個(gè)節(jié)點(diǎn)(Xeon 2.4 GHz)的Hadoop集群上時(shí)，可分解有10億個(gè)可觀測(cè)條目的5階張量.

本研究的目標(biāo)是在無需先驗(yàn)信息情況下高效分解大規(guī)模高維張量形式的時(shí)間序列大數(shù)據(jù).所提出方法的研究包括：1)運(yùn)行效率；2)分解的維度尺寸；3)可處理大型張量數(shù)據(jù)規(guī)模限制.

2 海量并行貝葉斯因子化分析方法(G -BF)

本節(jié)的主要內(nèi)容包括4個(gè)部分：1)介紹貝葉斯分解算法(Bayesian factorization， BF)的理論知識(shí)；2)對(duì)本文提出的面向時(shí)間序列大數(shù)據(jù)海量并行貝葉斯因子化分析方法(G -BF)的算法設(shè)計(jì)和并行加速設(shè)計(jì)的描述；3)簡要介紹G -BF和H-PARAFA框架的融合設(shè)計(jì)；4)對(duì)G -BF為滿足處理任意維度分解需求采用的高維數(shù)據(jù)映射方法進(jìn)行闡述.

2.1 貝葉斯分解算法

2.1.1 張量分解概率模型及貝葉斯先驗(yàn)

定義張量分解的線性模型表達(dá)式為

(1)

其中,Y表示大小為I1×I2×…×IN的N維張量，可以由因子張量X和噪聲張量組成，°表示外積，R表示因子矩陣秩，U(n)表示大小為In×R的因子矩陣集合，表示第n個(gè)因子矩陣的第r行.

(2)

其中，大小為R×R的矩陣σ和標(biāo)量λ為超參數(shù).

2.1.2 基于變分貝葉斯推斷的貝葉斯分解

概率模型在變分貝葉斯框架下采用確定性近似推理[20]，將自由能(Free)定義為p(Y)似然估計(jì)的邊界下限:

(3)

(4)

(5)

(6)

符號(hào)·表示內(nèi)積.

2.2 海量并行貝葉斯因子化分解方法

2.2.1 G -BF并行加速設(shè)計(jì)策略

G -BF是基于BF算法理論可以在GPU節(jié)點(diǎn)上運(yùn)行的大規(guī)模并行方法.在無需先驗(yàn)信息及精確因子初始化的情況下，G -BF算法可準(zhǔn)確快速迭代地分解大規(guī)模張量獲取因子矩陣.本方法涉及的大型矩陣運(yùn)算操作主要包括Khatri-Rao乘積、Hadamard乘積、逐元除法和矩陣乘法等.利用CUDA進(jìn)行密集矩陣運(yùn)算并行加速，實(shí)現(xiàn)因子矩陣集的快速更新迭代求解.G -BF將運(yùn)算過程作為設(shè)備并行模塊映射到圖1中CUDA核函數(shù).G -BF分解大規(guī)模張量獲取因子矩陣過程根據(jù)圖1內(nèi)容包含8個(gè)過程：

Fig. 1 G -BF method parallel acceleration design strategy diagram圖1 G -BF方法并行加速設(shè)計(jì)策略圖

1) G -BF的第1步發(fā)生在主機(jī)端.首先，讀取待處理張量Y，獲取Y的維度大小N、因子化分析過程中表示因子集大小的參數(shù)R以及初始化超參數(shù)λ.定義G -BF最大迭代次數(shù)ItemMax.當(dāng)?shù)螖?shù)達(dá)到ItemMax時(shí)，無論是否達(dá)到收斂條件都將停止因子矩陣更新.此參數(shù)是限制整個(gè)因子化過程的總體耗時(shí)的重要參數(shù)，當(dāng)不滿足收斂條件時(shí)，ItemMax設(shè)定值越大，通過迭代分解獲取的因子矩陣越接近最優(yōu)解，但耗時(shí)也將更大.

2) 該過程初始化G -BF中各項(xiàng)參數(shù)矩陣及參數(shù)張量，此過程發(fā)生在主機(jī)端.首先將因子矩陣U(n)隨機(jī)初始化大小為In×R的矩陣；S與G是3維張量，大小均為N×R×R,初始化方法為：當(dāng)?shù)?維和第3維的索引相同時(shí)，將值設(shè)為0.01，其余元素都設(shè)置為0.

3) 該過程在設(shè)備端負(fù)責(zé)計(jì)算初始化3維張量，其大小為N×R×R.E_Bar(n)的更新公式為

E_Bar(n)=U(n)TU(n)+InG(n),

這里涉及的運(yùn)算操作包括矩陣乘法和矩陣加法.根據(jù)E_Bar(n)的更新公式可知，各個(gè)E_Bar(n)的更新過程相互獨(dú)立.上角標(biāo)n代表更新的不同維度.這里采用并行方式更新，調(diào)用cudaStreams執(zhí)行多個(gè)維度的E_Bar(n)的并行更新，每1個(gè)cudaStream執(zhí)行任務(wù)相同.第n個(gè)cudaStream的執(zhí)行內(nèi)容為：將對(duì)應(yīng)第n維的因子矩陣U(n)和參數(shù)矩陣G(n)從主機(jī)端復(fù)制到設(shè)備端，主機(jī)端調(diào)用核函數(shù)在設(shè)備端計(jì)算E_Bar(n)，獲得結(jié)果后拷貝回主機(jī)端.在完成上述準(zhǔn)備后，開始迭代進(jìn)行因子矩陣U(n)的維度更新，設(shè)置item=1，n=1，item表示G -BF實(shí)際迭代周期數(shù)，n表示G -BF當(dāng)前更新維度.

4) 高維張量Y的存儲(chǔ)映射.G -BF算法設(shè)計(jì)要求可分解任意維度張量，但NVIDIA的CUDA計(jì)算架構(gòu)的數(shù)據(jù)管理固定為1維數(shù)據(jù).在G -BF中將張量Y通過張量矩陣化展開方式映射為1維向量進(jìn)行存儲(chǔ).當(dāng)維度更新時(shí)，張量元素坐標(biāo)根據(jù)維度間關(guān)系設(shè)定映射規(guī)則轉(zhuǎn)換.該過程也可以在迭代更新之外完成，如此可減少維度轉(zhuǎn)換的時(shí)間消耗，但這也將增加存儲(chǔ)空間的負(fù)擔(dān).詳細(xì)的維度映射方法將在2.2.4節(jié)給出.

5) G -BF的核心內(nèi)容，其目的是更新每個(gè)迭代循環(huán)中各維度的因子矩陣U(n).大致分為3個(gè)部分：計(jì)算過程參數(shù)矩陣T、更新當(dāng)前維度下參數(shù)矩陣G(n)、更新當(dāng)前維度下因子矩陣U(n).

① 計(jì)算過程參數(shù)矩陣的主要內(nèi)容為在主機(jī)端初始化過程參數(shù)矩陣T，其大小為R×R，數(shù)據(jù)值全部設(shè)為1λ.調(diào)用核函數(shù)更新設(shè)備端過程參數(shù)矩陣T，更新公式為T=T⊕U⊕-n,主要涉及的矩陣計(jì)算為Hadmard乘積.為減少核函數(shù)調(diào)用和數(shù)據(jù)傳輸帶來的時(shí)間開銷，將除當(dāng)前維度之外的因子矩陣合并成1維向量與過程參數(shù)矩陣T一同作為Hadmard乘運(yùn)算核函數(shù)的輸入，輸出結(jié)果為T的更新結(jié)果.

② 更新當(dāng)前維度下參數(shù)矩陣G(n)更新規(guī)則為G(n)=(T+S(n)-1)-1，需要通過CUDA進(jìn)行矩陣求逆和矩陣加法的GPU加速.

③ 更新當(dāng)前維度下因子矩陣U(n)其更新公式為U(n)={(Y(n)λ)U⊙-n}G(n)，這部分發(fā)生在設(shè)備端，調(diào)用Khatri-Rao乘積核函數(shù)直接計(jì)算U⊙-n，將N-1個(gè)因子矩陣拷貝到1個(gè)一維向量作為核函數(shù)的輸入?yún)?shù)以便一步直接計(jì)算U⊙-n，這樣的方法減少了核函數(shù)的調(diào)用.(Y(n)λ)U⊙-n由于矩陣乘法特性以及輸入矩陣規(guī)模過大導(dǎo)致本過程是整個(gè)G -BF中最耗時(shí)的步驟.調(diào)用設(shè)備端矩陣乘法核函數(shù)計(jì)算得到中間矩陣M，再次調(diào)用矩陣乘法核函數(shù)計(jì)算M×G(n)獲得當(dāng)前維度更新后的因子矩陣U(n).

6) 對(duì)當(dāng)前更新維度的因子矩陣U(n)進(jìn)行非負(fù)性和正則化處理.通過遍歷將因子矩陣中的負(fù)數(shù)設(shè)為0或某一極小值以及將因子矩陣U(n)以列優(yōu)先方式進(jìn)行歸一化.此操作目的是為減少計(jì)算復(fù)雜性以及滿足H-PARAFAC框架.通過對(duì)因子矩陣附加約束，保證了結(jié)合H-PARAFAC框架后結(jié)果的唯一性.

7) 參數(shù)矩陣E_Bar(n)和參數(shù)張量S的更新.該更新計(jì)算過程在設(shè)備端運(yùn)行.首先更新當(dāng)前維度中的參數(shù)矩陣E_Bar(n)，更新規(guī)則為E_Bar(n)=U(n)TU(n)+G(n)，使用CUDA調(diào)用主機(jī)端的核函數(shù)，在設(shè)備端進(jìn)行加速矩陣乘法和矩陣加法運(yùn)算操作.接下來，更新參數(shù)張量S.張量S實(shí)際上是1組數(shù)量為N且大小為R×R的矩陣集，其中單個(gè)矩陣的更新規(guī)則為S(n)=E_Bar(n)In+G(n)，任何2個(gè)S(n)更新過程彼此獨(dú)立，調(diào)用N個(gè)cudaStreams完成到設(shè)備端的并行傳輸并更新矩陣S(n).完成當(dāng)前維度中的所有更新操作后，n增加1，更新下一個(gè)維度.

當(dāng)G -BF中收斂到全局最小值或迭代次數(shù)達(dá)到預(yù)設(shè)最大值時(shí)，輸出的因子矩陣集為最優(yōu)因子.

2.2.2 并行加速線程映射策略

G -BF模型中在GPU設(shè)備端進(jìn)行的密集矩陣計(jì)算主要包括矩陣的Hadmard乘積、Khatri-Rao乘積、矩陣加減、逐元除法和矩陣乘法.G -BF中的矩陣運(yùn)算根據(jù)線程管理和調(diào)用策略的不同分為：

1) 直接映射.輸出矩陣的每個(gè)元素直接從輸入矩陣中的單個(gè)元素確定.例如，在矩陣Hadmard乘積運(yùn)算中，通過直接將2個(gè)輸入矩陣的第i個(gè)元素相乘獲得結(jié)果矩陣C的第i個(gè)元素，因此設(shè)備端線程直接映射到輸入元素.如果輸出矩陣有M個(gè)元素，并且核函數(shù)總共調(diào)用T個(gè)線程，那么每個(gè)線程最多負(fù)責(zé)計(jì)算(M+T-1)T元素.在G -BF中，Khatri-Rao乘積矩陣、Hadmard乘積矩陣、矩陣加減矩陣和逐元除法都屬于直接映射.算法1以Hadmard乘積為例給出直接映射的代碼設(shè)計(jì).

算法1.并行加速——矩陣Hadmard乘積.

輸入：矩陣A大小(AX,AY)、矩陣B大小(AX,AY);

輸出：矩陣C大小(AX,AY).

①row=blockIdx.y×blockDim.y+threadIdx.y；

②col=blockIdx.x×blockDim.x+threadIdx.x；

③offsetX=gridDim.x×blockDim.x；

④offsetY=gridDim.y×blockDim.y；

⑤ fori=coltoAXstepoffsetX

⑥ forj=rowtoAYstepoffsetY

⑦C[i×AY+j]=A[i×AY+j]×

B[i×AY+j];

⑧ end for

⑨ end for

2) 共享內(nèi)存.輸出矩陣的每個(gè)元素由輸入矩陣的1行或1列獲得.矩陣的乘法就是采用這種策略.矩陣乘法的輸出矩陣C的索引(i,j)是輸入矩陣A的第i行和輸入矩陣B的第j列的對(duì)應(yīng)元素乘積的和.若通過直接映射進(jìn)行矩陣乘法，則每個(gè)線程需要將1列和1行元素從全局內(nèi)存復(fù)制到局部內(nèi)存，這會(huì)帶來較高的時(shí)間開銷，大大降低計(jì)算性能.

本文在G -BF中使用共享內(nèi)存策略加快矩陣乘法的執(zhí)行效率.設(shè)備端Block中共享存儲(chǔ)器是線程共享的存儲(chǔ)單元，與設(shè)備端的全局內(nèi)存和局部內(nèi)存相比，具有更高的帶寬和更低的延遲.使用共享內(nèi)存能夠避免線程將數(shù)據(jù)從全局內(nèi)存復(fù)制到本地內(nèi)存的高額時(shí)間開銷.在矩陣的乘法中，許多不同的線程使用相同的數(shù)據(jù)部分，而共享內(nèi)存可以滿足同一個(gè)Block中的線程可以共享數(shù)據(jù).因此在對(duì)較大矩陣進(jìn)行乘法運(yùn)算時(shí)可以使用共享內(nèi)存上的“分塊矩陣”，每個(gè)塊物理上對(duì)應(yīng)于設(shè)備端的Block.假設(shè)2個(gè)輸入矩陣A(m×l)和B(l×n)，將每個(gè)矩陣劃分成一系列子矩陣(例如Aik和Bkj)，確保對(duì)應(yīng)的子矩陣塊的列等于行.

(7)

2.2.3 G -BF和H-PARAFAC

H-PARFAC框架是以大規(guī)模并行方式分解巨大稠密張量的一種解決方案，它使用現(xiàn)代并行和集群計(jì)算技術(shù)，實(shí)現(xiàn)了數(shù)據(jù)的近實(shí)時(shí)處理，并且使得提取任意大小和尺寸張量的全部因子成為可能.在大張量分解過程中，H-PARAFAC的分解和融合策略是非常有優(yōu)勢(shì)的.

G -BF在提取大規(guī)模高維數(shù)據(jù)低維信息時(shí)，具有優(yōu)良的迭代效率和性能.該方法能有效分解一個(gè)龐大的稠密張量，并且不需要精確的初始化過程即可快速迭代出因子的全局最優(yōu)解.將G -BF嵌入到H-PARAFAC框架中，它可以完成原來精確初始化模塊中并行直接三線性分解(direct trilinear decomposition，DTLD[6])和子張量分解2個(gè)部分的工作，進(jìn)而加速完成超大規(guī)模高維數(shù)據(jù)的低維因子化分解.圖2展示了G -BF結(jié)合H-PARAFAC框架的完整設(shè)計(jì).利用分解策略將原始高維張量變換為多個(gè)子張量；然后使用G -BF方法，將子張量因子化分解得到子因子矩陣；在子張量經(jīng)G -BF分解后，利用H-PARAFAC方法將各維度的子因子矩陣進(jìn)行融合，得到完整的原始張量的因子矩陣；最終完成超大規(guī)模高維時(shí)間序列張量大數(shù)據(jù)的分解.

Fig. 2 The design of H-PARAFAC framework with G -BF圖2 G -BF應(yīng)用H-PARAFAC架構(gòu)結(jié)合設(shè)計(jì)圖

2.2.4 高維張量的數(shù)據(jù)存儲(chǔ)映射

由于處理大規(guī)模時(shí)間序列大數(shù)據(jù)高維度分解設(shè)計(jì)目標(biāo)的限制，張量Y的數(shù)據(jù)存儲(chǔ)設(shè)計(jì)尤為重要.在分解過程中，幾乎所有的運(yùn)算都是矩陣運(yùn)算和向量運(yùn)算.在G -BF算法中，將張量Y擴(kuò)展為矩陣Y(n)進(jìn)行矩陣運(yùn)算，其映射策略為：假設(shè)N維張量Y中任意元素的坐標(biāo)為[i1,i2,…,iN]，張量每個(gè)維度的大小分別對(duì)應(yīng)為I1,I2,…,IN，則映射到一維向量的元素位置為

(8)

由式(8)可知，張量中各元素的N維索引與一維矢量Y中的坐標(biāo)是一一對(duì)應(yīng)的，然后可以直接利用相應(yīng)的元素進(jìn)行計(jì)算.基于此映射策略來完成4)中提到的張量維度(Y(1)→Y(n))的轉(zhuǎn)換.在嵌入G -BF的H-PARAFAC框架中子因子融合過程亦采用相似的維度轉(zhuǎn)換映射規(guī)則.

3 實(shí)驗(yàn)結(jié)果與分析

本次實(shí)驗(yàn)對(duì)本文提出方法的性能進(jìn)行了評(píng)估，并與G -HALS方法進(jìn)行了對(duì)比.實(shí)驗(yàn)包括：1)使用CUDA8.0重新編譯HALS算法G -HALS；2)使用CUDA8.0編譯本文的G -BF算法.所有的實(shí)驗(yàn)都采用了相同的編譯工具g++和nvcc.實(shí)驗(yàn)在1臺(tái)裝有NVIDIA GeForce TITAN X顯卡(Maxwell archite-cture)的工作站上進(jìn)行.具體的實(shí)驗(yàn)環(huán)境配置如表1所示:

Table 1 Configuration of Experimental Environment表1 實(shí)驗(yàn)環(huán)境配置

3.1 運(yùn)行效率

為了保持性能評(píng)估的一般性，本文實(shí)驗(yàn)數(shù)據(jù)均為由隨機(jī)數(shù)組成的多通道數(shù)據(jù)集.在使用G -BF方法進(jìn)行分析時(shí)，采用不同大小的數(shù)據(jù)集來測(cè)試該方法的適用性.隨機(jī)生成不同大小的4階張量，張量各個(gè)維度大小相同，包括40(404=2560000),50,60,70,80,90,100.

在所有的實(shí)驗(yàn)中，G -BF模型在每個(gè)張量上執(zhí)行20次，每次實(shí)驗(yàn)的迭代次數(shù)設(shè)置為100次.G -BF模型單次迭代的執(zhí)行時(shí)間如圖3所示，G -HALS和G -BF的單次迭代執(zhí)行時(shí)間隨著數(shù)據(jù)量增大的變化趨勢(shì)大致相同.

Fig. 3 The single execution time of G -HALS and G -BF processing the same tensor (order-4) with different sizes圖3 G -HALS和G -BF處理大小相同張量單次執(zhí)行時(shí)間(4-階)

G -BF的單次迭代時(shí)間范圍為[43.3 ms,645.4 ms]，分別對(duì)應(yīng)于最小體積和最大體積的張量.G -HALS對(duì)應(yīng)的單次迭代時(shí)間值為[46.4 ms，752.5 ms].結(jié)果表明，G -BF和G -HALS單次迭代的性能都有相同的趨勢(shì).在處理相對(duì)小的張量時(shí)，由于GPU的基本通信時(shí)間和計(jì)算能力沒有得到充分利用，G -BF的執(zhí)行時(shí)間接近G -HALS.隨著數(shù)據(jù)規(guī)模的增大，G -BF展現(xiàn)出性能上的優(yōu)勢(shì)，在數(shù)據(jù)規(guī)模非常大的情況下，G -BF在一定程度上能夠減少執(zhí)行時(shí)間.

使用相同的數(shù)據(jù)大小和截止條件來測(cè)試G -HALS和G -BF的迭代終止執(zhí)行時(shí)間.圖4顯示了用G -HALS和G -BF過程分解張量的執(zhí)行時(shí)間隨數(shù)據(jù)量增長的變化趨勢(shì).在大多數(shù)情況下，G -BF的性能優(yōu)于G -HALS，與上述單次迭代時(shí)的趨勢(shì)一致.同時(shí)與單次迭代時(shí)間相比，性能改善更為明顯，主要原因是G -BF在相同的條件下迭代次數(shù)較少.圖5顯示了張量分解中迭代次數(shù)的變化.相同條件下，在保證分解精度的同時(shí)G -BF算法的迭代次數(shù)更少，具有減少迭代次數(shù)的優(yōu)點(diǎn).

order-4:100×100×100×100 Fig. 4 The breakdown execution time of G -HALS andG -BF factoring analysis of tensor圖4 G -HALS和G -BF處理相同張量迭代終止時(shí)間

order-4:100×100×100×100 Fig. 5 The breakdown iterations of G -HALS andG -BF factoring analysis of tensor圖5 G -HALS和G -BF分解張量的迭代次數(shù)

3.2 可擴(kuò)展性

圖6顯示了相同數(shù)據(jù)規(guī)模在不同因子數(shù)量的分解下單次迭代的執(zhí)行時(shí)間.G -HALS算法的執(zhí)行時(shí)間范圍約為[725.3 ms，771.8 ms]，并且處理時(shí)間隨著因子數(shù)量的增大而有增加的趨勢(shì).隨著因子數(shù)量的增加，G -BF的執(zhí)行時(shí)間[620.5 ms，673.3 ms]具有相同的變化趨勢(shì).以上實(shí)驗(yàn)結(jié)果表明：G -BF和G -HALS都具有良好的可擴(kuò)展性.實(shí)驗(yàn)的數(shù)據(jù)集大小是固定的(1004=100 000 000)，因子數(shù)目作為唯一變量設(shè)置為4，8，12，16，20，24，28，32.

order-4:100×100×100×100 Fig. 6 The single execution time of G -HALS and G -BF processing the same tensor with different ranks圖6 不同秩下G -HALS和G -BF分解張量的單次執(zhí)行時(shí)間

Fig. 7 The single execution time of G -BF processing the same volume (107) of different dimensionality圖7 G -BF處理不同維度張量(107)的單次執(zhí)行時(shí)間

G -BF模型在數(shù)據(jù)維度上也表現(xiàn)出了良好的可擴(kuò)展性.圖7顯示了G -BF模型在處理數(shù)據(jù)規(guī)模相同、維度不同時(shí)的單次迭代執(zhí)行時(shí)間.隨著張量維數(shù)的增加，單次迭代的執(zhí)行時(shí)間也呈上升趨勢(shì).由于算法是基于維度數(shù)的迭代更新來求解最優(yōu)因子，雖然數(shù)據(jù)量是固定的，但是維度數(shù)增長帶來了更多的時(shí)間復(fù)雜度.實(shí)驗(yàn)表明：該方法適應(yīng)于處理不同維度的數(shù)據(jù)，使得任意維張量的分解成為可能.

3.3 處理超大規(guī)模張量的性能探究

基于已有工作提出的H-PARAFAC框架[13]，探究G -BF可以處理超大規(guī)模的高維張量數(shù)據(jù)性能.G -BF可同時(shí)在多個(gè)節(jié)點(diǎn)上運(yùn)行，圖8示出了基于H-PARAFAC在2個(gè)計(jì)算節(jié)點(diǎn)使用G -BF進(jìn)行張量分解的執(zhí)行時(shí)間.

Fig. 8 Total execution time of H-PARAFAC with G -BF deriving sub-factors and full factors of simulated augmenting tensor (order-4) when data scale>108圖8 當(dāng)數(shù)據(jù)規(guī)模>108時(shí)分解大張量總執(zhí)行時(shí)間

圖8上每個(gè)值表示處理張量的總時(shí)間.其中，數(shù)據(jù)規(guī)模從100×100×100×10增加到100×100×100×100 000,增加了10 000倍;而執(zhí)行時(shí)間從10.3 s增加到16 346.2 s，時(shí)間開銷增加了1 587倍.

GigaTensor[21]和CDTF[19]是與這項(xiàng)研究密切相關(guān)的2個(gè)優(yōu)秀的張量分解框架.GigaTensor使用100個(gè)計(jì)算節(jié)點(diǎn)可以分解3階稀疏張量(體積=109)，其時(shí)間開銷比數(shù)據(jù)量的增加快得多.CDTF使用運(yùn)行Hadoop的40個(gè)計(jì)算節(jié)點(diǎn)(Xeon 2.4 GHz)，可以分解具有109個(gè)數(shù)據(jù)元素的5階張量，并且其開銷隨著數(shù)據(jù)量線性增加.作為對(duì)比，G -BF使用H-PARAFAC框架處理109個(gè)數(shù)據(jù)元素大小的4階張量具有更高的效率.當(dāng)時(shí)間大小相同時(shí)，結(jié)合H-PARAFAC的G -BF能夠處理具有1011個(gè)數(shù)據(jù)元素的4階張量.更重要的是，在處理海量數(shù)據(jù)時(shí)，G -BF打破了數(shù)據(jù)規(guī)模的限制.

4 總 結(jié)

先驗(yàn)知識(shí)的缺乏已成為分析高維大規(guī)模時(shí)間序列的主要挑戰(zhàn)，大多數(shù)傳統(tǒng)的因子分解方法不能適應(yīng)數(shù)據(jù)的超高維和超大尺度.本文提出了一種海量并行貝葉斯分解方法(G -BF)來分析高維張量形式的時(shí)間序列大數(shù)據(jù).這種方法具有3個(gè)優(yōu)點(diǎn)：1)在沒有先驗(yàn)信息的情況下獲得因子矩陣；2)支持可變維度的張量分解；3)在處理大規(guī)模數(shù)據(jù)時(shí)保持性能和可擴(kuò)展性的優(yōu)勢(shì).

本研究使用GPU來加速密集矩陣的并行運(yùn)算，高效的存儲(chǔ)器配置方案和高維數(shù)據(jù)映射策略使得程序能夠滿足任意維度數(shù)據(jù)的大規(guī)模并行需求.同時(shí)，更新過程中的中間數(shù)據(jù)盡可能在設(shè)備存儲(chǔ)器中創(chuàng)建和操作，顯著降低了并行編程模型中主機(jī)與設(shè)備之間的通信消耗.

本文通過實(shí)驗(yàn)探究了G -BF方法在處理不同規(guī)模多維張量上的性能，與同一實(shí)驗(yàn)環(huán)境下HALS的并行優(yōu)化算法(G -HALS)實(shí)驗(yàn)作為對(duì)比.由實(shí)驗(yàn)結(jié)果可得出3個(gè)結(jié)論：

1) 與G -HALS(常規(guī)因子分解算法)相比，G -BF具有更好的性能，并且隨著數(shù)據(jù)量的增加，其優(yōu)越性更加明顯；

2) G -BF在數(shù)據(jù)規(guī)模、維度和因子數(shù)量等方面具有良好的可擴(kuò)展性；

3) 將G -BF方法與已有的H-PARAFAC方法相結(jié)合，可以對(duì)超大張量(2個(gè)節(jié)點(diǎn)上的體積可達(dá)1011)進(jìn)行因子分解，其性能相比傳統(tǒng)方法顯著提高.

總體而言，G -BF在對(duì)大規(guī)模多維張量的因子分解方面明顯優(yōu)于同類算法，針對(duì)時(shí)間序列大數(shù)據(jù)的因子降維分析方面具有很大的潛力.