Mixed-Model

當組在混合模型中被視為隨機與固定時,斜率估計存在很大差異

  • April 20, 2013

我了解,當我們認為某些模型參數在某些分組因素中隨機變化時,我們會使用隨機效應(或混合效應)模型。我希望擬合一個模型,其中響應已在分組因子中標準化和居中(不完美,但非常接近),但x沒有以任何方式調整自變量。這導致我進行了以下測試(使用捏造的數據),以確保我找到我正在尋找的效果,如果它確實存在的話。我運行了一個具有隨機截距(跨由 定義的組)的混合f效應模型和第二個以因子 f 作為固定效應預測因子的*固定效應模型。*我將 R 包lmer用於混合效果模型,以及基本功能lm()為固定效應模型。以下是數據和結果。

請注意y,無論組如何,都在 0 左右變化。並且在組內x始終y變化,但在組之間的變化比y

> data
     y   x f
1  -0.5   2 1
2   0.0   3 1
3   0.5   4 1
4  -0.6  -4 2
5   0.0  -3 2
6   0.6  -2 2
7  -0.2  13 3
8   0.1  14 3
9   0.4  15 3
10 -0.5 -15 4
11 -0.1 -14 4
12  0.4 -13 4

如果您對處理數據感興趣,這裡是dput()輸出:

data<-structure(list(y = c(-0.5, 0, 0.5, -0.6, 0, 0.6, -0.2, 0.1, 0.4, 
-0.5, -0.1, 0.4), x = c(2, 3, 4, -4, -3, -2, 13, 14, 15, -15, 
-14, -13), f = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 
4L, 4L, 4L), .Label = c("1", "2", "3", "4"), class = "factor")), 
.Names = c("y","x","f"), row.names = c(NA, -12L), class = "data.frame")

擬合混合效應模型:

> summary(lmer(y~ x + (1|f),data=data))
Linear mixed model fit by REML 
Formula: y ~ x + (1 | f) 
  Data: data 
  AIC   BIC logLik deviance REMLdev
28.59 30.53  -10.3       11   20.59
Random effects:
Groups   Name        Variance Std.Dev.
f        (Intercept) 0.00000  0.00000 
Residual             0.17567  0.41913 
Number of obs: 12, groups: f, 4

Fixed effects:
           Estimate Std. Error t value
(Intercept) 0.008333   0.120992   0.069
x           0.008643   0.011912   0.726

Correlation of Fixed Effects:
 (Intr)
x 0.000 

我注意到截距方差分量估計為 0,對我來說重要的是,x它不是y.

接下來,我將固定效應模型擬合f為預測變量,而不是隨機截距的分組因子:

> summary(lm(y~ x + f,data=data))

Call:
lm(formula = y ~ x + f, data = data)

Residuals:
    Min       1Q   Median       3Q      Max 
-0.16250 -0.03438  0.00000  0.03125  0.16250 

Coefficients:
           Estimate Std. Error t value Pr(>|t|)    
(Intercept) -1.38750    0.14099  -9.841 2.38e-05 ***
x            0.46250    0.04128  11.205 1.01e-05 ***
f2           2.77500    0.26538  10.457 1.59e-05 ***
f3          -4.98750    0.46396 -10.750 1.33e-05 ***
f4           7.79583    0.70817  11.008 1.13e-05 ***
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.1168 on 7 degrees of freedom
Multiple R-squared: 0.9484, Adjusted R-squared: 0.9189 
F-statistic: 32.16 on 4 and 7 DF,  p-value: 0.0001348 

現在我注意到,正如預期x的那樣,它是y.

我正在尋找的是對這種差異的直覺。我的想法在哪方面是錯誤的?為什麼我錯誤地期望x在這兩個模型中找到一個重要參數,但實際上只在固定效應模型中看到它?

這裡發生了幾件事。這些都是有趣的問題,但要解釋這一切需要相當多的時間/空間。

首先,如果我們繪製數據,這一切都會變得更容易理解。這是一個散點圖,其中數據點按組著色。此外,我們為每個組提供了單獨的組特定回歸線,以及粗體虛線的簡單回歸線(忽略組):

