基于偏相關(guān)性測試的遞歸式因果推斷算法

陳銘杰，張浩，彭昱忠，謝峰，龐悅

（1.東莞理工學院 計算機科學與技術(shù)學院，廣東 東莞 523808；2.廣東石油化工學院 計算機學院，廣東 茂名 525099；3.復旦大學 計算機科學技術(shù)學院，上海 200433；4.南寧師范大學 計算機與信息工程學院，南寧 530001；5.北京大學 數(shù)學科學學院，北京 100871；6.中國銀聯(lián)博士后科研工作站，上海 201201）

0 概述

在因果關(guān)系推斷問題中，研究者通常使用統(tǒng)計獨立性測試［1］和條件獨立性測試［2］對數(shù)據(jù)集間蘊含的因果關(guān)系進行推斷，因此，對于數(shù)據(jù)中因果信息的獲取尤為重要［3］。以PC 算法［4］為例，作為一個經(jīng)典的基于約束的因果推斷算法，其算法原理是根據(jù)數(shù)據(jù)集間變量的聯(lián)合概率分布與條件概率分布特點，通過條件獨立性（Conditional Independence，CI）測試判斷變量之間是否滿足獨立性或條件獨立性，在進行一系列基于方向規(guī)則的判斷后構(gòu)造出一個大致的因果網(wǎng)絡圖。此圖通常是一個部分有向無環(huán)圖（Partial Directed Acyclic Graph，PDAG），其中包含一組馬爾可夫（Markov）等價類，而不是一個完全準確的因果網(wǎng)絡圖。

因為CI 測試的條件集會較大程度影響到測試結(jié)果的準確率，所以CI 測試比統(tǒng)計獨立性測試的難度要大得多［5］。對于CI 測試，一般有三種使用場景：一種是數(shù)據(jù)離散分布的情況，可以根據(jù)對應的條件概率P(X，Y|Z)=P(X|Z)P(Y|Z)將條件集合Z變?yōu)閱蝹€條件變量，然后使用常用的條件互信息或卡方檢驗技術(shù)判斷條件獨立性關(guān)系是否成立［5］；另外兩種分別針對數(shù)據(jù)連續(xù)分布和混合連續(xù)離散分布的情況，會采用離散化的方法對數(shù)據(jù)進行離散化處理，然后通過傳統(tǒng)的CI 測試方法即第一種場景中的條件互信息和卡方檢驗等方法進行檢驗，但是存在的問題是在離散化的過程中容易丟失一些數(shù)據(jù)信息，導致最終推斷結(jié)果不夠準確。雖然可以通過計算條件密度估計值PX|YZ=PX|Z間的距離來檢驗CI 測試是否準確［6］，但是當數(shù)據(jù)量較大時，離散的區(qū)間比較難設定，過大容易導致數(shù)據(jù)丟失，過小則需要很多的樣本量，在Z較大的情況下，離散化數(shù)據(jù)和估算條件密度都是需要解決的重要問題。

近年來，研究者提出一系列基于核函數(shù)的CI 測試方法。這些方法可以利用高階矩信息將觀測變量映射到再生核希爾伯特空間（RKHSs），此時的變量具有幫助推斷高維變量分布的特點，如獨立性和同質(zhì)性［7］。文獻［8］提出一種利用希爾伯特-施密特范數(shù)的條件交叉協(xié)方差算子，它是可以對應到RKHSs函數(shù)下x和y條件協(xié)方差的一種度量。該方法證實了如果與之對應的RKHSs 是特征核［9-11］，則算子范數(shù)為0且x⊥y|Z。相關(guān)文獻以及深入研究表明，基于核方法的CI 測試可以從給定的變量中獲得比基于離散數(shù)據(jù)集的CI 測試更多且更完整的信息，同時也表明基于核方法的CI 測試在因果推斷場景中能夠發(fā)現(xiàn)更準確的因果關(guān)系。

