基于因果自回歸流模型的因果結(jié)構(gòu)學(xué)習(xí)算法

盧小金，陳薇，郝志峰，2，蔡瑞初*

（1.廣東工業(yè)大學(xué)計算機學(xué)院，廣東 廣州 510006；2.汕頭大學(xué)理學(xué)院，廣東 汕頭 515063）

0 引言

許多科學(xué)研究的目的為揭示某些事物的因果關(guān)系，進(jìn)而找到支配它們的規(guī)律［1］。干預(yù)或隨機實驗是發(fā)現(xiàn)因果關(guān)系的傳統(tǒng)方法［2-3］，然而通常存在高成本、高耗時甚至無法實現(xiàn)等問題。因此，通過分析觀察數(shù)據(jù)來揭示事物因果信息，也被稱為因果發(fā)現(xiàn)的相關(guān)研究近年來引起了很多關(guān)注［4］。

在過去的幾十年里，因果發(fā)現(xiàn)研究已經(jīng)取得了一系列跨學(xué)科的進(jìn)展，被廣泛應(yīng)用于生態(tài)學(xué)、生物基因?qū)W、流行病學(xué)、神經(jīng)科學(xué)等領(lǐng)域［5］。同時，一系列識別因果效應(yīng)的算法被提出［6-8］，典型的方法有基于約束［9-10］、基于評分［11］和基于梯度優(yōu)化［12-13］的因果學(xué)習(xí)算法，以及一些混合因果結(jié)構(gòu)學(xué)習(xí)算法［14-15］。其中：基于約束或評分的算法存在無法識別等價類的問題；基于梯度優(yōu)化的算法可解決等價類識別的問題，但通常依賴于加性噪聲項和原因變量互相獨立的因果函數(shù)模型假設(shè)。

KHEMAKHEM等［16］于2021 年提出了因果自回歸流模型（CAREFL），該模型基于仿射自回歸流拓展了現(xiàn)有的非線性函數(shù)因果模型，在結(jié)果變量的噪聲項與原因變量不獨立時仍然是可識別的，但只能識別兩個變量的因果對，無法學(xué)習(xí)高維變量下的因果網(wǎng)絡(luò)結(jié)構(gòu)。針對加性噪聲項受原因變量影響的多維變量因果網(wǎng)絡(luò)學(xué)習(xí)問題，本文結(jié)合基于約束的因果骨架學(xué)習(xí)方法和因果自回歸流函數(shù)模型，提出一種混合因果結(jié)構(gòu)學(xué)習(xí)算法。該算法從完全無向圖出發(fā)，基于條件獨立性得到典型因果骨架，進(jìn)而基于CAREFL 計算備選方向的邊緣似然度進(jìn)行因果方向推斷。

1 相關(guān)工作

非時序觀察數(shù)據(jù)的因果發(fā)現(xiàn)方法主要包括基于約束的因果發(fā)現(xiàn)方法、基于評分的因果發(fā)現(xiàn)方法和基于梯度優(yōu)化的因果發(fā)現(xiàn)算法。基于約束的方法包括PC 算法［9］和FCI 算法［10］。PC 和FCI 算法基于條件獨立性檢驗的方法，雖然能夠在滿足算法前提的情況下近似準(zhǔn)確地學(xué)習(xí)到因果結(jié)構(gòu)，但因為馬爾可夫等價類的存在（不同的因果結(jié)構(gòu)可以滿足同樣的條件獨立性），兩者都無法確保輸出全部的因果信息。針對沒有混淆因子的場景，基于評分的方法通過優(yōu)化合理定義的評分函數(shù)來找到具有完整因果信息的結(jié)構(gòu)，經(jīng)典的方法有GES［11］算法。基于梯度優(yōu)化的因果發(fā)現(xiàn)方法，如NOTEARS［12］和GOLEM 算法［13］，由于對數(shù)據(jù)分布做了額外的假設(shè)而不僅是基于條件獨立關(guān)系，因此能夠區(qū)分同一等價類中的不同因果結(jié)構(gòu)。基于梯度優(yōu)化的方法一般需要依賴某類函數(shù)因果模型（也被稱為結(jié)構(gòu)方程模型）來保證因果可識別性（最優(yōu)結(jié)構(gòu)存在且唯一）。

