NLME 函数中的错误：第 0 级、第 1 块的反求解中存在奇异性

Question

我正在尝试使用 nlme 函数将 ELISA 板数据拟合到非线性混合效应模型。之前，我询问过如何使语法适用于嵌套组，尽管答复让我觉得这是一个失败的原因。但是，当仅应用单个组时，我仍然难以生成模型对象。

示例数据可以在这里找到：https://pastebin.com/UzVa04TY

这是我运行 NLME 函数的最初设想，使用单独的变量来实现固定和随机效应（分别是大写字母和小写字母）。 x 是样品稀释倍数，y 是 ELISA 信号密度。

# Define the 3-parameter logistic function
logistic3PL <- function(x, A, B, C, a, b, c) {
  (A + a) / (1 + (x / (C + c))^ (B + b))
}

# Define the nonlinear mixed effects model
model <- nlme(y ~ logistic3PL(x, A, B, C, a, b, c),
              data = dataSet,
              fixed = A + B + C ~ 1,
              random = a + b + c ~ 1,
              groups =  ~ Plate,
              start = c(A = 1, B = 1, C = 1))

这产生了以下错误：

Error in nlme.formula(y ~ logistic3PL(x, A, B, C, a, b, c), data = dataSet,  : 
  Singularity in backsolve at level 0, block 1

我之前问题的评论者建议他实际上会用 3 参数函数来表达固定效应和随机效应，并用同一组字母表示。

# Define the 3-parameter logistic function
logistic3PL <- function(x, A, B, C) {
  A / (1 + (x / C) ^B)
}

# Define the nonlinear mixed effects model
model <- nlme(y ~ logistic3PL(x, A, B, C),
              data = dataSet,
              fixed = A + B + C ~ 1,
              random = A + B + C ~ 1,
              groups =  ~ Plate,
              start = c(A = 1, B = 1, C = 1))

这产生了相同的错误。

我应该注意到，我还尝试在这两种方法中将

random = A + B + C ~ 1,

替换为

... ~ Plate,

。每次，程序运行时都会思考大约 30 秒，然后才会给出完全相同的错误，而不是当该行以

~ 1

结尾时立即给出错误。

很乐意回答任何可能有助于解决问题的问题。预先感谢您的帮助！

Answer 1

我建议使用参考文献 [1] 中公式 1 给出的参数化。

#print the data
library(ggplot2)
p <- ggplot(dataSet, aes(x, y)) +
  geom_point(aes(color = factor(Plate)))
print(p)

#we use the other parametrization from the reference
#they say the base of the log can be arbitrary but
#if you know the dilution factor, you can use that for the base
foo <- function(x, A, B, C) A / (1 + exp((C - log10(x))/B))

#plot the function to see if the starting values are sensible 
curve(foo(x, 1, 10, 0))

library(nlme)
#first use nlsList
fit1 <- nlsList(y ~ foo(x, A, B, C) | Plate, data = dataSet, 
                start = list(A = 1, B = 10, C = 0))

#then use that fit as input for nlme
#specify no correlation between random effects with pdDiag
fit2 <- nlme(fit1, random = pdDiag(list(A ~ 1, B ~ 1, C ~ 1)), 
             groups = ~ Plate)
summary(fit2)

#plot predictions
p + lapply(1:3, 
           \(i) geom_function(fun = \(x) predict(fit2, newdata = data.frame(x = x, Plate = i)),
                              aes(color = factor(i))))

大多数人都会同意三个受试者不足以估计随机效应。您应该大幅增加盘子数量。

[1]：https://www.aphis.usda.gov/animal_health/vet_biologics/publications/STATWI0006.pdf