全文鏈接：http://tecdat.cn/?p=10809

本文用于比較六個不同統(tǒng)計軟件程序（SAS，Stata，HLM，R，SPSS和Mplus）的兩級分層線性模型的過程和輸出

下面介紹的六個模型都是兩級分層模型的變體，也稱為多級模型，這是混合模型的特殊情況。此比較僅對完全嵌套的數(shù)據(jù)有效（不適用于交叉或其他設(shè)計的數(shù)據(jù)，可以使用混合模型進行分析）。盡管HLM軟件的網(wǎng)站聲明可以用于交叉設(shè)計，但這尚未得到確認。下面的SAS，Stata，R，SPSS和Mplus中使用的過程是其多層次或混合模型過程的一部分，并且可以擴展為非嵌套數(shù)據(jù)。

但是出于比較的目的，我們將僅研究完全嵌套的數(shù)據(jù)集。除了HLM（完全由GUI運行）以外，所有程序的下面都包含用于每個模型的代碼/語法。我們提供了HLM和SPSS的屏幕截圖。此外，每個模型均以分層格式和混合格式指定。盡管模型的這兩個表達式是等效的，但一些研究領(lǐng)域更傾向于可視化層次結(jié)構(gòu)，因為它更容易看到層次之間的分離，而另一些研究領(lǐng)域則更喜歡混合格式，在其中容易區(qū)分固定效果和隨機效果。

模型注意事項將預(yù)測變量添加到本文檔討論的六個模型中時，我們選擇以均值居中為中心，這意味著我們從每個受試者的得分中減去了該變量的總體均值。

正如Enders＆Tofighi（2007）所詳細討論的那樣，以總體平均值為中心，而不是以組平均值（每個組的平均值均以該組中受試者的得分為準(zhǔn)）為中心，并不適合所有模型。。使用哪種居中方法的選擇應(yīng)由所詢問的具體研究問題決定。另一個考慮因素是這些程序使用的估計方法來產(chǎn)生參數(shù)估計，即最大似然（ML）或受限最大似然（REML）。每種都有自己的優(yōu)點和缺點。ML更適合不平衡的數(shù)據(jù)，但是會產(chǎn)生偏差的結(jié)果。REML是無偏的，但是在將兩個嵌套模型與似然比檢驗進行比較時，不能使用REML。

兩種方法將產(chǎn)生相同的固定效應(yīng)估計，但它們對隨機效應(yīng)的估計卻有所不同（Albright＆Marinova，2010）。正如我們將在下面討論的模型中看到的那樣，這兩種方法產(chǎn)生的結(jié)果非常相似，并且不會極大地影響隨機因素的p值。但是，重要的是要意識到，方法的選擇會影響隨機因素的估計，標(biāo)準(zhǔn)誤差和p值，并且可能會影響宣布隨機因素是否重要的決策。SAS，HLM，R和SPSS默認使用REML，而Stata和Mplus使用ML。在本文檔中的Stata示例中，我們告訴Stata使用REML以便將輸出與其他四個程序進行比較。

類內(nèi)相關(guān)系數(shù)

我們還報告了每種模型的類內(nèi)相關(guān)系數(shù)（ICC）ρ。ICC是結(jié)果變量中方差的比例，由分層模型的分組結(jié)構(gòu)解釋。它是根據(jù)組級別誤差方差與總誤差方差之比來計算的：

其中，是2級殘差的方差，是1級殘差的方差。換句話說，與總的無法解釋的方差（方差之內(nèi)和之間）相比，ICC報告了模型中任何可歸因于分組變量的預(yù)測變量無法解釋的變化量。

示例數(shù)據(jù)集

流行的數(shù)據(jù)集由來自不同班級的學(xué)生組成，并且由于每個學(xué)生都屬于一個唯一的班級，因此它是一個嵌套設(shè)計。因變量是“流行”，它是一個自評的流行度，范圍為0-10。預(yù)測指標(biāo)包括學(xué)生級別的性別（二分法）和Extrav（連續(xù)的自我評價的外向得分），以及班級的Texp（多年的老師經(jīng)驗， 是連續(xù)的）。

