optim {base}R Documentation

汎用的最適化

説明

Nelder–Mead 法、準ニュートン法と共役勾配法アルゴリズムに基づく 汎用的最適化プログラム。オプションとして矩形型制約下での 最適化も行う。

用法

optim(par, fn, gr = NULL,
      method = c("Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN"),
      lower = -Inf, upper = +Inf,
      control = list(), hessian = FALSE), ...)

引数

par 最適化されるべきパラメータの初期値。
fn 最初の引数として最適化されるパラメータベクトルを持つ、最小化されるべき関数。 スカラー値の返さなければならない。
gr グラディエントを返す関数。"Nelder-Mead" と "SANN" 法では不要。 もしこれが NULL であり、しかも必要とされるならば、グラディエントの 有限差分近似が使用される。 fn への呼び出しが行われると、同じ引数で gr が直ちに呼び出される ことが保証されている。
method 使われる方法。詳細 を見よ。
lower, upper "L-BFGS-B" 法に対する引数の限界。
control 制御パラメータのリスト。詳細 を見よ。
hessian 論理値。数値微分で求めたヘッセ行列を返すべきか？
... fn と gr に渡されるその他の引数。

詳細

既定ではこの関数は最小化を実行する。しかし、 もし control$fnscale が負ならば最大化を行う。

既定の方法は Nelder と Mead (1965) の方法の移植であり、関数値だけを 使い頑健であるが、相対的に遅い。微分不可能な関数に対してもそれなりに 使える。

"BFGS" 法は準ニュートン法 (別名 variable metric algorithm)であり、 特に 1970 年に Broyden, Fletcher, Goldfarb そして Shanno により同時に発表された ものである。これは関数値とグラディントを使い、最適化されるべき曲面の 図を作成する。

"CG" 法は Fletcher と Reeves (1964) に基づく共役勾配法である （但し、Polak–Ribiere もしくは Beale–Sorenson による更新法をオプションで 選べる）。共役勾配法は BFGS 法より一般に破綻しやすいが、行列を保持する必要が 無いため、かなり大規模な最適化問題でも成功することがある。

"L-BFGS-B" 法は Byrd et. al. (1994) の方法であり、 矩形型の拘束条件 、つまり各変数が上限もしくは下限を持つ、を許す。 初期値はこの拘束条件を満たす必要がある。 これは BFGS 準ニュートン法のメモリー節約型変形を用いる。 もし自明でない限界が与えられると、警告とともにこの方法が選択される。

"SANN" 法は Belisle (1992) によるシミュレーティッド・アニーリング法である。 シミュレーティッド・アニーリング法は確率的な大局的最適化法に属する。 これは関数値のみを使うが、かなり遅い。微分不可能な関数に対しても使える。 この移植は採択確率として Metropolis 関数を使う。次の候補点は現在の 温度パラメータに比例するスケールを持つガウス型マルコフ推移確率により 生成される。温度パラメータは Belisle (1992, p. 890) で与えられた対数的 冷却規則で減少する。"SANN" 法は制御パラメータの取り方に 本質的に依存することを注意しよう。これは汎用的な手法ではないが、非常に 荒い曲面上の望ましい点を得るのに非常に役にたつことがある。

関数 fn は、関数が与えられた値で評価不能なときは、NA もしくは Inf を返しても良い。しかし、初期値は fn が 計算可能な有限な値を持つものでなければならない。 （"L-BFGS-B" 法では関数値は常に有限でなければならない。）

optim は再帰的に使うことができ、一つでも複数のパラメータでも 使える。

引数 control はりストであり、次の成分のいずれかを提供する：