函數(shù)因果模型主要分為線性模型和非線性模型。經(jīng)典的線性模型為線性非高斯加性模型（LiNGAM）［7，17］。LiNGAM 對數(shù)據(jù)的生成方式做了線性和非高斯獨立噪聲的假設(shè)，并利用獨立成分分析進(jìn)行求解。經(jīng)典的非線性函數(shù)因果模型有加性噪聲模型（ANM）［6］和后非線性模型（PNL）［8］。該類模型對數(shù)據(jù)生成方式做了非線性和獨立噪聲的假設(shè)。在LiNGAM、ANM、PNL 等模型中，均要求加性噪聲項和原因變量獨立以保證模型的因果可識別性。CAREFL 基于噪聲項和父母結(jié)點獨立的假設(shè)，將加性噪聲模型進(jìn)行拓展，使得模型在加性噪聲項與原因變量不獨立時仍然能夠識別變量之間的因果方向，但只討論了二維變量的因果方向推斷場景，無法解決多維的因果結(jié)構(gòu)學(xué)習(xí)問題。

2 問題定義

本節(jié)對所研究問題進(jìn)行符號化定義和說明。

定義數(shù)據(jù)集D中包含n個可觀測變量的集合V={V1,V2,…,Vn}，樣本量大小為m。令數(shù)據(jù)集D蘊含的因果結(jié)構(gòu)為有向無環(huán)圖G={V,E}，其中，E={(Vi,Vj)|Vi→Vj,Vj?Vi}表示G中結(jié)點間邊的集合。Vpa（j）表示Vj的所有 的原因結(jié)點，即若Vi?Vpa（j），則Vi→Vj。同時假設(shè)滿足因果忠誠性假設(shè)：數(shù)據(jù)集D的分布P忠誠于因果圖G。這表明，對于任意的Vi,Vj?V(i≠j)和集合s?V，給定，如果Vi和Vj在集合s中的變量的條件下獨立，則Vi和Vj被集合s中的變量d-分離［18］。由于因果忠誠性假設(shè)保證了分布中的條件獨立性和因果圖的一致性，因此能夠利用基于條件獨立性的方法得到因果圖的骨架。

基于上述符號說明和假設(shè)，本文關(guān)注的是典型因果網(wǎng)絡(luò)學(xué)習(xí)問題，其定義如下：給定具有n個變量V={V1,V2,…,Vn}的數(shù)據(jù)集D，基于因果忠誠性假設(shè)，如何學(xué)習(xí)變量集V的因果結(jié)構(gòu)。

3 基于CAREFL 的因果結(jié)構(gòu)學(xué)習(xí)算法

針對因果網(wǎng)絡(luò)學(xué)習(xí)的目標(biāo)，本文提出一種基于因果自回歸流模型的混合因果結(jié)構(gòu)學(xué)習(xí)算法，簡稱SCARF 算法。圖1 給出了SCARF 算法的整體框架。該算法框架分為兩個階段：第一個階段基于無向完全圖，通過條件獨立性刪除不存在因果關(guān)系的變量之間的邊，得到因果骨架圖；第二個階段基于因果自回歸流模型，分別對骨架圖中的每一條無向邊進(jìn)行因果方向推斷。

圖1 SCARF 算法框架Fig.1 Framework of SCARF algorithm

3.1 因果骨架學(xué)習(xí)

受基于約束的方法啟發(fā)，本文考慮使用條件獨立性檢驗來判斷變量之間是否存在因果邊。首先從無向完全圖出發(fā)，對每個結(jié)點及其鄰接結(jié)點進(jìn)行獨立性檢驗，刪除不存在因果關(guān)系的結(jié)點之間的邊。因果骨架學(xué)習(xí)的算法偽代碼如算法1 所示。

算法1基于條件獨立性學(xué)習(xí)因果骨架

