使用带有长变量列表的 R 函数 nls答案

【问题标题】：using R function nls with long list of variables使用带有长变量列表的 R 函数 nls
【发布时间】：2018-06-12 08:35:17
【问题描述】：

我试图用 R 函数 nls 估计一个非线性模型。

我将通过 nls 的“帮助”工具中的一个示例来说明我的问题。

Treated <- Puromycin[Puromycin$state == "treated", ]
weighted.MM <- function(resp, conc, Vm, K)
{
    pred <- (Vm * conc)/(K + conc)
    (resp - pred) / sqrt(pred)
}

Pur.wt <- nls( ~ weighted.MM(rate, conc, Vm, K), data = Treated,
          start = list(Vm = 200, K = 0.1))

在这个例子中，weighted.MM 是一个参数数量非常有限的函数，对于我正在使用的模型类型实现类似的方法没有问题。

然而，我现在正试图转向一个更现实的问题，我实际上有几十个参数要传递给函数。我当然可以简单地列举它们，但我觉得这有点笨拙。

我考虑过先将它们放在一个单独的列表中。例如，使用上面的示例，我将首先定义：

MyArguments <- list(Vm, K)

然后将 MyArguments 作为参数传递给函数（并从函数内访问各个参数）。但是，这不起作用，因为我收到错误消息

Error: object 'Vm' not found

或者，

    MyArguments <- list("Vm", "K") 

weighted.MM1 <- function(resp, conc1, conc.1, thearguments)
{
  conc <- c(conc1, conc.1)
  pred <- (thearguments[[1]] * conc)/(thearguments[[2]] + conc)
  (resp - pred) / sqrt(pred)
}

Pur.wt1 <- nls( ~ weighted.MM1(rate, conc1, conc.1, MyArguments),
                data = lisTreat, start = list(Vm = 200, K = 0.1))

产量：

Error in thearguments[[1]] * conc : 
  non-numeric argument to binary operator
Called from: weighted.MM1(rate, conc1, conc.1, MyArguments)

有什么解决方法吗？

【问题讨论】：

标签： r arguments nls non-linear-regression

【解决方案1】：

由于您正在查找参数的值，因此无需定义它们：请看下面的代码：

weighted.MM1 <- function(resp, conc1, conc.1, thearguments)
{
  conc <- c(conc1, conc.1)
  pred <- (thearguments[1] * conc)/(thearguments[2] + conc)
  (resp - pred) / sqrt(pred)
}
nls( ~ weighted.MM1(rate, conc1, conc.1, MyArguments),
                data = lisTreat, start = list(MyArguments =c( 200, 0.1)))
Nonlinear regression model
  model: 0 ~ weighted.MM1(rate, conc1, conc.1, MyArguments)
   data: lisTreat
MyArguments1 MyArguments2 
   206.83468      0.05461 
 residual sum-of-squares: 14.6

Number of iterations to convergence: 5 
Achieved convergence tolerance: 3.858e-06

【讨论】：