基于核函數(shù)的方法，通過將變量集的條件概率分布空間映射到基于平方可積函數(shù)空間的高維的線性空間，找出一種單向或雙向映射關(guān)系，如KCIT［12］。這種方法準確率較高，但運算復雜度也非常高，一般來說，當變量個數(shù)超過10 或樣本量大于2 000 時，就很難在一般的計算環(huán)境中運行。因此，該方法難以在高維數(shù)據(jù)場景中高效地解決因果推斷問題。由此可見，因果推斷的一個關(guān)鍵難題就是設計一個合適的CI 測試方法，使之能夠在高維數(shù)據(jù)場景下也能高效運行。

本文提出一種基于偏相關(guān)性測試的因果推斷算法CIPCT，從遞歸分治的角度應對數(shù)據(jù)集和條件集過大的情況。采用變量/特征分割方法將變量集遞歸式地分割成多個子數(shù)據(jù)集，使用因果推斷算法對每個子數(shù)據(jù)集進行因果推斷，得到對應的子結(jié)果。在此基礎上，通過合并策略對所有的子結(jié)果進行合并整理，得到與原始數(shù)據(jù)集相對應的因果網(wǎng)絡。此外，在“分-治-合”過程中采用較為高效的偏相關(guān)性測試，避免計算復雜度高的核密度估算，進一步提高因果推斷的效率。

1 相關(guān)研究

本節(jié)將介紹一些關(guān)于提升因果推斷速度的工作。在研究因果關(guān)系相關(guān)問題時，函數(shù)模型性能也是其中一個決定性因素。本文采用加性噪聲模型（Additive Noise Model，ANM）作為因果函數(shù)模型，以此作為函數(shù)模型的優(yōu)勢在于：一是有獨特的模型特性；二是數(shù)據(jù)都有一個特殊的加噪聲項，可以借此獲取額外的信息［13］。由于基于最小殘差平方的回歸方法估算可得到關(guān)于Z的非線性平方可積函數(shù)，因此文獻［14］對（xi，Z）和（xj，Z）分別做非線性回歸得到函數(shù)f和g，通過檢驗f和g是否存在將CI 測試簡化為兩到三個統(tǒng)計獨立性測試［15］。當加性噪聲模型的特性被加以利用，即可避免估算條件概率密度，減少條件集Z所引起的高維難題，從而獲取完整的因果信息，加速CI 測試并打破Markov 等價類限制。

在打破Markov等價類限制的前提下，文獻［16-17］提出一種利用殘差獨立性發(fā)現(xiàn)條件獨立性的方法，其對于Z維度具有更好的魯棒性。Markov 等價類包含三個結(jié)構(gòu)：x1→x2→x3，x1←x2→x3，x1←x2←x3，這三個結(jié)構(gòu)有著相同的條件獨立性和統(tǒng)計獨立性，在傳統(tǒng)的因果推斷算法下判斷異常艱難，但是這三個結(jié)構(gòu)所具有的獨特的外加噪聲項是完全不同的。文獻［18］通過最小二乘回歸計算出殘差，利用外加噪聲項驗證條件獨立性等價于殘差獨立性，使得CI 測試進一步簡化成一個獨立性測試，從而進一步提升速度。

2 預備知識

2.1 因果網(wǎng)絡

因果網(wǎng)絡是一種概率推理模型，由一組隨機變量X=｛x1，x2，…，xn｝以及對應的頂點集V=｛v1，v2，…，vn｝和頂點間邊的集合E=｛e1，e2，…，em｝組成，記作G（X，V，E）?？梢钥闯?，因果網(wǎng)絡是由一組隨機變量X與其相對應的一個有向無環(huán)圖G（V，E）構(gòu)成的，其要求X的聯(lián)合概率分布P（X）與G（V，E）是對應的，其中，E=｛e（xi，xj）|xi，xj∈X｝表示有向無環(huán)圖G中兩個變量xi和xj之間的邊集合，e（xi，xj）表示變量xi和xj之間的關(guān)系，可以是xi-xj、xi→xj或xi←xj。需要注意的是，在因果推斷領域中，為了方便或提高閱讀性，在不影響理解的情況下，通常用X同時表示變量與節(jié)點，不分開論述。

2.2 馬爾可夫等價

