胡純嚴(yán) ,胡良平 ,2*
(1.軍事科學(xué)院研究生院,北京 100850;2.世界中醫(yī)藥學(xué)會(huì)聯(lián)合會(huì)臨床科研統(tǒng)計(jì)學(xué)專業(yè)委員會(huì),北京 100029*通信作者:胡良平,E-mail:lphu927@163.com)
基于因果圖中設(shè)定的變量之間的關(guān)系,不僅可以構(gòu)建和搜索調(diào)整集,也可以檢驗(yàn)調(diào)整集,還可以基于設(shè)定的多個(gè)因果圖模型尋找共用的調(diào)整集。本文先介紹因果圖的基礎(chǔ)知識(shí)和理論,再結(jié)合兩個(gè)實(shí)例并基于SAS軟件檢驗(yàn)調(diào)整集和尋找共用調(diào)整集。
在非隨機(jī)研究中,觀察單位不會(huì)被隨機(jī)分配到處理組和對(duì)照組。因此,一些協(xié)變量可以與組分配和結(jié)果變量相關(guān)聯(lián)。在這種情況下,結(jié)果變量的值由處理的因果效應(yīng)以及與協(xié)變量的虛假關(guān)聯(lián)決定。由于協(xié)變量混淆了因果關(guān)系,如果沒有某種形式的調(diào)整(例如匹配或分層)來消除虛假關(guān)聯(lián),就無法確定處理對(duì)結(jié)果的因果關(guān)系[1]。事實(shí)上,在不完善的隨機(jī)對(duì)照研究中,由混雜協(xié)變量引起的虛假關(guān)聯(lián)也有可能存在。
研究者可以使用因果圖過程來確定在何時(shí)、何種情況下估計(jì)因果效應(yīng)[2]。為此,需要以有向無環(huán)圖(Directed acyclic graph,DAG)的形式定義因果圖模型。因果圖過程接受輸入DAG,并輸出對(duì)因果分析和效應(yīng)估計(jì)有用的結(jié)果。
通常情況下,研究者可以運(yùn)用統(tǒng)計(jì)模型來描述和分析變量之間的依賴關(guān)系,例如,分析某個(gè)特定事件發(fā)生的概率有多大,以及隨著觀察到不同的變量,該概率將如何變化[3]。然而,統(tǒng)計(jì)模型中的某些變量之間可能確實(shí)存在一定程度的關(guān)聯(lián),而另一些變量之間也可能存在虛假的關(guān)聯(lián),此時(shí),僅憑統(tǒng)計(jì)方法不足以揭示變量之間的因果效應(yīng)[4]。這是因?yàn)樽兞恐g存在的虛假關(guān)聯(lián)會(huì)在相當(dāng)大的程度上產(chǎn)生誤導(dǎo)[5]。為了通過非隨機(jī)試驗(yàn)的數(shù)據(jù)來回答因果效應(yīng),必須用一組因果假設(shè)來補(bǔ)充聯(lián)合分布函數(shù),這些假設(shè)與專業(yè)上已知的信息共同構(gòu)成因果圖模型。
因果圖模型由DAG的形式表示,DAG由三部分組成[4],即節(jié)點(diǎn)、邊和缺失邊。
節(jié)點(diǎn):DAG中的每個(gè)節(jié)點(diǎn)代表一個(gè)變量,假設(shè)該變量在研究的過程中起因或果的作用;每個(gè)變量可以有任何分布;DAG中的每個(gè)變量可以是已觀測(cè)的,也可以是未觀測(cè)的,但未觀測(cè)的變量不宜過多;通常,誤差隨機(jī)變量(獨(dú)立誤差項(xiàng))不在DAG中表示。
邊:DAG中的所有邊都是定向的,也就是說,邊由從一個(gè)節(jié)點(diǎn)指向另一個(gè)節(jié)點(diǎn)的箭頭組成,且一對(duì)節(jié)點(diǎn)之間最多只能有一條邊;邊是因果假設(shè)的圖形表示,即DAG中的邊表示一個(gè)變量對(duì)另一個(gè)變量可能產(chǎn)生的直接因果效應(yīng);假設(shè)這些因果效應(yīng)是確定性的,但它們是完全非參數(shù)的,因?yàn)槊織l邊都可以有任何函數(shù)形式;因?yàn)镈AG中的每條邊都給出了因果解釋,所以每條邊都與一對(duì)節(jié)點(diǎn)的時(shí)間順序相關(guān)聯(lián),因此,DAG不能包含有向循環(huán)。CAUSALGRAPH過程對(duì)每個(gè)模型執(zhí)行語義驗(yàn)證,以驗(yàn)證它不包含有向循環(huán)。PROC CAUSALGRAPH允許指定雙向邊。雙向邊被解釋為兩個(gè)變量之間未測(cè)量的混雜,圖形仍然是DAG。例如,式(1)與式(2)中的符號(hào)“?”與“←L→”被解釋為一對(duì)邊:
其中,節(jié)點(diǎn)L代表一個(gè)或一些未測(cè)量的變量。
缺失邊:DAG中缺失邊(即兩個(gè)節(jié)點(diǎn)不通過邊直接連接)表示直接因果效應(yīng)為零的假設(shè)。因此,DAG中缺失的邊表示比邊具有更強(qiáng)的假設(shè)。這是圖形模型的強(qiáng)零假設(shè),它在計(jì)量經(jīng)濟(jì)學(xué)文獻(xiàn)中被稱為排除限制。缺失邊對(duì)因果模型所隱含的統(tǒng)計(jì)特性有影響。
DAG中的節(jié)點(diǎn)、邊和缺失邊一起形成一個(gè)因果圖模型,對(duì)研究者關(guān)于數(shù)據(jù)生成過程的假設(shè)進(jìn)行編碼。從這些生成數(shù)據(jù)的假設(shè)開始,研究者可以使用一組在DAG上操作的圖形規(guī)則來導(dǎo)出統(tǒng)計(jì)關(guān)聯(lián)關(guān)系。
因果圖中的基本術(shù)語:①相鄰,如果DAG中的兩個(gè)變量由一條邊直接連接,則它們是相鄰的。②路徑,路徑是一個(gè)有序的變量列表,其中沒有變量出現(xiàn)超過一次,列表中的連續(xù)變量在圖中相鄰,連接路徑中連續(xù)節(jié)點(diǎn)的邊可以指向任一方向。③因果關(guān)系與非因果關(guān)系,對(duì)于路徑上的每一對(duì)連續(xù)變量,如果連接兩個(gè)變量的箭頭指向后一個(gè)變量,那么路徑就是因果關(guān)系,否則,就稱為非因果關(guān)系。④正確路徑,如果路徑以處理變量(即研究者關(guān)注的重要影響因素)開始,且不包含任何其他處理變量,則該路徑是正確的[2]。⑤關(guān)系描述,DAG編碼因果圖模型中變量之間的特定關(guān)系的描述,對(duì)于一個(gè)邊,如P→Q,其中,P是Q的父項(xiàng),Q是P的子項(xiàng);如果變量S到變量T之間存在因果路徑,那么T是S的后代,S是T的祖先,因此,變量S的后代集是由S直接或間接引起的所有變量的集合,同樣,變量T的祖先集是T的所有直接或間接原因的集合。⑥碰撞器,對(duì)于路徑上的變量V,如果V有兩個(gè)指向它的箭頭(每側(cè)一個(gè)),則V是路徑上的碰撞器。否則,就是非碰撞器,碰撞器的定義是特定于路徑的,變量可以是一條路徑上的碰撞器,也可以是另一條路徑上的非碰撞器。⑦因果與非因果關(guān)聯(lián),變量集之間的統(tǒng)計(jì)關(guān)聯(lián)可分為因果部分和非因果部分或虛假部分,如果所有虛假的關(guān)聯(lián)都能被消除,那么因果關(guān)系就被認(rèn)為是確定的[2]。⑧調(diào)整或消除虛假關(guān)聯(lián),虛假關(guān)聯(lián)通常通過某種形式的統(tǒng)計(jì)調(diào)整或條件作用來消除,例如,研究者可以通過在回歸模型中包含一個(gè)變量作為回歸自變量,或者通過按變量的水平對(duì)資料分層來進(jìn)行調(diào)整計(jì)算,還可以使用因果圖模型來確定調(diào)整集中必須包含哪些變量。
3.1.1 實(shí)例及其背景信息
【例1】圖1所示的因果圖模型(改編自文獻(xiàn)[6])用于檢驗(yàn)受試者血清尿酸鹽(Urate)與心血管疾?。–VD)風(fēng)險(xiǎn)之間的關(guān)系。
圖1模型中各變量的含義如下:Urate,尿酸鹽(處理變量);CVD,心血管疾?。ńY(jié)果變量);Anti-HypertensiveUse,抗高血壓藥物使用;Creatinine,肌酐(測(cè)量的血清肌酐水平);Diabetes,糖尿??;Ethnicity,種族;Gender,性別;Gout,痛風(fēng);HbA1c,糖化血紅蛋白;MedicationPropensity,藥物傾向(反映個(gè)人服用處方藥傾向的潛在變量);Nutrition,營養(yǎng)(反映飲食或營養(yǎng)的潛在變量);Obesity,肥胖(體重指數(shù)指標(biāo));CurrentBP,當(dāng)前血壓;CurrentHDL,當(dāng)前高密度脂蛋白(高密度脂蛋白膽固醇指標(biāo));PreviousBP,先前血壓(研究前血壓指標(biāo));PreviousHDL,先前高密度脂蛋白(研究前高密度脂蛋白膽固醇指標(biāo));Smoking,吸煙(當(dāng)前吸煙狀態(tài)指標(biāo));StatinUse,他汀類藥物使用。
圖1 血清尿酸鹽對(duì)心血管疾病風(fēng)險(xiǎn)影響的因果圖模型Figure 1 Causal model of the effect of serum urate on risk of cardiovascular disease
變量MedicationPropensity和Nutrition對(duì)應(yīng)潛在變量,故無法觀察到。還假設(shè)未觀察到變量PreviousBP、PreviousHDL和Obesity。在考慮諸多協(xié)變量的前提下,試分析處理變量Urate對(duì)結(jié)果變量CVD是否具有因果效應(yīng)。
3.1.2 用SAS實(shí)現(xiàn)因果效應(yīng)分析
3.1.2.1 檢驗(yàn)調(diào)整集
根據(jù)圖1所顯示的因果圖模型,變量Urate和CVD之間的統(tǒng)計(jì)關(guān)聯(lián)反映了真實(shí)因果關(guān)聯(lián)和附加虛假或非因果關(guān)聯(lián)的組合。為了揭示這兩個(gè)變量之間的真正因果聯(lián)系,Thornley等[6]考慮對(duì)變量CurrentHDL、Ethnicity、Gender、HbA1c和 Smoking進(jìn)行調(diào)整。以下程序調(diào)用因果圖過程,以檢驗(yàn)該調(diào)整集是否可用于根據(jù)因果圖模型估計(jì)變量Urate對(duì)CVD的因果效應(yīng)。設(shè)所需要的SAS程序如下:
【SAS輸出結(jié)果及解釋】
檢驗(yàn)調(diào)整集是否為有效的因果效應(yīng)的輸出結(jié)果,見表1。
表1 檢驗(yàn)調(diào)整集是否為有效的因果效應(yīng)的輸出結(jié)果Table 1 Output results for testing whether the adjustment set is a valid causal effect
【表1中有關(guān)內(nèi)容的說明】第2列的“大小”指協(xié)變量的個(gè)數(shù)(本例有5個(gè));第3列的“有效”指協(xié)變量的集合對(duì)研究變量Urate對(duì)CVD因果效應(yīng)是否有效(本例經(jīng)檢驗(yàn),其結(jié)果為“無效”,輸出中用“否”表示);第4列的“最小”指所找到的調(diào)整集是否為最小的調(diào)整集(本例的調(diào)整集包含5個(gè)協(xié)變量,故它不是最小的調(diào)整集);表1中5個(gè)協(xié)變量對(duì)結(jié)果變量的影響具有統(tǒng)計(jì)學(xué)意義,即在研究變量Urate對(duì)CVD因果效應(yīng)中是不可忽視的,但由第3列上的“否”可知,包含5個(gè)協(xié)變量的調(diào)整集不足以估計(jì)變量Urate對(duì)CVD的因果效應(yīng)。
3.1.2.2 枚舉路徑和搜索有效調(diào)整集
要了解建議的調(diào)整集無效的原因,可以請(qǐng)求枚舉將處理變量與模型中的結(jié)果變量聯(lián)系起來的正確路徑,也可以使用該過程搜索有效的調(diào)整集。以下程序調(diào)用CAUSALGRAPH過程來執(zhí)行這兩項(xiàng)任務(wù),SAS程序與上文“3.1.2.1節(jié)”的程序基本相同,僅以下兩句稍有改變:
proc causalgraph compact list;
testid Gender HbA1c Ethnicity Smoking CurrentHDL/paths=(noncausal nonblocked);
【SAS輸出結(jié)果及解釋】
輸出了兩條路徑,第1條路徑為:
Urate<==Nutrition==>PreviousHDL==>StatinUse<==MedicationPropensity==> AntiHypertensiveUse==>CurrentBP==>CVD
第2條路徑為:
Urate<==Nutrition==>Obesity==>PreviousBP==>AntiHypertensiveUse==>CurrentBP==>CVD
注:沒有滿足“Thor12”的指定準(zhǔn)則的調(diào)整集。
【關(guān)于上述兩條路徑的說明】第一條路徑?jīng)]有被阻塞,因?yàn)樽兞縎tatinUse是路徑上的碰撞器,并且它的一個(gè)子變量CurrentHDL出現(xiàn)在調(diào)整集中。第二條路徑不包含任何碰撞器,但不會(huì)被阻塞,因?yàn)樗话ㄗh調(diào)整集中的任何變量。這兩條路徑都是非因果關(guān)系的路徑。如果使用建議調(diào)整集,變量Urate和CVD之間的某些關(guān)聯(lián)可能歸因于這兩條非因果路徑,并且因果效應(yīng)無法正確估計(jì)。
3.1.2.3 減少不可測(cè)量變量
為了獲得可識(shí)別性,研究者可以考慮收集附加數(shù)據(jù)。例如,如果研究者要收集變量Obesity和PreviousHDL的數(shù)據(jù),以便這兩個(gè)變量成為可測(cè)量的變量,這等價(jià)于在上面的輸出結(jié)果中阻塞兩條非因果路徑。下面的程序可實(shí)現(xiàn)這一目標(biāo):
SAS程序與前面“3.1.2.1節(jié)”的程序基本相同,僅以下兩句稍有改變:
unmeasured Nutrition PreviousBP Medication Propensity;
testid Gender HbA1c Ethnicity Smoking Current HDL PreviousHDL Obesity;
【SAS輸出結(jié)果及解釋】
TESTID語句中包含7個(gè)變量時(shí)的輸出結(jié)果,見表2。
表2 TESTID語句中包含7個(gè)變量時(shí)輸出的結(jié)果Table 2 Output results when the testid statement contains 7 variables
由表2可知,建議的調(diào)整集(含7個(gè)協(xié)變量)在“有效”列中標(biāo)記為“是”,說明該調(diào)整集足以估計(jì)變量Urate對(duì)CVD的因果效應(yīng)。
3.2.1 因果圖模型確切結(jié)構(gòu)未知時(shí)的解決方案
當(dāng)研究者不確定因果圖模型的確切結(jié)構(gòu)或數(shù)據(jù)可能存在其他因果圖模型時(shí),如何使用因果圖過程找到可用于估計(jì)因果效應(yīng)的調(diào)整集?例如,研究者可能不確定是否要在因果圖模型中包含哪些邊,邊具有什么方向,或者模型包含哪些協(xié)變量,則可以在過程中使用多個(gè)MODEL語句來指定每個(gè)合理的因果圖模型,然后使用PROC CAUSALGRAPH語句中的公共選項(xiàng)來搜索對(duì)所有模型有效的調(diào)整集。
對(duì)于單個(gè)因果圖模型,如果存在調(diào)整集,則可以使用該調(diào)整集來估計(jì)因果效應(yīng)。如果多個(gè)因果圖模型共享一個(gè)調(diào)整集,那么無論哪個(gè)因果圖模型反映了真實(shí)的數(shù)據(jù)生成過程,都可以使用該調(diào)整集來估計(jì)因果效應(yīng)。
3.2.2 實(shí)例及其背景信息
【例2】沿用例1的實(shí)例和背景信息,在圖1的基礎(chǔ)上,添加一些新的信息。見圖2。
圖2 血清尿酸鹽對(duì)心血管疾病風(fēng)險(xiǎn)影響的兩種可能因果模型Figure 2 Two possible causal models of the effect of serum urateon risk of cardiovascular disease
圖2中的模型源自Thornley等[6]開發(fā)的更大模型。這些模型研究了受試者血清Urate與CVD風(fēng)險(xiǎn)之間的關(guān)系。有關(guān)模型中變量的含義,請(qǐng)參見本文例1。圖2描述了變量StatinUse對(duì)CVD效應(yīng)的兩個(gè)合理因果假設(shè)。一個(gè)模型假設(shè)變量StatinUse對(duì)CVD的所有因果效應(yīng)都是由變量CurrentHDL介導(dǎo)的;另一個(gè)模型,除介導(dǎo)效應(yīng)外,還包括變量StatinUse對(duì)CVD的直接影響。這兩個(gè)模型的不同之處在于一條邊,如圖2中的虛線箭頭所示。在本例的兩個(gè)模型中,假設(shè)變量Nutrition和PreviousHDL是不可測(cè)量的,故由虛線輪廓的節(jié)點(diǎn)表示。
3.2.3 用SAS實(shí)現(xiàn)因果效應(yīng)分析
調(diào)用PROC CAUSALGRAPH為圖2中的兩個(gè)模型構(gòu)造調(diào)整集。所需要的SAS程序與例1中的程序基本相同,此處從略。
【SAS輸出結(jié)果及解釋】
因篇幅所限,具體輸出結(jié)果從略?,F(xiàn)將輸出結(jié)果分3部分概述如下。
第1部分輸出結(jié)果:圖2中定義的第1個(gè)因果圖模型(Thor12SimpleHDL)的計(jì)算結(jié)果,由6個(gè)可觀測(cè)的協(xié)變量(CurrentHDL、Ethnicity、Gender、Gout、Smoking及StatinUse)構(gòu)造出18個(gè)調(diào)整集,其中,有4個(gè)調(diào)整集屬于最小調(diào)整集(每個(gè)最小調(diào)整集僅包含2個(gè)協(xié)變量)。第2部分輸出結(jié)果:圖2中定義的第2個(gè)因果圖模型(Thor12AltHDL)的計(jì)算結(jié)果,由6個(gè)可觀測(cè)的協(xié)變量(CurrentHDL、Ethnicity、Gender、Gout、Smoking及StatinUse)構(gòu)造出12個(gè)調(diào)整集,其中,有2個(gè)調(diào)整集屬于最小調(diào)整集(每個(gè)最小調(diào)整集僅包含2個(gè)協(xié)變量)。第3部分輸出結(jié)果:2種因果圖模型共用的調(diào)整集。巧合的是,在兩個(gè)模型中有效的調(diào)整集與在具有額外邊的模型中有效的調(diào)整集相同。
由于至少存在一個(gè)共用調(diào)整集,該分析表明,通過“第3部分輸出結(jié)果”中的任何一個(gè)調(diào)整集,無論變量StatinUse是否對(duì)CVD有直接效應(yīng),都可以估計(jì)變量Urate對(duì)CVD的因果效應(yīng)。
基于因果圖模型所做的統(tǒng)計(jì)分析,其結(jié)論正確與否,主要取決于以下關(guān)鍵點(diǎn):圖中所設(shè)定的變量之間的關(guān)系(包括變量的先后順序,發(fā)出箭頭和接收箭頭的變量,變量是否為觀測(cè)變量,哪些變量是處理變量、協(xié)變量和結(jié)果變量);是否有專業(yè)知識(shí)為依據(jù);是否符合基本常識(shí)和常規(guī)邏輯。這幾個(gè)方面在多重線性回歸分析[7]、多重Logistic回歸分析[8]、多重 Cox模型回歸分析[9]和結(jié)構(gòu)方程模型分析中[10],也都是極為重要的。
本文介紹了因果圖模型的基礎(chǔ)知識(shí)和理論,涉及因果圖的組成和術(shù)語。針對(duì)兩個(gè)實(shí)例并借助SAS軟件,實(shí)現(xiàn)了兩項(xiàng)因果圖分析任務(wù),即檢驗(yàn)調(diào)整集和尋找多個(gè)模型的共用調(diào)整集。