Re: Problem behavior with FindMaximum
- To: mathgroup at smc.vnet.net
- Subject: [mg59702] Re: Problem behavior with FindMaximum
- From: "Carl K. Woll" <carlw at u.washington.edu>
- Date: Thu, 18 Aug 2005 00:16:34 -0400 (EDT)
- Organization: University of Washington
- References: <ddur1e$oc4$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
"James H. Steiger" <jsteiger at bellsouth.net> wrote in message news:ddur1e$oc4$1 at smc.vnet.net... > Hello all: > > I wonder if you could give me some advice about behavior of FindMaximum[] > that I cannot seem to decipher. > > There is a broad class of problems in statistics that involves finding the > maximum > of a function of several parameters, all but one (call it "a") of which > are > *nuisance parameters*. > The function is evaluated at any value of "a" by maximizing it w.r.t. all > the nuisance parameters. > > A simple example (constructed just for Mathgroup -- the actual > functions I work with are messier) should make this clear. > > > myFunc1[a_, b_] := 1/((a - 16)^2 + 1) - (a - b)^2 > myFunc2[a_] := FindMaximum[myFunc1[a, b], {b, 0}][[1]] > > FindMaximum returns a list, the first element of which is the maximized > value of the function, > the second of which is a replacement rule specifying the value of b at > which > the maximum occurs. > > As you can quickly verify, myFunc2 is well behaved, and you can plot > myFunc2 > without incident. > > Plot[myFunc2[a],{a,14,18}] produces a nice plot with no error messages > > Here is where the problem arises. Suppose you want to use FindMaximum[] to > obtain the maximum of myFunc2 which clearly occurs at > a=16. > The simplest answer is to add ?NumericQ to your function definition: myFunc2[a_?NumericQ] := FindMaximum[myFunc1[a, b], {b, 0}][[1]] This prevents your myFunc2 code from executing when its argument is nonnumeric, avoiding the error messages you list below. When FindMaximum is given a function it first takes a look at the symbolic form to see if it can find the derivative. The ?NumericQ pattern prevents errors during this phase. Once FindMaximum determines that it can't evaluate derivatives, it will use a nonderivative routine to compute the maximum. Note that nonderivative methods require many more function evaluations than derivative methods: In[21]:= FindMaximum[myFunc2[x],{x,10},EvaluationMonitor:>Sow[x]]//Reap Out[21]= {{1.,{x->16.}},{{10.,10.,10.,10.0438,10.0438,10.0438,10.2191,10.2191, 10.2191,10.9204,10.9204,10.9204,13.7253,13.7253,13.7253,24.9452,24.9452, 24.9452,16.8133,16.8133,16.8133,14.2454,14.2454,14.2454,16.142,16.142, 16.142,15.999,15.999,15.999,17.756,17.756,17.756,16.0006,16.0006,16.0006, 16.,16.,16.,16.,16.,16.}}} The derivative of myFunc2 is simply the derivative of myFunc1 with respect to its first argument evaluated at the appropriate a and b: myFunc2p[a_?NumericQ] := myFunc1p[a, b] /. FindMaximum[myFunc1[a, b], {b, 0}][[2]] myFunc1p[a_, b_] = D[myFunc1[a, b], a]; myFunc2' = myFunc2p; myFunc1p is the appropriate derivative of myFunc1, and the line myFunc2'=myFunc2p teaches Mathematica how to take derivatives of myFunc2. Now we can use derivative methods: In[25]:= FindMaximum[myFunc2[x],{x,10},EvaluationMonitor:>Sow[x]]//Reap Out[25]= {{1.,{x->16.}},{{10.,10.0438,10.2191,10.9204,13.7253,24.9452,16.8133, 14.2454,16.142,15.999,17.7557,16.0006,16.}}} It takes about 1/3 as many function evaluations (and 1/3 as long) now. Of course, myFunc2p is coded inefficiently since we make the identical call to FindMaximum as is done in myFunc2. We should instead change the myFunc2 code to compute both myFunc2 and myFunc2p so only one call to FindMaximum is necessary to compute both myFunc2 and myFunc2p: myFunc2[a_?NumericQ] := Module[{fm}, fm = FindMaximum[myFunc1[a, b], {b, 0}]; myFunc2p[a] = myFunc1p[a, b] /. fm[[2]]; fm[[1]]] myFunc1p[a_, b_] = D[myFunc1[a, b], a]; myFunc2' = myFunc2p; Good luck! > If you input the command > > FindMaximum[myFunc2[x],{x,14}] > > you obtain a pair of error messages (can anyone tell me how to copy these > in > Mathematica as text?) > > ----------------- > > FindMaximum::nnum: The function value 1/((1+<<1>>)^2) - (0.+a)^2 is not a > number at {b}={0.} > > FindMaximum::nnum: The function value 0.2 - (14. -b)^2 is not a number at > {a}={14.} > > ----------------- > > Is there some problem of "scope" of these variables that I am not aware > of? > Or is there some bug in FindMaximum[]? > Is there a fix? > > Thanks to all, > > Jim > > > James H. Steiger, Professor and Director > Quantitative Methods and Evaluation > Dept. of Psychology and Human Development > Vanderbilt University > Peabody College #512 > Nashville, TN, 37203 > > Phone: 615-322-7060 > email: james.h.steiger at vanderbilt.edu > Carl Woll Wolfram Research