R语言nlme函数groups参数嵌套分组语法问题求助

问题背景 #

使用R语言的nlme函数拟合ELISA板数据的非线性混合效应模型时，遇到groups参数嵌套分组的语法错误，现有示例仅覆盖单分组场景，缺乏嵌套分组参考。

模型设定 #

数据中每个板（j）包含Ref和Test两种样本类型（i）的Sigmoid曲线，采用三参数逻辑斯蒂（3PL）模型：
$$(A+a)/(1+(k/(C+c))^{(B+b)})$$
其中：

• 大写字母A、B、C为固定效应
• 小写字母a、b、c为基于板的随机效应
• k为稀释因子

报错情况 #

1. 使用groups = ~i/j时，出现以下错误：

Error in names(reSt) <- namGrp : 
  'names' attribute [2] must be the same length as the vector [1]

2. 使用groups = ~i + j时，出现以下错误：

Error in getGroups.data.frame(dataMix, eval(parse(text = paste("~1", deparse(groups[[2]]),  : 
  invalid formula for groups

样本数据结构 #

structure(list(i = c("Ref", "Ref", "Ref", "Ref", "Ref", "Ref", 
"Ref", "Ref", "Ref", "Ref", "Ref", "Test", "Test", "Test", "Test", 
"Test", "Test", "Test", "Test", "Test", "Test", "Test", "Ref", 
"Ref", "Ref", "Ref", "Ref", "Ref", "Ref", "Ref", "Ref", "Ref", 
"Ref", "Test", "Test", "Test", "Test", "Test", "Test", "Test", 
"Test", "Test", "Test", "Test"), j = c("1", "1", "1", "1", "1", 
"1", "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", "1", 
"1", "1", "1", "1", "2", "2", "2", "2", "2", "2", "2", "2", "2", 
"2", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2", "2"
), k = c("5", "10", "20", "40", "80", "160", "320", "640", "1280", 
"2560", "5120", "5", "10", "20", "40", "80", "160", "320", "640", 
"1280", "2560", "5120", "5", "10", "20", "40", "80", "160", "320", 
"640", "1280", "2560", "5120", "5", "10", "20", "40", "80", "160", 
"320", "640", "1280", "2560", "5120"), y = c("0.0128", "0.0269", 
"0.0611", "0.1783", "0.4333", "0.6237", "0.846", "0.8952", "0.9145", 
"0.9443", "0.9831", "0.00470000000000001", "0.0087", "0.0223", 
"0.0557", "0.137", "0.3419", "0.555", "0.6829", "0.7723", "0.8426", 
"0.8433", "0.0332", "0.042", "0.0755", "0.172", "0.4082", "0.6905", 
"0.8783", "0.9687", "1.0235", "1.0798", "1.0744", "0.0114", "0.0165", 
"0.0368", "0.1016", "0.1442", "0.3734", "0.6292", "0.7851", "0.8644", 
"0.8932", "0.9376")), row.names = c(NA, 44L), class = "data.frame")

解决方案 #

1. 修正数据类型 #

首先将字符型的k和y转换为数值型，否则无法正常拟合：

dataMix$k <- as.numeric(dataMix$k)
dataMix$y <- as.numeric(dataMix$y)

2. 正确的嵌套分组语法 #

根据数据结构，板（j）是外层分组，样本类型（i）是板内的内层分组，嵌套分组的正确语法为：

groups = ~j/i

该公式表示i嵌套在j中，即每个板内的不同样本类型作为子分组。

3. 匹配模型与分组逻辑 #

若随机效应a、b、c仅基于板（j），而非板内的样本类型，无需嵌套分组，直接指定groups = ~j即可。此时需确保模型公式中样本类型i作为固定效应纳入（例如让A、B、C随i变化），示例模型公式定义如下：

logistic3PL <- function(k, A_Ref, A_Test, B_Ref, B_Test, C_Ref, C_Test, a, b, c) {
  A <- ifelse(dataMix$i == "Ref", A_Ref, A_Test)
  B <- ifelse(dataMix$i == "Ref", B_Ref, B_Test)
  C <- ifelse(dataMix$i == "Ref", C_Ref, C_Test)
  (A + a)/(1 + (k/(C + c))^(B + b))
}

调用nlme时的参数设置：

library(nlme)
fit <- nlme(y ~ logistic3PL(k, A_Ref, A_Test, B_Ref, B_Test, C_Ref, C_Test, a, b, c),
            data = dataMix,
            fixed = A_Ref + A_Test + B_Ref + B_Test + C_Ref + C_Test ~ 1,
            random = a + b + c ~ 1 | j,
            groups = ~j,
            start = c(A_Ref = 1, A_Test = 0.9, B_Ref = 1.5, B_Test = 1.2, C_Ref = 80, C_Test = 160))

4. 错误原因说明 #

• 使用groups = ~i/j时，分组顺序颠倒（将样本类型作为外层、板作为内层），不符合数据结构，导致随机效应分组名称与数量不匹配。
• 使用groups = ~i + j时，nlme的groups参数不支持用+表示分组关系，交叉分组需用~i*j，嵌套分组需用外层/内层的语法。

内容的提问来源于stack exchange，提问作者AndyH