僅截距模型（無條件模型）

無條件混合模型規(guī)范類似于單因素方差分析，其總體均值和類效應(yīng)。但是，我們將其視為隨機效應(yīng)（均值為零的正態(tài)分布變量），而不是像方差分析中那樣的固定因子效應(yīng)。因此，我們將估計值解釋為每個類別的平均數(shù)在總體平均人氣得分附近的方差。

估算值是每個班級的“大眾”平均值的平均值，而不是研究中所有學(xué)生的平均值。如果數(shù)據(jù)完全平衡（即每個班級的學(xué)生人數(shù)相同），則無條件模型的結(jié)果將與方差分析程序的結(jié)果相同。

SAS結(jié)果

?需要“ covtest”選項來報告方差分量估計的標(biāo)準(zhǔn)誤差。另外， 需要指定非結(jié)構(gòu)化協(xié)方差矩陣類型，這是HLM和R默認情況下使用的類型，我們在這里使用它進行比較。SAS的輸出等于Hox的書表2.1中的結(jié)果。我們可以得出結(jié)論，各類別之間的平均人氣得分為5.078，并且各類別之間的差異（1.221）比不同類別之間的差異（0.702）多。當(dāng)我們?yōu)樵撃Ｐ陀嬎鉏CC時，將對此進行進一步討論。

Stata結(jié)果

?Stata的xtmixed命令需要因變量，后跟“ ||”? 指定固定變量和隨機變量之間的分隔。我們必須包括方差選項以查看輸出中方差分量的估計值，以及reml選項以使用受限的最大似然估計。還要注意，Stata不會輸出隨機分量估計的p值，但是可以通過置信區(qū)間中是否包含零來確定有效值。這些結(jié)果與SAS的結(jié)果完全匹配

點擊標(biāo)題查閱往期內(nèi)容

用SPSS估計HLM多層（層次）線性模型模型

左右滑動查看更多

01

02

03

04

HLM結(jié)果

HLM報告方差組件的標(biāo)準(zhǔn)偏差，而不是標(biāo)準(zhǔn)誤差。同樣，對于隨機效應(yīng)，他 僅報告截距的卡方統(tǒng)計量和p值。這些結(jié)果與其他程序的結(jié)果相同。

R結(jié)果

R報告方差成分（例如HLM）的標(biāo)準(zhǔn)偏差，而lme4軟件包報告固定效應(yīng)的t統(tǒng)計量。?

SPSS結(jié)果

?屏幕截圖：

?

?需要在“隨機”窗口中指定非結(jié)構(gòu)化協(xié)方差類型。這些結(jié)果與其他程序和本文得出的結(jié)果相同。請注意，像SAS和Mplus一樣，SPSS報告方差分量的標(biāo)準(zhǔn)誤差，而HLM和R報告標(biāo)準(zhǔn)差。我們無法得出結(jié)論，哪個更適合報告，但是差異不會影響這些參數(shù)的p值。

?因為這是一個無條件模型，所以我們不需要指定任何WITHIN或BETWEEN變量。下面列出了在MODEL語句中列出變量的標(biāo)準(zhǔn)。在以下各節(jié)中，我們將看到前三個示例：

1.％WITHIN％– 1級固定因子（非隨機斜率）2.具有潛在斜率變量的％WITHIN％– 1級隨機因子3.％BETWEEN％– 2級固定因子4.在任一個陳述–在學(xué)生水平上測得的變量，但具有1級和2級方差估計 。

上表顯示了Mplus輸出底部的“模型結(jié)果”部分的結(jié)果。Mplus確實會報告每個估計的p值，并且所有估計都與其他程序的p值匹配，但隨機截距的方差估計相差約0.007。這種差異是由于Mplus使用ML估計這一事實造成的。盡管存在這種差異，但我們看不到任何變量的重要性發(fā)生變化。

匯總