trace
論理値。 もし真なら、最適化の途中過程に関する追跡情報が生成される。
fnscale
最適化の間 fn と gr の値に 常に適用される尺度スケール。もし負ならば、問題は最大化問題になる。 最適化は値 fn(par)/fnscale に対して実行される。
parscale
パラメータに対するスケール値のベクトル。 最適化は par/parscale に対して行われ、一つの要素の基本単位分の 変更が、スケール化された値の基本単位分の変更をもたらすという意味で、 相互に比較可能でなければならない。
ndeps
par/parscale 尺度で与えられた、 グラディエントの有限差分近似で使うステップサイズのベクトル。 既定値は 1e-3 。
maxit
繰り返しの最大値。既定値は、微分値を使う 方法では 100 回、"Nelder-Mead" 法では 500 回。 "SANN" では maxit は関数評価の総数を与える。 他の停止条件は無い。既定値は 10000。
abstol
収束の許容絶対値。非負値関数に対してのみ、ゼロへの到達の許容値として役にたつ。
reltol
収束の許容相対値。アルゴリズムは一回のステップで値を reltol * (abs(val) + reltol) 以上減少できなければ 停止する。 既定値は sqrt(.Machine$double.eps) で、典型的には約 1e-8。
alpha, beta, gamma
"Nelder-Mead" 法に対するスケール化パラメータ。 alpha は reflection 因子（既定値は 1.0）、 beta は contraction 因子 (0.5) そして gamma は expansion 因子 (2.0)。
REPORT
control$trace が真のときの "BFGS" 法の レポートの頻度。既定値は 10 回の繰り返し毎。
type
グラディエント共役法に対応。 1 なら Fletcher–Reeves の更新法、2 なら Polak–Ribiere の、そして 3 なら Beale–Sorenson の更新方法。
lmm
整数値で "L-BFGS-B" 法で保持される BFGS 更新の数を与える。 retained in the "L-BFGS-B" method, 既定値は 5。
factr
"L-BFGS-B" 法の収束を制御する。 収束は目的値の減少が機械許容度のこの数値倍以内であるとき起こる。 既定値は約 1e-8。
pgtol
"L-BFGS-B" 法の収束の制御を助ける。 これは現在の探索方向へのグラディエントの射影に対する許容値である。 既定値はゼロで、チェックが行われない。
temp
"SANN" 法に対する制御値。冷却計画の初期温度である。 既定値は 10。
tmax
"SANN" 法における各温度での関数評価回数。 既定値は 10。

返り値

次の成分を持つリスト：
par 得られた最良のパラメータ値。
value par に対する fn の値。
counts 二つの成分を持つ整数値ベクトルで、それぞれ fn と gr の呼び出し回数。 これは、ヘッセ行列を要求したときのヘッセ行列の計算に使われた回数、そして グラディエントを有限差分近似するための fn のすべての呼び出し を除く。
convergence 整数値コード。0 は収束が成功理に終ったことを意味する。 エラーコードは
1
は繰り返しの最大値 maxit に達したことを意味する。
10
は Nelder–Mead 単体が退化したことを意味する。
51
は "L-BFGS-B" 法からの警告。 詳細は成分 message を見よ。
52
は "L-BFGS-B" 法からの警告。 詳細は成分 message を見よ。
message 最適化プログラムからの任意の追加情報を与える文字列、もしくは NULL。
hessian 引数 hessian が真のときのみ。 解が見付かったときのヘッセ行列の推定値を与える対称行列。 これは矩形型制約が有効なときでも、非拘束条件下でのヘッセ行列を 与えることを注意せよ。

注意

"Nelder-Mead" 法、"BFGS"法、そして "CG" 法に対するオリジナルのコードは Nash (1990) によるパスカルコードに 基づき、p2c で変換後手作業で最適化された。Dr Nash はコードが自由に 配布されることに同意している。

"L-BFGS-B" 法に対するコードは Netlib に含まれる Zhu, Byrd, Lu-Chen and Nocedal の フォートランコード(ファイル opt/lbfgs_bcm.shar: 他のバージョンが toms/778 にある)に基づく。

"SANN" 法に対するコードは A. Trapletti が提供した。

