Chapter 7

Optimization

Scientists and engineers often have to solve for the maximum or minimum of a multi-

dimensional function, sometimes with constraints on the values of some or all of the

variables. This is known as optimization, and is a rich, highly developed, and often

difﬁcult problem. Generally the problem is phrased as a minimization, which shows

its kinship to the least-squares data ﬁtting procedures discussed in a subsequent chap-

ter. If a maximum is desired, one simply solves for the minimum of the negative of

the function. The greatest difﬁculties typically arise if the multi-dimensional surface

has local minima in addition to the global minimum, because there is no way to show

that the minimum is local except by trial and error.

R has three functions in the base installation: optimize() for one-dimensional

problems, optim() for multi-dimensional problems, and constrOptim() for opti-

mization with linear constraints. (Note that even “unconstrained” optimization nor-

mally is constrained by the limits on the search range.) We shall consider each of

these with suitable examples, and introduce several add-on packages that expand the

power of the basic functions. Optimization is a sufﬁciently large and important topic

to deserve its own task view in R, at http://cran.r-project.org/web/views/

Optimization.html.

In addition to the packages considered in this chapter, the interested reader

should become acquainted with the nloptr package, which is considered one of

the strongest and most comprehensive optimization packages in R. According to its

synopsis,“nloptr is an R interface to NLopt. NLopt is a free/open-source library for

nonlinear optimization, providing a common interface for a number of different free

optimization routines available online as well as original implementations of various

other algorithms.”

7.1 One-dimensional optimization

The base R function for ﬁnding minima (the default) or maxima of functions of a

single variable is optimize(). As a concrete example, given a 16 ×20 sheet of

cardboard, ﬁnd the size x of the squares to be cut from the corners that maximizes

the volume of the open box formed when the sides are folded up. x must be in the

range (0,8).

> optimize(function(x) x*(20-2*x)*(16-2*x), c(0,8), maximum=T)

159

160 OPTIMIZATION

0.0 1.0 2.0 3.0

-2 -1 0 1 2

x

f(x)

Figure 7.1: Plot of function f (x) = x sin(4x) showing several maxima and minima.

$maximum

[1] 2.944935

$objective

[1] 420.1104

Consider next the use of optimize with the function

> f = function(x) x*sin(4*x)

which plotted looks like this (Figure 7.1):

> curve(f,0,3)

It has two minima in the x = 0 −3 range, with the global minimum near 2.8,

and two maxima, with the global maximum near 2.0. Applying optimize() in the

simplest way yields

> optimize(f,c(0,3))

$minimum

[1] 1.228297

$objective

[1] -1.203617

which gives the local minimum because it is the ﬁrst minimum encountered by the

search algorithm (Brent’s method, which combines root bracketing, bisection, and

inverse quadratic interpolation). Because we have a plot of the function, we can see

that we must exclude the local minimum from the lower and upper endpoints of the

search interval.

> optimize(f,c(1.5,3))

$minimum

[1] 2.771403

$objective

[1] -2.760177

To ﬁnd the global maximum we enter

ONE-DIMENSIONAL OPTIMIZATION 161

> optimize(f,c(1,3),maximum=TRUE)

$maximum

[1] 1.994684

$objective

[1] 1.979182

We could have obtained the same result by minimizing the negative of the function

> optimize(function(x) -f(x),c(1,3))

$minimum

[1] 1.994684

$objective

[1] -1.979182

which ﬁnds the maximum in the right place but, of course, yields the negative of the

function value at the maximum.

If necessary, the desired accuracy can be adjusted with the tol option in the

function call.

The pracma package contains the function findmins(), which ﬁnds the posi-

tions of all the minima in the search interval by dividing it n times (default n = 100)

and applying optimize in each interval. To ﬁnd the values at those minima, evaluate

the function.

> require(pracma)

> f.mins = findmins(f,0,3)

> f.mins # x values at the minima

[1] 1.228312 2.771382

> f(f.mins[1:2]) # function evaluated at the minima

[1] -1.203617 -2.760177

The Examples section of the help page for optimize() shows how to include

a parameter in the function call. It also shows how, for a function with a very ﬂat

minimum, the wrong solution can be obtained if the search interval is not properly

chosen. Unfortunately, there is no clear way to choose the search interval in such

cases, so if the results are not as expected from inspecting the graph of the function,

other intervals should be explored.

The function f (x) = |x

2

−8| yields some interesting behavior (Figure 7.2).

> f = function(x) abs(x^2-8)

> curve(f,-4,4)

Straightforward solving for the maximum over the entire range yields the result at

x = 0.

> optimize(f,c(-4,4),maximum=T)

$maximum

[1] -1.110223e-16

$objective

[1] 8

Excluding the middle from the search interval ﬁnds the maxima at the extremes.

Get *Using R for Numerical Analysis in Science and Engineering* now with O’Reilly online learning.

O’Reilly members experience live online training, plus books, videos, and digital content from 200+ publishers.