總體而言，這六個程序?qū)τ趦H截取模型產(chǎn)生了非常相似的結(jié)果（唯一的差異發(fā)生在隨機效應(yīng)的Mplus估計中）。唯一的區(qū)別是他們?nèi)绾螆蟾骐S機方差估計的精度。此模型的ICC等于：

這告訴我們，“流行”課程總變化的大約三分之一可以由每個學(xué)生所在的班級解釋。

固定Level-1因子的隨機截距（非隨機斜率）模型

該模型增加了一個學(xué)生級別的固定因子Extrav，即自我報告的外向得分?；旌夏Ｐ涂雌饋硐袷腔趲в袇f(xié)變量Extrav的類的ANCOVA，但請記住，我們?nèi)匀徽J為這是隨機效應(yīng)，而不是固定效應(yīng)。因此， 估算值與ANCOVA程序所得出的估算值不同。?

?在此數(shù)據(jù)的實際應(yīng)用中，Extrav應(yīng)該具有固定的效果而不是隨機的效果是沒有意義的，因為學(xué)生外向性的水平應(yīng)隨班級而變化。但是，出于比較這四個程序的目的，我們?nèi)匀幌Ｍ{(diào)查一個具有一個學(xué)生級別固定因子的案例。

SAS結(jié)果

現(xiàn)在，我們對Extrav的固定效果進行了估算。學(xué)生報告的外向得分每增加一個單位，他們的受歡迎度得分就會增加0.486。這些結(jié)果等于使用REML的其他程序的結(jié)果。

Stata結(jié)果

當(dāng)我們向Stata中的模型添加預(yù)測變量時，我們添加了cov（un）選項，指定了非結(jié)構(gòu)化協(xié)方差矩陣。我們將 Extraversion變量放在“ ||”之前 表示它是一個固定因子（具有非隨機斜率）。這些結(jié)果與其他程序的結(jié)果相同。

?HLM結(jié)果

?

R結(jié)果

SPSS結(jié)果

Mplus結(jié)果

?現(xiàn)在，我們在VARIABLE語句的WITHIN選項中包括居中的Extrav變量。對于內(nèi)部MODEL規(guī)范，我們必須使用“ ON”選項，以告知Mplus Extrav是固定的1級因子。? 可以看到由于使用ML估計而不是REML，許多估計和估計的標(biāo)準(zhǔn)誤差（以及t統(tǒng)計量）存在細微差異。由于方差的估計值與其他程序不同，因此Mplus報告的ICC與下面報告的有所不同。

** **

匯總

對于此模型，前五個程序的結(jié)果完全相同，而Mplus的估算值相差很小。此模型的ICC大于無條件模型的ICC（正如預(yù)期的那樣，因為我們通過添加固定因子來控制某些學(xué)生水平的變化）：

?

使用一個學(xué)生水平的固定因子，“流行”總變化的幾乎一半可以由該學(xué)生的班級和學(xué)生水平的固定因子“外向”解釋

一級因子的隨機截距和斜率模型

?

該模型包含Extrav的隨機斜率，這意味著我們允許回歸方程的斜率隨類而變化。該模型比以前的模型更適合于所使用的變量，因為可以直觀地假設(shè)外向因類而異。

SAS結(jié)果

隨機Extrav斜率的估計值很重要（p值為0.003），因此我們可以說學(xué)生的外向性得分隨班級的不同而變化。這些結(jié)果與其他程序的結(jié)果完全匹配，除了固定效果的t統(tǒng)計量存在一些細微差異。

Stata結(jié)果

HLM結(jié)果

SPSS結(jié)果

Mplus結(jié)果

?這次我們在WITHIN語句中包括一個潛在的斜率變量，以將Extrav指定為隨機因子，該變量告訴Mplus不要在數(shù)據(jù)集中尋找“ randoms1”，因為沒有觀察到它。可以將此變量的輸出解釋為Extrav的隨機斜率分量。我們必須這樣做，因為Mplus是為結(jié)構(gòu)方程模型設(shè)計的，其多級模型功能是其潛在潛伏分析程序的改編。

匯總

總體而言，前五個程序?qū)υ撃Ｐ彤a(chǎn)生相同的結(jié)果，而Mplus再次由于ML估計而相差很小。此模型的ICC為：

