|
Classroom Tips and Techniques: Nonlinear Fit, Optimization, and the DirectSearch Package
Robert J. Lopez
Emeritus Professor of Mathematics and Maple Fellow
Maplesoft
|
Introduction
|
| |
In this month's article, I revisit a nonlinear curve-fitting problem that appears in my Advanced Engineering Mathematics ebook, examine the role of Maple's Optimization package in that problem, and then explore the DirectSearch package from Dr. Sergey N. Moiseev.
The nonlinear curve-fitting problem was brought to me by one of the statisticians at the Rose-Hulman Institute of Technology well before my 2003 retirement from the classroom. Originating in the newly established biomed faculty at RHIT, its experimentally obtained data points were to be fit by the Michaelis-Menton function  .
When I initially examined the problem in Maple V Release 4, I found that the fit command in the (now deprecated) stats package only applied to linear models. Hence, when this problem became a section in the 2001 Addison-Wesley version of my Advanced Engineering Math text, it was solved from first principles by forming and minimizing  , the sum of squares of deviations. Since the text predated the Optimization package,  was minimized by techniques of the calculus. (The essence of the discussion in the text is a comparison between the nonlinear fit and two different linearizations. The linearizations prove to be poor fits!)
By the time the AEM ebook was prepared in Maple 10, the Optimization package had been added to Maple, and its use for this problem was then included in the ebook. The results for the minimization of  , as one would expect, agree for the two approaches. Moreover, the NonlinearFit command in the Statistics package (new in Maple 10) provided an immediate fit of the Michaelis-Menton model to the given data.
The recent announcement in MaplePrimes of version 2 of the DirectSearch package caught my attention. I was curious to see to what extent this package would have been of use for the curve-fitting problem in the AEM ebook. I was surprised (and delighted) to discover how robust this package actually is.
|
|
Michaelis-Menton Fit to Biomed Data
|
| |
Figure 1 is a graph of the 46 data points that are to be fit with the Michaelis-Menton function  . The plot command and the actual data are "behind" the table containing the graph. To see them, select the "Show input" option in the Table Properties dialog. The data,  , are given in floating-point form, to six places.
|
|
|
Figure 1 Forty-six data points  to be fit with 
|
The following definition of the Michaelis-Menton function also takes place in a table for which the input is not shown.
Prior to Maple 10, there was no direct functionality for a nonlinear fit in Maple. Hence, from first principles, we let  be the sum of squares of deviations for  and the given data:
Figure 2 is a graph of a portion of the surface  . It shows the existence of a minimum point, and also shows that this minimum occurs in a steep, but narrow, "valley." Again, the graph appears in a table that conceals input.
 
|
|
|
Figure 2 The sum of squares of the deviations graphed as a function of the parameters  and 
|
The ordinary tools of the calculus provide the minimizing values of  and  ,
|
|
(1) |
and hence of the fitting function itself
 =  
and the minimum value of  , namely,
 =  
Applying the Optimization package's Minimize command to  gives
This extreme value is not as accurate as we obtained in (1). An improvement can be obtained with the Optimization package's NLPSolve command, which became available in Maple 9.5.
 
From Maple 10 onwards, the NonlinearFit command from the Statistics package immediately gives
 