plot(y ~ x, data=dat, col=f, pch=19)
abline(coef(lm(y ~ x, data=dat)), lwd=3, lty=2)
by(dat, dat$f, function(i) abline(coef(lm(y ~ x, data=i)), col=i$f))

數據

固定效應模型

固定效應模型將如何處理這些數據相當簡單。的效果估計“控制”組。換句話說,首先相對於組啞元正交化,然後這個新的正交化的斜率是估計的。在這種情況下,這種正交化將消除(具體來說,集群間的變異性),因為組假人與. (為了直觀地認識到這一點,想想如果我們回歸會發生什麼僅在一組虛擬模型上,離開出等式。從上面的情節來看,我們當然希望有一些高-此回歸中每個虛擬係數的統計數據!)

所以基本上這對我們來說意味著只有**集群內的可變性用於估計效果. 簇間變異性(正如我們在上面看到的那樣,這是實質性的),被“控制”在分析之外。所以我們得到的斜率lm()是 4 條類內回歸線的平均值,在這種情況下,所有這些線都相對陡峭。

混合模型

混合模型的作用稍微複雜一些。混合模型試圖同時使用集群內和集群間的可變性估計效果. 順便說一句,這確實是該模型的賣點之一,因為它整合這些額外信息的能力/意願意味著它通常可以產生更有效的估計。但不幸的是,當和平均群內效應不太同意,就像這裡的情況一樣。注意:這種情況是面板數據的“豪斯曼檢驗”試圖診斷的!

具體來說,混合模型在這裡嘗試做的是估計平均簇內斜率之間的某種折衷。以及忽略集群的簡單回歸線(粗體虛線)。混合模型所依賴的這個折衷範圍內的確切點取決於隨機截距方差與總方差的比率(也稱為類內相關性)。當該比率接近 0 時,混合模型估計值接近簡單回歸線的估計值。當比率接近 1 時,混合模型估計接近平均集群內斜率估計。

以下是簡單回歸模型的係數(圖中的粗體虛線):

> lm(y ~ x, data=dat)

Call:
lm(formula = y ~ x, data = dat)

Coefficients:
(Intercept)            x  
  0.008333     0.008643  

如您所見,此處的係數與我們在混合模型中獲得的係數相同。這正是我們期望找到的,因為正如您已經指出的,我們對隨機截距的方差估計為 0,使前面提到的比率/類內相關性為 0。因此,在這種情況下,混合模型估計只是簡單的線性回歸估計,正如我們在圖中看到的那樣,這裡的斜率遠沒有集群內的斜率那麼明顯。

這給我們帶來了最後一個概念問題……

為什麼隨機截距的方差估計為 0?

這個問題的答案有可能變得有點技術性和困難性,但我會盡量保持簡單和非技術性(為了我們倆的緣故!)。但它可能還會有點囉嗦。

我之前提到了類內相關性的概念。這是另一種思考依賴的方式(或者更準確地說,是模型的錯誤)由聚類結構引起。類內相關性告訴我們,相對於從數據集中任何位置(即,可能或可能不在同一個簇中)抽取的兩個錯誤的平均相似度,從同一個簇中抽取的兩個錯誤的平均相似度有多大。一個正的類內相關性告訴我們,來自同一個集群的錯誤往往相對更相似;如果我從一個集群中得出一個錯誤並且它具有很高的價值,那麼我可以預期我從同一個集群中得出的下一個錯誤也將具有很高的價值。雖然不太常見,但類內相關性也可能是負的;從同一個集群中提取的兩個錯誤與整個數據集通常預期的相似度較低(即,值相差更遠)。

我們正在考慮的混合模型沒有使用表示數據依賴關係的類內相關方法。相反,它描述了方差分量方面的依賴性。只要類內相關性是正的,這一切都很好。在這些情況下,類內相關性可以很容易地用方差分量來寫,特別是前面提到的隨機截距方差與總方差的比率。(參見關於類內相關的 wiki 頁面有關這方面的更多信息。)但不幸的是,方差分量模型很難處理我們具有負類內相關性的情況。畢竟,將類內相關性寫成方差分量涉及將其寫成方差的比例,而比例不能為負。

