王 威 綜述 楊 帆 審校

在基于臨床數(shù)據(jù)的研究中,無(wú)論是回顧性還是前瞻性的研究,常常會(huì)出現(xiàn)數(shù)據(jù)缺失的情況。當(dāng)出現(xiàn)了缺失值,我們應(yīng)給予妥當(dāng)?shù)奶幚?以便做出可靠的統(tǒng)計(jì)推斷。目前,應(yīng)對(duì)缺失值的方法是多種多樣的,主要有刪除法和填補(bǔ)法。刪除法又稱完整案例分析(complete case analysis)或列表刪除(listwise deletion),即刪除存在缺失值的所有觀測(cè),對(duì)保留下來(lái)的無(wú)缺失值的觀測(cè)進(jìn)行統(tǒng)計(jì)分析。此方法的操作是最簡(jiǎn)單易行的,因此為大多數(shù)研究者所采用,也是多數(shù)統(tǒng)計(jì)分析軟件所默認(rèn)采用的方法。然而,此方法在實(shí)際運(yùn)用中可能會(huì)存在兩點(diǎn)不足:一是缺失值的存在可能不是隨機(jī)的,即存在缺失值的觀測(cè)與完整觀測(cè)之間存在某些方面的差異,且這些差異會(huì)給后續(xù)的參數(shù)估計(jì)帶來(lái)偏倚;二是在研究變量數(shù)較多的情況下,所有研究變量均無(wú)缺失的觀測(cè)可能會(huì)比較少,運(yùn)用此法會(huì)舍棄過多觀測(cè),不僅會(huì)嚴(yán)重降低投入產(chǎn)出比,而且會(huì)增大參數(shù)估計(jì)的標(biāo)準(zhǔn)差和置信區(qū)間,降低統(tǒng)計(jì)功效。只有在存在缺失值的觀測(cè)占比很小的時(shí)候(比如<5%),所造成的參數(shù)估計(jì)的偏倚和統(tǒng)計(jì)功效的降低才可以近似忽略[1]。因此,多建議對(duì)缺失值進(jìn)行填補(bǔ)后再進(jìn)行分析。填補(bǔ)又稱插補(bǔ),可分為單一填補(bǔ)法和多重填補(bǔ)法。單一填補(bǔ)包括均值填補(bǔ),啞變量填補(bǔ)和基于回歸模型的單一填補(bǔ)等。單一填補(bǔ)因?yàn)闀?huì)降低被填補(bǔ)變量的不確定性,縮小參數(shù)估計(jì)的標(biāo)準(zhǔn)差,已逐漸被淘汰,本文不作詳細(xì)介紹[2]。多重填補(bǔ)法因其考慮了缺失的不確定性等優(yōu)點(diǎn)越來(lái)越受到大家的推崇[3]。另外,廣義線性模型(generalized linear models, GLM),囊括線性回歸模型、二元Logistic回歸模型、Poisson回歸模型等,在臨床數(shù)據(jù)分析中占據(jù)著重要地位,其中二元Logistic回歸常常作為主要的多因素模型,用以產(chǎn)出因果推斷的統(tǒng)計(jì)學(xué)依據(jù)[4-5]。目前,雖然二元Logistic回歸與多重填補(bǔ)技術(shù)的結(jié)合在危險(xiǎn)因素分析的醫(yī)學(xué)研究中的應(yīng)用越來(lái)越多,如Foerster等[6]通過二元logistic回歸和多重填補(bǔ)方法建立了一個(gè)風(fēng)險(xiǎn)分層模型以更好的識(shí)別適合上尿路上皮癌內(nèi)鏡下保留腎臟手術(shù)的患者,但以實(shí)例分析介紹兩者如何結(jié)合應(yīng)用的文獻(xiàn)較少。本文就數(shù)據(jù)缺失的模式和比例、多重填補(bǔ)的流程進(jìn)行簡(jiǎn)單的梳理,并使用急性腎損傷的數(shù)據(jù)演示多重填補(bǔ)(mice程序包)與廣義線性模型的結(jié)合應(yīng)用,以期為含缺失值的臨床數(shù)據(jù)的統(tǒng)計(jì)分析提供參考依據(jù)。

數(shù)據(jù)缺失的模式與比例