當有兩個有向無環(huán)圖是馬爾可夫等價類時，即兩個有向無環(huán)圖滿足包含同樣節(jié)點、邊（不考慮方向）和V結(jié)構(gòu)的條件時，它們具有完全一樣的獨立性與條件獨立性。以圖1、圖2 為例，這兩個網(wǎng)絡結(jié)構(gòu)是屬于馬爾可夫等價的，他們的獨立性和條件獨立性一致：1）獨立性：X不獨立Y，X不獨立Z和Y不獨立Z；2）條件獨立性：X不獨立Y|Z，X⊥Z|Y，Y不獨立Z|X。因此，對于馬爾可夫等價類，由于它們都具有相同的獨立性與條件獨立性，因此CI 測試無法正確判斷出真實結(jié)構(gòu)。

圖1 順連結(jié)構(gòu)Fig.1 Consequent structure

圖2 分連結(jié)構(gòu)Fig.2 Disjuncition structure

2.3 忠實性假設

已知一個因果網(wǎng)絡G，其中有三個變量（集）X、Y、Z，他們的聯(lián)合概率分布為P（XYZ），如果滿足從X⊥Y|Z可以推導出X、Y被Z集合D分離的條件，那么就稱這個聯(lián)合概率分布P（XYZ）滿足關(guān)于因果網(wǎng)絡G的因果忠實性假設。

3 偏相關(guān)性測試

本節(jié)將討論如何在因果推斷中應用偏相關(guān)性測試檢測條件獨立性是否成立，這里先給出DAUDIN的關(guān)于CI 的一個經(jīng)典理論［19-20］，其中描述了如何在平方可積函數(shù)空間中把相關(guān)性與獨立性相關(guān)聯(lián)，將以該理論支撐本文后續(xù)的理論分析。

3.1 條件獨立性特性

條件獨立性特性（Characterization of Conditional Independence，CCI）的相關(guān)定義如下：

令X、Y、Z為3 個隨機變量或隨機變量集，定義：

CCI 描述的是一種在平方可積函數(shù)空間中相關(guān)性與獨立性的等價關(guān)系。后續(xù)將給出一些相關(guān)的理論分析與結(jié)果，即分別在高斯和非高斯情況下偏相關(guān)性與條件（不）獨立性的關(guān)系。

3.2 偏相關(guān)性與條件獨立性的關(guān)系

定理1給定m+2個隨機變量：x、y和Z=｛z1，z2，…，zm｝，它們都是基于獨立的隨機變量si（i=1，2，…，l）的線性組合，如果所有的si都符合聯(lián)合高斯分布，那么x⊥y|Z成立，當且僅當σxy.Z=0。

證明考慮到給定變量集Z后x和y的偏相關(guān)性系數(shù)為ρxy.Z=，其中σ**.Z（給定Z的協(xié)方差）可以認為是投影到Z上的x與y殘差的協(xié)方差。因此，由獨立與相關(guān)關(guān)系可得出σxy.Z=cov(x-E(x|Z)，y-E(y|Z))=0。在聯(lián)合高斯分布下，σxy.Z=0 ?x⊥y|Z。

另一方面，如果x⊥y|Z，那么根據(jù)式（7）進行推斷，因為E(x-E(x|Z)|Z)=0和E(y-E(y|Z)|Z)=0，所以有x-E(x|Z)∈E1和y-E(y|Z)∈E2。此外，考慮到給定Z關(guān)于x和y的偏相關(guān)性，由于給定Z后的偏協(xié)方差σxy.Z等價于x和y在Z形成的線性空間上投影的協(xié)方差，因此有：

定理2偏相關(guān)性不成立與獨立成立在聯(lián)合高斯分布條件下是等價的，即獨立性測試可以由偏相關(guān)性測試替換。

推論1給定m+2 個隨機變量：x、y以及Z=｛z1，z2，…，zm｝，它們都是基于獨立的隨機變量si（i=1，2，…，l）的線性組合，如果所有的si符合聯(lián)合非高斯分布，那么σxy.Z≠0 ?x不獨立y|Z。

推論1 可以由定理1 直接推導出，表明在非高斯場景下，偏相關(guān)性不成立是條件不獨立的充分條件，即可以通過偏相關(guān)性測試判定獨立性不成立。所以，在非?？量痰臈l件——非高斯且偏相關(guān)性成立的情況下，無法判斷條件獨立的情況。因此，在絕大多數(shù)情況下，都可以用偏相關(guān)性對獨立性進行判斷，從而大幅縮短檢測獨立性是否成立的時間。

