Services & ResourcesWolfram ForumsMathGroup Archive

[Date Index] [Thread Index] [Author Index]

Re: FindRoot within a FindRoot Problem

• To: mathgroup at smc.vnet.net
• Subject: [mg84234] Re: [mg84212] FindRoot within a FindRoot Problem
• From: Daniel Lichtblau <danl at wolfram.com>
• Date: Fri, 14 Dec 2007 15:00:16 -0500 (EST)
• References: <200712140856.DAA07770@smc.vnet.net>

stephen.a.mccormack at baesystems.com wrote:
> I need to find mu1 and sigma1 such that the two errors are zero (or at
> least minimized). The following code works but is manual and iterative
> in that I need to manually guess the values for mu1 and sigma1, check
> the errors and guess again, and so on.
> 
> \[Mu]1 = 80.(* a guessed value *);
> \[Sigma]1 = 107.0     (* a guessed value *);
> \!\(TraditionalForm\`\(\(solution1 = FindRoot[{\[Mu]1 ==
> \*SuperscriptBox[\(E\), \(M1 +
> \*FractionBox[
> SuperscriptBox[\(S1\), \(2\)], \(2\)]\)], \[Sigma]1 ==
> \*SqrtBox[\(
> \*SuperscriptBox[\(E\), \(
> \*SuperscriptBox[\(S1\), \(2\)] + 2*M1\)]*\((
> \*SuperscriptBox[\(E\),
> SuperscriptBox[\(S1\), \(2\)]] -
>             1)\)\)]}, {M1,  .1}, {S1,  .1}];\)\[IndentingNewLine]
> \(m1 = ReplaceAll[M1,
>       solution1[\([\)\(1\)\(]\)]];\)\[IndentingNewLine]
> s1 = ReplaceAll[S1, solution1[\([\)\(2\)\(]\)]]\)\);
> dist1 = LogNormalDistribution[m1, s1];
> dist1pdf = PDF[dist1, x];
> ivalue2a = Integrate[dist1pdf, {x, 0, 10}]*100
> ivalue13a = 100 - Integrate[dist1pdf, {x, 0, 1000}]*100
> error1 = (ivalue2a - 2.78*10^-16)/(2.78*10^-16)*100
> 
> error2 = (ivalue13a - .2)/.2*100
> {m1, s1}
> 
> I tried wrapping the code within FindRoot but get error messages (see
> below). Does anyone know how to force Mathematica to solve for mu1 and
> sigma1 while setting both errors to zero ?
> 
> Thanks
> 
> \!\(TraditionalForm\`\(\(solution1 = FindRoot[{\[Mu]1 ==
> \*SuperscriptBox[\(E\), \(M1 +
> \*FractionBox[
> SuperscriptBox[\(S1\), \(2\)], \(2\)]\)], \[Sigma]1 ==
> \*SqrtBox[\(
> \*SuperscriptBox[\(E\), \(
> \*SuperscriptBox[\(S1\), \(2\)] + 2*M1\)]*\((
> \*SuperscriptBox[\(E\),
> SuperscriptBox[\(S1\), \(2\)]] -
>             1)\)\)]}, {M1,  .1}, {S1,  .1}];\)\[IndentingNewLine]
> \(m1 = ReplaceAll[M1,
>       solution1[\([\)\(1\)\(]\)]];\)\[IndentingNewLine]
> s1 = ReplaceAll[S1, solution1[\([\)\(2\)\(]\)]]\)\);
> dist1 = LogNormalDistribution[m1, s1];
> dist1pdf = PDF[dist1, x];
> ivalue2a = Integrate[dist1pdf, {x, 0, 10}]*100
> ivalue13a = 100 - Integrate[dist1pdf, {x, 0, 1000}]*100
> solution2 =
>  FindRoot[{0. == (ivalue2a - 2.78*10^-16)/(2.78*10^-16)*100.,
>    0. == (ivalue13a - .2)/.2*100. == 0.}, {\[Mu]1, 95.}, {\[Sigma]1,
>    107.}]

I find this quite hard to follow. In general it is vastly easier on the 
eyes if you use plain text ascii instead of the human-unreadable mess 
above. To do this, you might do Cell > Convert To > InputForm before you 
copy/paste into email.

What you want to do is define these various entities as functions of mu 
and sigma. Start with solution1[mu,sigma]. It will do a FindRoot, 
returning the {m1,s1} pair appropriate for the given values of mu and 
sigma. We will use this in dist1, etc.

Here is the code I used, and the resulting {mu,sigma} pair.

solution1[(mu_)?NumberQ, (sigma_)?NumberQ] :=
    {m1, s1} /. FindRoot[{mu == E^(m1 + s1^2/2),
          sigma == Sqrt[E^(s1^2 + 2*m1)*(E^s1^2 - 1)]}, {m1, 0.1},
        {s1, 0.1}]

dist1[(mu_)?NumberQ, (sigma_)?NumberQ] :=
      Module[{m, s}, {m, s} = solution1[mu, sigma];
         LogNormalDistribution[m, s]];

dist1pdf[(mu_)?NumberQ, (sigma_)?NumberQ] :=
      PDF[dist1[mu, sigma], x];

ivalue2a[(mu_)?NumberQ, (sigma_)?NumberQ] :=
    NIntegrate[dist1pdf[mu, sigma], {x, 0, 10}]*100
ivalue13a[(mu_)?NumberQ, (sigma_)?NumberQ] :=
    100 - NIntegrate[dist1pdf[mu, sigma], {x, 0, 1000}]*100

error1[(mu_)?NumberQ, (sigma_)?NumberQ] :=
    ((ivalue2a[mu, sigma] - 2.78/10^16)/(2.78/10^16))*100
error2[(mu_)?NumberQ, (sigma_)?NumberQ] :=
    ((ivalue13a[mu, sigma] - 0.2)/0.2)*100

FindMinimum[error1[mu, sigma]^2 + error2[mu, sigma]^2,
    {mu, 80}, {sigma, 107}]

Out[16]= {3.24884*10^-23, {mu -> 342.766, sigma -> 142.691}}


Daniel Lichtblau
Wolfram Research

• References:
  • FindRoot within a FindRoot Problem
    • From: stephen.a.mccormack@baesystems.com