我正在尝试使用R来估计一个具有手动规范的多项logit模型。我找到了一些包,可以用来估计MNL模型here或here。
我找到了其他一些关于“滚动”你自己的MLE函数here的文章。然而,根据我的深入研究,所有这些函数和包都依赖于内部的optim
函数。
在我的基准测试中,optim
是瓶颈。使用一个包含大约16000个观察值和7个参数的模拟数据集,R在我的机器上大约需要90秒。Biogeme中的等效模型大约需要10秒。一位用Ox编写自己的代码的同事报告说,对于相同的模型,大约有4秒。
有没有人有编写自己的最大似然函数的经验,或者可以给我指出一些在默认optim
函数之外进行优化的方向(没有双关语)?
如果任何人想要R代码重新创建模型,请让我知道-我很乐意提供它。我没有提供它,因为它与优化optim
函数和节省空间的问题没有直接关系……
编辑:感谢大家的想法。基于下面无数的评论,对于更复杂的模型,我们能够获得与Biogeme相同的R,并且对于我们运行的几个较小/更简单的模型,R实际上更快。我认为这个问题的长期解决方案将涉及到编写一个独立的最大化函数,该函数依赖于fortran或C库,但当然也可以采用其他方法。
发布于 2010-09-21 19:21:59
已经尝试过nlm()函数了吗?不知道它是否快很多,但它确实提高了速度。还要检查选项。optim默认使用较慢的算法。通过使用准牛顿算法(method="BFGS")而不是默认算法,可以获得大于5倍的加速比。如果您不太关心最后一位数,您还可以将nlm()的容差级别设置得更高,以获得额外的速度。
f <- function(x) sum((x-1:length(x))^2)
a <- 1:5
system.time(replicate(500,
optim(a,f)
))
user system elapsed
0.78 0.00 0.79
system.time(replicate(500,
optim(a,f,method="BFGS")
))
user system elapsed
0.11 0.00 0.11
system.time(replicate(500,
nlm(f,a)
))
user system elapsed
0.10 0.00 0.09
system.time(replicate(500,
nlm(f,a,steptol=1e-4,gradtol=1e-4)
))
user system elapsed
0.03 0.00 0.03
发布于 2010-09-21 20:17:02
你考虑过CRAN Task View for Optimization上的材料了吗?
发布于 2010-09-21 20:46:11
顺便说一句,我已经用C-ish用OPTIF9完成了这项工作。你很难再跑得更快了。有很多方法可以让事情变得更慢,比如运行R这样的解释器。
补充:从评论中可以清楚地看到,OPTIF9被用作优化引擎。这意味着很可能大部分时间都花在计算R中的目标函数上。虽然可能在底层使用C函数来执行一些操作,但仍然存在解释器开销。有一种快速的方法可以确定R中的哪些代码行和函数调用负责大部分时间,也就是用key键暂停并检查堆栈。如果一条语句花费了X%的时间,那么它在堆栈上的时间就是X%。您可能会发现,有些操作不会转到C,而应该转到C。当你找到一种方法来并行化R执行时,你以这种方式获得的任何加速因子都将被保留下来。
https://stackoverflow.com/questions/3757321
复制相似问题