從圖中可以看出,這些數據中的類內相關性看起來會略微為負。(我在得出這個結論時所看到的事實是,在每個集群內,但集群中相對較小的差異意味著,因此從同一個簇中提取的兩個錯誤往往具有幾乎跨越,而從不同集群中得出的誤差往往會有更溫和的差異。)因此,您的混合模型在實踐中正在做混合模型在這種情況下經常做的事情:它給出的估計值與負類內相關性一致因為它可以集合,但它停止在 0 的下限(這個約束通常被編程到模型擬合算法中)。所以我們最終得到一個估計的隨機截距方差為 0,這仍然不是一個很好的估計,但它與我們可以使用這種方差分量類型的模型得到的一樣接近。

所以,我們能做些什麼?

一種選擇是只使用固定效應模型。這在這裡是合理的,因為這些數據有兩個獨立的特徵,這對於混合模型來說很棘手(隨機組效應與和負類內相關)。

另一種選擇是使用混合模型,但設置它的方式是我們分別估計集群間和集群內的斜率而不是笨拙地試圖將它們匯集在一起。在這個答案的底部,我參考了兩篇討論這種策略的論文;我遵循 Bell & Jones 第一篇論文中提倡的方法。

為此,我們採取我們的預測器並將其分成兩個預測器,這將只包含集群間的變化, 和這將只包含集群內的變化. 這是這樣的:

> dat <- within(dat, x_b <- tapply(x, f, mean)[paste(f)])
> dat <- within(dat, x_w <- x - x_b)
> dat
     y   x f x_b x_w
1  -0.5   2 1   3  -1
2   0.0   3 1   3   0
3   0.5   4 1   3   1
4  -0.6  -4 2  -3  -1
5   0.0  -3 2  -3   0
6   0.6  -2 2  -3   1
7  -0.2  13 3  14  -1
8   0.1  14 3  14   0
9   0.4  15 3  14   1
10 -0.5 -15 4 -14  -1
11 -0.1 -14 4 -14   0
12  0.4 -13 4 -14   1
> 
> mod <- lmer(y ~ x_b + x_w + (1|f), data=dat)
> mod
Linear mixed model fit by REML 
Formula: y ~ x_b + x_w + (1 | f) 
  Data: dat 
  AIC   BIC logLik deviance REMLdev
6.547 8.972  1.726   -23.63  -3.453
Random effects:
Groups   Name        Variance Std.Dev.
f        (Intercept) 0.000000 0.00000 
Residual             0.010898 0.10439 
Number of obs: 12, groups: f, 4

Fixed effects:
           Estimate Std. Error t value
(Intercept) 0.008333   0.030135   0.277
x_b         0.005691   0.002977   1.912
x_w         0.462500   0.036908  12.531

Correlation of Fixed Effects:
   (Intr) x_b  
x_b 0.000       
x_w 0.000  0.000

這裡有幾點需要注意。首先,係數為與我們在固定效應模型中得到的完全一樣。到現在為止還挺好。二、係數為是我們從回歸中得到的回歸斜率僅在集群手段的向量上. 因此,它並不完全等同於我們第一個圖中的粗虛線,它使用了總方差,但它很接近。第三,雖然係數為比簡單回歸模型的係數小,標準誤差也小很多,因此-統計量更大。這也不足為奇,因為在這個混合模型中,殘餘方差要小得多,因為隨機組效應會消耗掉簡單回歸模型必須處理的大量方差。

最後,由於我在上一節中闡述的原因,我們仍然對隨機截距的方差進行了估計為 0。我不確定至少在不切換到除lmer(). 也許其他用戶可以對這個問題提出一些想法。

參考

  • Bell, A. 和 Jones, K. (2014)。解釋固定效應:時間序列橫截面和麵板數據的隨機效應建模。政治科學研究​​和方法。 PDF格式
  • Bafumi, J. 和 Gelman, AE (2006)。當預測變量和組效應相關時擬合多級模型。PDF格式

引用自:https://stats.stackexchange.com/questions/56695

comments powered by Disqus