4 因果推斷算法

本節(jié)提出基于偏相關(guān)性測試的遞歸式因果推斷（CIPCT）算法，對算法基本框架、因果分割、因果方向?qū)W習與子圖合并4 個部分分別進行介紹。

4.1 算法基本框架

CIPCT 算法是一種遞推式算法，其采用傳統(tǒng)的“分-治-合”方法［21-22］進行因果學習，并在具體的算法步驟中融合偏相關(guān)性測試。CIPCT 算法流程如圖3所示，具體步驟如下：

圖3 CIPCT 算法流程Fig.3 Procedure of CIPCT algorithm

1）采用分治策略將龐大的原始數(shù)據(jù)集進行遞歸，分解成更小的子數(shù)據(jù)集，使用高效的偏相關(guān)性測試提高效率。

2）在分出若干子數(shù)據(jù)集的基礎上，對每一個子數(shù)據(jù)集進行因果方向?qū)W習，形成子結(jié)果。

3）以比較顯著性值的方式作為合并策略，對所有的子結(jié)果進行合并整理后，得到原始數(shù)據(jù)集對應的完整因果網(wǎng)絡圖。

4.2 因果分割

由于傳統(tǒng)的因果推斷算法在面對高維數(shù)據(jù)時運行時間較長，因此本文在進行因果分割時采用分治策略，并融合偏相關(guān)性測試提升效率。因果分割階段的目的就是盡可能地將數(shù)據(jù)集分割成子數(shù)據(jù)集，子數(shù)據(jù)集會通過與父數(shù)據(jù)集相同或更高階的偏相關(guān)性測試繼續(xù)進行因果分割，從而減少偏相關(guān)性測試的次數(shù)，縮短算法在高維環(huán)境下的耗時。在基于CI 測試的因果推斷中，一般將CI 測試的執(zhí)行次數(shù)作為算法時間復雜度，因為與CI測試相比，其他操作在執(zhí)行時間上可以忽略不計。假設原始變量集V=｛v1，v2，…，vn｝被遞歸分割成m個子集｛V1，V2，…，Vm｝，其中，對于所有m，有|Vm|≤n。當假設CIPCT 使用PC 算法進行因果推斷時，則求解子問題的算法時間復雜度為，其中kmax=max（|V1|，|V2|，…，|Vm|）。另一方面，需要在每次迭代中都構(gòu)造一個CI 測試表，在最壞的情況下，需要針對原始變量集V計算一個σ階CI 測試表（σ是預先設定的參數(shù)）。因此，遞歸式推斷算法CIPCT 的時間復雜度為如果只使用PC 算法解決因果推斷問題，則算法復雜度為O（n22n-2）。由于kmax通常要比n小得多，因此CIPCT 算法的算法時間復雜度較傳統(tǒng)基于CI 測試的算法大幅降低。

因果分割的具體實現(xiàn)如算法1 所示，主要步驟如下：

1）構(gòu)造一個0 階偏相關(guān)性測試表M，其通過對數(shù)據(jù)集V中v1～vn進行偏相關(guān)性測試后得到。Mij=1 表示在0 階條件下vi和vj相關(guān)系數(shù)為0，即corr（vi，vj）=0；Mij=0 表示在0 階條件下vi和vj不相關(guān)。

2）從Mij=1 中取vi∈A，vj∈B，把數(shù)據(jù)集V分為子數(shù)據(jù)集ABC，其中C中斷了AB之間所有聯(lián)系，要求C必須是滿足條件的最小集，但不一定為是A或者B的D 分離集。令D=V，從D中去除掉所有與C不相關(guān)的變量，令V1=A∪D，V2=B∪D。

3）若在上一步中無法對V進行因果分割，則構(gòu)造更高一階的偏相關(guān)性測試表，再執(zhí)行步驟2。通過以上操作，最后可得一個因果分割V=｛V1，V2｝。

算法1因果分割

4.3 因果方向?qū)W習

