#LyX 1.3 created this file. For more info see http://www.lyx.org/
\lyxformat 221
\textclass article
\begin_preamble
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\newtheorem{theorem}{Theorem}
\newtheorem{acknowledgement}[theorem]{Acknowledgement}
\newtheorem{axiom}[theorem]{Axiom}
\newtheorem{case}[theorem]{Case}
\newtheorem{claim}[theorem]{Claim}
\newtheorem{conclusion}[theorem]{Conclusion}
\newtheorem{condition}[theorem]{Condition}
\newtheorem{conjecture}[theorem]{Conjecture}
\newtheorem{corollary}[theorem]{Corollary}
\newtheorem{criterion}[theorem]{Criterion}
\newtheorem{definition}[theorem]{Definition}
\newtheorem{example}[theorem]{Example}
\newtheorem{exercise}[theorem]{Exercise}
\newtheorem{lemma}[theorem]{Lemma}
\newtheorem{notation}[theorem]{Notation}
\newtheorem{problem}[theorem]{Problem}
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{remark}[theorem]{Remark}
\newtheorem{solution}[theorem]{Solution}
\newtheorem{summary}[theorem]{Summary}
\usepackage{graphicx}
\usepackage{amsmath}
\usepackage{html}
\usepackage{dcolumn}
\usepackage{mathpazo}
\usepackage[ps2pdf,pdftitle={Econometrics},urlcolor=blue,linktocpage,a4paper,colorlinks=true]{hyperref}
\end_preamble
\language english
\inputencoding default
\fontscheme palatino
\graphics default
\paperfontsize 11
\spacing onehalf
\papersize a4paper
\paperpackage a4wide
\use_geometry 0
\use_amsmath 1
\use_natbib 0
\use_numerical_citations 0
\paperorientation portrait
\secnumdepth 3
\tocdepth 3
\paragraph_separation indent
\defskip medskip
\quotes_language swedish
\quotes_times 2
\papercolumns 1
\papersides 1
\paperpagestyle default
\layout Title
Unconstrained Optimization with MINTOOLKIT for GNU Octave
\layout Author
Michael Creel
\begin_inset Foot
collapsed false
\layout Standard
Dept.
of Economics and Economic History, Universitat Autònoma de Barcelona.
michael.creel@uab.es
\end_inset
\layout Date
March 2004
\layout Abstract
The paper documents MINTOOLKIT for GNU Octave.
MINTOOLKIT provides functions for minimization and numeric differentiation.
The main algorithms are BFGS, LBFGS, and simulated annealing.
Examples are given.
\layout Section
Introduction
\layout Standard
Unconstrained optimization is a basic tool in many disciplines, and there
is no need to discuss its importance.
This paper discusses how unconstrained optimization may be done using GNU
Octave (Eaton,
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{www.octave.org}{http://www.octave.org}
\end_inset
) using the package MINTOOLKIT (Creel,
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{http://pareto.uab.es/mcreel/MINTOOLKIT}{http://pareto.uab.es/mcreel/MINTOOLKIT}
\end_inset
).
If you would just like to see some examples of how to use the algorithms,
skip to section
\begin_inset LatexCommand \ref{sec:Examples}
\end_inset
.
Otherwise, here's some introductory information that explains how algorithms
we selected for inclusion into MINTOOLKIT.
\layout Subsection
Types of problems
\layout Standard
We first briefly discuss types of optimization problems, in order to identify
the cases where GNU Octave will be a good platform for analysis, and which
of the algorithms in MINTOOLKIT will likely work well for given cases.
\layout Subsubsection
Large/small
\layout Standard
What is the number of parameters to be minimized? Let
\begin_inset Formula $k$
\end_inset
be the number of parameters.
Memory usage of an algorithm will be some function
\begin_inset Formula $f(k).$
\end_inset
If this function is growing rapidly for a certain algorithm, that algorithm
will cease to be useful for
\begin_inset Quotes sld
\end_inset
large
\begin_inset Quotes srd
\end_inset
problems, since memory will be exhausted.
Of course,
\begin_inset Quotes sld
\end_inset
large
\begin_inset Quotes srd
\end_inset
is a relative term that increases over time as memory becomes cheaper.
For
\begin_inset Quotes sld
\end_inset
small
\begin_inset Quotes srd
\end_inset
problems, the speed of convergence of the algorithm will be of primary
importance, since memory resources will not be a bottleneck.
\layout Subsubsection
Continuous/discontinuous
\layout Standard
Gradient-based methods such as Newton's algorithm or quasi-Newton methods
rely on the function being differentiable.
This will not hold if the function is not continuous.
Search-type algorithms will be appropriate here.
\layout Subsubsection
Convex/nonconvex
\layout Standard
A convex objective function will have a single global minimizer, whereas
nonconvex functions may have additional local minima.
Quasi-Newton methods only use local information in their updates, so they
may well converge to a non-global minimum, depending upon starting values.
A possible solution is to try a number of starting values.
This is likely to work well if the nonconvexity problem is not too severe.
When there are
\emph on
many
\emph default
local minima, a search-type algorithm may become more efficient, since the
problem of local minima is dealt with automatically and doesn't require
the analysts' intervention.
After all, who's time is more important, yours, or your computer's?
\layout Subsubsection
Costly/cheap
\layout Standard
Can the objective function be evaluated quickly, or is it time-consuming?
Octave is an interpreted language, and is in general slower than FORTRAN,
C, or similar.
So expensive objective functions are best not implemented in pure Octave.
One might go to a different environment for analysis, but in fact it is
relatively easy to convert objective functions written in Octave to C++,
and call them dynamically from Octave scripts.
In this way, the expensive calculations are done using a fast language,
while the user deals with the convenient, friendly Octave environment.
In fact, C++ functions may make use of the Octave classes, so converting
an Octave function to C++ is not very difficult.
The algorithms in MINTOOLKIT serve as examples of how this may be done.
\layout Section
Algorithms
\layout Standard
The algorithms in MINTOOLKIT were chosen based upon many sources of information,
two of which are
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{Nocedal (1992)}{http://www.ece.northwestern.edu/~nocedal/PDFfiles/acta.pdf}
\end_inset
, for continuous, convex problems, and
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{Mittelmann}{http://plato.la.asu.edu/topics/problems/global.html}
\end_inset
for global minimization.
The goal of MINTOOLKIT is to be able to solve well-posed problems quickly
and robustly, using the smallest set of algorithms possible.
\begin_inset Quotes sld
\end_inset
Well-posed
\begin_inset Quotes srd
\end_inset
is an important adjective here - MINTOOLKIT is not meant to be able to
solve poorly conditioned problems or problems that can easily fall into
numeric precision traps.
Please pay attention to how your data is scaled, try to
\begin_inset Quotes sld
\end_inset
bullet-proof
\begin_inset Quotes srd
\end_inset
your objective function to avoid divisions by zero, etc.
Nevertheless, if you find bugs, or have suggestions or comments, please
contact the author.
\layout Subsection
BFGSMIN
\layout Standard
The BFGS algorithm is probably the most widely-used quasi-Newton method
for moderately-sized continuous problems that are not extremely nonconvex.
It is more robust than other quasi-Newton methods such as DFP (Nocedal,
1992), and it is faster than Newton's method, since the Hessian matrix
need not be calculated.
One could easily create a Newton algorithm using the source code for
\family typewriter
bfgsmin
\family default
.
\layout Subsection
LBFGSMIN
\layout Standard
For large problems, the BFGS algorithm may not be feasible, since it requires
storing a
\begin_inset Formula $k\times k$
\end_inset
matrix.
The LBFGS (
\begin_inset Quotes sld
\end_inset
L
\begin_inset Quotes srd
\end_inset
is for
\begin_inset Quotes sld
\end_inset
limited memory
\begin_inset Quotes srd
\end_inset
) method is able to store all information for updates in vectors, which
substantially reduced memory requirements.
It may also be faster than the BFGS algorithm in some circumstances, since
it may use fewer floating point operations per iteration, depending upon
the size of the problem.
While it will usually require more iterations that BFGS, it may be faster
if each iteration is faster.
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{Liu and Nocecal
\backslash
(1989
\backslash
)}{http://www.ece.northwestern.edu/~nocedal/PDFfiles/limited-memory.pdf}
\end_inset
is a reference.
\layout Subsection
SAMIN
\layout Standard
When problems are mildly nonconvex, the quasi-Newton methods above, combined
with a number of trial start values, may be the fastest way to find the
global minimum.
This solution is implemented in
\family typewriter
battery.m.
\family default
But when the problem becomes less smooth, with many local minima, this solution
may fail.
Simulated annealing is one of the algorithms that works well in this case.
The implementation by
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{Goffe}{http://www.netlib.no/netlib/opt/simann.f}
\end_inset
has been used widely, and is the basis for the version included in MINTOOLKIT.
An additional reference is Goffe (1996).
\layout Standard
\family typewriter
samin
\family default
differs from the Goffe code in two important ways:
\layout Enumerate
The
\begin_inset Quotes sld
\end_inset
temperature
\begin_inset Quotes srd
\end_inset
is determined automatically.
In a first stage, the temperature is increased until the active search
region covers the entire parameter space defined as the
\begin_inset Formula $k$
\end_inset
-dimensional rectangle
\begin_inset Formula $\times_{j=1}^{k}(lb_{j},ub_{j})$
\end_inset
where
\begin_inset Formula $lb_{j}$
\end_inset
and
\begin_inset Formula $ub_{j}$
\end_inset
are the
\begin_inset Formula $j$
\end_inset
th elements of the LB and UB vectors that are the lower and upper bounds
for the parameters (these are user-specified arguments to samin.
Once this is achieved, the temperature decreases as usual.
\layout Enumerate
Convergence is defined as
\emph on
two
\emph default
conditions holding simultaneously.
\begin_deeper
\layout Enumerate
The last NEPS best function values cannot differ by more than FUNCTOL.
This is as in Goffe's code.
\layout Enumerate
The width of the search interval must be less than PARAMTOL for each parameter.
This allows to avoid accepting points on a flat plateau.
\end_deeper
\layout Section
Obtaining the code
\layout Itemize
MINTOOLKIT is available directly from the author, at
\begin_inset ERT
status Open
\layout Standard
\backslash
htmladdnormallink{http://pareto.uab.es/mcreel/MINTOOLKIT}{http://pareto.uab.es/mcreel/MINTOOLKIT}
\end_inset
.
If you get it this way, uncompress the file where you like, change to the
MINTOOLKIT directory, and compile by typing
\begin_inset Quotes sld
\end_inset
make all
\begin_inset Quotes srd
\end_inset
(this supposes that you have Octave installed already).
Make sure that Octave knows where the MINTOOLKIT directory is.
This option guarantees that you have the most recent version.
\layout Itemize
Otherwise, you can obtain MINTOOLKIT as part of the
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{octave-forge}{http://sourceforge.net/project/showfiles.php?group_id=2888}
\end_inset
package.
\layout Itemize
If you happen to be running
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{Debian}{http://www.debian.org}
\end_inset
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{Linux}{http://www.linux.org}
\end_inset
, you can install a pre-compiled version of octave-forge and all required
files by typing
\begin_inset Quotes sld
\end_inset
apt-get install octave-forge
\begin_inset Quotes srd
\end_inset
.
This is the easiest (and recommended) option.
\layout Section
\begin_inset LatexCommand \label{sec:Examples}
\end_inset
Examples
\layout Standard
\added_space_bottom medskip
MINTOOLKIT contains some functions for use by users, and some other functions
that users can ignore.
The functions for users are
\layout Standard
\added_space_bottom medskip
\begin_inset Tabular
|
\begin_inset Text
\layout Standard
Function
\end_inset
|
\begin_inset Text
\layout Standard
Purpose
\end_inset
|
|
\begin_inset Text
\layout Standard
bfgsmin
\end_inset
|
\begin_inset Text
\layout Standard
Ordinary BFGS algorithm
\end_inset
|
|
\begin_inset Text
\layout Standard
battery
\end_inset
|
\begin_inset Text
\layout Standard
Calls bfgsmin with a set of starting values
\end_inset
|
|
\begin_inset Text
\layout Standard
lbfgsmin
\end_inset
|
\begin_inset Text
\layout Standard
Limited-memory BFGS, for large problems
\end_inset
|
|
\begin_inset Text
\layout Standard
samin
\end_inset
|
\begin_inset Text
\layout Standard
Simulated annealing, for global minimization
\end_inset
|
|
\begin_inset Text
\layout Standard
numgradient
\end_inset
|
\begin_inset Text
\layout Standard
numeric first derivative of vector-valued function
\end_inset
|
|
\begin_inset Text
\layout Standard
numhessian
\end_inset
|
\begin_inset Text
\layout Standard
numeric second derivative matrix
\end_inset
|
\end_inset
\layout Standard
This section gives some very simple examples of the use of the algorithms
and functions in MINTOOLKIT.
The first examples are intended to clearly illustrate how to use the algorithms.
Realism is not important.
Then some more difficult problems are considered.
\layout Standard
The functions in MINTOOLKIT allow minimization or differentiation with respect
to any of the arguments of a function, holding the other arguments fixed.
The other arguments can include data or fixed parameters of the function,
for example.
The argument with respect to which minimization or differentiation is done
is denoted by
\begin_inset Formula $minarg$
\end_inset
, which by default is equal to 1.
Any function to be minimized or differentiated by algorithms in MINTOOLKIT
must follow one of the forms
\layout Standard
\begin_inset Formula \begin{eqnarray*}
value & = & \, f(arg_{1},\, arg_{2},\,...,\, arg_{p})\\
{}[value,\, return_{2},\,...,\, return_{n}] & = & \, f(arg_{1},\, arg_{2},\,...,\, arg_{p})\end{eqnarray*}
\end_inset
\emph on
Special case
\emph default
: If the second form is used
\emph on
and
\emph default
\begin_inset Formula $return_{2}$
\end_inset
is a
\begin_inset Formula $k\times1$
\end_inset
vector, where
\begin_inset Formula $k$
\end_inset
is the dimension of
\begin_inset Formula $minarg$
\end_inset
, then is assumed to be the gradient of
\begin_inset Formula $f$
\end_inset
with respect to
\begin_inset Formula $minarg$
\end_inset
, if the algorithm called uses the gradient.
Otherwise, it (and any other returns from
\begin_inset Formula $f$
\end_inset
) are ignored by MINTOOLKIT.
\layout Subsection
Minimization
\layout Subsubsection
\family typewriter
bfgsmin
\layout Standard
\family typewriter
bfgsmin
\family default
is called as
\begin_inset Formula \[
[theta,\, value,\, convergence]=bfgsmin("f",\,\{ args\},\,\{ control\})\]
\end_inset
The first argument
\begin_inset Formula $"f"$
\end_inset
is a string variable that holds the name of the function to be minimized.
The second argument,
\begin_inset Formula $args$
\end_inset
, is a cell array that hold the arguments of
\begin_inset Formula $f$
\end_inset
.
The third argument
\begin_inset Formula $control$
\end_inset
is an optional cell array of
\begin_inset Formula $4$
\end_inset
elements.
The elements of
\begin_inset Formula $control$
\end_inset
are described in Table
\begin_inset LatexCommand \ref{cap:Controls-for-bfgsmin}
\end_inset
.
\begin_inset Float table
wide false
collapsed false
\layout Caption
\begin_inset LatexCommand \label{cap:Controls-for-bfgsmin}
\end_inset
Controls for
\family typewriter
bfgsmin
\layout Standard
\begin_inset Tabular
|
\begin_inset Text
\layout Standard
\size footnotesize
Element
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
Purpose
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
Default Value
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
Other possible values
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
1
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
maxiters
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
-1 (infinity)
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
any positive integer
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
2
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
verbosity
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
0
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
1: summary every iteration; 2: only final summary
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
3
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
criterion
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
1: strict convergence (f, g,
\begin_inset Formula $\Delta p$
\end_inset
)
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
2: weak convergence (only f)
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
4
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
minarg
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
1: first argument
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
int:
\begin_inset Formula $1\le minarg\le k$
\end_inset
,
\begin_inset Formula $k=\# args$
\end_inset
\end_inset
|
\end_inset
\end_inset
The outputs of
\family typewriter
bfgsmin
\family default
are obvious, except the code values that
\begin_inset Formula $convergence$
\end_inset
can take on.
These are -1 for no convergence, maxiters exceeded; 1: convergence according
to the specified strong or weak criterion; 2: no convergence due to failure
of the algorithm (
\emph on
e.g.,
\emph default
the gradient calculation fails, or a stepsize cannot be found).
\layout Standard
Consider a simple example - minimizing a quadratic function.
The program
\begin_inset ERT
status Open
\layout Standard
\backslash
htmladdnormallink{bfgsmin-example.m}{http://pareto.uab.es/mcreel/MINTOOLKIT/bfgsmin-example.m}
\end_inset
follows:
\begin_inset Include \verbatiminput{bfgsmin-example.m}
preview false
\end_inset
\layout Itemize
The first example uses numeric derivatives, and minimizes with respect to
\begin_inset Formula $x,$
\end_inset
the first argument of the objective function.
The second argument,
\begin_inset Formula $y,$
\end_inset
is treated as fixed.
\layout Itemize
The second example uses analytic derivatives, since it calls objective2,
and minimizes with respect to
\begin_inset Formula $x,$
\end_inset
the first argument of the objective function.
The second argument,
\begin_inset Formula $y,$
\end_inset
is treated as fixed.
\layout Itemize
The third example uses numeric derivatives, and minimizes with respect to
\begin_inset Formula $y,$
\end_inset
the second argument of the objective function, since the 4th element of
\family typewriter
control
\family default
,
\family typewriter
minarg
\family default
, is 2 .
The first argument,
\begin_inset Formula $x,$
\end_inset
is treated as fixed.
\layout Standard
The output of running this example is
\layout Standard
\begin_inset Include \verbatiminput{bfsg-example.out}
preview false
\end_inset
Notice that analytic gradients lead to faster convergence that do numeric
gradients.
Also note in the third example, where
\family typewriter
minarg=2
\family default
, that minimization can be with respect to any of the arguments of the objective
function.
\layout Subsubsection
\family typewriter
lbfgsmin
\layout Standard
When the problem is very large, a limited-memory bfgs algorithm may be needed,
if
\family typewriter
bfgsmin
\family default
is not feasible due to memory limitations.
\family typewriter
lbfgsmin
\family default
is called as
\begin_inset Formula \[
[theta,\, value,\, convergence]=lbfgsmin("f",\,\{ args\},\,\{ control\})\]
\end_inset
The first argument
\begin_inset Formula $"f"$
\end_inset
is a string variable that holds the name of the function to be minimized.
The second argument,
\begin_inset Formula $args$
\end_inset
, is a cell array that hold the arguments of
\begin_inset Formula $f$
\end_inset
.
The third argument
\begin_inset Formula $control$
\end_inset
is an optional cell array of
\begin_inset Formula $5$
\end_inset
elements.
The elements of
\begin_inset Formula $control$
\end_inset
are the same as for
\family typewriter
bfgsmin,
\family default
except that there is one more element that controls how many iterations
are used to form the quasi-Hessian matrix (this is the
\emph on
memory
\emph default
of the method).
The control vector is fully described in Table
\begin_inset LatexCommand \ref{cap:Controls-for-lbfgsmin}
\end_inset
.
You can easily modify the above example to use the
\family typewriter
lbfgsmin
\family default
method.
\begin_inset Float table
wide false
collapsed false
\layout Caption
\begin_inset LatexCommand \label{cap:Controls-for-lbfgsmin}
\end_inset
Controls for
\family typewriter
lbfgsmin
\layout Standard
\begin_inset Tabular
|
\begin_inset Text
\layout Standard
\size footnotesize
Element
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
Purpose
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
Default Value
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
Other possible values
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
1
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
maxiters
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
-1 (infinity)
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
any positive integer
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
2
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
verbosity
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
0
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
1: summary every iteration; 2: only final summary
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
3
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
criterion
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
1: strict convergence (f, g,
\begin_inset Formula $\Delta p$
\end_inset
)
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
2: weak convergence (only f)
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
4
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
minarg
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
1: first argument
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
int:
\begin_inset Formula $1\le minarg\le k$
\end_inset
,
\begin_inset Formula $k=\# args$
\end_inset
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
5
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
memory
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
5
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
any positive integer
\end_inset
|
\end_inset
\end_inset
\layout Standard
It is possible that
\family typewriter
lbfgsmin
\family default
can outperform
\family typewriter
bfgsmin
\family default
even when memory is not an issue.
Remember that both of these algorithms are
\emph on
approximating
\emph default
the Hessian matrix using previous gradient evaluations.
If the true Hessian is changing rapidly, then a limited memory approximation
may be better than a long memory approximation.
The Rosenbrock function is such a case.
The program
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{lbfgsmin-example.m}{http://pareto.uab.es/mcreel/MINTOOLKIT/lbfgsmin-example.m}
\end_inset
minimizes a 200-dimensional Rosenbrock function using both algorithms.
The output
\begin_inset Include \verbatiminput{lbfgsmin-example.out}
preview false
\end_inset
shows that the limited memory algorithm uses significantly more iterations
that the ordinary BFGS algorithm, but it is almost 4 times as fast.
In general, though, the ordinary BFGS algorithm is recommended when memory
limitations are not a problem.
\layout Subsubsection
\family typewriter
samin
\layout Standard
For discontinuous and/or seriously nonconvex problems, the quasi-Newton
methods are not likely to work well.
\family typewriter
samin
\family default
is called as
\begin_inset Formula \[
[theta,value,convergence]=samin("f",args,control)\]
\end_inset
\layout Standard
The controls for
\family typewriter
samin
\family default
are summarized in Table
\begin_inset Float table
wide false
collapsed false
\layout Caption
\begin_inset LatexCommand \label{cap:Controls-for-lbfgsmin}
\end_inset
Controls for
\family typewriter
samin
\layout Standard
\begin_inset Tabular
|
\begin_inset Text
\layout Standard
Element
\end_inset
|
\begin_inset Text
\layout Standard
Name
\end_inset
|
\begin_inset Text
\layout Standard
Purpose
\end_inset
|
\begin_inset Text
\layout Standard
Description
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
1
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
lb
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
lower bounds
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
vector of lower bounds for parameters
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
2
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
ub
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
upper bounds
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
vector of upper bounds for parameters
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
3
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
nt
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
control looping
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
loops per temp.
reduction, e.g., nt=20
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
4
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
ns
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
control looping
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
loops per stepsize adjustment, e.g., ns=5
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
5
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
rt
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
reduce temp.
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
0 < rt < 1, e.g., rt=0.75
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
6
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
maxevals
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
limit evaluations
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
usually, a large number, unless just exploratory, e.g., 1e10
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
7
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
neps
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
convergence
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
positive integer.
Higher is stricter criterion, e.g., neps=5
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
8
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
functol
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
convergence
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
last neps function values must be this close to eachother
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
9
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
paramtol
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
convergence
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
width of search interval must be less than this value
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
10
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
verbosity
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
output
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
0: no outout; 1: intermediate; 2: only final
\end_inset
|
|
\begin_inset Text
\layout Standard
\size footnotesize
11
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
minarg
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
arg.
for min.
\end_inset
|
\begin_inset Text
\layout Standard
\size footnotesize
which arg to min.
w.r.t., usually = 1
\end_inset
|
\end_inset
\end_inset
\layout Standard
The example program
\begin_inset ERT
status Open
\layout Standard
\backslash
htmladdnormallink{sa-example.m}{http://pareto.uab.es/mcreel/MINTOOLKIT/sa-example.m}
\end_inset
is listed here:
\begin_inset Include \verbatiminput{sa-example.m}
preview false
\end_inset
The objective function is the sum of
\begin_inset Formula $k$
\end_inset
exponentiated cosine waves, each shifted down so the minimum is zero, with
some curvature added in to create a global minimum of
\begin_inset Formula $f(x)=0$
\end_inset
at
\begin_inset Formula $x=(0,0,...,0).$
\end_inset
The (edited to shorten) output of the example is here:
\begin_inset Include \verbatiminput{sa-example.out}
preview false
\end_inset
\layout Standard
You can see that the minimum was found correctly.
\layout Subsubsection
\begin_inset LatexCommand \label{sub:A-more-difficult}
\end_inset
A more difficult problem
\layout Standard
The Moré-Garbow-Hillstrom test suite contains some relatively difficult
minimization problems.
\family typewriter
bfgsmin
\family default
by itself can solve some of these problems, but not all of them, since some
have multiple local minima, or completely flat regions where a gradient-based
method will not be able to find a decreasing direction of search.
The
\begin_inset Quotes sld
\end_inset
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{Biggs EXP6}{http://www.uni-graz.at/imawww/kuntsevich/solvopt/results/moreset.html}
\end_inset
\begin_inset Quotes sld
\end_inset
problem #18 is one for which
\family typewriter
bfgsmin
\family default
fails to find the global minimum.
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{This program}{http://pareto.uab.es/mcreel/MINTOOLKIT/solvebiggs.m}
\end_inset
shows how the global minimum may be found by combining an initial search
that uses
\family typewriter
samin
\family default
to find good starting values with refinement using
\family typewriter
bfgsmin
\family default
to sharpen up the final results.
The
\family typewriter
samin
\family default
results from running this program, which use a fast temperature reduction
and a fairly low limit on function evaluations are:
\layout Standard
\paragraph_spacing single
\family typewriter
NO CONVERGENCE: MAXEVALS exceeded
\layout Standard
\paragraph_spacing single
\family typewriter
Stage 2, decreasing temperature
\layout Standard
\paragraph_spacing single
\family typewriter
Obj.
fn.
value 0.000006
\layout Standard
\paragraph_spacing single
\family typewriter
parameter search width
\layout Standard
\paragraph_spacing single
\family typewriter
9.844228 0.000000
\layout Standard
\paragraph_spacing single
\family typewriter
4.294696 0.000000
\layout Standard
\paragraph_spacing single
\family typewriter
-5.231675 0.000000
\layout Standard
\paragraph_spacing single
\family typewriter
-3.114330 0.000000
\layout Standard
\paragraph_spacing single
\family typewriter
1.076916 0.000000
\layout Standard
\paragraph_spacing single
\family typewriter
1.118343 0.000000
\layout Standard
Then come the BFGS iterations to sharpen up the results.
The final BFGS results are:
\layout Standard
\paragraph_spacing single
\family typewriter
BFGSMIN intermediate results: Iteration 33
\layout Standard
\paragraph_spacing single
\family typewriter
Stepsize 0.0000000
\layout Standard
\paragraph_spacing single
\family typewriter
Using analytic gradient
\layout Standard
\paragraph_spacing single
\family typewriter
Objective function value 0.0000000000
\layout Standard
\paragraph_spacing single
\family typewriter
Function conv 1 Param conv 1 Gradient conv 1
\layout Standard
\paragraph_spacing single
\family typewriter
params gradient change
\layout Standard
\paragraph_spacing single
\family typewriter
10.0000 0.0000 0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
4.0000 -0.0000 -0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
-5.0000 0.0000 0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
-3.0000 -0.0000 0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
1.0000 -0.0000 0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
1.0000 0.0000 0.0000
\layout Itemize
The minimum is found, but note that the solution values are in a different
order than those given on the SolvOpt web page, with some negative signs.
This problem suffers from a lack of identification - there are multiple
values that give exactly the same value of zero.
\layout Standard
An alternative which will often be faster, but is less sure to find the
global minimum, is to call
\family typewriter
bfgsmin
\family default
with many random starting values and a limited number of iterations.
This is implemented in
\family typewriter
battery.m
\family default
.
You can see an example in
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{This program}{http://pareto.uab.es/mcreel/MINTOOLKIT/solvebiggs2.m}
\end_inset
.
This leads to the results
\layout Standard
\paragraph_spacing single
\family typewriter
BFGSMIN intermediate results: Iteration 130
\layout Standard
\paragraph_spacing single
\family typewriter
Stepsize 0.0000000
\layout Standard
\paragraph_spacing single
\family typewriter
Using analytic gradient
\layout Standard
\paragraph_spacing single
\family typewriter
Objective function value 0.0000000000
\layout Standard
\paragraph_spacing single
\family typewriter
Function conv 1 Param conv 1 Gradient conv 1
\layout Standard
\paragraph_spacing single
\family typewriter
params gradient change
\layout Standard
\paragraph_spacing single
\family typewriter
4.0000 -0.0000 0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
10.0000 0.0000 -0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
3.0000 0.0000 0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
5.0000 -0.0000 0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
1.0000 -0.0000 0.0000
\layout Standard
\paragraph_spacing single
\family typewriter
1.0000 0.0000 0.0000
\layout Standard
The minimum is found correctly, and you can see that the problem is not
identified.
\layout Subsubsection
Tips for successful minimization
\layout Paragraph
Scaling
\layout Standard
Scaling the data and other constant parameters of the objective function
so that the elements of the gradient are of approximately the same order
of magnitude will help improve accuracy of the Hessian approximation.
This can help a lot in obtaining convergence, and the results will have
higher accuracy.
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{This program}{http://pareto.uab.es/mcreel/MINTOOLKIT/tips.m}
\end_inset
illustrates.
The output is
\begin_inset Include \verbatiminput{tip1.out}
preview false
\end_inset
You can see that the scaled data gives a more accurate solution, using less
than half the iterations.
\layout Paragraph
Bullet-proofing
\layout Standard
Writing your objective function so that it cannot return NaN or otherwise
crash can save a lot of grief.
Insert something like
\family typewriter
if (((abs(obj_value) == Inf)) || (isnan(obj_value)))
\layout Standard
\family typewriter
obj_value = realmax;
\layout Standard
\noindent
\family typewriter
endif
\layout Standard
\noindent
at the end of the objective function, and then return
\family typewriter
obj_value
\family default
.
This way, parameter values that lead to crashes are penalized, and will
be avoided automatically.
\layout Subsection
Numeric differentiation
\layout Standard
\family typewriter
numgradient
\family default
and
\family typewriter
numhessian
\family default
can be used for numeric differentiation.
\family typewriter
numgradient
\family default
returns the derivative of an
\begin_inset Formula $n\times1$
\end_inset
vector-valued function with respect to a
\begin_inset Formula $k\times1$
\end_inset
vector in a
\begin_inset Formula $n\times k$
\end_inset
matrix.
\family typewriter
numhessian
\family default
returns the derivative of a real-valued function with respect to a
\begin_inset Formula $k\times1$
\end_inset
vector in a
\begin_inset Formula $k\times k$
\end_inset
matrix.
Both functions are quite accurate.
\begin_inset ERT
status Open
\layout Standard
\backslash
htmladdnormallink{numderivatives.m}{http://pareto.uab.es/mcreel/MINTOOLKIT/numderivatives.m}
\end_inset
, which follows, shows how it can be done.
\begin_inset Include \verbatiminput{numderivatives.m}
preview false
\end_inset
\layout Standard
The results are:
\begin_inset Include \verbatiminput{numderivative.out}
preview false
\end_inset
\layout Section
Testing the code
\layout Standard
The program
\begin_inset ERT
status Open
\layout Standard
\backslash
htmladdnormallink{mgh-test.m}{http://pareto.uab.es/mcreel/MINTOOLKIT/mgh-test.m}
\end_inset
allows testing the algorithms using the Moré-Garbow-Hillstrom test suite,
obtained from the
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{SolvOpt}{http://www.uni-graz.at/imawww/kuntsevich/solvopt/}
\end_inset
source code.
You can compare the output with
\begin_inset ERT
status Collapsed
\layout Standard
\backslash
htmladdnormallink{these results}{http://www.uni-graz.at/imawww/kuntsevich/solvopt/results/table1.html}
\end_inset
, if you like.
Note that simply applying BFGS with a single start value will sometimes
lead to a failure of convergence, or convergence to a non-global minimum.
This is expected, considering the nature of the problems.
See section
\begin_inset LatexCommand \ref{sub:A-more-difficult}
\end_inset
for an appropriate means of proceeding with these problems.
\layout Standard
If you find any bugs in the code, please contact me.
\layout Bibliography
\bibitem {key-5}
Eaton, J.W.,
\begin_inset LatexCommand \url{http://www.octave.org/}
\end_inset
\layout Bibliography
\bibitem {key-7}
Liu and Nocedal,
\begin_inset LatexCommand \url{http://www.ece.northwestern.edu/~nocedal/PDFfiles/limited-memory.pdf}
\end_inset
\layout Bibliography
\bibitem {key-4}
Mittelmann,
\begin_inset LatexCommand \url{http://plato.asu.edu/topics/problems/global.html}
\end_inset
\layout Bibliography
\bibitem {key-3}
Nocedal (1992),
\begin_inset LatexCommand \url{http://www.ece.northwestern.edu/~nocedal/PDFfiles/acta.pdf}
\end_inset
\layout Bibliography
\bibitem {key-2}
Goffe,
\begin_inset LatexCommand \url{http://www.netlib.no/netlib/opt/simann.f}
\end_inset
\layout Bibliography
\bibitem {key-6}
Goffe, "SIMANN: A Global Optimization Algorithm using Simulated Annealing
" Studies in Nonlinear Dynamics & Econometrics, Oct96, Vol.
1 Issue 3.
\the_end