哦，太棒了！谢谢你。
@LaurentFranckx 如果这回答了您的问题，您可以接受它作为解决方案
我认为问题已经解决了。我首先认为这是由于我弄错了，但示例代码实际上也出现了同样的问题。要看到这一点，请使用 nlsLM 而不是 nls 运行回归。如果这样做，您将收到错误消息：“rownames<-(*tmp*, value = "MyArguments") 中的错误：'dimnames' [1] 的长度不等于数组范围”。如果你深入 nls.lm 的代码，你会发现只有列表的第一个元素会随着迭代而变化。
对不起。这是否意味着这不能解决[问题？
恐怕不行。我现在已经详细阐述了我在对原始问题的“答案”中遇到的新问题。

【解决方案2】：

我已尝试根据@Onyambu 提出的建议实施新的解决方案。

不过，这也带来了新的问题。

首先，我尝试使用 nls 实现一个解决方案。这是我使用的实际代码：

DiffModel <- nls(COPERTFreq ~ CalculateProbaVehChoiceDiffusion(MyArguments, years = RegistYear),
                 data = DataSetForModel , start = list(MyArguments = c(ASC_Mat, ASC_No_Size)))

其中 CalculateProbaVehChoiceDiffusion() 是在别处定义的非线性函数，RegistYear 是常数，而 MyArguments 是要估计的系数列表，其中 c(ASC_Mat, ASC_No_Size) 作为初始值。

这会导致以下错误消息：

Error in numericDeriv(form[[3L]], names(ind), env) : 
  Missing value or an infinity produced when evaluating the model

现在，我在其他地方读到，这个问题可以通过使用 nlsLM 来解决。这会导致一条新的错误消息：

Error in `rownames<-`(`*tmp*`, value = "MyArguments") : 
  length of 'dimnames' [1] not equal to array extent

好的，所以我在调试模式下使用 nls.lm 再次运行了模型。这表明错误消息来自以下代码行：

names(out$par) <- rownames(out$hessian) <- colnames(out$hessian) <- names(out$diag) <- names(par)

然而，正是通过检查“out”对象才能清楚问题所在。首先， out$hessian 只是一个标量，而我希望它的行数和列数等于参数的数量。其次，out$par$MyArguments 表明，除了第一个元素之外，MyArguments 的值不会从一个迭代到另一个迭代。

这是一个已知的错误还是我必须修改将 MyArguments 传递给函数调用的方式？

请注意，据我所知，当我将 nlsLm 应用于@Onyambu 提供的示例时也会出现此问题：

> undebug(nls.lm)
> Treated <- Puromycin[Puromycin$state == "treated", ]
> 
> lisTreat <- with(Treated,
+                  list(conc1 = conc[1], conc.1 = conc[-1], rate = rate))
> 
> weighted.MM1 <- function(resp, conc1, conc.1, thearguments)
+ {
+   conc <- c(conc1, conc.1)
+   pred <- (thearguments[1] * conc)/(thearguments[2] + conc)
+   (resp - pred) / sqrt(pred)
+ }
> nls( ~ weighted.MM1(rate, conc1, conc.1, MyArguments),
+      data = lisTreat, start = list(MyArguments =c( 200, 0.1)))
Nonlinear regression model
  model: 0 ~ weighted.MM1(rate, conc1, conc.1, MyArguments)
   data: lisTreat
MyArguments1 MyArguments2 
   206.83468      0.05461 
 residual sum-of-squares: 14.6

Number of iterations to convergence: 5 
Achieved convergence tolerance: 3.858e-06
> 
> nlsLM( ~ weighted.MM1(rate, conc1, conc.1, MyArguments),
+        data = lisTreat, start = list(MyArguments =c( 200, 0.1)))
Error in `rownames<-`(`*tmp*`, value = "MyArguments") : 
  length of 'dimnames' [1] not equal to array extent

因此，虽然 nls 可用于此示例，但 nlsLM 不能，并且会产生与我的代码相同的错误消息。然后我在调试模式下使用 nls.lm 再次运行 nlsLM。在以下行之后，

out <- .Call("nls_lm", par, lower, upper, fn1, jac1, ctrl, 
    new.env(), PACKAGE = "minpack.lm")

我检查了输出对象并看到：

$par
$par$MyArguments
[1] 244.5117   0.1000


$hessian
[1] 0.02859739

$fvec
 [1] -5.5215487 -0.9787468 -0.5543382 -1.5986605  0.4486608 -0.9651245  0.7020058  1.2419040  1.1430780  0.4488084  1.1445818  1.6121474

$info
[1] 1

$message
[1] "Relative error in the sum of squares is at most `ftol'."

$diag
$diag[[1]]
[1] 0.2077949


$niter
[1] 4

$rsstrace
[1] 112.59784  43.41211  42.89350  42.89349  42.89349

因此，第二个参数的值在 4 次迭代后没有改变。这可能当然是正确的解决方案。但我发现我的模型也发生了同样的情况，这是一个有趣的巧合。

最终编辑：我最终决定依靠对所有论点的完整枚举。正如我在问题的陈述中所写的，它不是很优雅，但是，至少它可以与 nls.lm 一起使用（尽管仍然不能与 nls 一起使用）。

【讨论】：