参考文献

Belisle, C. J. P. (1992) Convergence theorems for a class of simulated annealing algorithms on Rd. J Applied Probability, 29, 885–895.

Byrd, R. H., Lu, P., Nocedal, J. and Zhu, C. (1995) A limited memory algorithm for bound constrained optimization. SIAM J. Scientific Computing, 16, 1190–1208.

Fletcher, R. and Reeves, C. M. (1964) Function minimization by conjugate gradients. Computer Journal 7, 148–154.

Nash, J. C. (1990) Compact Numerical Methods for Computers. Linear Algebra and Function Minimisation. Adam Hilger.

Nelder, J. A. and Mead, R. (1965) A simplex algorithm for function minimization. Computer Journal 7, 308–313.

関連事項

nlm, optimize

例

fr <- function(x) {   ## Rosenbrock のバナナ関数
    x1 <- x[1]
    x2 <- x[2]
    100 * (x2 - x1 * x1)^2 + (1 - x1)^2
}
grr <- function(x) { ## `fr' のグラディエント
    x1 <- x[1]
    x2 <- x[2]
    c(-400 * x1 * (x2 - x1 * x1) - 2 * (1 - x1),
       200 *      (x2 - x1 * x1))
}
optim(c(-1.2,1), fr)
optim(c(-1.2,1), fr, grr, method = "BFGS")
optim(c(-1.2,1), fr, NULL, method = "BFGS", hessian = TRUE)
optim(c(-1.2,1), fr, grr, method = "CG")
optim(c(-1.2,1), fr, grr, method = "CG", control=list(type=2))
optim(c(-1.2,1), fr, grr, method = "L-BFGS-B")

flb <- function(x)
    { p <- length(x); sum(c(1, rep(4, p-1)) * (x - c(1, x[-p])^2)^2) }
## 25-次元の矩形型制約
optim(rep(3, 25), flb, NULL, "L-BFGS-B",
      lower=rep(2, 25), upper=rep(4, 25)) # par[24] is *not* at boundary

## "wild" な関数、大局的最小値は約 -15.81515
fw <- function (x)
    10*sin(0.3*x)*sin(1.3*x^2) + 0.00001*x^4 + 0.2*x+80
plot(fw, -50, 50, n=1000, main = "optim() minimising `wild function'")

res <- optim(50, fw, method="SANN",
             control=list(maxit=20000, temp=20, parscale=20))
res
## ここで局所的に改良する
(r2 <- optim(res$par, fw, method="BFGS"))
points(r2$par, r2$val, pch = 8, col = "red", cex = 2)
[Package Contents]
R: General-purpose Optimization
optim {base}R Documentation

General-purpose Optimization

Description

General-purpose optimization based on Nelder–Mead, quasi-Newton and conjugate-gradient algorithms. It includes an option for box-constrained optimization.

Usage

optim(par, fn, gr = NULL,
      method = c("Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN"),
      lower = -Inf, upper = +Inf,
      control = list(), hessian = FALSE), ...)

Arguments

par Initial values for the parameters to be optimized over.
fn A function to be minimized, with first argument the vector of parameters over which minimization is to take place. It should return a scalar result.
gr A function to return the gradient. Not needed for the "Nelder-Mead" and "SANN" method. If it is NULL and it is needed, a finite-difference approximation will be used. It is guaranteed that gr will be called immediately after a call to fn at the same parameter values.
method The method to be used. See Details.
lower, upper Bounds on the variables for the "L-BFGS-B" method.
control A list of control parameters. See Details.
hessian Logical. Should a numerically differentiated Hessian matrix be returned?
... Further arguments to be passed to fn and gr.

Details

By default this function performs minimization, but it will maximize if control$fnscale is negative.

The default method is an implementation of that of Nelder and Mead (1965), that uses only function values and is robust but relatively slow. It will work reasonably well for non-differentiable functions.

