復(fù)雜抽樣數(shù)據(jù)多水平模型分析方法及其應(yīng)用

于石成廖加強(qiáng)于 妺郭 瑩肖革新△金承剛馮國(guó)雙胡躍華馬林茂

復(fù)雜抽樣數(shù)據(jù)多水平模型分析方法及其應(yīng)用

于石成1廖加強(qiáng)2于 妺1郭 瑩1肖革新1△金承剛3馮國(guó)雙1胡躍華1馬林茂1

目的本文通過抽樣調(diào)查實(shí)例，闡述多階段抽樣、不等抽樣概率和事后分層特性不同產(chǎn)生的復(fù)雜抽樣數(shù)據(jù)，其應(yīng)用多水平模型分析的原理和方法。方法對(duì)我國(guó)某省行為危險(xiǎn)因素抽樣調(diào)查的數(shù)據(jù)，應(yīng)用未加權(quán)和加權(quán)的隨機(jī)截距l(xiāng)ogistic回歸模型分析了某些因素與跌倒性傷害的關(guān)系。結(jié)果實(shí)際分析包括50個(gè)區(qū)縣（PSU），250個(gè)鄉(xiāng)鎮(zhèn)街道（2水平），12086個(gè)體（1水平）。未加權(quán)估計(jì)結(jié)果顯示：對(duì)跌倒性傷害有統(tǒng)計(jì)學(xué)影響的變量是健康狀況中等和差、未被雇傭和未婚，年齡為負(fù)相關(guān)，即年齡越大，發(fā)生跌倒性傷害的危險(xiǎn)性越小；復(fù)雜抽樣2水平logistic回歸分析顯示：對(duì)跌倒性傷害有統(tǒng)計(jì)學(xué)影響的變量與未加權(quán)的結(jié)果基本一致，但未婚失去了統(tǒng)計(jì)學(xué)意義。體重指數(shù)、性別和受教育程度與跌倒性傷害的發(fā)生沒有統(tǒng)計(jì)學(xué)聯(lián)系。結(jié)論與未加權(quán)的結(jié)果比，加權(quán)分析對(duì)跌倒性傷害有統(tǒng)計(jì)學(xué)影響的變量基本一致，但加權(quán)復(fù)雜抽樣PMLE估計(jì)的標(biāo)準(zhǔn)誤偏大，結(jié)果更保守；對(duì)性別的分析發(fā)現(xiàn)，加權(quán)后的結(jié)果符合目前對(duì)跌倒性傷害發(fā)生機(jī)制的認(rèn)識(shí)，因此納入權(quán)重的多水平分析方法對(duì)該資料可能更合理。

復(fù)雜抽樣 多水平模型 多階段抽樣 隨機(jī)效應(yīng)logistic回歸

目前主流統(tǒng)計(jì)分析軟件，如SAS，SPSS，Stata，MPlus和SUDAN，都將抽樣權(quán)重納入統(tǒng)計(jì)分析過程，除可進(jìn)行復(fù)雜抽樣數(shù)據(jù)的描述性統(tǒng)計(jì)分析外，還可進(jìn)行復(fù)雜抽樣數(shù)據(jù)的多元線性回歸、logistic回歸、Poisson回歸和Cox回歸等，使得復(fù)雜抽樣數(shù)據(jù)的統(tǒng)計(jì)推斷方法越來越多地在數(shù)據(jù)分析中得到應(yīng)用［1-4］。復(fù)雜抽樣數(shù)據(jù)大多具有層次結(jié)構(gòu)即多水平，其特點(diǎn)是反應(yīng)變量的分布在個(gè)體間不具獨(dú)立性，存在地理距離內(nèi)、行政區(qū)劃內(nèi)或特定空間范圍內(nèi)的聚集性［5］。多水平模型在醫(yī)學(xué)領(lǐng)域已有多年的應(yīng)用，在處理層次結(jié)構(gòu)數(shù)據(jù)上已發(fā)揮了重要的作用［6-8］；復(fù)雜抽樣數(shù)據(jù)的分析既要考慮抽樣權(quán)重，又要兼顧數(shù)據(jù)的層次結(jié)構(gòu)，一般是將權(quán)重納入廣義線性混合模型（generalized linearm ixed models，GLMMs）來處理這類數(shù)據(jù)［9-11］。當(dāng)GLMMs將抽樣權(quán)重納入模型后，使GLMMs能處理復(fù)雜抽樣數(shù)據(jù)，解決了復(fù)雜抽樣數(shù)據(jù)多水平模型統(tǒng)計(jì)分析和計(jì)算問題。鑒于目前我國(guó)全國(guó)性的流行病學(xué)抽樣調(diào)查多采用多階段抽樣設(shè)計(jì)，并且數(shù)據(jù)具有層次結(jié)構(gòu)，因此復(fù)雜抽樣數(shù)據(jù)的多水平模型分析已有明顯的應(yīng)用價(jià)值。

