0 
obj1<-function(monthly.savings,
               success,
               start.capital,
               target.savings,
               monthly.mean.return,
               monthly.ret.std.dev,
               monthly.inflation,
               monthly.inf.std.dev,
               n.obs,
               n.sim=1000){

  req = matrix(start.capital, n.obs+1, n.sim) #matrix for storing target weight

  monthly.invest.returns = matrix(0, n.obs, n.sim)
  monthly.inflation.returns = matrix(0, n.obs, n.sim)

  monthly.invest.returns[] = rnorm(n.obs * n.sim, mean = monthly.mean.return, sd = monthly.ret.std.dev)
  monthly.inflation.returns[] = rnorm(n.obs * n.sim, mean = monthly.inflation, sd = monthly.inf.std.dev)

  #for loop to be 
  for (a in 1:n.obs){
    req[a + 1, ] = req[a, ] * (1 + monthly.invest.returns[a,] - monthly.inflation.returns[a,]) + monthly.savings
  }

  ending.values=req[nrow(req),]
  suc<-sum(ending.values>target.savings)/n.sim
  value<-success-suc

  return(abs(value))
}

我有上面想要最小化的目标函数。它试图解决给定成功概率所需的每月节省。给定以下输入假设

success<-0.9
start.capital<-1000000
target.savings<-1749665
monthly.savings=10000
monthly.mean.return<-(5/100)/12
monthly.ret.std.dev<-(3/100)/sqrt(12)
monthly.inflation<-(5/100)/12
monthly.inf.std.dev<-(1.5/100)/sqrt(12)
monthly.withdrawals<-10000
n.obs<-10*12  #years * 12 months in a year
n.sim=1000

我使用了以下符号：

optimize(f=obj1,
        success=success,
        start.capital=start.capital,
        target.savings=target.savings,
        monthly.mean.return=monthly.mean.return,
        monthly.ret.std.dev=monthly.ret.std.dev,
        monthly.inflation=monthly.inflation,
        monthly.inf.std.dev=monthly.inf.std.dev,
        n.obs = n.obs,
        n.sim = n.sim,
        lower = 0,
        upper = 10000,
        tol = 0.000000001,maximum=F)

我得到 7875.03

由于我是从正态分布中采样的，因此每次输出都会有所不同，但它们应该大致相同，或者取几个百分点。我遇到的问题是我不能任意指定上限。上述示例的上限（10000）是经过多次试验后采摘的樱桃。如果说我输入了 100000 的上限（我知道这不合理），它将返回该数字，而不是找到全局最低储蓄。我在错误地构建目标函数的任何想法？

谢谢，

4

1 回答 1

4

对于给定的输入，您的函数并不总是返回相同的输出这一事实可能会带来一些问题（它会产生很多虚假的局部最小值）：您可以通过在函数（例如，set.seed(1)），或者通过存储随机数并每次重用它们，或者通过使用低差异序列（例如，randtoolbox::sobol）。

由于它是一个变量的函数，您可以简单地绘制它以查看发生了什么：它在 10,000 之后有一个平台——优化算法无法区分平台和局部最优。

f <- function(x) {
  set.seed(1)
  obj1(x,
      success             = success,
      start.capital       = start.capital,
      target.savings      = target.savings,
      monthly.mean.return = monthly.mean.return,
      monthly.ret.std.dev = monthly.ret.std.dev,
      monthly.inflation   = monthly.inflation,
      monthly.inf.std.dev = monthly.inf.std.dev,
      n.obs               = n.obs,
      n.sim               = n.sim 
  )
}
g <- Vectorize(f)
curve(g(x), xlim=c(0, 20000))

您最初的问题实际上不是最小化问题，而是求根问题，这要容易得多。

obj2 <- function(monthly.savings) {
  set.seed(1)
  req = matrix(start.capital, n.obs+1, n.sim)
  monthly.invest.returns <- matrix(0, n.obs, n.sim)
  monthly.inflation.returns <- matrix(0, n.obs, n.sim)
  monthly.invest.returns[] <- rnorm(n.obs * n.sim, mean = monthly.mean.return, sd = monthly.ret.std.dev)
  monthly.inflation.returns[] <- rnorm(n.obs * n.sim, mean = monthly.inflation, sd = monthly.inf.std.dev)
  for (a in 1:n.obs)
    req[a + 1, ] <- req[a, ] * (1 + monthly.invest.returns[a,] - monthly.inflation.returns[a,]) + monthly.savings
  ending.values <- req[nrow(req),]
  suc <- sum(ending.values>target.savings)/n.sim
  success - suc
}
uniroot( obj2, c(0, 1e6) )
# [1] 7891.187
于 2013-06-26T22:21:07.693 回答