Services & Resources / Wolfram Forums / MathGroup Archive
-----

MathGroup Archive 2012

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

Search the Archive

Re: Fast vs. Slow NonlinearModelFit models

  • To: mathgroup at smc.vnet.net
  • Subject: [mg123990] Re: Fast vs. Slow NonlinearModelFit models
  • From: "Dan O'Brien" <danobrie at gmail.com>
  • Date: Wed, 4 Jan 2012 05:02:48 -0500 (EST)
  • Delivered-to: l-mathgroup@mail-archive0.wolfram.com
  • References: <201112211152.GAA02211@smc.vnet.net> <201112220926.EAA13844@smc.vnet.net>

Ok this idea occurred to me, but I can't expand it so that it's a good 
approximation over a large enough domain.  Any other thoughts?

On 12/22/2011 3:26 AM, DrMajorBob wrote:
> You might fit using a short Series expansion in place of Erfc, then use
> the fitted parameters as initial guess in a second fit for the exact model.
>
> Bobby
>
> On Wed, 21 Dec 2011 05:52:34 -0600, Dan O'Brien<danobrie at gmail.com>  wrote:
>
>> Hello,
>>
>> I'm wondering if anyone has a suggestion for speeding up the
>> NonlinearModelFit computation for model2 below.  The model used to fit
>> my experimental data can be a simple complex function who's magnitude
>> squared is a lorentzian (model1), or, some may argue a more complete
>> model is the same function convolved with a gaussian which gives a
>> function that involves the complementary error function (model2), the
>> Voigt profile.  As indicated below with AbsoluteTiming using these two
>> different models show VERY different computational times.  I imagine it
>> has something to do with mathematica's internal implementation of Erfc
>>
>> In[1]:= $Version
>> (*Lorentzian Lineshape*)
>> \[Chi]R[\[Omega]_, {A_, \[Phi]A_, \[CapitalGamma]_, \[Omega]v_}] := (
>>     A E^(I \[Phi]A))/(-(\[Omega] - \[Omega]v) - I \[CapitalGamma]);
>>
>> (*Voigt Lineshape*)
>> z[\[Sigma]_, \[CapitalGamma]_, \[Omega]v_, \[Omega]_] := (\[Omega] - \
>> \[Omega]v + I \[CapitalGamma])/(Sqrt[2] \[Sigma])
>> FadeevaF[z_] := Exp[-z^2] Erfc[-I z]
>> IH\[Chi]R[\[Omega]_, {\[Sigma]_,
>>      A_, \[Phi]A_, \[CapitalGamma]_, \[Omega]v_}] :=
>>    I Sqrt[\[Pi]/2] (
>>     A E^(I \[Phi]A))/\[Sigma]  FadeevaF[
>>      z[\[Sigma], \[CapitalGamma], \[Omega]v, \[Omega]]]
>> (*Redefine functions so A is peak amplitude of the function magnitude \
>> squared*)
>> A\[Chi]R[\[Omega]_, {\[Sigma]_,
>>      A_, \[Phi]A_, \[CapitalGamma]_, \[Omega]v_}] := \[Chi]R[ \[Omega], \
>> {Sqrt[Abs[A]] Abs[\[CapitalGamma]], \[Phi]A,
>>      Abs[\[CapitalGamma]], \[Omega]v}]
>> AIH\[Chi]R[\[Omega]_, {\[Sigma]_,
>>      A_, \[Phi]A_, \[CapitalGamma]_, \[Omega]v_}] :=
>>    IH\[Chi]R[\[Omega], {\[Sigma], (
>>      Sqrt[Abs[A]] E^(-(1/2) \[CapitalGamma]^2/\[Sigma]^2) Sqrt[2/\[Pi]]
>>        Abs[\[Sigma]])/
>>      Erfc[Abs[\[CapitalGamma]]/(
>>       Sqrt[2] Abs[\[Sigma]])], \[Phi]A, \[CapitalGamma], \[Omega]v}]
>> (*Generate data with noise*)
>> dataplot =
>>    ListPlot[data =
>>      Table[{\[Omega],
>>        Abs[AIH\[Chi]R[\[Omega], {4, 1.1, 0, 2.1, 0}] +
>>           AIH\[Chi]R[\[Omega], {4, .6, \[Pi], 2.1, 18}]]^2 +
>>         RandomReal[{-0.05, 0.05}]}, {\[Omega], -30, 50, .5}],
>>     Joined ->  True]
>> (*Models*)
>> model1 = Abs[
>>       A\[Chi]R[#, {4, a1, 0, \[CapitalGamma]1, \[Omega]vo1}] +
>>        A\[Chi]R[#, {4, a2, \[Pi], \[CapitalGamma]2, \[Omega]vo2}]]^2&;
>> model2 =
>>     Abs[AIH\[Chi]R[#, {4, a1, 0, \[CapitalGamma]1, \[Omega]vo1}] +
>>        AIH\[Chi]R[#, {4, a2, \[Pi], \[CapitalGamma]2, \[Omega]vo2}]]^2&;
>>
>> (*Fit*)
>> AbsoluteTiming[
>>     mod1 = NonlinearModelFit[data,
>>       model1[\[Omega]], {{a1, 1}, {\[CapitalGamma]1, 2}, {\[Omega]vo1,
>>         1}, {a2, .8}, {\[CapitalGamma]2, 2}, {\[Omega]vo2,
>>         15}}, \[Omega]]][[1]]
>> AbsoluteTiming[
>>     mod2 = NonlinearModelFit[data,
>>       model2[\[Omega]], {{a1, 1}, {\[CapitalGamma]1, 2}, {\[Omega]vo1,
>>         1}, {a2, .8}, {\[CapitalGamma]2, 2}, {\[Omega]vo2,
>>         15}}, \[Omega]]][[1]]
>>
>> Show[{dataplot,
>>     Plot[{Normal[mod1[\[Omega]]],
>>       Normal[mod2[\[Omega]]]}, {\[Omega], -30, 50}]}]
>>
>> Out[1]= "8.0 for Microsoft Windows (32-bit) (October 6, 2011)"
>>
>> Out[8]=****Plot Snipped
>>
>> Out[11]= 0.0390625 *Here is AbsoluteTiming for model1*
>>
>> Out[12]= 75.0078125 *Here is AbsoluteTiming for model2*
>>
>> Out[13]=****Plot Snipped
>



  • Prev by Date: Re: Mathematica on tablet computer
  • Next by Date: Re: Mathematica on tablet computer
  • Previous by thread: Re: Why does the Front End freeze when using a TCPIP-mode connection with the kernel?
  • Next by thread: How to plot divergence of gradient as contour plot