多水平復(fù)雜抽樣數(shù)據(jù)的線性模型，Pfeffermann應(yīng)用偽最大似然估計(jì)（pseudo-maximum-likelihood estimation）算法來估計(jì)模型的參數(shù)，效果很好。但對(duì)廣義線性混合模型（generalized linearm ixed models），認(rèn)為較好的模型估計(jì)方法是全偽最大似然估計(jì)（full pseudo-maximum-likelihood estimation），它應(yīng)用自適應(yīng)積分法（adaptive quadrature）估計(jì)模型參數(shù)，標(biāo)準(zhǔn)誤的估計(jì)采用泰勒線性化三明治估計(jì)量［10］。目前國(guó)外復(fù)雜抽樣數(shù)據(jù)多水平模型理論和應(yīng)用研究已有二十幾年的歷史，我國(guó)在這方面的理論研究鮮有報(bào)告，可能由于數(shù)據(jù)本身和應(yīng)用的限制，在醫(yī)學(xué)研究領(lǐng)域國(guó)內(nèi)還沒有應(yīng)用復(fù)雜抽樣多水平模型來處理流行病學(xué)抽樣調(diào)查數(shù)據(jù)。本文介紹了復(fù)雜抽樣數(shù)據(jù)多水平模型分析的原理和方法，并用STATA軟件對(duì)我國(guó)某省行為危險(xiǎn)因素調(diào)查的數(shù)據(jù)擬合復(fù)雜抽樣數(shù)據(jù)多水平模型，并與未考慮復(fù)雜抽樣的多水平模型結(jié)果進(jìn)行了比較和解釋，以闡述具有層次結(jié)構(gòu)的復(fù)雜抽樣數(shù)據(jù)應(yīng)用復(fù)雜抽樣數(shù)據(jù)多水平模型的合理性。

模型基本原理

廣義線性混合模型參數(shù)估計(jì)構(gòu)造一個(gè)常邊際對(duì)數(shù)似然函數(shù)（usualmarginal log likelihood function），其公式如下：

復(fù)雜抽樣數(shù)據(jù)的廣義線性混合模型參數(shù)估計(jì)構(gòu)造一個(gè)加權(quán)對(duì)數(shù)偽似然函數(shù)（log pseudo-likelihood function w ith weights），其公式如下：

這里，wj＝1/πj，wi｜j＝1/πi｜j，i和j分別代表1水平個(gè)體和2水平的層。一般情況下，復(fù)雜抽樣數(shù)據(jù)多水平模型采用牛頓-拉夫遜最大算法（New ton-Raphson maximum algorithm）使加權(quán)對(duì)數(shù)偽似然函數(shù)達(dá)到最大值，即PMLE（pseudo-maximum-likelihood estimation）參數(shù)估計(jì)值；其標(biāo)準(zhǔn)誤估計(jì)采用泰勒線性化三明治估計(jì)量（sandw ich estimator）。

