1
$\begingroup$

If I have the following data (I have updated the data):

  data={{{60, 1852.94}, {65, 178.035}, {70, 7.97143}, {75, 48.9479}, {80, 
  133.561}, {85, 8.65079}, {87, 1.78915}}

How can I fit it to the following equation?

eq = 1/((x*c1)/((Exp[B1/(x - T0)])^0.75)*(1 - Exp[((Tm - x)*c2)/Tm^2])) /. {T0 -> 259.246,B1 -> 2595.89, Tm -> 88.2 + 273.15}

I am trying the following:

  fun[x_] = NonlinearModelFit[data,eq, {c1, c2}, x] // Normal (*c1=732975, c2=2.65721*10^6*)

Which gives me an erroneous c1 and c2 of c1=732975 & c2=2.65721*10^6. When I do it in excel it gives me a c1 of -65514626.34 and a c2 of 68.4 which fits the data relative well. Why mathematica is not fitting the data or finding the correct c1 and c2 values?

$\endgroup$
6
  • 1
    $\begingroup$ Try to plot this: Show[{ ListPlot[data], Plot[(E^(2595.89/(-259.246 + x)))^0.75/(c1 (1 - E^(5.092118394960718*^-6 c2 (443.15 - x))) x) /. {c1 ->313634, c2 -> -8434.47}, {x, First[data][[1]], Last[data][[1]]}, PlotRange -> All]. It shows that the values of c1 and c2 you obtained from excel yield a very bad fit, if any. $\endgroup$ Commented May 21, 2021 at 18:13
  • $\begingroup$ @AlexeiBoulbitch thank you for you comment. Then what would the best fit be? and why Mathematica cannot find it? $\endgroup$
    – John
    Commented May 21, 2021 at 19:05
  • $\begingroup$ Often such things indicate that the model is not optimal. I would try to think about another to fit the data. $\endgroup$ Commented May 21, 2021 at 19:23
  • $\begingroup$ @AlexeiBoulbitch I have updated the data. I think the reason it was not working is because the y axis of the data was in log scale. I have updated the data now. I think now the c1 and c2 values that I found in excel are the ones that fit the data but still Mathematica does not find it $\endgroup$
    – John
    Commented May 21, 2021 at 19:29
  • $\begingroup$ I do not see that. Try to look at the plot (the code from my first comment) with your new data. $\endgroup$ Commented May 21, 2021 at 21:04

1 Answer 1

3
$\begingroup$

Given that the response and predictor variables are both positive and structure of the model, one of the coefficients (c1 and c2) needs to be negative and the other needs to be positive. Excel's answer at least has that property but I don't see why you say that the Excel answer fits the data relatively well.

data = {{60, 1852.94}, {65, 178.035}, {70, 7.97143}, {75, 48.9479}, {80, 133.561}, 
  {85, 8.65079}, {87, 1.78915}};

eq = 1/((x*c1)/((Exp[B1/(x - T0)])^0.75)*(1 - Exp[((Tm - x)*c2)/Tm^2])) /. 
  {T0 -> 259.246, B1 -> 2595.89, Tm -> 88.2 + 273.15};

excel = {c1 -> -65514626.34, c2 -> 68.4}; 
Show[ListPlot[data, PlotRange -> All], 
 Plot[eq /. excel, {x, 60, 87}, PlotStyle -> Red]]

Data and Excel fit

If you start off with estimates of the coefficients with opposite signs, then NonlinearModelFit finds a better fit.

nlm = NonlinearModelFit[data, eq, {{c1, -1}, {c2, 1}}, x, MaxIterations -> 10000];
nlm["BestFitParameters"]
(* {c1 -> 0.000589194, c2 -> -0.000841265} *)
nlm["EstimatedVariance"]^0.5
(* 543.51 *)
Show[ListPlot[data, PlotRange -> All], Plot[nlm[x], {x, 60, 87}]]

Data and better fit

Because the estimators of the coefficients given the data and model are nearly perfectly correlated (look at nlm["CorrelationMatrix"]), there is a great deal of numerical instability. In fact, just setting c2 to some arbitrary number and fitting with just c1 gives an even better fit:

nlm = NonlinearModelFit[data, eq /. c2 -> 1, {{c1, -1}}, x, MaxIterations -> 10000];
nlm["BestFitParameters"]
(* {c1 -> -4.95197*10^-7} *)
nlm["EstimatedVariance"]^0.5
(* 496.162 *)
Show[ListPlot[data, PlotRange -> All], Plot[nlm[x], {x, 60, 87}]]

Even better fit

But all fits are horrible. To paraphrase Bullwinkle the Moose: "You need a better 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.