问R具有潜在警告/错误的全局优化问题和tryCatch的使用
EN

Stack Overflow用户

提问于 2021-08-31 21:55:42

回答 1查看 54关注 0票数 1

简单地说，我有一个函数f(x, t2)，我想找出x的值，使f(x, t2)相对于t2的积分最大化。选择粒子群算法进行优化。可执行代码如下所示

library(pso)
xl=0; xu=2000; n=1; t2l=100; t2u=2000; t1=1
g<-function(x, t2) t1*x/(t2+x)
h<-function(z) 1/z^n
gdot<-function(x, t2){
  c(x/(t2+x),-t1*x/(t2+x)^2)
}
logdetHinv<-function(dp, dw, t2){
  gmat=mapply(function(x) gdot(x,t2),dp)
  D0=gmat%*%diag(dw)%*%t(gmat)
  D1=gmat%*%diag(1/h(g(dp,t2)))%*%diag(dw)%*%t(gmat)
  2*log(det(D1))-log(det(D0))
}
obj<-function(x){
  dp=x[1:2]; dw=c(x[3],1-x[3])
  fitness_value=-integrate(Vectorize(function(t2) logdetHinv(dp, dw, t2)*1/(t2u-t2l)), t2l, t2u)$value
  return(ifelse(dw[2]>0, fitness_value, fitness_value+1e3))
}

x <- psoptim(rep(1,3), fn = obj, lower = c(rep(xl,2),0.1), upper = c(rep(xu,2), 0.9))$par
x

因为全局优化涉及一些随机过程，所以它有时会报告正确的结果

> x
[1] 2000.0000  754.4146    0.5000

其他时候它会报告错误

Error in integrate(Vectorize(function(t2) logdetHinv(dp, dw, t2) * 1/(t2u -  : 
  non-finite function value
In addition: There were 11 warnings (use warnings() to see them)
> warnings()
Warning messages:
1: In log(det(D1)) : NaNs produced
2: In log(det(D0)) : NaNs produced
3: In log(det(D1)) : NaNs produced
4: In log(det(D0)) : NaNs produced

我假设该算法尝试对logdetHinv中的一些负值进行记录，这将返回NaN并显示一条警告消息，而不是错误，最终导致integrate中出现错误。

我想要避免这样的值，可能是使用tryCatch，比如如果在函数logdetHinv中有警告，它返回一个非常小的值，但不是NaN，所以它不会在integrate中导致错误，并且psoptim在最大化目标函数(最小化-integrate(logdetHinv))时不太可能选择这样的值。在如此复杂的情况下，我对tryCatch并不熟悉。我应该把tryCatch放在哪里？谢谢。

此外，我想知道R中是否有一些调试技术可以让我知道在这种情况下是什么随机值(D0/D1)导致了错误。我猜它在log中是一些负值，但它不应该是负值，因为log中是一个正定矩阵的行列式。在回溯模式下，在browse中，如果我键入D0，将找不到对象'D0‘。

debugging

optimization

try-catch

腾讯特效SDK 2.5折起

美颜基础/原子能力套餐低至1500元/月，提供丰富的美化能力，支持全平台集成

回答 1

Stack Overflow用户

回答已采纳

发布于 2021-08-31 22:33:19

在这种情况下，我不会使用tryCatch，因为它在测试中通常比在主代码中更合适。为什么不简单地测试函数中的行列式呢？像这样的东西应该是有效的：

logdetHinv<-function(dp, dw, t2){
    gmat=mapply(function(x) gdot(x,t2),dp)
    D0=gmat%*%diag(dw)%*%t(gmat)
    D1=gmat%*%diag(1/h(g(dp,t2)))%*%diag(dw)%*%t(gmat)
    detD1 <- max(0.01, det(D1))
    detD0 <- max(0.01, det(D0)) 
    2*log(detD1)-log(detD0)
}