算法1包含了3 層關(guān)鍵的循環(huán)：第1 層循環(huán)的作用是從條件集的空集開始遍歷，逐個增加條件集個數(shù)，從而找到使得變量獨立的最小條件集；第2 層循環(huán)遍歷所有滿足條件|adj(C,Vi){Vj}|≥Llength的無向邊(Vi,Vj)，尋找是否有使得Vi和Vj獨立的條件集，因為對于(Vi,Vj)而言，其條件集只會出現(xiàn)在Vi的鄰居中，因此此判斷條件可加快遍歷條件集的速度；第3 層循環(huán)的作用是對于Vi的鄰接結(jié)點遍歷長度為Llength的所有組合，進(jìn)而將每個組合作為條件集，檢驗Vi與Vj的獨立性，刪除檢驗結(jié)果為條件獨立的邊。

在本文算法中，需要對有限的樣本評估條件獨立性，做法是首先計算樣本偏相關(guān)系數(shù)，利用Fisher轉(zhuǎn)換得到偏相關(guān)系數(shù)的概率分布，進(jìn)而基于假設(shè)檢驗進(jìn)行獨立性測試。文獻(xiàn)［10，19］為該獨立性檢驗方法的可靠性提供了理論證明。

3.2 因果方向推斷

在得到由無向邊構(gòu)成的因果骨架后，需要對無向邊進(jìn)行因果定向。在基于結(jié)構(gòu)方程模型的因果發(fā)現(xiàn)算法中，由于需要具備因果關(guān)系可識別性和計算可實踐性，結(jié)構(gòu)方程模型通常具有特定的形式，如被廣泛使用的加性噪聲模型：

其中：fj是關(guān)于Vj的原因變量的非線性函數(shù)；nj是Vj對應(yīng)的加性噪聲項，各個變量對應(yīng)的噪聲獨立同分布。ANM 模型將關(guān)于原因變量的非線性函數(shù)與其噪聲項相加，簡潔高效，但在加性噪聲項與父母結(jié)點不獨立的場景下，無法保證因果可識別性。KHEMAKHEM等［16］提出的因果自回歸流模型則弱化了加性噪聲模型的假設(shè)，使其在因果方向判斷上也具備可識別性。

3.2.1 自回歸標(biāo)準(zhǔn)化流

在介紹因果自回歸流模型前，先介紹標(biāo)準(zhǔn)化流。標(biāo)準(zhǔn)化流基于隱變量z?Rn，通過一系列可微可逆的變換T表達(dá)觀測數(shù)據(jù)V?Rn的分布。z通常是簡單的基分布pz(z)，從而可以基于變量轉(zhuǎn)換獲得V的分布：

其中：T或者T-1通常使用神經(jīng)網(wǎng)絡(luò)實現(xiàn)。為了提高網(wǎng)絡(luò)擬合能力和保持變換的可逆可微性質(zhì)，標(biāo)準(zhǔn)化流模型通常將同族的一系列映射T1,T2,…,Tk鏈接起來，從而使得T=T1?T2?…?Tk。同時，T的雅克比矩陣行列式|detJT-1(V)|可以從子變換Tl的行列式計算中得到。因此，在設(shè)計標(biāo)準(zhǔn)化流模型時，研究者需要考慮的重要因素就是Tl雅克比矩陣行列式的計算復(fù)雜度。

自回歸流模型［20-22］固定變量的輸入順序，以及設(shè)計特定的模型結(jié)構(gòu)將雅克比矩陣限制為下三角矩陣以達(dá)到簡化行列式計算的目的。自回歸流模型結(jié)構(gòu)包括簡單的加法和仿射變換［23］，以及更復(fù)雜的神經(jīng)樣條變換［24］。在文獻(xiàn)［22］中，自回歸流模型的轉(zhuǎn)換函數(shù)T表示如下：

其中：π是關(guān)于V中變量自回歸結(jié)構(gòu)的排列；V＜π（j）表示在π中排在Vj之前的變量；函數(shù)τ（j也被稱為轉(zhuǎn)換器）相對于 第一個參數(shù)zj是可逆 的并且由V＜π（j）參數(shù)化。

