R中的哪个优化算法/程序包可以解决此问题？

问题描述

我正在尝试求解随机系数logit模型以估计需求参数。我的教授建议在Matlab中执行此操作，但也可以使用其他语言，因此我选择R是因为我对此更加熟悉。

我无需过多介绍细节，我想我需要为这段代码和数据添加一些背景知识。基本思想是，数据集中产品的市场份额取决于对所有消费者都相同的均值效用（取决于产品特征），以及与该均值的偏差（取决于从商品中得出的消费者特征）。正态分布）。然后将市场份额大致计算为

marketshare = exp（mean + dev实用程序）/ [1 +该市场中所有产品exp（mean + dev实用程序）的总和。]

然后的问题是找到最适合数据的产品特征和消费者特征的系数（还有其他一些计算，但这在这里并不重要）。

数据的每一行对应一个产品，并且产品被分组到彼此不交互的市场中。

下面的代码解决了他的问题。您可以下载data.Rdata here。如果要运行它，可以通过将tol_inner增加到1E-8之类的方法来显着减少所需的时间，并且仍然可以保持合理的距离。

library(nloptr)
library(tictoc)

# Import data
load(".../data.Rdata")

# Simulated market share function
shares <- function(meanutil,devutil) {
  num <- exp(devutil + meanutil)
  mktsum <- rowsum(num,group = mktids)
  denom <- 1 + unname(mktsum[rep(1:nrow(mktsum),times = mktprods),])
  s_jm <- rowMeans(num/denom)
  return(s_jm)
}

# Inner loop settings/initialization
tol_inner <- 1E-14
max_iter <- 10000
delta_guess <- rep(0,nrow(x)) # Guess for mean utilities for first outer loop iteration

# Objective function
obj_fun <- function(theta){
  theta1 <- theta[1:Kbeta]
  theta2 <- theta[(Kbeta+1):(Kbeta+Ktheta)]
  
  mu <- x[,-1] %*% (theta2 * v) # Deviations from mean utility
  delta_h <- delta_guess        # Start inner loop with guess or previous converged delta
  
  iter <- 1                     # Initialize counter
  max_diff <- 1
  
  # Inner loop
  while (max_diff > tol_inner & iter < max_iter){
    delta_h1 <- delta_h + logdatashare - log(shares(meanutil = delta_h,devutil = mu)) # Update mean utility
    max_diff <- max(abs(delta_h1 - delta_h))                                           # Check for difference with previous iteration
    delta_h <- delta_h1                                                                # Update delta for next iteration
    iter  = iter + 1                                                                   # Increment counter
  }
  
  cat("Inner loop iterations: ",iter-1)
  
  delta_guess <<- delta_h      # Store converged mean utility to global env. for next outer loop iteration
  xi <- delta_h - x %*% theta1 # Compute unobserved prod. characteristics
  g <- t(IV) %*% xi            # Compute moments
  f <- t(g) %*% W %*% g        # Objective function value
  return(f)
}

# Outer loop optimization
tic()
res <- nloptr(
  x0 = rep(1,9),eval_f = obj_fun,lb = c(-Inf,-Inf,0),ub = c(Inf,Inf,Inf),opts = list("algorithm" = "NLOPT_LN_BOBYQA","xtol_rel" = 1E-10,"ftol_rel" = 1E-10,"print_level" = 3,"maxeval" = 0)
)
toc()

# Some other optimization algorithms I tried: NLOPT_LN_COBYLA,NLOPT_LN_BOBYQA,NLOPT_LN_SBPLX,NLOPT_LN_NELDERMEAD,NLOPT_LN_NEWUOA_BOUND,NLOPT_GN_CRS2_LM,NLOPT_GN_MLSL_LDS,NLOPT_GN_DIRECT_L


# This should be close to the true parameters:
solution <- c(3.19298639055865,3.3072606294676,0.546084673417042,0.576567495578966,2.02309552738769,1.0427987901479,0.622031348646723,0.52060518402574,0.605067406591273)
obj_fun(solution)

但是，我有两个问题：

此代码在Intel Xeon E5-1620 v2 @ 3.5GHz上运行大约需要20分钟，而我的同学们可以在不到4分钟的时间里在笨拙的大学计算机上运行等效的Matlab。我认为我已经通过消除所有不必要的循环，应用调用和向量化所有东西来优化了一切。如果可能的话，我想进一步加快速度，但是我认为这取决于优化算法。我还注意到，运行此代码时，我的CPU利用率仅为20％。
结果似乎非常取决于算法和起始值。我有95％的把握从理论上讲应该只有一个最小值，但是当从全零或runif(9)开始时，到目前为止，我尝试过的所有算法都截断了目标函数的错误但看似平坦的部分。我的Matlab朋友没有这个问题。

我尝试了许多不同的算法和软件包，但是到目前为止，nloptr和NLOPT_LN_BOBYQA似乎是最好的。是否有其他算法或软件包更适合此类问题？我可以以某种方式使用并行化来加快速度吗？ Matlab如何更快地解决这个问题？

解决方法

暂无找到可以解决该程序问题的有效方法，小编努力寻找整理中！

如果你已经找到好的解决方法，欢迎将解决方案带上本链接一起发送给小编。

小编邮箱:dio#foxmail.com (将#修改为@）

nonlinear-optimization r

R中的哪个优化算法/程序包可以解决此问题？

问题描述

解决方法

相关问答