經(jīng)過因果分割后將會得到若干個子數(shù)據(jù)集，此處將使用V結(jié)構(gòu)性質(zhì)進行因果方向?qū)W習。首先從局部結(jié)構(gòu)x-z-y入手，如果在給定某個條件集Z之后，能令x、y滿足偏相關(guān)性，即pcorr（x，y|Z）=0，根據(jù)D 分離準則，此時若z不屬于Z，就可以判斷出x-z-y結(jié)構(gòu)為x→z←y；然后進行一致性傳播［3］，推導出剩下邊的方向。在因果分割過程中會將這幾個相似的局部結(jié)構(gòu)保存，以此來完成下一步子圖合并。

4.4 子圖合并

本文在文獻［21］研究的基礎上，將子圖合并中的CI 測試替換為偏相關(guān)測試，根據(jù)邊與邊之間的顯著性值sv 的大小決定連接方向。顯著性值sv 表示了方向被接受的概率，給出相關(guān)的理論分析如下：

定理3對于任意一個無向圖或局部無向圖，假設給定一個局部結(jié)構(gòu)v1-v2-v3，在對該結(jié)構(gòu)進行基于V結(jié)構(gòu)的方向判斷時，接受兩個方向v1→v2和v3→v2的概率相同，均為svvstructure：=P(v1⊥v3|S)×P(v1⊥v3|(S∪v2))，其中S為滿足CI 的條件集。

證明根據(jù)V結(jié)構(gòu)定義，v1與v3在給定某一個條件集S滿足CI，而條件集增加中間節(jié)點集后不滿足CI 這兩個特征即可直接推斷得到。首先假設給定一個局部V結(jié)構(gòu)v1-v2-v3，當v1→v2和v3→v2時，它們的顯著性值分別定義為：

其中：S是對應的D 分離集。然后對比顯著性值，取顯著性值更大的方向進行方向連接。最后通過合并所有的子因果網(wǎng)絡，得到關(guān)于原始數(shù)據(jù)集的完整因果網(wǎng)絡結(jié)構(gòu)圖。在子圖合并過程中，顯著性值的計算是在偏相關(guān)測試的基礎上進行的，通過檢測上一步因果分割中的偏相關(guān)測試得到顯著性值再保存結(jié)果，因此，不會增加額外的計算時間。子圖合并的具體實現(xiàn)如算法2 所示。

算法2子圖合并

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

本文所有的實驗均在MATLAB R2016a 中進行，運行環(huán)境是一臺CPU 為i5-7200U，主頻為2.50 GHz、內(nèi)存12 GB 的64 位操作系統(tǒng)的筆記本，將本文提出的CIPCT 算法和目前具有代表性的遞歸式因果推斷算法CAPA［21］進行對比。

實驗所用數(shù)據(jù)為10 組不同的真實因果網(wǎng)絡結(jié)構(gòu)［22-23］，涵蓋了多個領域，表1 中列出了這些數(shù)據(jù)的網(wǎng)絡結(jié)構(gòu)信息，其中，節(jié)點數(shù)代表維度。

表1 10 個因果網(wǎng)絡結(jié)構(gòu)的統(tǒng)計特性Table 1 Statistical characteristics of ten causal network structures

實驗分為兩種方式進行：首先在不同數(shù)據(jù)場景下固定相同樣本量，將樣本量控制在2|V|=2n，這是一種經(jīng)典的評估樣本量、節(jié)點個數(shù)與準確率關(guān)系的方法［14］，并同時對比兩種推斷算法的耗時；然后在相同數(shù)據(jù)環(huán)境下，分別比較不同樣本量下的實驗結(jié)果，評估樣本量對準確率、耗時的影響。本文實驗使用召回率R、精確率P和F1值F1對算法的綜合性能進行評估。其中：召回率反映算法發(fā)現(xiàn)因果關(guān)系的查全率；精確率反映算法發(fā)現(xiàn)因果關(guān)系的查準率；F1 值反映算法的綜合的準確率。這三個指標的計算公式如下：

5.1 不同數(shù)據(jù)相同樣本量的實驗