Method "BFGS" is a quasi-Newton method (also known as a variable metric algorithm), specifically that published simultaneously in 1970 by Broyden, Fletcher, Goldfarb and Shanno. This uses function values and gradients to build up a picture of the surface to be optimized.

Method "CG" is a conjugate gradients method based on that by Fletcher and Reeves (1964) (but with the option of Polak–Ribiere or Beale–Sorenson updates). Conjugate gradient methods will generally be more fragile that the BFGS method, but as they do not store a matrix they may be successful in much larger optimization problems.

Method "L-BFGS-B" is that of Byrd et. al. (1994) which allows box constraints, that is each variable can be given a lower and/or upper bound. The initial value must satisfy the constraints. This uses a limited-memory modification of the BFGS quasi-Newton method. If non-trivial bounds are supplied, this method will be selected, with a warning.

Method "SANN" is a variant of simulated annealing given in Belisle (1992). Simulated-annealing belongs to the class of stochastic global optimization methods. It uses only function values but is relatively slow. It will also work for non-differentiable functions. This implementation uses the Metropolis function for the acceptance probability. The next candidate point is generated from a Gaussian Markov kernel with scale proportional to the actual temperature. Temperatures are decreased according to the logarithmic cooling schedule as given in Belisle (1992, p. 890). Note that the "SANN" method depends critically on the settings of the control parameters. It is not a general-purpose method but can be very useful in getting to a good value on a very rough surface.

Function fn can return NA or Inf if the function cannot be evaluated at the supplied value, but the initial value must have a computable finite value of fn. (Except for method "L-BFGS-B" where the values should always be finite.)

optim can be used recursively, and for a single parameter as well as many.

The control argument is a list that can supply any of the following components:

trace
Logical. If true, tracing information on the progress of the optimization is produced.
fnscale
An overall scaling to be applied to the value of fn and gr during optimization. If negative, turns the problem into a maximization problem. Optimization is performed on fn(par)/fnscale.
parscale
A vector of scaling values for the parameters. Optimization is performed on par/parscale and these should be comparable in the sense that a unit change in any element produces about a unit change in the scaled value.
ndeps
A vector of step sizes for the finite-difference approximation to the gradient, on par/parscale scale. Defaults to 1e-3.
maxit
The maximum number of iterations. Defaults to 100 for the derivative-based methods, and 500 for "Nelder-Mead". For "SANN" maxit gives the total number of function evaluations. There is no other stopping criteria. Defaults to 10000.
abstol
The absolute convergence tolerance. Only useful for non-negative functions, as a tolerance for reaching zero.
reltol
Relative convergence tolerance. The algorithm stops if it is unable to reduce the value by a factor of reltol * (abs(val) + reltol) at a step. Defaults to sqrt(.Machine$double.eps), typically about 1e-8.
alpha, beta, gamma
Scaling parameters for the "Nelder-Mead" method. alpha is the reflection factor (default 1.0), beta the contraction factor (0.5) and gamma the expansion factor (2.0).
REPORT
The frequency of reports for the "BFGS" method in control$trace is positive. Defaults to every 10 iterations.
type
for the conjugate-gradients method. Takes value 1 for the Fletcher–Reeves update, 2 for Polak–Ribiere and 3 for Beale–Sorenson.
lmm
is an integer giving the number of BFGS updates retained in the "L-BFGS-B" method, It defaults to 5.
factr
controls the convergence of the "L-BFGS-B" method. Convergence occurs when the reduction in the objective is within this factor of the machine tolerance. Default is 1e7, that is a tolerance of about 1e-8.
pgtol
helps controls the convergence of the "L-BFGS-B" method. It is a tolerance on the projected gradient in the current search direction. This defaults to zero, when the check is suppressed.
temp
controls the "SANN" method. It is the starting temperature for the cooling schedule. Defaults to 10.
tmax
is the number of function evaluations at each temperature for the "SANN" method. Defaults to 10.

