尝试在 R 中运行 lmer() 时出错

Question

所以这是我的问题。我在 R 中有一个数据集，需要在其上运行混合效果模型。代码如下：

data <- read.csv("D:/blahblah.csv")
analysis.data <- lmer(intdiff ~ stress_limit * word_position * follows + (1|speaker), data)
summary(analysis.data)

当我尝试运行脚本时，它返回以下错误：

 Error in mer_finalize(ans) : Downdated X'X is not positive definite, 15.

我已将错误跟踪到“follows”参数，因为当我只使用stress_limit 和word_position 时，它运行良好。如果有帮助，“follows”中的数据只有 3 个字符串：n 或 l、辅音、元音。我尝试用 _ 替换空格，但没有成功。在这种情况下， lmer() 函数的内部工作是否阻止了“跟随”的使用？任何帮助都会很棒！

更多信息：intdiff 包含数值，stress_limit 是字符串（重读或非重读），单词位置也是字符串（Word Medial 或 Word Initial）。

编辑：这是重现错误的数据示例：

structure(list(intdiff = c(11.45007951, 12.40144758, 13.47898367, 
6.279497762, 18.19461897, 16.15539707), word_position = structure(c(2L, 
2L, 2L, 1L, 1L, 1L), .Label = c("Word Initial", "Word Medial"
), class = "factor"), follows = structure(c(4L, 4L, 4L, 1L, 2L, 
4L), .Label = c("Consonant", "n or l", "Pause", "Vowel"), class = "factor"), 
stress_limit = structure(c(2L, 1L, 1L, 2L, 2L, 2L), .Label = c("Stressed", 
"Unstressed"), class = "factor"), speaker = structure(c(2L, 
2L, 2L, 2L, 2L, 2L), .Label = c("f11r", "f13r", "f15a", "f16a", 
"m09a", "m10a", "m12r", "m14r"), class = "factor")), .Names = c("intdiff", 
"word_position", "follows", "stress_limit", "speaker"), row.names = c(NA, 
6L), class = "data.frame")

我也尝试了 lme() 函数，但它返回了这个错误：

Error in MEEM(object, conLin, control$niterEM) : 
Singularity in backsolve at level 0, block 1

我原始帖子中的代码是我正在使用的确切代码，减去 library(lme4) 调用，所以我不会留下任何我能想到的信息。

我的 R 版本是 2.15.2

Answer 1

如果没有可重复的示例，很难确定：How to make a great R reproducible example?

但是，猜测：这类问题通常是由于设计矩阵中的共线性造成的。将您的连续预测器 (intdiff) 居中可能会有所帮助。也可以直接探索设计矩阵

X <- model.matrix( ~ stress_limit * word_position * follows, data)

对之间的共线性：cor(X)。不幸的是，尽管您可以查看计算工具方差膨胀因子（例如library("sos"); findFn("VIF")）。

作为交叉检查，lme 也应该能够处理您的模型：

library(nlme)
lme(intdiff ~ stress_limit * word_position * follows, 
   random=~1|speaker, data=data)

当我在 lme4 (available on github) 的开发版本中运行您的测试数据时，我得到了Error in lmer(intdiff ~ stress_limit * word_position * follows + (1 | : rank of X = 5 < ncol(X) = 12。另一方面，对于这么小的输入数据集（6 个观察值），不可能拟合 12 个参数。确切地说出问题出在哪里有点困难。您的 3 个变量的所有 12 种组合是否都实际出现在您的数据中？如果缺少一些，那么您需要遵循开发版本帮助中给出的建议：

不像一些更简单的建模框架，如“lm”和“glm”，它们 自动检测完全共线的预测变量， “[gn]lmer”无法处理小于满秩的设计矩阵。为了 例如，在模型具有未观察到的交互的情况下 级别的组合，由用户定义新变量 （例如，根据结果在数据中创建“ab” 'droplevels(interaction(a,b))').

特别是，您可以按如下方式拟合此模型：

data <- transform(data,
       allcomb=interaction(stress_limit,word_position,follow,drop=TRUE))
lme(intdiff ~ allcomb, random=~1|speaker, data=data)

这将为您提供单向方差分析，将数据中实际存在的独特级别组合视为类别。您必须自己弄清楚它们的含义。

另一种方法是减少模型中的交互次数，直到你得到一个没有任何缺失组合的集合；如果你很幸运(stress_limit+word_position+follow)^2（所有双向交互）会起作用，但你可能需要进一步缩小模型（例如stress_limit + word_position*follow）。

另一种测试方法是在您提出的模型上使用lm()，并检查估计系数中是否没有NA 值。

您将在这些方式中失去的主要是方便/易于解释，因为无论如何都无法从数据中估计缺失组合的参数......