1
$\begingroup$

I asked this question. Answers and comments were appreciated.

But (some) of the comments did not make satisfied. [Maybe due to my bad English I did not describe my problem very well].

Let me try again with another way of explanation, and I hope you can help me:

  • Start with the following example:

When we have these data points: $(0,5),(1,7),(2,9),(3,11)$ then Mr.A asked us to find the best fit of these points that is of the form $y=ax+b$. Then we can tell Mr.A: "the equation that you look for is $y=2x+5$".

Another person, Mr.B, asked us to find the best fit (of the same data points) that is of the form $y=ax^2+bx$. Then we should not say: "it is not good fit because the points represent a straight line, but your equation form is quadratic, and therefore we can not find because it is bad".

Yes, we can suggest Mr.B to go for straight line, but if he was strict, then we can tell him: "the equation that you look for is $y=-1.31579x^2+7.52632x$"

Again Mr.C asked us to find the best fit (of the same data points) that is of the form $y=ax^2+\frac{b}{x+c}$, then we can tell him: "the equation that you look for is $y=0.153282x^2+\frac{-38.869}{x+(-6.93418)}$".

What I want to say is; we should not force Mr.B and Mr.C to go for straight line.

Come to my original problem;

data = {{595098, 335}, {595158, 336}, {595218, 336}, {595338, 
    344}, {595458, 347}, {595638, 351}, {595818, 352}, {596082, 
    356}, {596322, 357}, {596922, 362}, {597521, 363}, {598481, 
    363}, {599322, 371}, {600523, 376}, {601723, 378}, {603523, 
    380}, {605323, 384}, {608924, 390}, {612523, 392}, {619724, 
    399}, {626926, 404}, {648527, 413}, {670129, 415}, {691731, 
    419}, {712906, 424}, {734504, 425}, {756104, 424}, {776690, 
    426}, {798291, 426}, {819890, 429}, {841490, 431}, {863090, 
    432}, {884692, 433}, {906290, 434}, {927892, 437}, {949492, 
    438}, {971090, 437}, {992691, 437}, {1014291, 438}, {1035891, 
    442}, {1039491, 441}, {1041236, 442}};
model = a*x + b;
fit = FindFit[data, model, {a, b}, x];
Show[Plot[Evaluate[model /. fit], {x, 595070, 1041236}], 
 ListPlot[data, PlotStyle -> Red]]

is fine

data = {{595098, 335}, {595158, 336}, {595218, 336}, {595338, 
    344}, {595458, 347}, {595638, 351}, {595818, 352}, {596082, 
    356}, {596322, 357}, {596922, 362}, {597521, 363}, {598481, 
    363}, {599322, 371}, {600523, 376}, {601723, 378}, {603523, 
    380}, {605323, 384}, {608924, 390}, {612523, 392}, {619724, 
    399}, {626926, 404}, {648527, 413}, {670129, 415}, {691731, 
    419}, {712906, 424}, {734504, 425}, {756104, 424}, {776690, 
    426}, {798291, 426}, {819890, 429}, {841490, 431}, {863090, 
    432}, {884692, 433}, {906290, 434}, {927892, 437}, {949492, 
    438}, {971090, 437}, {992691, 437}, {1014291, 438}, {1035891, 
    442}, {1039491, 441}, {1041236, 442}};
model = a*x^2 + b/x;
fit = FindFit[data, model, {a, b}, x];
Show[Plot[Evaluate[model /. fit], {x, 595070, 1041236}], 
 ListPlot[data, PlotStyle -> Red]]

is fine

data = {{595098, 335}, {595158, 336}, {595218, 336}, {595338, 
    344}, {595458, 347}, {595638, 351}, {595818, 352}, {596082, 
    356}, {596322, 357}, {596922, 362}, {597521, 363}, {598481, 
    363}, {599322, 371}, {600523, 376}, {601723, 378}, {603523, 
    380}, {605323, 384}, {608924, 390}, {612523, 392}, {619724, 
    399}, {626926, 404}, {648527, 413}, {670129, 415}, {691731, 
    419}, {712906, 424}, {734504, 425}, {756104, 424}, {776690, 
    426}, {798291, 426}, {819890, 429}, {841490, 431}, {863090, 
    432}, {884692, 433}, {906290, 434}, {927892, 437}, {949492, 
    438}, {971090, 437}, {992691, 437}, {1014291, 438}, {1035891, 
    442}, {1039491, 441}, {1041236, 442}};
model = a + b*x + c*x^2 + d*x^3 + e/x^4;
fit = FindFit[data, model, {a, b, c, d, e}, x];
Show[Plot[Evaluate[model /. fit], {x, 595070, 1041236}], 
 ListPlot[data, PlotStyle -> Red]]

is fine too.

But this is not fine:

data = {{595098, 335}, {595158, 336}, {595218, 336}, {595338, 
    344}, {595458, 347}, {595638, 351}, {595818, 352}, {596082, 
    356}, {596322, 357}, {596922, 362}, {597521, 363}, {598481, 
    363}, {599322, 371}, {600523, 376}, {601723, 378}, {603523, 
    380}, {605323, 384}, {608924, 390}, {612523, 392}, {619724, 
    399}, {626926, 404}, {648527, 413}, {670129, 415}, {691731, 
    419}, {712906, 424}, {734504, 425}, {756104, 424}, {776690, 
    426}, {798291, 426}, {819890, 429}, {841490, 431}, {863090, 
    432}, {884692, 433}, {906290, 434}, {927892, 437}, {949492, 
    438}, {971090, 437}, {992691, 437}, {1014291, 438}, {1035891, 
    442}, {1039491, 441}, {1041236, 442}};
model = a + b*Exp[c*(x^d + e)^f];
fit = FindFit[data, model, {a, b, c, d, e, f}, x];
Show[Plot[Evaluate[model /. fit], {x, 595070, 1041236}], 
 ListPlot[data, PlotStyle -> Red]]

I am sure that my model is fine to represent these data points. When I re-scaled the data (dividing x by 1000000 , and dividing y by 100), DESMOS found the parameters, this means my model is representative for my original data! https://i.imgur.com/klJatY4.png

@JimB said one of my parameters is redundant, and he was right. I want to stick with

model = a + b*Exp[c*(x^d + e)^f]

I am a new user of Mathematica, please just suggest me how to rescale my data if that will work, or suggest me to use good "starting values", or suggest me to use other codes. Do not suggest me to use other models like suggesting Mr.B and Mr.C.

Edit: Desmos result of the original (not scaled data):

enter image description here

Your help would be really appreciated. Thanks!

Improve this question
$\endgroup$
7
  • 1
    $\begingroup$ I took the freedom to edit the post, mostly to convert words and sentences in upper case to lower case. Nowadays whole sentences in upper case are interpreted as shouting and thus as offending. $\endgroup$ – Henrik Schumacher Jul 2 '20 at 18:19
  • 1
    $\begingroup$ @HenrikSchumacher Thanks a lot sir for editing the post, sorry for using upper case. I do not know if that is interpreted as shouting (I am Arabic). Also thanks for this information. I will use lower case letters next time. $\endgroup$ – Hussain-Alqatari Jul 2 '20 at 18:23
  • $\begingroup$ You may have a look at my answer [here] (mathematica.stackexchange.com/questions/26336/…) $\endgroup$ – Romke Bontekoe Jul 2 '20 at 19:04
  • $\begingroup$ @RomkeBontekoe Can you please tell me how is that related to my problem? $\endgroup$ – Hussain-Alqatari Jul 2 '20 at 19:35
  • 2
    $\begingroup$ I understand that you want to stick with your model, but please take this as constructive criticism: I don't think the images you show really support the idea that this model fits your data. The line vaguely captures the points, but the trend of the data is really not captured very well. The line is too flat on the right side, for example. I honestly think you need to reconsider this model. $\endgroup$ – Sjoerd Smit Jul 3 '20 at 13:39
3
$\begingroup$

With

xmin = Min@data[[All, 1]];
xmax = Max@data[[All, 1]];
ymin = Min@data[[All, 2]];
ymax = Max@data[[All, 2]];

After scaling

datanew = Table[{(data[[i, 1]] - xmin)/(xmax - xmin), (data[[i, 2]] - ymin)/(ymax - ymin)}, {i, 1, Length[data]}]

model = a + b*Exp[c*(Abs[x^d + e])^f];
modelk[k_] := a + b*Exp[c*(Abs[datanew[[k, 1]]^d + e])^f] - datanew[[k, 2]]
res = Sum[Abs[modelk[k]], {k, 2, Length[datanew]}];
fit = NMinimize[res, {a, b, c, d, e, f},  Method -> "DifferentialEvolution"]
model0 = model /. fit[[2]];
plot1 = ListPlot[datanew, PlotStyle -> Red];
plot2 = Plot[model0, {x, 0.001, 1}, PlotStyle -> {Thick, Blue}, PlotRange -> All];
Show[plot1, plot2]

enter image description here

NOTE

Considering the data set

data = {{1, 5}, {3, 6}, {5, 8}, {7, 12}, {9, 20}, {11, 28}, {13, 32}, {15, 34}, {17, 35}, {25, 36}};

and after scaling as we can't use the first scaled point, we include an additional near point as follows:

data = {{1, 5}, {1.1, 5}, {3, 6}, {5, 8}, {7, 12}, {9, 20}, {11, 28}, {13, 32}, {15, 34}, {17, 35}, {25, 36}};

and then we proceed without additional restrictions