3.2.2 因果自回歸流模型

結(jié)合上文可以看到式（1）和式（3）的高度相似性，這體現(xiàn)在兩個模型都顯式定義了基于變量的某種結(jié)構(gòu)或順序。鑒于這個相似性，基于自回歸流模型構(gòu)造的函數(shù)因果模型CAREFL 被提出：

其中：zj是噪聲變量，各結(jié)點的噪聲變量是獨立同分布 的；Vpa（j）是在圖G中的Vj的父母結(jié)點；tj和uj是 任意函數(shù)。在CAREFL 模型中，噪聲項不是簡單地與關(guān)于原因變量的非線性函數(shù)相加，而是可以受到原因變量的調(diào)控，這更符合現(xiàn)實場景。特別地，在各變量的t均為0 的情況下，式（2）是加性噪聲模型［見式（1）］的一部分。

假設(shè)因果順序為Vi→Vj，結(jié)合式（2）和式（4），可以得到因果自回歸流模型的邊緣對數(shù)似然度：

3.2.3 基于流模型的似然度比較

對因果骨架中的某條無向邊定向，實際上是在Vi→Vj和Vj→Vi兩個備選模型中選擇其一。本文利用標(biāo)準(zhǔn)化流分別擬合備選模型。對于某個備選方向，使用極大似然度的方法去優(yōu)化每一層流中神經(jīng)網(wǎng)絡(luò)的參數(shù)。為了避免過擬合，將樣本切分為測試集和訓(xùn)練集，用訓(xùn)練后的模型評估測試集的似然度。綜上所述，因果方向判斷的數(shù)學(xué)形式被定義為：

4 SCARF 算法描述

結(jié)合算法1 中基于約束的方法和標(biāo)準(zhǔn)化流技術(shù)，本文提出SCARF 算法。SCARF 算法偽代碼詳見算法2。SCARF 算法分為兩個階段：第一個階段是通過算法1 得到均為無向邊的因果骨架（第1 行）；第二個階段是遍歷因果骨架中的所有邊（第3～14 行），基于標(biāo)準(zhǔn)化流技術(shù)，通過極大似然法優(yōu)化模型中神經(jīng)網(wǎng)絡(luò)的參數(shù)，通過比較邊緣似然度的經(jīng)驗期望（第8～13 行）進(jìn)行因果方向推斷。最終，算法輸出的是表征變量間因果信息的有向無環(huán)圖。

算法2SCARF 算法

5 實驗與結(jié)果分析

為了驗證SCARF 算法在因果結(jié)構(gòu)學(xué)習(xí)上的有效性，本節(jié)對該算法在仿真因果結(jié)構(gòu)數(shù)據(jù)集和真實因果結(jié)構(gòu)數(shù)據(jù)集上的實驗效果進(jìn)行分析和評估。在評估實驗中，本文選取的對比方法有基于條件獨立性的約束算法PC、基于函數(shù)因果模型的ICALINGAM 和GOLEM。對于隨機生成和真實的因果結(jié)構(gòu)，實驗數(shù)據(jù)的生成機制服從如下非線性結(jié)構(gòu)方程：

其中：ui～N(0,1) 是隨機生成的噪聲變量；ai～U(-2,-0.5)∪(0.5,2)與bi～U(-2,-0.5)∪(0.5,2)是每個父母結(jié)點對Vi的隨機因果權(quán)重。每組實驗的運行次數(shù)在40 次以上。采用F1 值作為因果結(jié)構(gòu)學(xué)習(xí)的評價指標(biāo)，計算公式如下：

其中：TP表示因果結(jié)構(gòu)圖G中方向預(yù)測正確的邊數(shù)量；FP表示將反向預(yù)測成正向的數(shù)量；FN表示將正向預(yù)測成反向或無向的數(shù)量。

5.1 仿真因果結(jié)構(gòu)數(shù)據(jù)集實驗