復(fù)雜抽樣數(shù)據(jù)廣義線性混合模型用上式（2）做PMLE估計(jì)時(shí)，須考慮2水平權(quán)重；且不能直接使用1水平的個(gè)體權(quán)重，這樣可能影響參數(shù)估計(jì)的準(zhǔn)確性［11］，解決方法是對(duì)1水平個(gè)體權(quán)重進(jìn)行權(quán)重的縮放（scaling of weights）。Longford等［12］1995年提出的權(quán)重縮放方法1如下：

1998年P(guān)feffermann等［13］提出了權(quán)重縮放方法2如下：

應(yīng)用實(shí)例及結(jié)果解釋

我國(guó)某省在2010年進(jìn)行了行為危險(xiǎn)因素的調(diào)查，該省共有150個(gè)區(qū)縣，1900個(gè)鄉(xiāng)鎮(zhèn)，人口7000多萬(wàn)。根據(jù)計(jì)算本調(diào)查所需樣本量為15000人，使用了三階段整群抽樣設(shè)計(jì)，第一階段隨機(jī)抽取區(qū)縣，第二階段抽取鄉(xiāng)、鎮(zhèn)和街道，第三階段抽取個(gè)體。第一階段隨機(jī)抽取產(chǎn)生了50個(gè)區(qū)縣（primary sampling units，PSU）；第二階段，在每個(gè)區(qū)縣內(nèi)按鄉(xiāng)、鎮(zhèn)、街道分層，在每層內(nèi)按概率比例規(guī)模抽樣（PPS）方法，整群抽取2個(gè)鄉(xiāng)、1個(gè)鎮(zhèn)和2個(gè)街道，共產(chǎn)生了250個(gè)鄉(xiāng)鎮(zhèn)街道；在抽取的每個(gè)鄉(xiāng)、鎮(zhèn)和街道中，隨機(jī)抽取60名年齡45歲及以上的成年人，共15000個(gè)體。該調(diào)查將在過去的3個(gè)月內(nèi)發(fā)生1次或2次及以上跌倒性傷害為反應(yīng)變量，記為1，未發(fā)生任何跌倒性傷害為0，反應(yīng)變量為二分類變量，可用復(fù)雜抽樣logistic回歸分析該數(shù)據(jù)；但該數(shù)據(jù)顯然具有地區(qū)和個(gè)體的層次結(jié)構(gòu)，因此，合適的統(tǒng)計(jì)分析方法應(yīng)是考慮復(fù)雜抽樣2水平logistic回歸模型。這里個(gè)體為1水平，250個(gè)鄉(xiāng)、鎮(zhèn)、街道為2水平，最高水平區(qū)縣（PSU）在這里不作為水平考慮，但在分析中可作為層變量進(jìn)行調(diào)整處理。

本文分析了跌倒后造成的傷害與體重指數(shù)（BM I）及一些感興趣變量的關(guān)系，因變量為在過去三個(gè)月內(nèi)是否發(fā)生過跌倒性傷害，自變量有體重指數(shù)、性別、年齡、健康狀況、受教育程度、雇傭狀況和婚姻狀況等7個(gè)變量。該數(shù)據(jù)擬合復(fù)雜抽樣2水平隨機(jī)截距l(xiāng)ogistic回歸模型，分析了250個(gè)鄉(xiāng)、鎮(zhèn)和街道，由于有缺失數(shù)據(jù)，實(shí)際分析的樣本量為12086例。

調(diào)查采用了不等概率抽樣，需要根據(jù)抽樣設(shè)計(jì)對(duì)樣本進(jìn)行抽樣加權(quán)，加權(quán)過程涉及到個(gè)體抽樣權(quán)重，未包括不應(yīng)答權(quán)重和事后分層加權(quán)。三階段抽樣，用w代表權(quán)重，腳注1、2和3分別代表縣、鄉(xiāng)鎮(zhèn)街道和個(gè)體各階段的抽樣權(quán)重，用i表示某一樣本個(gè)體，s代表某一樣本個(gè)體所在的層。樣本區(qū)縣的抽樣權(quán)重（wsi1）值為分層簡(jiǎn)單隨機(jī)抽樣下樣本區(qū)縣抽樣概率的倒數(shù)，其計(jì)算公式如下：

