R pso全局优化NaN报错的tryCatch处理及调试方法咨询

一、tryCatch的插入位置 #

建议在两个层级插入tryCatch做双重保障，避免异常值传递导致报错：

1. 第一层：在logdetHinv函数内拦截行列式计算警告 #

这里是NaN产生的源头，直接拦截警告返回极小值，避免异常值进入后续积分步骤：

logdetHinv<-function(dp, dw, t2){
  tryCatch({
    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)
    return(2*log(det(D1))-log(det(D0)))
  }, warning = function(w) {
    # 出现警告时返回极小值，对应积分后obj会生成大惩罚值，PSO会自动避开该组参数
    return(-1e9)
  }, error = function(e) {
    return(-1e9)
  })
}

补充：你遇到的行列式为负属于数值计算浮点误差，理论上正定矩阵行列式为正，但当矩阵接近奇异时，浮点运算的舍入误差会让原本接近0的正行列式变成绝对值极小的负数。除了用tryCatch，也可以直接对det结果加兜底逻辑：2*log(max(det(D1), 1e-30)) - log(max(det(D0), 1e-30))，不需要捕获警告也能避免NaN问题。

2. 第二层：在obj函数内拦截积分报错 #

防止个别异常t2取值导致integrate直接崩溃，给积分逻辑也加一层拦截：

obj<-function(x){
  dp=x[1:2]; dw=c(x[3],1-x[3])
  # 提前给dw加边界判断，避免dw[2]接近0导致的矩阵奇异问题
  if(dw[2] < 1e-6) return(1e9)
  tryCatch({
    fitness_value=-integrate(Vectorize(function(t2) logdetHinv(dp, dw, t2)*1/(t2u-t2l)), t2l, t2u)$value
    return(fitness_value)
  }, error = function(e) {
    # 积分出错时返回大惩罚值，引导PSO避开该参数区间
    return(1e9)
  })
}

二、定位报错的调试技巧 #

你在traceback里找不到D0是因为报错触发时，logdetHinv的执行环境已经被销毁，可以用下面几种方法捕获异常状态：

• 方法1：将警告转为错误，触发时自动进入调试环境
  先执行options(warn = 2)，此时所有警告都会被当作错误抛出，再运行代码，报错后执行debug(logdetHinv)标记该函数，重新运行后每次进入函数都会暂停，你可以逐行运行检查D0、D1的取值。
• 方法2：在函数内增加条件断点
  修改logdetHinv，在det为负时自动进入浏览模式：

  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)
    d0 = det(D0)
    d1 = det(D1)
    if(d0 <=0 || d1 <=0) {
      # 触发异常时自动进入交互浏览，此时可以直接打印dp、dw、t2、D0、D1的值
      browser()
    }
    return(2*log(d1)-log(d0))
  }
  

• 方法3：保存全局运行环境
  在代码开头执行options(error = dump.frames)，报错后当前所有函数的运行环境都会保存到last.dump对象里，执行debugger(last.dump)就可以选择进入logdetHinv的执行环境，查看所有局部变量的取值。

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