通過4 個方面的控制變量實驗評估SCARF 算法的性能。每組實驗都根據(jù)實驗參數(shù)隨機生成因果結(jié)構(gòu)和數(shù)據(jù)，實驗參數(shù)設(shè)置如下：在數(shù)據(jù)生成機制實驗組中，控制變量是加性噪聲項是否與原因變量獨立（默認(rèn)加性噪聲項不與原因變量獨立），結(jié)點維度有10、16、20、30，因果圖的平均入度有0.5、1、1.5、2，樣本大小有500、1 000、2 000、3 000、4 000、5 000，其中，加粗字體為默認(rèn)實驗參數(shù)。

為了對比SCARF 在噪聲項和原因結(jié)點獨立時的結(jié)構(gòu)學(xué)習(xí)性能，本文引入了噪聲項不受父母結(jié)點調(diào)控的實驗效果，其數(shù)據(jù)生成機制服從以下結(jié)構(gòu)方程：

其中：wjl～U(-2,-0.5)∪(0.5,2),l=1,2,3，是每個父母結(jié)點對Vi的隨機因果權(quán)重，此時加性噪聲項ui不受父母調(diào)控。同時，對照組根據(jù)式（7）生成噪聲項受父母結(jié)點調(diào)控的數(shù)據(jù)。

圖2 給出了噪聲項是否受父母結(jié)點調(diào)控的仿真實驗結(jié)果。從圖2 的結(jié)果可以看到，無論加性噪聲項是否受到原因變量的影響，SCARF 算法的F1 值都顯著高于對比方法，保持在0.74 以上。同時，在加性噪聲項受到原因變量影響的場景中，基于條件獨立性的PC 算法指標(biāo)沒有受到顯著影響?；诤瘮?shù)因果模型的GOLEM 和ICALINGAM，F(xiàn)1 值較低并且指標(biāo)下降趨勢更明顯。這是因為GOLEM 基于ANM 模型，ICALINGAM 基于LiNGAM 模型，兩者都屬于加性噪聲模型，在噪聲項與原因變量不獨立時會破壞兩者的因果可識別性。實驗結(jié)果表明，基于因果自回歸流模型的因果推斷方法不僅適用于現(xiàn)有的加性噪聲模型，在噪聲項受父母結(jié)點影響時依然能夠識別因果方向，拓展了加性噪聲模型的應(yīng)用場景。

圖2 不同因果機制仿真數(shù)據(jù)的實驗結(jié)果Fig.2 Experimental results of different causal mechanisms simulation data

在圖3 中，實驗的控制變量是結(jié)點的維度?？梢钥吹?，隨著結(jié)點維度的增加，所有算法的指標(biāo)都有不同程度的下降，絕大部分因果算法對結(jié)構(gòu)結(jié)點都有一定敏感度。但本文算法在這10、15、20、30 這4 個結(jié)點維度下，結(jié)構(gòu)學(xué)習(xí)的F1 值始終維持在0.73 以上，顯著高于對比算法，這驗證了本文算法的魯棒性。

圖3 不同的結(jié)點維度下的實驗結(jié)果Fig.3 Experimental results at different node dimensions

圖4 給出了控制變量為平均入度的實驗結(jié)果。隨著結(jié)構(gòu)的平均入度增大，因果圖更加稠密，結(jié)構(gòu)學(xué)習(xí)難度會增大。實驗結(jié)果也表明隨著平均入度增大，所有算法的F1 值都有不同程度的下降。本文算法的F1 值維持在0.71 以上，高于所有對比方法的最高值。這說明本文算法對結(jié)構(gòu)平均入度有一定敏感度，但是依然能夠保持較強的結(jié)構(gòu)學(xué)習(xí)能力。

圖4 不同平均入度下的實驗結(jié)果Fig.4 Experimental results at different average in-degrees

圖5 給出了不同樣本數(shù)量下的實驗結(jié)果?？梢钥吹?，大部分算法對樣本數(shù)量較為敏感，但當(dāng)樣本規(guī)模提高到一定數(shù)量后，F(xiàn)1 值保持平穩(wěn)波動。SCARF 在樣本數(shù)量達(dá)到1 000 個后，F(xiàn)1 值在0.74 以上保持波動，高于其他對比方法的最高值，這是因為標(biāo)準(zhǔn)化流對分布具有強大的擬合能力，對樣本數(shù)量敏感度較低。

