為什麼這個線性混合模型是奇異的？
我試圖理解為什麼當線性混合效應模型擬合到下面的數據時我得到一個奇異的擬合。
我使用了 R
lme4::lmer，模型非常簡單，只有截距作為固定效應，因子變量作為隨機。這是數據集（可以復制並粘貼到 R）
data <- read.table(text= " group_id y 1 6.38 1 10.83 1 13.25 1 2.96 1 11.29 1 11.52 1 8.28 1 8.36 1 8.31 1 7.33 2 8.57 2 7.00 2 7.67 2 10.19 2 12.88 2 9.67 2 8.47 2 7.27 2 7.49 2 17.25 3 10.40 3 8.53 3 8.68 3 11.38 3 7.92 3 5.66 3 11.72 3 6.93 3 9.95 3 7.19 4 13.31 4 8.57 4 7.87 4 8.50 4 5.11 4 6.50 4 3.46 4 5.98 4 9.12 4 8.60 5 14.35 5 6.79 5 7.43 5 9.16 5 7.02 5 7.09 5 6.68 5 6.24 5 8.43 5 8.51", header= TRUE, colClasses= c('factor', 'numeric'))這是擬合模型：
library(lme4) fit <- lmer(data= data, y ~ 1 + (1|group_id)) boundary (singular) fit: see ?isSingular <<<<<< summary(fit) Linear mixed model fit by REML ['lmerMod'] Formula: y ~ 1 + (1 | group_id) Data: data REML criterion at convergence: 239 Scaled residuals: Min 1Q Median 3Q Max -2.139 -0.604 -0.093 0.467 3.242 Random effects: Groups Name Variance Std.Dev. group_id (Intercept) 0.00 0.00 Residual 7.05 2.66 Number of obs: 50, groups: group_id, 5 Fixed effects: Estimate Std. Error t value (Intercept) 8.641 0.376 23 optimizer (nloptwrap) convergence code: 0 (OK) boundary (singular) fit: see ?isSingular幫助
isSingular說，方差 - 協方差矩陣的一些“維度”已被估計為零，我想在數據中看到為什麼會發生這種情況。
正如您所發現的，當方差分量之一被估計為零時，就會發生這種情況。這通常有以下兩種解釋之一：
- 隨機效應結構過擬合——通常是因為隨機斜率太多
- 多個方差分量之一實際上非常接近於零，並且沒有足夠的數據來估計它高於零。
顯然，第一種情況與您的數據不同，因為您只有隨機截取。
因此，隨機截距的實際變化很可能
group_id非常接近於零。如果是，那麼只有 5 個組，軟件可能無法估計零以上的方差。如果通過繪製數據，這是一個很好的起點：
我們已經可以看到，與組內的變化相比，組均值的變化很小。
我們可以通過三種方式（至少）更正式地對此進行調查：
首先，讓我們看一下每組數據的均值：
library(tidyverse) data %>% group_by(group_id) %>% summarize(mean = mean(y)) ## 1 8.85 ## 2 9.65 ## 3 8.84 ## 4 7.70 ## 5 8.17請注意，所有組之間的差異很小，但請注意第 1 組和第 3 組的平均值幾乎相同。讓我們刪除第 1 組，看看會發生什麼：
data %>% filter(group_id != 1) %>% lmer(y ~ 1 + (1|group_id), data = .) %>% summary() ## Random effects: ## Groups Name Variance Std.Dev. ## group_id (Intercept) 0.03789 0.1947 ## Residual 6.75636 2.5993 ## Number of obs: 40, groups: group_id, 4 ## Fixed effects: ## Estimate Std. Error t value ## (Intercept) 8.5885 0.4224 20.34所以模型收斂沒有奇異性，但方差分量
group_id非常小，正如我們所懷疑的那樣。接下來，我們可以為
group_id組件添加一些額外的方差。這樣做的問題是，只有 5 個組，如果我們要從例如rnorm(5, 0, 1)（標準差為 1）中抽取 5 個觀測值，則樣本標準差可能不會接近 1，而平均值可能不會接近接近於零。解決這個問題的一個好方法是使用蒙特卡羅模擬（基本上只需多次進行並取平均值）。在這裡，我們將進行 100 次模擬：
n.sim <- 100 simvec_rint <- numeric(n.sim) # vector to hold the random intercepts variances simvec_fint <- numeric(n.sim) # vector to hold the fixed intercepts for (i in 1:n.sim) { set.seed(i) data$y1 = data$y + rep(rnorm(5, 0, 1), each = 10) m0 <- lmer(y1 ~ 1 + (1|group_id), data = data) if (!isSingular(m0)) { # If the model is not singular then extract the random and fixed effects VarCorr(m0) %>% as.data.frame() %>% pull(vcov) %>% nth(1) -> simvec_rint[i] summary(m0) %>% coef() %>% as.vector() %>% nth(1) -> simvec_fint[i] } else { simvec_rint[i] <- simvec_fint[i] <- NA } }因此，我們向方差為 1 的組添加了隨機噪聲。蒙特卡羅估計值為：
> mean(simvec_rint, na.rm = TRUE) [1] 1.116416 > mean(simvec_fint, na.rm = TRUE) [1] 8.637063注意：
- 隨機截距的方差為 1.12。然而，我們已經向等於 1 的組添加了方差，因此這意味著原始數據中隨機截距的方差接近於零，正如我們所懷疑的那樣。
- 固定截距為 8.64，與原始數據擬合的模型基本相同。
最後，讓我們看一個沒有隨機效應的模型，它顯然只是一個 ANOVA：
> aov(y ~ group_id, data = data) %>% summary() Df Sum Sq Mean Sq F value Pr(>F) group_id 4 22.0 5.489 0.763 0.555 Residuals 45 323.6 7.190因此，幾乎沒有證據表明這 5 個組的均值彼此不同。另一種看待這個問題的方法是：
> lm(y ~ group_id, data = data) %>% summary() Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 8.8510 0.8479 10.438 1.33e-13 *** group_id2 0.7950 1.1992 0.663 0.511 group_id3 -0.0150 1.1992 -0.013 0.990 group_id4 -1.1490 1.1992 -0.958 0.343 group_id5 -0.6810 1.1992 -0.568 0.573因此，也很少有證據表明第 2、3、4 和 5 組與第 1 組具有不同的均值。這兩個模型都與混合模型中隨機截距的非常小的變化一致。
因此，綜上所述，在這種情況下，我們可以得出結論，由於組數較少且組間估計的變化較小，軟件無法估計零以上的隨機截距變化，因此模型是奇異，儘管模型估計似乎是可靠的。