多重填補(bǔ)的創(chuàng)始人Little和Rubin教授將數(shù)據(jù)缺失的原因分為以下三種模式:完全隨機(jī)缺失(missing completely at random, MCAR)、隨機(jī)缺失(missing at random, MAR)和非隨機(jī)缺失(missing not at random, MNAR)。完全隨機(jī)缺失指的是數(shù)據(jù)的缺失與否既不能歸因于已觀測(cè)的變量,也不能歸因于未觀測(cè)的變量,即存在缺失的觀測(cè)與完整觀測(cè)來(lái)自于同一個(gè)分布的總體,完整案例分析方法僅適用于這種缺失模式[7]。但有學(xué)者認(rèn)為這種缺失模式只是理論上存在的,是隨機(jī)缺失的一種特殊情況[8]。隨機(jī)缺失指的是數(shù)據(jù)的缺失與否取決于已觀測(cè)的變量,即存在缺失的觀測(cè)和完整觀測(cè)的差別可以被其他觀測(cè)變量解釋。有研究者認(rèn)為,“隨機(jī)缺失”一詞的命名有歧義,其 “隨機(jī)”一詞與“完全隨機(jī)缺失”中的“隨機(jī)”意義是不同的,前者是一種有條件的、可控制的缺失,更精確的表達(dá)方式為隨機(jī)條件缺失(conditionally missing at random)[7-9]。非隨機(jī)缺失指的是數(shù)據(jù)的缺失與否取決于未觀測(cè)到的變量,即存在缺失的觀測(cè)和完整觀測(cè)的差別無(wú)法用所有已知變量解釋。將多重填補(bǔ)技術(shù)應(yīng)用于含有前兩種缺失(即完全隨機(jī)缺失和隨機(jī)缺失)的數(shù)據(jù),可得到近似無(wú)偏的參數(shù)估計(jì)結(jié)果,而非隨機(jī)缺失則不在多重填補(bǔ)技術(shù)的應(yīng)用范圍。在實(shí)際應(yīng)用中,事先幾乎不可能知道數(shù)據(jù)的缺失究竟屬于哪種類型,那么在應(yīng)用多重填補(bǔ)法對(duì)缺失數(shù)據(jù)進(jìn)行處理之時(shí),需要先假定數(shù)據(jù)缺失符合隨機(jī)缺失。只要構(gòu)建填補(bǔ)模型時(shí)納入足夠多的變量,數(shù)據(jù)缺失的模式就會(huì)非常接近隨機(jī)缺失,可以使用多重填補(bǔ)法進(jìn)行處理[10-11]?；蛘呖梢赃x擇多個(gè)可能適用的多重填補(bǔ)模型分別填補(bǔ),對(duì)各個(gè)模型所得完整數(shù)據(jù)集進(jìn)行相同的統(tǒng)計(jì)分析,將各自的結(jié)果匯總,進(jìn)行敏感性分析[12]。

在應(yīng)用多重填補(bǔ)法處理缺失數(shù)據(jù)時(shí),數(shù)據(jù)缺失的比例對(duì)多重填補(bǔ)的影響目前尚存在爭(zhēng)論[13]。人們所推薦的多重填補(bǔ)最多可接受的數(shù)據(jù)缺失的比例在5%～50%[14-16]。近年有學(xué)者對(duì)模擬數(shù)據(jù)進(jìn)行了大量的填補(bǔ)運(yùn)算,發(fā)現(xiàn)相比于數(shù)據(jù)缺失的比例,納入填補(bǔ)模型的輔助變量所起的作用更加明顯。即便是在缺失數(shù)據(jù)比例高達(dá)90%的情況下,只要在填補(bǔ)模型中納入合適的輔助變量,填補(bǔ)后的數(shù)據(jù)依然可以進(jìn)行無(wú)偏倚的統(tǒng)計(jì)推斷[13]。然而此研究所填補(bǔ)變量均為符合正態(tài)分布的連續(xù)變量,且為模擬數(shù)據(jù),其結(jié)論可能不適用于分類變量和其他分布類型的變量。

因此,在進(jìn)行多重填補(bǔ)之前,應(yīng)對(duì)所需填補(bǔ)數(shù)據(jù)的缺失值進(jìn)行多方面的評(píng)估,必要時(shí)需了解各個(gè)變量收集、整理、清洗的全過程,盡可能納入更多輔助變量進(jìn)入填補(bǔ)模型,并使用多個(gè)填補(bǔ)模型,進(jìn)行敏感性分析,方能盡可能減少偏倚。

多重填補(bǔ)的流程

以R語(yǔ)言程序包mice(multivariate imputation by chained equation)為例,由該程序包及相應(yīng)函數(shù)實(shí)現(xiàn)數(shù)據(jù)多重填補(bǔ)的原理和過程可知 (圖1),一個(gè)完整的數(shù)據(jù)多重填補(bǔ)過程應(yīng)包括缺失值的填補(bǔ)和使用所得數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析兩大部分[12,17]。因此,當(dāng)我們使用多重填補(bǔ)法對(duì)含缺失值的數(shù)據(jù)進(jìn)行回歸分析時(shí),建議按照如下流程進(jìn)行數(shù)據(jù)填補(bǔ)和回歸分析:

圖1 基于多重填補(bǔ)的廣義線性模型的實(shí)現(xiàn)流程