在不同數(shù)據(jù)場景下固定相同樣本量進行對比實驗，選取全部10 個因果網(wǎng)絡結(jié)構(gòu)并將樣本量控制在2|V|=2n，將提出的CIPCT 算法與CAPA 算法［21］進行對比。對比兩種推斷算法的準確率，并列出兩種算法的運行時間進行分析比較，實驗結(jié)果如表2 所示，其中加粗數(shù)據(jù)表示最優(yōu)值。

表2 CIPIT 和CAPA 算法在不同數(shù)據(jù)集中的性能對比Table 2 Performance comparison of CIPIT algorithm and CAPA algorithm on different datasets

分析表2 結(jié)果可以發(fā)現(xiàn)，CIPCT 算法的耗時相對較少，但是精確率相對略低，這是因為偏相關(guān)性測試對于高維數(shù)據(jù)時存在錯誤添加節(jié)點的現(xiàn)象。但是由于其采用了計算偏相關(guān)性系數(shù)的方式，因此在維度增加時可以篩選出更合適的節(jié)點，從而加快構(gòu)造網(wǎng)絡結(jié)構(gòu)的速度，進而加快對高維數(shù)據(jù)的處理速度。

5.2 相同數(shù)據(jù)不同樣本量的實驗

設計在相同數(shù)據(jù)場景下不同樣本量的對比實驗，挑選Barley（兩者準確率相近）和Andes（CAPA 優(yōu)于CIPIT）分別在不同樣本量（n、2n、3n、5n、7n）的情況下進行對比實驗，實驗結(jié)果如表3 所示，其中加粗數(shù)據(jù)表示最優(yōu)值。

表3 CIPIT 和CAPA 算法在不同樣本量下的性能對比Table 3 Performance comparison of CIPIT algorithm and CAPA algorithm under different sample sizes

從表3 可以看出，兩種算法在Barley 數(shù)據(jù)集上準確率幾乎持平，而在Andes 數(shù)據(jù)集上CAPA 優(yōu)于CIPIT，因為偏相關(guān)性測試是一種近似的CI 測試，在特定的網(wǎng)絡結(jié)構(gòu)下，如Andes，這個網(wǎng)絡結(jié)構(gòu)要求嚴格的CI，則偏相關(guān)性測試效果會一定程度下降。此外，可以看到當維度為n時，在Andes 結(jié)構(gòu)中CI 測試的速度較快，原因是CAPA 采用的是啟發(fā)式因果分割，在低樣本量情況下，加上CI 測試的耗時也可能要比CIPIT 快。

通過對比表3 中的F1 值和運行時間可以看出，雖然本文算法在Andes 網(wǎng)絡結(jié)構(gòu)中的F1值略低于CAPA 算法，但對比兩種算法的處理時間可知，CIPIT可以在高維網(wǎng)絡結(jié)構(gòu)中進行更快速的有效構(gòu)建，并且在維度越高的數(shù)據(jù)環(huán)境下這種時間差距越明顯。由于目前其他的因果推斷算法都難以在高維數(shù)據(jù)環(huán)境下進行快速有效的處理，使用偏相關(guān)性測試后耗時減少明顯，說明本文算法在運算時間上具有顯著優(yōu)勢。

6 結(jié)束語

本文提出的CIPCT 算法是一種基于偏相關(guān)性測試的快速因果推斷算法，其在傳統(tǒng)的“分-治-合”框架中融入偏相關(guān)性測試，從而保證準確率同時提高運算效率。在對因果結(jié)構(gòu)進行降維處理“分”的過程中，將整體的因果網(wǎng)絡快速拆分成若干個子因果網(wǎng)絡，對每一個子因果網(wǎng)絡都進行因果推斷，使得因果推斷能夠在高維數(shù)據(jù)環(huán)境下運行，避免傳統(tǒng)CI 測試運算復雜度和時間復雜度高的缺點。實驗結(jié)果表明，在高維數(shù)據(jù)環(huán)境下，本文算法與目前具有代表性的遞歸式因果推斷算法CAPA 在準確率幾乎持平的情況下，速度可提升2～10 倍。后續(xù)將改進推斷過程及算法框架，將本文算法與神經(jīng)網(wǎng)絡相結(jié)合，進一步提升準確率和降低算法復雜度。