通過將Extrav的影響從固定變?yōu)殡S機，由于我們正在考慮在學(xué)生級別上更多隨機變化，因此ICC會略有增加。

兩個1級因子的隨機斜率模型

?

?對于此模型，我們包括第二個學(xué)生級別的變量Sex，該變量也具有隨機斜率。這意味著我們既要考慮學(xué)生的性別，又要考慮他們的外向得分，并且允許這兩個因素的斜率隨班級而變化。

SAS結(jié)果

?

?在此輸出中，我們可以看到性別確實對學(xué)生自我報告的知名度有重大影響（p值<0.001）。對Sex的固定估計意味著，在Extrav不變的情況下，女學(xué)生（Sex = 1）的普遍得分比男學(xué)生（基線組，Sex = 0）高1.244。

SAS不喜歡在該模型中，Sex的估計方差非常接近零，因此沒有輸出標(biāo)準(zhǔn)誤差或p值。因為非常接近于零，所以我們可以得出結(jié)論，性別不會因類別而顯著變化。?

Stata結(jié)果

Stata在運行該模型時引用了一個錯誤：標(biāo)準(zhǔn)誤差計算失敗，這意味著未計算隨機效應(yīng)的標(biāo)準(zhǔn)誤差。我們發(fā)現(xiàn)通過刪除cov（un）選項，不會出現(xiàn)此錯誤。但是，該輸出中的所有估計均與其他程序不同，因此我們選擇使用非結(jié)構(gòu)化協(xié)方差矩陣規(guī)范報告輸出。我們不確定這是否是在Stata中運行此類模型的常見問題，但重要的是要意識到它會發(fā)生。

HLM結(jié)果

這些估計大致等于其他計劃的結(jié)果，但隨機性別影響的估計除外。由于這種影響非常接近于零，因此程序不會報告完全相同的值，但是所有結(jié)果都表明該值遠非重要。

R結(jié)果

SPSS結(jié)果

Mplus結(jié)果

這次，我們在WITHIN語句中包括了兩個潛在的斜率變量，以將Extrav和Sex指定為隨機因子。我們可以將“ randoms1”的輸出解釋為Extrav的估計，將“ randoms2”的輸出解釋為Sex的估計。

Mplus針對此模型的輸出所得出的估計值與先前模型中的其他程序相距甚遠。我們看到，由于模型必須估計更多隨機參數(shù)，因此估計程序（ML與REML）之間的差異變得更加明顯。但是，Mplus同意其他程序的觀點，即“性別”的隨機方差部分以外的所有估計值都非常重要。

匯總

對于方差非常接近零的隨機效應(yīng)，六個程序以不同的方式處理估計值。SAS和Stata無法報告隨機效應(yīng)的標(biāo)準(zhǔn)誤差或p值，而其他變量的估計值和標(biāo)準(zhǔn)誤差均具有相當(dāng)大的差異。Mplus結(jié)果也顯示出比以前的模型更大的差異。此模型的ICC為：

同樣，當(dāng)我們在模型中添加另一個學(xué)生級別的效果（包括隨機斜率）時，ICC略有增加。?

2級因子和兩個隨機1級因子（無交互）

?這是我們看到的第一個具有2級（班級）變量的模型：教師的多年經(jīng)驗（Texp），也是以均值為中心的。正如Enders和Tofighi（2007）指出的那樣，級別2變量的唯一居中選項是均值居中。? 無法對均值中心Texp進行分組，因為它已經(jīng)在班級水平上進行了度量，這意味著“分組均值”將等于原始值。

在分層格式中， 可以看到它具有固定的斜率系數(shù)，并且對于每個類j都是唯一的。該模型在教師的經(jīng)驗和學(xué)生水平的變量之間沒有任何相互作用。如果我們有理由相信Texp不會緩和Sex和Extrav對Popular的影響，那么我們將使用此模型，這意味著我們的學(xué)生水平變量的斜率是相同的，無論學(xué)生是否有新教師或新教師。一位擁有多年經(jīng)驗的人。

