获取跟踪以在 R 中存储一个值

Question

我是运行R中的nloptr包，这是一个非线性优化包。在 nloptr 包中，我使用 cobyla() 命令来执行我的优化。现在，我可以使用 maxeval 参数指定算法的最大迭代次数，但是当 cobyla() 命令执行时，它只允许我查看最终评估的输出。但是，我希望能够在该方法的每次迭代中看到所有输出。

有人向我推荐我使用 trace() 命令来查看所有中间输出，但我不知道如何获取 R 来存储此信息。任何人都可以向我建议如何使用 trace() （或任何其他方式）来完成此操作。

这是我正在使用的代码：

library(nloptr)


#########################################################################################
### 
#########################################################################################
f = function(x){
    ans = cos(pi*x/5)+.2*sin(4*pi*x/5)

    return(ans)
}

const = function(x){
    ans = numeric(1)
    ans[1] = -(sin(pi*x/5)+.2*cos(4*pi*x/5))
    return(ans)
}

#########################################################################################
###
#########################################################################################

x0 = runif(1,0,10)

COB = cobyla(x0,f,hin=const,lower=0,upper=10,nl.info = TRUE,
      control = list(xtol_rel = 1e-16, maxeval = 2000))

所以我希望能够存储运行 cobyla() 迭代 1,2,...,maxeval

的输出

所以我认为我需要做的是使用

trace(f)

或

trace(cobyla)

在执行代码之前，但我不确定如何存储运行和 trace() 命令的值

Answer 1

假设 Josh 的评论是正确的，那么您需要自己进行循环。它看起来像这样...

maxeval<-2000
x0<-c()
x0[1] = runif(1,0,10)
for (cureval in 1:maxeval) {
    COB = cobyla(x0[cureval],f,hin=const,lower=0,upper=10,nl.info = TRUE,
      control = list(xtol_rel = 1e-16, maxeval = 1))
    x0[cureval+1]  <-COB$par

  }

Answer 2

在不调整 nloptr 代码的情况下，最好的办法是将 print_level 选项设置为 3，将输出发送到 sink 并将其抓取到取回你需要的东西。

请注意，要使用 print_level，您不能使用包装函数 cobyla——而是必须下拉到 nloptr 函数。

最后，请注意 cobyla 和具有 algorithm = NLOPT_LN_COBYLA 的 nloptr 函数表现不同——前者将不等式约束设为 >= 0，而后者将它们设为<=0.

在下面的代码中，我做的第一件事是使用核心 nloptr 函数并显示它与使用 cobyla 包装器代码得出相同的答案。

然后，我设置 print_level = 3 并将其写入文件。优化完成后，我读入该文件，并使用正则表达式从输出中提取解决方案路径。

这是它的样子：

library(tidyr)
library(nloptr)
library(assertthat)
library(dplyr)
library(ggplot2)

# starting parameter values
x0.hs100 = c(1, 2, 0, 4, 0, 1, 1)

# objective function
fn.hs100 = function(x) {
  (x[1]-10)^2 + 5*(x[2]-12)^2 + x[3]^4 + 3*(x[4]-11)^2 + 10*x[5]^6 +
    7*x[6]^2 + x[7]^4 - 4*x[6]*x[7] - 10*x[6] - 8*x[7]
}

# inequality constraints
# NOTE: nloptr with algorithm = NLOPT_LN_COBYLA takes g(x) <= 0
hin.hs100 = function(x) {
  h = numeric(4)
  h[1] = 127 - 2*x[1]^2 - 3*x[2]^4 - x[3] - 4*x[4]^2 - 5*x[5]
  h[2] = 282 - 7*x[1] - 3*x[2] - 10*x[3]^2 - x[4] + x[5]
  h[3] = 196 - 23*x[1] - x[2]^2 - 6*x[6]^2 + 8*x[7]
  h[4] = -4*x[1]^2 - x[2]^2 + 3*x[1]*x[2] -2*x[3]^2 - 5*x[6]    +11*x[7]
  return(-h)
}

# compute the solution
sink(file = "Output/cobyla_trace.txt", type = "output")
S1 = nloptr(x0 = x0.hs100, 
            eval_f = fn.hs100,
            eval_g_ineq = hin.hs100,
            opts = list(algorithm = "NLOPT_LN_COBYLA", 
                        xtol_rel = 1e-8, maxeval = 2000, print_level = 3))
sink()

# inequality constraints
# NOTE: cobyla takes g(x) >= 0
hin.hs100_minus = function(x) {
  h = numeric(4)
  h[1] = 127 - 2*x[1]^2 - 3*x[2]^4 - x[3] - 4*x[4]^2 - 5*x[5]
  h[2] = 282 - 7*x[1] - 3*x[2] - 10*x[3]^2 - x[4] + x[5]
  h[3] = 196 - 23*x[1] - x[2]^2 - 6*x[6]^2 + 8*x[7]
  h[4] = -4*x[1]^2 - x[2]^2 + 3*x[1]*x[2] -2*x[3]^2 - 5*x[6]    +11*x[7]
  return(h)
}

# compute the solution
S2 = cobyla(x0.hs100, fn.hs100, hin = hin.hs100_minus,
            nl.info = TRUE, control = list(xtol_rel = 1e-8, maxeval = 2000))

# check that both return the same solution
assertthat::are_equal(S1$solution, S2$par)

# get the solutions from the output file
solutionPath = readLines(con = file("Output/cobyla_trace.txt"))

# extract the solution path data out of the raw output
solutionPathParamRaw = solutionPath[grepl("^\tx", solutionPath)]
solutionPathParamMatch = gregexpr("(-)?[0-9]+(\.[0-9]+)?", solutionPathParamRaw, perl = TRUE)
solutionPathParam = as.data.frame(
  t(
    sapply(
      regmatches(
        solutionPathParamRaw, solutionPathParamMatch
      ), 
      as.numeric, simplify = TRUE
    )
  )
)

# give the columns some names
names(solutionPathParam) = paste("x", seq(1, 7), sep = "")
solutionPathParam$IterationNum = seq(1, dim(solutionPathParam)[1])

# plot the solutions
solutionPathParam %>%
  gather(Parameter, Solution, -IterationNum) %>%
  ggplot(aes(x = IterationNum, y = Solution, group = Parameter, color = Parameter)) +
  geom_line() + 
  theme_bw()

ggsave("Output/SolutionPath.png")

获取跟踪以在 R 中存储一个值

Getting trace to store a value in R

r

mathematical-optimization