樣本鄉(xiāng)鎮(zhèn)街道的抽樣權(quán)重（wsi2）值為與人口數(shù)成比例的PPS抽樣下樣本鄉(xiāng)鎮(zhèn)街道抽樣概率的倒數(shù)，用下式計(jì)算：

樣本個(gè)體的抽樣權(quán)重（wsi3）值為樣本個(gè)體抽樣概率的倒數(shù)?？捎?jì)算如下：

依據(jù)上述各階段抽樣權(quán)重，最終樣本個(gè)體的抽樣權(quán)重為下式（5）：

wsi＝wsi1×wsi2×wsi3＝樣本個(gè)體i所在區(qū)縣分層抽樣比倒數(shù)×

樣本個(gè)體i所在鄉(xiāng)鎮(zhèn)街道45歲及以上人口百分比構(gòu)成，如無(wú)法得到實(shí)際數(shù)據(jù)，可用該縣或省的構(gòu)成（Pop45）計(jì)算。

計(jì)算鄉(xiāng)鎮(zhèn)街道（2水平單位）πj，個(gè)體被抽中的概率πij及它們的權(quán)重，即概率的倒數(shù)，其計(jì)算公式為：

ws2＝wsi1×wsi2＝樣本個(gè)體i所在區(qū)縣分層抽樣比倒數(shù)×

Stata提供了分析復(fù)雜抽樣2水平隨機(jī)截距l(xiāng)ogistic回歸模型的分析模塊gllamm和說明書，其未加權(quán)最大似然估計(jì)（MLE），stata分析語(yǔ)句為：gllamm injury bmi age gender health marriage education employed，i（ID_level_2）link（logit）fam ily（binom）nip（12）adapt

穩(wěn)健標(biāo)準(zhǔn)誤估計(jì)（Robust standard errors）：gllamm，robust

其加權(quán)偽最大似然估計(jì)（PMLE）stata分析語(yǔ)句為：

gllamm injury bm i age gender health marriage education employed，i（ID_level_2）cluster（w t2）link（logit）family（binom）pweight（w t1_sw）nip（12）adapt

調(diào)整PSU后的穩(wěn)健標(biāo)準(zhǔn)誤估計(jì)：gllamm，robust cluster（psu）

其中，因變量為二分類injury，自變量有：體重指數(shù)BM I（啞變量）、健康狀況health（啞變量）、婚姻狀況marriage（啞變量）、受教育程度education（啞變量）和雇傭狀況employed（啞變量）。i指出2水平變量為ID_level_2；cluster定義2水平變量的權(quán)重為w t2；link指出連接函數(shù)為logit；family指出擬合二項(xiàng)分布（binom），若擬合Poisson分布，則family（Poisson）；pweight定義1水平個(gè)體權(quán)重為w t1_sw，其為經(jīng)縮放后的值。前面提到1水平權(quán)重不能直接應(yīng)用，要進(jìn)行權(quán)重的縮放，下面是應(yīng)用式（3）和（4）兩種權(quán)重縮放方法的stata語(yǔ)句。

方法1：

方法2：