SAS結(jié)果

?現(xiàn)在，我們在固定效果表中看到了Texp，估計值為0.089，p值很大。這意味著，在使學(xué)生的性別和性取向得分保持不變的情況下，每增加一年的教師經(jīng)驗，該學(xué)生的熱門得分就會增加0.089。

同樣，我們看到SAS無法處理隨機性別效應(yīng)的很小變化。因此，沒有報告標(biāo)準(zhǔn)誤差，z統(tǒng)計量或p值。

Stata結(jié)果

**
**

?與以前的模型一樣，我們收到一個錯誤，告訴我們Stata無法計算方差分量的標(biāo)準(zhǔn)誤差。但是，這些估計值與其他程序的估計值大致相同。

HLM結(jié)果

?

?

這些估計值與其他程序的結(jié)果略有不同。

** **

R結(jié)果

SPSS結(jié)果

?

Mplus結(jié)果

對于Level-2因子，我們在BETWEEN語句中包括Ctexp。我們再次看到這些估計數(shù)與其他五個計劃的輸出有微小出入

匯總

使用REML的五個程序的輸出實質(zhì)上是相等的，僅相差幾千個單位。與以前的模型一樣，最大的差異出現(xiàn)在隨機性別效應(yīng)的方差估計中，因為它非常接近零。

請注意，此模型的ICC比以前的模型有所降低（= 0.542）：

請記住，ICC是衡量 所在的班級可以解釋多少無法解釋的變化的方法。通過添加班級級別的預(yù)測變量，我們可以解釋不同班級中較大比例的變化。因此，與沒有任何2級預(yù)測變量的模型相比，該模型的隨機截距存在較少的變異，因此ICC也較低。

具有相互作用的一個2級因子和兩個隨機1級因子

?

?這是我們在班級變量Texp與學(xué)生級變量Sex和Extrav之間進行跨級交互的唯一模型。例如，如果我們想找出具有更多經(jīng)驗的教師是否比新教師對學(xué)生的外向性或性別與他們自我報告的知名度之間的關(guān)系有不同的影響，則可以使用此模型。換句話說，教師的經(jīng)歷是否適度了性格外向或性別對受歡迎程度的影響？

?可以看到，在分層格式中，Texp在三個方程式的每個方程式內(nèi)都有一個斜率系數(shù)。這與混合模型中的交互項有關(guān)，即通過外向的教師體驗和按性別的教師體驗。

SAS結(jié)果

?在固定效果表中，有兩個交互作用項，其中一個（）遠不重要，p值> 0.5。?

在隨機方差分量表中，我們看到外向隨機斜率的估計值和性別隨機斜率的估計值與零沒有顯著差異。這意味著沒有證據(jù)表明這兩個因素實際上在該模型中因類別而異。?

Stata結(jié)果

Stata無法自動識別變量之間的交互項，因此我們必須為兩個跨級別的交互手動創(chuàng)建變量（請參見上面的代碼中的gen語句）。當(dāng)我們使用帶有非結(jié)構(gòu)化協(xié)方差矩陣選項的xtmixed命令運行時，Stata給出了一個錯誤，指出Hessian不是負半定性，一致性錯誤，并且沒有產(chǎn)生任何輸出。?

HLM結(jié)果

這些估計值大致等于其他程序的結(jié)果。

R結(jié)果

SPSS結(jié)果

?對于SPSS 19而言，此模型實在太多了。對于具有非結(jié)構(gòu)化協(xié)方差矩陣的更復(fù)雜的模型，其他程序可能會運行更有效的算法，因此優(yōu)于SPSS。

****

Mplus結(jié)果

?現(xiàn)在，我們在BETWEEN模型部分中包含兩個ON語句，以指示與教師經(jīng)驗的跨層次交互。同樣，我們發(fā)現(xiàn)與其他輸出之間的細微差異，但Mplus同意Texp和Sex之間的固定交互作用不顯著，而Extrav和Sex的隨機組成也并不重要。

