R语言拟合Dose-response curve时遇Convergence failed错误求助

问题与错误信息 #

拟合剂量反应曲线时，所有模型均出现以下收敛错误：

Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt, :
initial value in 'vmmin' is not finite
Error in drmOpt(opfct, opdfct1, startVecSc, optMethod, constrained, warnVal, :
Convergence failed

用户使用的R代码如下：

install.packages("drc")
library(drc)
library(ggplot2)

Sulf_Conc <- c(-3, -2.958607315, -2.698970004, -1.958607315, -0.995678626, 
  -0.300162274, 0.699056855)
Enzyme_Activities <- c(100, 55.30724303, 69.24958323, 82.88291229, 100.627504, 
  104.0535174, 67.16165961)
missing_enzyme <- sum(is.na(Enzyme_Activities))
cat("Missing values in Enzyme_Activities:", missing_enzyme, "\n")
missing_sulf <- sum(is.na(Sulf_Conc))
cat("Missing values in Sulf_Conc:", missing_sulf, "\n")
initial_params <- c(A = min(Enzyme_Activities),  # Lower asymptote
                    D = max(Enzyme_Activities),  # Upper asymptote
                    B = 1,  # Reasonable slope value (can be adjusted)
                    C = median(Sulf_Conc))  # EC50 (e.g., use median of concentrations)

ll.4 <- drm(
  Enzyme_Activities ~ Sulf_Conc, 
  data = DHA_Sulphadiazine, 
  fct = LL.4(), 
  start = initial_params
)
ll.4 <- drm(
  Enzyme_Activities ~ Sulf_Conc, 
  data = DHA_Sulphadiazine, 
  fct = LL.4()
)

错误原因与解决方案 #

1. 修复未定义的数据框 #

代码中引用了DHA_Sulphadiazine数据框，但未提前创建。需先将变量组合为数据框：

DHA_Sulphadiazine <- data.frame(Sulf_Conc, Enzyme_Activities)

2. 数据趋势不匹配LL.4模型 #

LL.4模型是典型的S型曲线模型，要求剂量与响应呈单调相关（要么随剂量增加持续上升，要么持续下降）。但你的数据：

• 浓度从-3升至-0.3时，酶活性从100升至104
• 浓度进一步升至0.699时，酶活性骤降至67
  这种非单调的波动是拟合失败的核心原因，优化算法无法找到符合S型假设的参数。

3. 调整模型或数据 #

• 改用无单调假设的模型：使用GN()广义非线性模型，它不强制趋势单调性，更适配当前数据：

gn_model <- drm(Enzyme_Activities ~ Sulf_Conc, data = DHA_Sulphadiazine, fct = GN())
summary(gn_model)

• 验证数据有效性：检查浓度和酶活性的测量是否存在误差，或是否需要对数据做转换（如浓度取反、调整坐标轴方向）。
• 补充样本量：仅7个数据点的样本量过小，拟合稳定性差，建议补充更多数据点。

4. 强制拟合LL.4的优化尝试（谨慎使用） #

如果必须使用LL.4模型，需调整初始参数并增加迭代次数，但拟合结果可能无生物学意义：

initial_params <- c(A = 50, D = 105, B = -2, C = 0)  # B设为负，匹配高浓度活性下降趋势

ll.4 <- drm(
  Enzyme_Activities ~ Sulf_Conc, 
  data = DHA_Sulphadiazine, 
  fct = LL.4(), 
  start = initial_params,
  control = list(maxit = 1000)  # 增加迭代次数
)

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