應(yīng)用2水平隨機(jī)截距l(xiāng)ogistic回歸模型分析了鄉(xiāng)鎮(zhèn)街道（2水平）的個(gè)體（1水平）跌倒性傷害與感興趣的因素關(guān)系，用區(qū)縣變量（PSU）分層調(diào)整。分析時(shí)沒有將區(qū)縣作為一個(gè)水平來對(duì)待，主要考慮區(qū)縣的變異可能不大，但在標(biāo)準(zhǔn)誤估計(jì)時(shí)，調(diào)整了區(qū)縣（PSU）的聚集效應(yīng)。由于缺失數(shù)據(jù)的存在，最終的分析在50個(gè)PSU，250個(gè)鄉(xiāng)鎮(zhèn)（2水平）和12086個(gè)體（1水平）中進(jìn)行；最終結(jié)果報(bào)告了未加權(quán)2水平logistic回歸的參數(shù)估計(jì)-基于模型標(biāo)準(zhǔn)誤估計(jì)和穩(wěn)健標(biāo)準(zhǔn)誤估計(jì)，以及應(yīng)用權(quán)重縮放方法的PMLE參數(shù)估計(jì)，結(jié)果見表1。從表1可見，未加權(quán)2水平logistic回歸的參數(shù)估計(jì)，基于模型和穩(wěn)健標(biāo)準(zhǔn)誤估計(jì)結(jié)果很接近，對(duì)跌倒性傷害有統(tǒng)計(jì)學(xué)顯著影響的變量：健康狀況中等和差、未被雇傭和未婚，年齡為負(fù)相關(guān)，即年齡越大，發(fā)生跌倒性傷害的危險(xiǎn)性越??；復(fù)雜抽樣2水平logistic回歸參數(shù)PMLE估計(jì)，權(quán)重縮放方法1和2估計(jì)結(jié)果非常接近；與未加權(quán)的方法比，對(duì)跌倒性傷害有統(tǒng)計(jì)學(xué)顯著影響的變量與未加權(quán)的結(jié)果基本一致，但未婚失去了統(tǒng)計(jì)顯著性，兩種縮放方法估計(jì)的P值在0.05～0.10之間。其它分析變量，如體重指數(shù)、性別和受教育程度與跌倒性傷害的發(fā)生沒有統(tǒng)計(jì)學(xué)聯(lián)系。加權(quán)復(fù)雜抽樣PMLE估計(jì)，大多數(shù)參數(shù)估計(jì)值比未加權(quán)估計(jì)增加，但其標(biāo)準(zhǔn)誤比未加權(quán)估計(jì)值增加明顯，導(dǎo)致大部分變量P值增大，估計(jì)保守。

表1 2水平隨機(jī)截距l(xiāng)ogistic回歸和復(fù)雜抽樣2水平隨機(jī)截距l(xiāng)ogistic回歸擬合

討 論

流行病學(xué)抽樣調(diào)查一般應(yīng)用多階段抽樣，第一階段先抽取地區(qū)或群（clusters），第二階段抽取群下面的亞群（subclusters）；最后階段抽取基本抽樣單位，如個(gè)體。這種抽樣方法產(chǎn)生了多水平數(shù)據(jù)，基本抽樣單位為1水平或最低水平，最高水平為PSU。對(duì)這類資料的分析不但考慮復(fù)雜抽樣的權(quán)重（包括抽樣權(quán)重、不應(yīng)答權(quán)重和事后分層權(quán)重），而且考慮數(shù)據(jù)的多水平結(jié)構(gòu)，即應(yīng)用復(fù)雜抽樣數(shù)據(jù)多水平模型分析方法。本文應(yīng)用我國(guó)某省行為危險(xiǎn)因素調(diào)查的數(shù)據(jù)，應(yīng)用Rabe-Hesketh描述的考慮權(quán)重的全PMLE參數(shù)估計(jì)方法，通過適應(yīng)積分法估計(jì)任何水平模型的參數(shù)和標(biāo)準(zhǔn)誤。Stata軟件gllamm程序分析發(fā)現(xiàn)：健康狀況中等和差、未被雇傭和未婚與跌倒性傷害的關(guān)系與以前的文獻(xiàn)或研究一致［14］；但一般認(rèn)為年齡大發(fā)生跌倒性傷害的可能性大［15］，但本結(jié)果與之相反，這也許與本研究選取的個(gè)體年齡大有關(guān)。文獻(xiàn)認(rèn)為男性較女性更易發(fā)生跌倒性傷害［16］，本研究未加權(quán)的分析結(jié)果是男性發(fā)生跌倒性傷害的危險(xiǎn)性小；加權(quán)的2水平隨機(jī)截距l(xiāng)ogistic回歸模型顯示男性與跌倒性傷害是正相關(guān)；雖然加權(quán)和未加權(quán)的結(jié)果均沒有統(tǒng)計(jì)顯著性，但加權(quán)后的結(jié)果合理，符合目前對(duì)跌倒性傷害發(fā)生機(jī)制的認(rèn)識(shí)。這也部分說明了復(fù)雜抽樣的層次結(jié)構(gòu)數(shù)據(jù)，分析時(shí)納入權(quán)重的多水平分析方法的合理性。