查看數(shù)據(jù)缺失情況當(dāng)我們拿到一個(gè)數(shù)據(jù)集后,首先應(yīng)當(dāng)明確數(shù)據(jù)的缺失模式,即前文中提到的MCAR、MAR和 MNAR,但這一點(diǎn)往往難以做到。不過,我們能做的是查看數(shù)據(jù)集的缺失情況,包括明確含缺失值的變量、變量類型和缺失比例,這對(duì)是否進(jìn)行填補(bǔ)和多重填補(bǔ)中一些命令行參數(shù)的設(shè)定很重要。比如,若因變量有缺失,一般我們會(huì)選擇刪除相應(yīng)的觀測(cè)行,而不對(duì)因變量進(jìn)行填補(bǔ);若缺失的變量既有連續(xù)變量,也有分類變量,那我們可以通過函數(shù)factor()先對(duì)分類變量進(jìn)行定義,或在多重填補(bǔ)時(shí)為每個(gè)需要填補(bǔ)的變量設(shè)定正確的參數(shù)defaultMethod;若某個(gè)或某些變量的缺失比例較大(>50%),我們一般不建議進(jìn)行數(shù)據(jù)填補(bǔ),而建議放棄使用該變量。

缺失值的多重填補(bǔ)當(dāng)缺失比例較低時(shí),我們便可以采用函數(shù)mice()嘗試缺失值的多重填補(bǔ)。在這個(gè)過程中,有3個(gè)重要的參數(shù)需要設(shè)定。第一個(gè)為參數(shù)m,用來(lái)設(shè)定多重填補(bǔ)生成的完整數(shù)據(jù)集的數(shù)量,默認(rèn)為5。第二個(gè)為method或defaultMethod,用來(lái)指定變量的填補(bǔ)方法,比如pmm對(duì)應(yīng)定量變量、logreg對(duì)應(yīng)二分類變量、polyreg對(duì)應(yīng)無(wú)序多分類變量、polr對(duì)應(yīng)有序多分類變量。第三個(gè)為seed,設(shè)置種子數(shù),使得下次再運(yùn)行代碼時(shí)所得到的數(shù)據(jù)集與前次一致,保證可重復(fù)性。

完整數(shù)據(jù)集的統(tǒng)計(jì)分析在完成數(shù)據(jù)的填補(bǔ)之后,我們會(huì)得到一個(gè)mids對(duì)象,該對(duì)象包含了通過鏈?zhǔn)蕉嘣钛a(bǔ)得到的多個(gè)完整數(shù)據(jù)集。由于函數(shù)glm()中的參數(shù)data所識(shí)別的數(shù)據(jù)集類型為data frame (數(shù)據(jù)框)、list(列表)或環(huán)境 (environment),那么我們就無(wú)法直接使用函數(shù)glm()對(duì)mids對(duì)象進(jìn)行回歸分析。解決的方法是,聯(lián)用函數(shù)with()和glm(),對(duì)mids對(duì)象中的每個(gè)數(shù)據(jù)集均進(jìn)行回歸分析。

統(tǒng)計(jì)分析結(jié)果的匯總在使用with()完成了每個(gè)數(shù)據(jù)集的回歸分析后,我們會(huì)得到每個(gè)數(shù)據(jù)集對(duì)應(yīng)的回歸分析結(jié)果,即多套回歸分析結(jié)果。接下來(lái),我們需要使用函數(shù)pool()對(duì)多個(gè)回歸分析結(jié)果進(jìn)行匯總,綜合為包含回歸系數(shù)、標(biāo)準(zhǔn)誤、P值等在內(nèi)的一套回歸分析結(jié)果。

在廣義線性模型中的應(yīng)用

案例情況如下:某腎內(nèi)科醫(yī)師回顧性收集了109例在院內(nèi)發(fā)生急性腎損傷和231例未發(fā)生急性腎損傷患者的相關(guān)信息,包括一般特征(性別、年齡、BMI等級(jí))和入院時(shí)的實(shí)驗(yàn)室檢查結(jié)果(中性粒細(xì)胞和淋巴細(xì)胞比值、尿素氮、尿肌酐、尿酸、鈉離子、鉀離子、超敏反應(yīng)蛋白、血清乳酸),擬通過二元logistic回歸分析急性腎損傷的影響因素。

