喬 杰,蔡瑞初,郝志峰
(1.廣東工業(yè)大學計算機學院,廣州 510006;2.佛山科學技術學院數學與大數據學院,廣東佛山 528000)
因果網絡是一種描述多維數據間因果關系的網絡[1-3],在地球系統(tǒng)學[4-5]、生物基因學[6]等領域均具有廣泛應用。隨機試驗是因果關系推斷的黃金標準[7],然而通常存在實現成本高、實際操作困難等問題。因此,如何從觀測數據中推斷出因果結構是一個重要且具有挑戰(zhàn)性的課題。目前,研究人員已提出一系列從數據中發(fā)現因果關系的方法[8],但僅能觀測到整個因果結構的部分變量,還有研究人員針對存在隱變量的因果結構學習提出FCI[1]、RFCI[9]、GFCI[10]等算法,然而這類算法具有大量不可識別的因果邊。對于某個不存在v-結構[6]的因果對,如果其結果變量無法觀測,只能觀測到結果變量的子代,那么就存在一條級聯式傳遞的因果關系,而這樣的因果關系通常會破壞非線性函數式因果模型的可識別性。CAI 等[11]在2019 年提出級聯加性噪聲模型(Cascade Additive Noise Model,CANM),該模型在數據服從非線性加性噪聲假設下級聯結構仍是可識別的,但只能識別兩個結點的因果對,無法學習因果結構。本文針對包含隱變量的級聯加性噪聲模型,結合基于約束的因果骨架學習算法以及級聯函數式因果模型,提出一種混合因果結構學習算法。
從觀測數據中學習因果網絡主要包括面向時序序列(如格蘭杰因果[12]、基于神經網絡的聯合識別方法[13]等)以及非時序序列兩類算法。本文主要針對非時序數據,從非時序觀測數據中學習因果網絡分為以下兩類方法:第一類是基于獨立性的方法,包括基于約束[1]、基于評分[14]等方法,此類方法存在無法有效區(qū)分馬爾科夫等價類的問題[15],導致無法推斷部分因果結構;第二類是基于函數的方法,此類方法通常假設原因變量與結果變量間存在某種特定的函數映射關系,使得結果變量是由相互獨立的原因變量與噪聲變量經過映射得到的。通過判斷該噪聲與原因變量的獨立性,在某些條件下能夠推斷出因果關系的方向,將這類可以推斷因果方向的模型稱為可識別模型,典型的函數式因果模型一般只考慮兩個結點的因果對關系,主要包括線性非高斯無環(huán)模型
(Linear Non-Gaussian Acyclic Model,LiNGAM)[16-17]、
非線性加性噪聲模型(Additive Noise Model,ANM)[18]、后非線性(Post-Nonlinear,PNL)因果模型[19]。
定義因果結構為圖G={X,E},其中,X={X1,X2,…,Xn} 表示n個因果變量結點,E={(i,j)|Xi→Xj}表示結點間的因果邊集合,數據集為,樣本量大小為m。假設Xi是Xj的父母,若滿足直接因果關系則Xi→Xj,假設Xi是Xj的祖先,若存在間接因果邊則Xi→Xk→Xj或Xi←Xk←Xj。定義Xi⊥G Xj|S為給定條件集S,Xi與Xj在圖G上條件獨立。同時,給出因果忠誠性假設[1]、因果充分性假設等常見的基本假設。因果忠誠性假設保證了真實數據分布與因果結構獨立的一致性,這一假設使得通過條件獨立性進行因果骨架學習成為可能。因果充分性假設保證了任意兩個可觀測的因果結點不存在一個共同的隱變量父親根結點,這一假設使得在隱變量下利用級聯因果模型進行方向推斷成為可能。
本文使用一種典型因果結構[20]來表達隱因果結構網絡。簡而言之,對于因果路徑,若Xi→Xk→Xj且Xk是隱變量,則在典型因果結構上認為Xi→Xj,對于因果路徑Xi←Xk→Xj且Xk是隱變量,同時Xk沒有任何可觀測的父親或祖先結點,則認為Xi?Xj,反復迭代直到結構不再變化。在因果充分性假設下,典型因果結構并不存在形如Xi?Xj的雙向邊。
本文的目標是學習典型因果結構網絡。圖1 給出了基于級聯非線性加性噪聲模型的因果結構學習算法SCANM 框架。該框架分為兩個階段:第一階段是將輸入的觀測數據通過因果骨架來學習因果變量間的骨架;第二階段利用級聯非線性加性噪聲模型,針對存在隱變量的數據進行因果方向推斷。
圖1 因果結構學習框架Fig.1 Causal structure learning framework
因果骨架學習的基本流程為:首先初始化一個完全圖作為因果結構;其次使用條件獨立性檢驗對獨立邊進行刪邊,直到無法再刪除邊為止。
算法1因果骨架學習
在算法1中有3層關鍵的循環(huán):第1層循環(huán)(第3行)的目的是從小到大遍歷不同長度的條件集以測試條件獨立性,條件集長度從0 開始,即初始條件集為空集。這么做是因為如果條件越多,那么條件獨立性的判斷會越不準確,所以傾向于從小到大遍歷。第2 層循環(huán)(第5 行)的目的是找到每條滿足給定條件集的邊,此外,對于某條邊(Xi,Xj)而言,其條件集應該只會出現在Xi的鄰居中,即屬于集合Y?adj(G,Xi){Xj},通過該方法可以加快遍歷條件集的速度。第3 層循環(huán)(第7 行)的目的是遍歷所有滿足集合長度等于length 的條件集,即|Y|=length,并使用遍歷的條件對結點Xi和Xj進行獨立性檢驗,如果發(fā)現存在一個條件使得它們獨立,則刪除這條邊。通過以上步驟,可以學習到典型因果骨架圖G。
通過因果骨架學習找到典型因果骨架后,需要對每條邊進行因果定向。由于隱變量的存在且在因果充分性假設下,因此會遇到X→Z→Y和X→Y←Z兩類存在隱變量的結構,其中Z均是不可觀測的。
3.2.1 結構級聯非線性加性噪聲模型
對于第1 類結構已在非線性級聯加性噪聲模型中得到解決。對于第2 類結構,CANM 是沒有考慮的,即X→Y←Z*,在此使用Z*加以區(qū)分。該結構滿足ANM 模型,即Y=f(X,Z*)+ε,且Z*是隱變量。因此,本文認為非線性級聯加性噪聲模型可以應用于該結構,通過對第1 類結構的邊緣分布進行變換使其得到等價于第2 類結構的形式:
其中:nz是在第1 類結構分布下X→Z的噪聲,在第2 個等式中使用一般的f表達了中間傳遞過程,即Y=fy(fz(X)+Nz)=f(X,Nz)。關鍵在于第3 個等式,如果令噪聲分布與Z*的分布相等,則可以得到等價于第2 類結構的邊緣分布。換而言之,可以使用CANM框架對第2 類結構進行建模。在下文中將給出CANM 的變分下界以及基于變分自編碼機的優(yōu)化方案。
3.2.2 變分下界
根據式(1)可進一步將單個隱變量推廣到多個,并且它們的噪聲用向量n表示,因此非線性級聯加性噪聲模型的邊緣似然度可表示如下:
基于式(2)可以給出在樣本點(x(i),y(j))下的邊緣似然度的變分下界:
其中:qφ(n|x(i),y(i))是變分分布的,其作用是近似后驗分布pθ(n|x(i),y(i)),當且僅當KL(qφ(n|x(i),y(i))||pθ(n|x(i),y(i)))=0 時,上述下界的等號成立,即只要能夠找到一個足夠近似的變分分布,就能近似最大化模型的邊緣似然度。
3.2.3 變分自編碼機
為更好地優(yōu)化邊緣似然度的同時近似后驗分布pθ(n|x(i),y(i)),使用基于變分自編碼機(Variational Auto-Encoder,VAE)[21]的求解方案。具體地,使用多層感知機(Multi-Layer Perceptron,MLP)來分別作為編碼機qφ(n|x(i),y(i))和解碼器pθ(n|x(i),y(i))的全局近似器[22]。同時,針對期望項中的梯度求解問題,使用一種重參數化技術,假設先驗分布p(n)~N(0,I)服從高斯分布,使得式(3)改寫如下:
在得到各種隱變量結構下仍然適用的因果方向推斷方法后,將這種結構級聯非線性加性噪聲模型應用到典型因果骨架上,得到一種基于級聯非線性加性噪聲模型的混合因果結構學習算法。
算法2混合因果結構學習
算法2 分為兩個階段:第一階段是使用算法1 學習出典型因果骨架(第1 行);第二階段是使用級聯非線性加性噪聲模型對每個因果邊進行方向推斷(第2~12 行)。最終算法的輸出是典型因果結構(第13 行),該結構揭示了變量間存在隱變量下的因果結構關系。
使用虛擬因果結構與真實結構數據集對模型進行評估,將本文SCANM 算法與以下主流因果結構學習算法進行對比:1)基于盲源分離方法在線性非高斯數據下進行因果結構學習的LiNGAM 算法[16];2)基于爬山法與貝葉斯評分進行因果結構搜索的HC 算法[23];3)將約束與貝葉斯評分兩者集成進行因果結構搜索的MMHC 算法[24]。為驗證隱變量建模的有效性,在因果結構上的每一個因果對之間均會引入額外的隱變量,該隱變量滿足級聯加性噪聲的形式。在實驗中,每組參數都至少運行80 次以上,并采用F1 值(F)作為評價指標,計算公式如下:
在虛擬因果結構數據集中設計4 個變量控制實驗,每個實驗都會隨機生成不同的虛擬因果結構,具體設置為:隱變量數量為0、1、2、3、4,結構維度為6、10、15、20,平均入度為1.0、1.5、2.0,樣本數量為250、500、1 000、2 000、3 000、4 000、5 000,其中加粗數據表示在各控制實驗中的默認設置。
圖2 給出了不同隱變量數量下的實驗結果。由圖2 可以看出:一方面,SCANM 算法在0 個和1 個隱變量下F1 值都在0.87 附近,驗證了隱變量實驗的有效性;另一方面,隨著隱變量數量的增加,在隱變量數量為4 時,SCANM 算法的F1 值下降至0.73,這是因為隱變量的增加會導致信噪比的降低,從而影響模型學習效果。
圖2 不同隱變量數量下的實驗結果Fig.2 Experimental results under different number of hidden variables
圖3 給出了不同結構維度下的實驗結果,可以看出在不同的結構維度下SCANM 算法的F1 值在0.83 附近波動,這意味著級聯非線性加性噪聲模型在不同的結構維度下具有較強的魯棒性。
圖3 不同結構維度下的實驗結果Fig.3 Experimental results under different structural dimensions
圖4 給出了不同平均入度下的實驗結果。由圖4 可以看出,HC 與MMHC 算法對平均入度較為敏感,而且隨著平均入度的增加F1 值下降比較明顯,尤其在平均入度在1.0~1.5 時,F1 值下降了近0.1,其原因是當邊數較多時,若算法無法識別存在的隱變量,則出錯的邊數會增多,從而使得F1 值下降,而SCANM 算法由于對于隱變量的方向也可識別,因此對平均入度不敏感。
圖4 不同平均入度下的實驗結果Fig.4 Experimental results under different average in-degree
圖5 給出了不同樣本數量下的實驗結果。由圖5可以看出,所有算法對于樣本數量都較為敏感,樣本數量是一個較為重要的變量。盡管在樣本數量較少時,所有算法的F1 值均不太理想,但SCANM 算法的F1 值仍至少比其他算法高0.15,體現出其性能的優(yōu)越性。
圖5 不同樣本數量下的實驗結果Fig.5 Experimental results under different number of samples
選擇4 種真實結構進行測試,真實結構數據來自https://www.bnlearn.com/bnrepository/,具體信息如表1所示。
表1 真實因果結構數據集信息Table 1 Real causal structure dataset information
在真實因果結構數據集實驗中使用F1 值、準確率、召回率作為評價指標,計算公式如下:
在真實因果結構數據集中,使用與虛擬因果結構數據集實驗相同的設置,實驗結果如圖6~圖8 所示,總體而言,SCANM 算法在不同真實結構的平均F1 值為0.82,相比其他算法提升了51%。特別地,可以看出:在不同的真實結構中SCANM 算法均取得了最好的效果,尤其在sachs 結構中,其原因是sachs 結構的平均入度更大,級聯非線性加性噪聲模型在該情況下具有更好的效果;對比算法在不同的真實結構中召回率偏高而準確率偏低,尤其在child 結構中最為明顯,主要原因為沒有考慮隱變量的存在,容易學習出冗余邊。上述結果驗證了SCANM 算法的有效性。
圖6 真實因果結構數據集中的算法F1 值比較Fig.6 Comparison of F1 values for algorithms in real causal structure dataset
圖7 真實因果結構數據集中的算法準確率比較Fig.7 Comparison of precision for algorithms in real causal structure dataset
圖8 真實因果結構數據集中的算法召回率比較Fig.8 Comparison of recall for algorithms in real causal structure dataset
本文提出一種基于結構級聯非線性加性噪聲模型的因果結構學習算法,通過聯合傳統(tǒng)基于獨立性的因果骨架學習算法以及級聯加性噪聲模型,解決了存在隱變量的因果結構學習問題。實驗結果表明,在虛擬結構數據集和真實結構數據集下,該算法相比主流因果結構學習算法準確率更高、魯棒性更強。后續(xù)將研究不滿足因果充分性假設下的因果結構學習算法,進一步擴展級聯非線性加性噪聲模型的適用范圍。