Value

A list with components:
par The best set of parameters found.
value The value of fn corresponding to par.
counts A two-element integer vector giving the number of calls to fn and gr respectively. This excludes those calls needed to compute the Hessian, if requested, and any calls to fn to compute a finite-difference approximation to the gradient.
convergence An integer code. 0 indicates successful convergence. Error codes are
1
indicates that the iteration limit maxit had been reached.
10
indicates degeneracy of the Nelder–Mead simplex.
51
indicates a warning from the "L-BFGS-B" method; see component message for further details.
52
indicates an error from the "L-BFGS-B" method; see component message for further details.
message A character string giving any additional information returned by the optimizer, or NULL.
hessian Only if argument hessian is true. A symmetric matrix giving an estimate of the Hessian at the solution found. Note that this is the Hessian of the unconstrained problem even if the box constraints are active.

Note

The code for methods "Nelder-Mead", "BFGS" and "CG" was based originally on Pascal code in Nash (1990) that was translated by p2c and then hand-optimized. Dr Nash has agreed that the code can be made freely available.

The code for method "L-BFGS-B" is based on Fortran code by Zhu, Byrd, Lu-Chen and Nocedal obtained from Netlib (file opt/lbfgs_bcm.shar: another version is in toms/778).

The code for method "SANN" was contributed by A. Trapletti.

References

Belisle, C. J. P. (1992) Convergence theorems for a class of simulated annealing algorithms on Rd. J Applied Probability, 29, 885–895.

Byrd, R. H., Lu, P., Nocedal, J. and Zhu, C. (1995) A limited memory algorithm for bound constrained optimization. SIAM J. Scientific Computing, 16, 1190–1208.

Fletcher, R. and Reeves, C. M. (1964) Function minimization by conjugate gradients. Computer Journal 7, 148–154.

Nash, J. C. (1990) Compact Numerical Methods for Computers. Linear Algebra and Function Minimisation. Adam Hilger.

Nelder, J. A. and Mead, R. (1965) A simplex algorithm for function minimization. Computer Journal 7, 308–313.

See Also

nlm, optimize

Examples

fr <- function(x) {   ## Rosenbrock Banana function
    x1 <- x[1]
    x2 <- x[2]
    100 * (x2 - x1 * x1)^2 + (1 - x1)^2
}
grr <- function(x) { ## Gradient of `fr'
    x1 <- x[1]
    x2 <- x[2]
    c(-400 * x1 * (x2 - x1 * x1) - 2 * (1 - x1),
       200 *      (x2 - x1 * x1))
}
optim(c(-1.2,1), fr)
optim(c(-1.2,1), fr, grr, method = "BFGS")
optim(c(-1.2,1), fr, NULL, method = "BFGS", hessian = TRUE)
optim(c(-1.2,1), fr, grr, method = "CG")
optim(c(-1.2,1), fr, grr, method = "CG", control=list(type=2))
optim(c(-1.2,1), fr, grr, method = "L-BFGS-B")

flb <- function(x)
    { p <- length(x); sum(c(1, rep(4, p-1)) * (x - c(1, x[-p])^2)^2) }
## 25-dimensional box constrained
optim(rep(3, 25), flb, NULL, "L-BFGS-B",
      lower=rep(2, 25), upper=rep(4, 25)) # par[24] is *not* at boundary

## "wild" function , global minimum at about -15.81515
fw <- function (x)
    10*sin(0.3*x)*sin(1.3*x^2) + 0.00001*x^4 + 0.2*x+80
plot(fw, -50, 50, n=1000, main = "optim() minimising `wild function'")

res <- optim(50, fw, method="SANN",
             control=list(maxit=20000, temp=20, parscale=20))
res
## Now improve locally
(r2 <- optim(res$par, fw, method="BFGS"))
points(r2$par, r2$val, pch = 8, col = "red", cex = 2)
[Package Contents]