0
$\begingroup$

I'm trying to fit the array data as follow:

NonlinearModelFit[
  data, 
  s*t^r + (A/(w*(Pi/2)^0.5))*Exp[-2*((t - u)/w)^2], 
  {s, r, A, w, u}, 
  t];

But when I run it, it gives the following problem:

Exp[] is too small to represent as a normalized machine number; precision may be lost.

Do you know how to solve it?

Thank you!

Here the code:

data={{9331200, 2.25014*10^-14}, {9590400, 2.94221*10^-14}, {12441600, 
  2.25158*10^-14}, {13219200, 1.88385*10^-14}, {13564800, 
  8.74689*10^-15}, {13824000, 1.65873*10^-14}, {13910400, 
  1.91969*10^-14}, {14083200, 1.85371*10^-14}, {14256000, 
  1.21808*10^-14}, {22464000, 1.55479*10^-14}};

errors={2.35061*10^-15, 4.03579*10^-15, 2.57128*10^-15, 3.17395*10^-15, 
 4.44869*10^-15, 3.86327*10^-15, 3.73115*10^-15, 1.55858*10^-15, 
 4.62527*10^-15, 1.50732*10^-15};

dataPlot = ListLogLogPlot[data, PlotRange -> Full];

(*First Fit with Powerlaw, it works*)
Pow = NonlinearModelFit[data, k*t^h, {k,h}, t, 
   MaxIterations -> 1000, Weights -> 1/errors^2];
PowPlot = 
  LogLogPlot[Pow[t], {t, 100*24*3600, 260*24*3600}, PlotStyle -> Red, 
   PlotRange -> Full];
"Plot powerlaw:"
Show[{PowPlot, dataPlot}, PlotRange -> Full]

Pow["ParameterTable"]

chi1r = Pow["EstimatedVariance"]


(*Second Fit with Powerlaw+Gaussian, where it gives that problem*)
Pgauss = NonlinearModelFit[data, 
   s*t^r + (A/(w*(Pi/2)^(1/2)))*Exp[-2*((t - u)/w)^2], {s, r, A, w, 
    u}, t, Weights -> 1/errors^2, MaxIterations -> 1000, 
   WorkingPrecision -> 50, Method -> NMinimize];

PgaussPlot = 
  LogLogPlot[Pgauss[t], {t, 100*24*3600, 260*24*3600}, PlotStyle -> Red];

"Plot powerlaw+gaussiana:"
Show[{PgaussPlot, dataPlot}]

Pgauss["ParameterTable"]
chi2r = Pgauss["EstimatedVariance"]

Thank you very much!

$\endgroup$
4
  • 3
    $\begingroup$ Please provide Mathematica-code and your data! $\endgroup$ Commented May 26, 2020 at 18:43
  • 1
    $\begingroup$ Giving reasonable starting values for the parameters, e.g. NonlinearModelFit[data, model, {{s, 1.}, {r, 2.}, etc}}] or choosing Method->"NMinimize" might help. $\endgroup$
    – aooiiii
    Commented May 26, 2020 at 19:36
  • 1
    $\begingroup$ Rationalize your data, change 0.5 to 1/2, and specify a WorkingPrecision to use arbitrary-precision: NonlinearModelFit[data, st^r + (A/(w (Pi/2)^(12))) Exp[-2 ((t - u)/w)^2], {s, r, A, w, u}, t, WorkingPrecision -> 25] Even higher precision may be needed. $\endgroup$
    – Bob Hanlon
    Commented May 26, 2020 at 21:17
  • $\begingroup$ Some good tips here: mathematica.stackexchange.com/questions/139038/… $\endgroup$
    – JimB
    Commented May 27, 2020 at 16:59

1 Answer 1

2
$\begingroup$

Higher precision is not needed. One just needs to scale the data so that numerical instabilities don't occur. This is especially the case here where the scales of the predictor and response variables are so far apart.

For both models multiplying the predictor and response by appropriate constants don't change the underlying model. (Sometimes subtracting a constant is what you want which doesn't change the underlying model but that's not the case here.)

(* Original data *)
data = {{9331200, 2.25014*10^-14}, {9590400, 2.94221*10^-14}, {12441600, 2.25158*10^-14},
   {13219200, 1.88385*10^-14}, {13564800, 8.74689*10^-15}, {13824000,  1.65873*10^-14},
   {13910400, 1.91969*10^-14}, {14083200, 1.85371*10^-14}, {14256000,  1.21808*10^-14},
   {22464000, 1.55479*10^-14}};
errors = {2.35061*10^-15, 4.03579*10^-15, 2.57128*10^-15, 3.17395*10^-15, 4.44869*10^-15,
   3.86327*10^-15, 3.73115*10^-15, 1.55858*10^-15, 4.62527*10^-15, 1.50732*10^-15};

(* Scale data *)
errors2 = errors 10^15;
data2 = data;
data2[[All, 1]] = data[[All, 1]]/10^7;
data2[[All, 2]] = data[[All, 2]] 10^14;

(* Power law fit *)
Pow = NonlinearModelFit[data2, k*t^h, {k, h}, t, Weights -> 1/errors2^2];
Pow["ParameterTable"]

Parameter table for power law fit

(* More complicated model *)
Pgauss = NonlinearModelFit[data2, s*t^r + (A/(w (Pi/2)^(1/2)))*Exp[-2*((t - u)/w)^2],
   {s, r, A, w, u}, t, Weights -> 1/errors2^2];

(* Results *)
Pgauss["ParameterTable"]

Parameter table for more complicated model

 Show[ListPlot[data2],
   Plot[{Pow[t], Pgauss[t]}, {t, Min[data2[[All, 1]]], Max[data2[[All, 1]]]}, PlotRange -> All,
   PlotLegends -> {"Pow model", "Pgauss model"}], PlotRange -> All]

Data and fits from Pow and Pgauss models

The Pgauss fit looks pretty implausible as an explanation. (The fit for the Pgauss model is fine. It's the consideration of using it to describe the relationship is what is implausible.) The $AIC_c$ statistics bear that out:

Pow["AICc"]
(* 15.9845 *)
Pgauss["AICc"]
(* 41.1028 *)

as the $AIC_c$ statistic for the Pow model is far lower than that of the Pgauss model.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.