復(fù)雜抽樣數(shù)據(jù)的統(tǒng)計(jì)分析也有百年的歷史［17］，1990年代以后更有新統(tǒng)計(jì)理論和方法發(fā)展，復(fù)雜抽樣的統(tǒng)計(jì)方法已包括：列聯(lián)表資料的對(duì)數(shù)線性模型和相關(guān)的方法、廣義線性模型（logistic回歸和Poisson回歸）、生存分析、一般線性混合模型、結(jié)構(gòu)方程模型、隱變量模型，這些方法充分考慮了抽樣權(quán)重，對(duì)參數(shù)標(biāo)準(zhǔn)誤和可信區(qū)間的估計(jì)更準(zhǔn)確。目前SAS、Stata、SPSS和SUDAN統(tǒng)計(jì)軟件包都包含了復(fù)雜抽樣數(shù)據(jù)的統(tǒng)計(jì)描述、一般線性回歸、logistic回歸（二分類、有序和無(wú)序）、Poisson回歸（零膨脹、負(fù)二項(xiàng)）、生存分析等，可以滿足大部分復(fù)雜抽樣數(shù)據(jù)的統(tǒng)計(jì)分析。

但是目前流行的統(tǒng)計(jì)軟件還都沒有包括復(fù)雜抽樣多水平模型擬合程序。當(dāng)抽樣權(quán)重納入模型時(shí)，一般構(gòu)造偽似然函數(shù)（pseudo-maximum-likelihood）來估計(jì)模型參數(shù)。在實(shí)踐中的一個(gè)主要問題是大多數(shù)流行病學(xué)抽樣調(diào)查僅給出了基本抽樣單位或1水平的抽樣權(quán)重，而沒有高水平的權(quán)重；但在進(jìn)行復(fù)雜抽樣數(shù)據(jù)多水平模型擬合時(shí)，需要高水平的抽樣權(quán)重。因此，在今后的流行病學(xué)調(diào)查設(shè)計(jì)中和實(shí)施時(shí)，應(yīng)注意收集計(jì)算各水平權(quán)重指標(biāo)數(shù)據(jù)，以便在數(shù)據(jù)分析時(shí)做加權(quán)處理。

1.Rao JNK.Interplay between sample survey theory and practice；anappraisal.Survey Methodology，2005，31：117-138.

2.呂筠，何平平，李立明.復(fù)雜抽樣調(diào)查數(shù)據(jù)實(shí)例分析.中華流行病學(xué)雜志，2008，29（8）：832.

3.繆凡，童峰.復(fù)雜抽樣數(shù)據(jù)的logistic回歸分析方法及其應(yīng)用.中國(guó)衛(wèi)生統(tǒng)計(jì)，2008，25（6）：577-579.

4.胡躍華，匡翔宇，金承剛，等.復(fù)雜抽樣Poisson回歸分析方法及應(yīng)用.中國(guó)衛(wèi)生統(tǒng)計(jì)，2012，29（5）：650-653.

5.楊珉，李曉松主編.醫(yī)學(xué)和公共衛(wèi)生研究常用多水平統(tǒng)計(jì)模型.北京，北京大學(xué)醫(yī)學(xué)出版社，2007.

6.賈改珍，閆陽(yáng)，徐天和，等.多水平模型在大學(xué)生預(yù)防艾滋病健康教育影響因素分析中的應(yīng)用.中國(guó)衛(wèi)生統(tǒng)計(jì)，2013，30（1）：37-39.

7.Gebremariam MK，Andersen LF，Bielland M，et al.Does the school food environment influence the dietary behaviours of Norwegian 11-yearolds.The HEIA study.Scand JPublic Health，2012，40（5）：491-497.

8.Nansel TR，Lipsky LM，Lannotti RJ.Cross-sectional and longitudinal relationships of bodymass index with glycemic control in children andadolescentswith type 1 diabetesmellitus，2013，100（1）：126-132.