圖5 不同樣本數(shù)量下的實驗結(jié)果Fig.5 Experimental results at different sample sizes

5.2 真實因果結(jié)構(gòu)數(shù)據(jù)集實驗

從https：//www.bnlearn.com/bnrepository/中選取4 個真實因果結(jié)構(gòu)進(jìn)行實驗，驗證本文算法在真實因果結(jié)構(gòu)下的學(xué)習(xí)能力，這4 個真實因果結(jié)構(gòu)的信息如表1 所示。

表1 真實因果結(jié)構(gòu)信息Table 1 True causal structure information

在真實結(jié)構(gòu)的數(shù)據(jù)集中，根據(jù)式（7）生成樣本數(shù)量為3 000 的數(shù)據(jù)集。在圖6～圖8 中分別比較了各算法在4 個真實因果結(jié)構(gòu)上的F1 值、召回率和準(zhǔn)確率。圖6 表明，SCARF 在真實因果結(jié)構(gòu)上的F1 值顯著優(yōu)于對比方法，F(xiàn)1 平均值比ICALINGAM 高出48%。SCARF 對于不同的因果結(jié)構(gòu)的F1 值在0.7 以上，相對穩(wěn)定，這驗證了SCARF 在真實因果結(jié)構(gòu)中的魯棒性。在圖7 和圖8 中，對比方法中基于函數(shù)因果模型的ICALINGAM 和GOLEM 的準(zhǔn)確率很高但召回率卻很低，這是因為在噪聲項受原因變量影響的場景下，不符合這兩個算法的可識別性前提，導(dǎo)致結(jié)構(gòu)中許多父母結(jié)點的信息無法識別。

圖6 真實因果結(jié)構(gòu)數(shù)據(jù)集中不同算法的F1值Fig.6 F1 scores of different algorithms in true causal structure dataset

圖7 真實因果結(jié)構(gòu)數(shù)據(jù)集中不同算法的召回率Fig.7 Recall of different algorithms in true causal structure dataset

圖8 真實因果結(jié)構(gòu)數(shù)據(jù)集中不同算法的準(zhǔn)確率Fig.8 Accuracy of different algorithms in true causal structure dataset

6 結(jié)束語

本文提出一種基于結(jié)構(gòu)因果自回歸流模型的因果網(wǎng)絡(luò)學(xué)習(xí)算法SCARF，通過結(jié)合傳統(tǒng)的基于獨立性檢驗的因果骨架學(xué)習(xí)方法以及因果自回歸流模型，解決噪聲項受父母結(jié)點影響場景下的因果網(wǎng)絡(luò)學(xué)習(xí)問題。在虛擬因果結(jié)構(gòu)數(shù)據(jù)集和真實因果結(jié)構(gòu)數(shù)據(jù)集上的實驗結(jié)果表明，SCARF 算法較對比算法召回率更高、魯棒性更強。實驗中的真實因果結(jié)構(gòu)來自于生物、醫(yī)學(xué)等領(lǐng)域的經(jīng)典數(shù)據(jù)集，這表明SCARF 具備在實際應(yīng)用場景中進(jìn)行因果發(fā)現(xiàn)任務(wù)的能力，可應(yīng)用于神經(jīng)科學(xué)、生物信息學(xué)等領(lǐng)域。本文算法中對于網(wǎng)絡(luò)的全局搜索依賴于傳統(tǒng)的條件獨立性方法，隨著結(jié)點維度增大，算法的時間效率會受到約束，下一步將研究如何利用標(biāo)準(zhǔn)化流中神經(jīng)網(wǎng)絡(luò)可使用反向梯度傳播的屬性，基于梯度優(yōu)化學(xué)習(xí)全局網(wǎng)絡(luò)的算法，達(dá)到提升算法的使用范圍和降低時間復(fù)雜度的目的。