|
|
The DirectSearch Package
|
| |
|
Description
|
| |
Version 2 of the DirectSearch package is available from the Maple Application Center. Installation instructions are provided in the README.txt file that accompanies the code. Although written for Maple 13, I've used it in Maple 14 and Maple 15 with no ill effects.
The code for the commands in the package are in the archive DirectSearch.mla; the archive DirectSearch.hdb contains the files for the help pages, which are very good. There are three different optimizers: Search, GlobalOptima, and Global Search, all of which use derivative-free direct search methods so all the computations in the package apply to functions that may not be differentiable or continuous. In contrast, the only derivative-free algorithm in Maple's Optimization package is the sequential simplex technique of Nelder and Mead, implemented by the nonlinearsimplex option in the NLPSolve command.
For multiobjective optimization, the DirectSearch package has the seven commands Minimax, CompromiseProgramming, WeightedSum, BoundedObjective, ModifiedTchebycheff, WeightedProduct, and ExponentialWeightedSum. In addition, there is the SolveEquations command for solving systems of equations, and the DataFit command for nonlinear fits to data. The DataFit command supports nine different fitting methods, best described via Table 1, extracted from the help page for this command.
|
Fit method
|
Method description
|
Objective function
|
|
lms
|
Least mean squares method. This method is effective but not robust to outliers.
|
|
|
lad
|
Least absolute deviations method. This method is more robust to outliers than LMS method.
|
|
|
lts
|
Least trimmed squares method. Robustness of this method depends on upper percentile value p ( by default).
|
where TrimmedMean is trimmed mean
|
|
lws
|
Least Winsorized squares method. Robustness of this method depends on upper percentile value p ( by default).
|
where WinsorizedMean is Winsorized mean
|
|
minimax
|
Minimax method. This method is very sensitive to outliers.
|
|
|
median
|
Least median squares method. This method is most robust to outliers.
|
where is sample median
|
|
quantile
|
Least quantile squares method. Robustness of this method depends on quantile order p ( by default).
|
where is sample quantile of order p
|
|
trimean
|
Trimean method. This method is robust to outliers.
|
where is sample median and quantile of order p
|
|
lfad
|
Least function of absolute deviations method. Robustness of this method depends on the function  specified (by default,  ).
|
|
|
Table 1 Options to the DataFit command in the DirectSearch package
|
|
|
Application 1
|
| |
Our first application of DirectSearch is to the problem iconized by Figure 1. According to Table 1, the "ordinary" least squares method is implemented in DirectSearch as
The first number is the minimum value of the objective function  defined in the top row of Table 1. Note that this function differs from our  by the factor  , where, for this calculation,  . Hence, we compute
 = 
results compatible with those obtained earlier.
|
|
Application 2
|
| |
Exercise 3 in Section 47.3 of my Advanced Engineering Mathematics ebook is the isoperimetric problem asking for the extremal  that renders a functional  stationary, satisfies  , has arc length  , and intersects the transversal  , where  is given by
An extremal  , an arc of the circle  , is determined by a triple  that is a solution of the following three equations.
Maple's fsolve command (numeric solver) finds one solution, namely,
Unfortunately, it fails to find any other solution. There is a second solution, found by solving the middle equation for a, and so eliminating a from the first and third equations. This leads to two equations of the form  , whose graph reveals that there are two intersections of the two implicit functions so defined. The second solution is then found to be
The SolveEquations command in DirectSearch seeks minima of  , and returns an m × 4 matrix, each row of which details an extreme point: the value of R, a vector of residuals for each  , the coordinates of the point, and the number of evaluations of R it took to find that point. The minima are arranged in order of increasing value of R.
With n = 3, our equations yield the matrix
 
which has some 83 rows. The first two solutions are
 
and
 
for which the R-values are on the order of  and  , respectively. The value of R in the third row is 1.7 ×  ; in the last row it is greater than 7. Hence, we conclude that there are just two solutions of our three equations and these solutions correspond to the two found in Maple.
|
|
Application 3
|
| |
Some time ago I was asked if Maple could solve the following set of five equations for  .
Maple's fsolve command does not provide a solution, and its LSSolve command from the Optimization package gives different solutions (depending on the initial point), generally with sum of squares approximately 3 ×  . I never could bring my self to believe that any of these minima for the sum of squares was actually a solution of the equations.
Of course, when I started exploring DirectSearch, I tried
which results in a 22 × 4 matrix of extrema. The first such extremum, namely,
has a residual of 9.8 ×  , but no solution produced by SolveEquations matches any returned by LSSolve. The existence of a solution now becomes an urgent issue.
To this end, consider the difference
|
|
(2) |
It would appear from the difference in signs on the left and right sides of (2) that no solution exists. Hence, we must be careful with numeric "solutions" obtained by minimizing  . Just because a powerful minimization tool is used does not necessarily mean that a solution has been found for the set of equations  .
|
|
Legal Notice: © Maplesoft, a division of Waterloo Maple Inc. 2011. Maplesoft and Maple are trademarks of Waterloo Maple Inc. This application may contain errors and Maplesoft is not liable for any damages resulting from the use of this material. This application is intended for non-commercial, non-profit use only. Contact Maplesoft for permission if you wish to use this application in for-profit activities.
|