9.Steven G.Heeringa，Wagner J，Torres M，et al.Sample designs and samplingmethods for the Collaborative Psychiatric Epidem iology Studies（CPES）.Int.J.Methods Psychiatr.Res.，2004，13（4）：221-240.

10.Sophia Rabe-Hesketh.Multilevelmodeling of complex survey data.J. R.Statist.Soc.，2006，169：805-827.

11.Moshe Feder，Gad Nathan，Danny Pfeffermann.Survey Methodology，2000，26（1）：53-65.

12.Longford NT.Model-basedmethods for analysis of data from 1990 NAEP Trial State Assessment.Research and Development Report NCES 95-696.Washington DC：National Center for Education Statistics.

13.Pfeffermann D，Holmes CJ，Goldstein DJ，et al.Weighting for unequal selection probabilities in multilevelmodels.J.R.Statist.Soc.B，1998，60：23-40.

14.Roe B，Howell F，Riniotis K，et al.Older people and falls：health status，quality of life，lifestyle，care networks，prevention and views on service use follow ing a recent fall.JClin Nurs，2009，18：2261-2272.

15.Hausdorff JM，Rios DA，Edelberg HK.Gait variability and fall risk in community-living older adults：A 1-year prospective study.A rch Phys Med Rehabil，2001，82：1050-1056.

16.Stevens JA，Sogolow ED.Gender differences for non-fatal unintentional fall related injuries among older adults.Injury Prevention，2005，11：115-119.

17.Fisher RA.StatisticalMethods for Research Work.Oliver and Boyd，Edinburgh，1925.

（責(zé)任編輯：劉 壯）

App lication of M ultilevel M odeling to Com plex Sam ple Survey Data

Yu Shicheng，Liao Jiaqiang，Yu Mo，et al（Chinese Center for Disease Control and Prevention（102206），Beijing）

ObjectiveTo illustrate the principal and application ofmultilevelmodeling of complex survey data thatwere derived from multistage sampling，unequal sampling probabilities and different features of post-stratification.MethodsWeighted and un-weighted random intercept logistic regressionmodelswere applied to complex survey data of behavioral risk factors in a province to look at the association of fall injuries w ith some factors of interest.ResultsThere were 12086 subjects（level 1）aged 45 years or above nested w ithin 250 villages，towns and sub-districts（level 2）from 50 counties/districts（PSU）.Un-weighted results showed that variables significantly and positively associated w ith the risk of fall injurieswere fair or poor health，unemployed situation，unmarried；age was significantly and negatively associated w ith the risk of fall injuries，or one less likely got injured when getting older.The results from 2-level random intercept logisticmodel demonstrated that the variables associated w ith the risk of fall injuries were sim ilar to those from un-weighted models，but the variable of unmarried m itigated its significance to be insignificant.Body mass index，beingmale，educational levelwere notassociated w ith the risk of fall injuries from the analyses.ConclusionIn contrast to the results from un-weighted methods，statistically significant variables from weightedmethodswere analogous to those from weighted ones；however，estimates using full pseudo-maximum-likelihood estimation（PMLE）weremore conservative as opposed to un-weighted ones.As for gender，weighted result was in consistent w ith the currentunderstanding of themechanism for the developmentof fall injuries，therefore，it soundedmore reasonable to employ multilevel modeling for the complex survey data.

Complex survey data；Multilevelmodel；Multistage sampling；Random intercept logistic regression

*：淮河流域癌癥綜合防治項(xiàng)目（1310800003）

1.中國(guó)疾病預(yù)防控制中心公共衛(wèi)生監(jiān)測(cè)與信息服務(wù)中心（102206）

2.四川大學(xué)華西公共衛(wèi)生學(xué)院衛(wèi)生統(tǒng)計(jì)教研室（610041）

3.北京師范大學(xué)社會(huì)發(fā)展與公共政策學(xué)院（100875）

△通信作者：肖革新，E-mail：biocomputer＠126.com