xmin = Min@data[[All, 1]];
xmax = Max@data[[All, 1]];
ymin = Min@data[[All, 2]];
ymax = Max@data[[All, 2]];
datanew = Table[{(data[[i, 1]] - xmin)/(xmax - xmin), (data[[i, 2]] - ymin)/(ymax - ymin)}, {i, 1, Length[data]}];
modelf = a + b*Exp[c*(x^d + e)^f];
modelkf[k_] := a + b*Exp[c*(datanew[[k, 1]]^d + e)^f] - datanew[[k, 2]]
res = Sum[Abs[modelkf[k]]^4 , {k, 2, Length[datanew]}];
fit = NMinimize[res Length[datanew], {a, b, c, d, e, f}, Method -> "DifferentialEvolution"]
modelf0 = modelf /. fit[[2]];
plot1 = ListPlot[datanew, PlotStyle -> Red];
plot2 = Plot[modelf0, {x, 0, 1}, PlotStyle -> {Thick, Blue}, PlotRange -> All];
Show[plot1, plot2]

enter image description here

Improve this answer
$\endgroup$
11
  • $\begingroup$ Yes, better (+) $\endgroup$ – yarchik Jul 3 '20 at 15:33
  • 2
    $\begingroup$ This is certainly a better fit and the OP ought to consider it. But a problem is that you no longer have the same model. In this case standardizing the predictor variable results in a different model. You can't reconstruct the parameters of the original model from the fitted model. (Of course, maybe I'm wrong.) $\endgroup$ – JimB Jul 3 '20 at 16:03
  • $\begingroup$ The introduction of Abs[] into the exponentiation was done in order to avoid complex solutions. After that, the model without Abs[] behaves well as expected. $\endgroup$ – Cesareo Jul 3 '20 at 16:32
  • $\begingroup$ Also the first element modelk[[1]] is not considered due to problems handling the origin $(0,0)$ from datanew. $\endgroup$ – Cesareo Jul 3 '20 at 16:35
  • $\begingroup$ Your model is different from mine. See this example, desmos shows i.imgur.com/IQTfeJS.png ... but your way shows i.imgur.com/1rfVhSZ.png .. I believe you can correct what is wrong. Thanks in advance. $\endgroup$ – Hussain-Alqatari Jul 3 '20 at 16:44
2
$\begingroup$

I think your model is much too complex for the data. But with the right starting values and a specific restriction on one of the parameters, Mathematica provides a much better fit than DESMOS. And to my surprise the fit seems fine in terms of predicting but I have no confidence in the estimate of any specific parameter (or even their signs) because of the high correlations in the parameter correlation matrix (with all correlations being very close to -1 or +1).

data = {{595098, 335}, {595158, 336}, {595218, 336}, {595338, 
    344}, {595458, 347}, {595638, 351}, {595818, 352}, {596082, 
    356}, {596322, 357}, {596922, 362}, {597521, 363}, {598481, 
    363}, {599322, 371}, {600523, 376}, {601723, 378}, {603523, 
    380}, {605323, 384}, {608924, 390}, {612523, 392}, {619724, 
    399}, {626926, 404}, {648527, 413}, {670129, 415}, {691731, 
    419}, {712906, 424}, {734504, 425}, {756104, 424}, {776690, 
    426}, {798291, 426}, {819890, 429}, {841490, 431}, {863090, 
    432}, {884692, 433}, {906290, 434}, {927892, 437}, {949492, 
    438}, {971090, 437}, {992691, 437}, {1014291, 438}, {1035891, 
    442}, {1039491, 441}, {1041236, 442}};

(* Using approximately your parameter estimates as starting values *)
{a0, b0, c0, d0, f0, g0} = {495, -54, 0.2, 0.21, -15.83, -1.3}

nlm = NonlinearModelFit[data, {a + b Exp[c (x^d + f)^g], f > -16}, 
  {{a, a0}, {b, b0}, {c, c0}, {d, d0}, {f, f0}, {g, g0}}, x, MaxIterations -> 10000];

Warning message

Even after 10,000 iterations there is no convergence. (And, there's still no convergence after 50,000 iterations.) However, current estimates are provided.

nlm["BestFitParameters"]
(* {a -> 1170.12, b -> -376.038, c -> 0.65887, d -> 0.156187, f -> -7.97678, g -> -0.0292045} *)

nlm["RSquared"]
(* 0.999973 *)

Show[ListPlot[data], 
  Plot[nlm[x], {x, Min[data[[All, 1]]], Max[data[[All, 1]]]}, PlotRange -> All]]

data and fit

The correlation matrix spells trouble:

nlm["CorrelationMatrix"] // MatrixForm

Parameter correlation matrix

We have (estimated) perfect correlation among several coefficients and nearly perfect correlation (close to -1 or close to +1) for all of the others. This can be a symptom of an overparameterized model given the available data.

I short the combination of the available data and the particular model results in an unstable estimation process for any piece of software. If DESMOS allows starting values, you might want to consider trying the "better" estimates from Mathematica.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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