** **

匯總

加上兩個跨層交互項，Stata和SPSS無法使用非結(jié)構(gòu)化協(xié)方差選項運行模型。這并不是說不應(yīng)該將它們用于這種類型的分析，但是在向具有非結(jié)構(gòu)化協(xié)方差矩陣的模型中添加更復(fù)雜的參數(shù)時，應(yīng)謹慎使用。

與以前的模型一樣，SAS，HLM和R的結(jié)果相對接近相等，而Mplus的估計略有不同。另外，ICC與模型5幾乎完全相同，這意味著交互作用項不會改變按類別說明的差異比例：

總結(jié)

進行比較的目的是調(diào)查來自六個不同統(tǒng)計軟件程序的嵌套兩級層次模型的過程和結(jié)果可能存在的差異?？傮w而言，我們發(fā)現(xiàn)SAS，Stata（帶有reml選項），HLM，R和SPSS產(chǎn)生的實際估計值之間沒有太大差異。Mplus使用另一種估算方法ML，這導(dǎo)致其估算值與其他估算值有所不同。另外，重要的是要注意以下幾點：

1.對于方差估計非常接近零的隨機效應(yīng)，SAS無法產(chǎn)生標(biāo)準(zhǔn)誤差或p值。其他三個程序在估計這些參數(shù)方面的差異與其他效果相比更大。

1. Stata和SPSS無法處理最復(fù)雜的模型，該模型包含兩個跨級別的交互項。建議使用其他程序來分析復(fù)雜模型并指定非結(jié)構(gòu)化協(xié)方差矩陣。

此外，我們研究了每種模型中類內(nèi)相關(guān)系數(shù)的值。通過添加1級預(yù)測因子，ICC有所增加。但是，當(dāng)我們添加2級預(yù)測變量時，ICC會大大降低，甚至比無條件模型更低。這是由于在類級別添加了預(yù)測變量時，無法解釋的Level-2變異（隨機截距項）減少了。

盡管本文檔可以用作為嵌套數(shù)據(jù)集運行各種兩級分層模型的指南，但我們強烈建議讀者僅在適合回答您的特定研究問題時使用這些模型。在確定固定因素和隨機因素之間，以及對于中心平均值為1的總體平均值或組平均值時，必須謹慎使用。?

** **

參考文獻

Enders, Craig K. and Tofighi, Davood (2007). “Centering Predictor Variables in Cross-Sectional Multilevel Models: A New Look at an Old Issue.” Psychological Methods, vol. 12, pg. 121-138.

Hox, Joop J. (2010). Multilevel Analysis (2nd ed.). New York: Routledge

本文摘選?《?使用SAS，Stata，HLM，R，SPSS和Mplus的分層線性模型HLM?》?，點擊“閱讀原文”獲取全文完整代碼和數(shù)據(jù)資料。

點擊標(biāo)題查閱往期內(nèi)容

R語言用WinBUGS 軟件對學(xué)術(shù)能力測驗建立層次（分層）貝葉斯模型
WinBUGS對多元隨機波動率模型：貝葉斯估計與模型比較
R語言BUGS/JAGS貝葉斯分析: 馬爾科夫鏈蒙特卡洛方法（MCMC）采樣
WINBUGS對隨機波動率模型進行貝葉斯估計與比較
R語言通過WinBUGS對MGARCH和MSV模型進行貝葉斯估計和比較
R語言用WinBUGS 軟件對學(xué)術(shù)能力測驗（SAT）建立分層模型
R語言中RStan貝葉斯層次模型分析示例
用SPSS估計HLM多層（層次）線性模型模型
R語言使用貝葉斯層次模型進行空間數(shù)據(jù)分析
R語言基于copula的貝葉斯分層混合模型的診斷準(zhǔn)確性研究
R語言分層線性模型案例
R語言用WinBUGS 軟件對學(xué)術(shù)能力測驗（SAT）建立分層模型
使用SAS，Stata，HLM，R，SPSS和Mplus的分層線性模型HLM