在廣義線性模型中,將glm()函數(shù)的famliy參數(shù)設(shè)置為binomial時(shí),所擬合的模型便是二元logistic回歸模型。接下來(lái),我們以此類回歸模型為例,展示基于多重填補(bǔ)數(shù)據(jù)集、缺失數(shù)據(jù)集和原始數(shù)據(jù)集的廣義線性模型的結(jié)果,并將基于多重填補(bǔ)的回歸分析結(jié)果與基于無(wú)缺失值原始數(shù)據(jù)的回歸分析結(jié)果進(jìn)行比較(本文所用的R代碼均可通過github網(wǎng)站查看https://github.com/hamody316/mice/blob/main/README.md#mice)。

由于本次所使用的的數(shù)據(jù)集并無(wú)缺失值,因此我們使用R程序包simFrame中的函數(shù)setNA()為性別、年齡、體質(zhì)量指數(shù)(BMI)等級(jí)、尿肌酐、鈉離子和血清乳酸構(gòu)造了缺失值,使得他們均有5%的隨機(jī)缺失。另外,在該數(shù)據(jù)集中,性別為二分類變量,BMI等級(jí)為有序多分類變量(等級(jí)變量),其他均為連續(xù)變量。

假定我們已確定了納入到多因素回歸分中的變量為性別、年齡、BMI等級(jí)、尿肌酐、鈉離子和血清乳酸。在設(shè)定了填補(bǔ)數(shù)據(jù)集數(shù)量(m=10)、變量的填補(bǔ)方法和種子數(shù)后,我們得到了10個(gè)完整的數(shù)據(jù)集,使用這10個(gè)數(shù)據(jù)集進(jìn)行回歸分析的匯總結(jié)果如表1。

表1 基于多重填補(bǔ)的回歸分析結(jié)果

為了能夠說明多重填補(bǔ)的效果,我們用含缺失值的數(shù)據(jù)集和無(wú)缺失值的原始數(shù)據(jù)集分別進(jìn)行了回歸分析,所得結(jié)果見表2和表3;同時(shí),我們還對(duì)表3和表1的回歸系數(shù)進(jìn)行了差值計(jì)算(表4)。

表2 基于含缺失值數(shù)據(jù)集的回歸分析結(jié)果

表3 基于原始數(shù)據(jù)的回歸分析結(jié)果

表4 回歸系數(shù)的差值

由表3的結(jié)果可知,所納入的多個(gè)變量在基于多重填補(bǔ)的回歸分析中具有統(tǒng)計(jì)學(xué)意義,與急性腎損傷有統(tǒng)計(jì)學(xué)關(guān)聯(lián),這與表1中的結(jié)論一致。由表2的結(jié)果可知,在基于含缺失值的回歸分析結(jié)果中,年齡和BMI等級(jí)與急性腎損傷無(wú)統(tǒng)計(jì)學(xué)關(guān)聯(lián),這與基于原始數(shù)據(jù)的回歸分析的結(jié)果不一致,說明表2中的結(jié)果發(fā)生了偏倚。另外,由表4可知,基于多重填補(bǔ)數(shù)據(jù)得到的變量系數(shù)與基于原始數(shù)據(jù)得到的變量系數(shù)的差值很小,進(jìn)一步說明多重填補(bǔ)得到結(jié)果的穩(wěn)健性和可靠性。當(dāng)然,本文多重填補(bǔ)的結(jié)果并不意味著多重填補(bǔ)對(duì)所有含缺失值的數(shù)據(jù)集而言都是最佳的填補(bǔ)方法,我們應(yīng)該根據(jù)數(shù)據(jù)類型、缺失情況等選擇合適的填補(bǔ)方法,比如隨機(jī)森林填補(bǔ)法、K鄰近值法等[18-19]。

小結(jié):本文對(duì)缺失值的常見處理方法、缺失值模式、缺失值比例、多重填補(bǔ)的流程進(jìn)行了簡(jiǎn)單總結(jié),并通過急性腎損傷的影響因素分析的實(shí)例,展示了基于多重填補(bǔ)的廣義線性模型的分析過程及結(jié)果。對(duì)比基于多重填補(bǔ)的回歸分析結(jié)果、基于原始數(shù)據(jù)的回歸分析結(jié)果和基于缺失數(shù)據(jù)的回歸分析結(jié)果可知,在缺失比例較低時(shí),前兩者的回歸分析結(jié)果(如各變量的回歸系數(shù))雖然有一定的差值,但整體結(jié)論并未出現(xiàn)偏倚,即在回歸方程中各變量的統(tǒng)計(jì)學(xué)意義是一致的。這也表明,在條件適當(dāng)?shù)那闆r下,可以在腎病相關(guān)數(shù)據(jù)的統(tǒng)計(jì)分析中使用多重填補(bǔ)。本文結(jié)合腎病相關(guān)案例數(shù)據(jù)展示了如何使用R語(yǔ)言實(shí)現(xiàn)該方法,并對(duì)該方法所得結(jié)果的穩(wěn)健性進(jìn)行了驗(yàn)證,希望能為廣大醫(yī)護(hù)或科研人員在處理缺失數(shù)據(jù)的思路和實(shí)踐上提供參考。