6

I asked a question in a previous post that was closed because "The title of the question is not correct and the issue here is too trivial to help anyone later".

https://mathematica.stackexchange.com/questions/19297/nonlinearmodelfit

It's possible because I'm not so smart in Mathematica. But in this period I've tried a lot of alternative solutions.

Recall my problem

 data={{220., 542.821}, {225., 531.898}, {230., 521.391}, {235., 
 510.596}, {240., 499.603}, {245., 488.585}, {250., 479.002}, {255., 
 468.574}}
 func= a0*Exp[a1 t + a2 t^2]
 nlm = NonlinearModelFit[data, func , {a0, a1, a2}, t]

which gives the result: 

FittelModel = 0.

Some answers told me that it is not an exponential, but linear problem. But this problem is a physical problem that comes out from a lot of papers that describes this phenomenum as having exponential behavior. So also if my question was trivial, I've continued to solve it.

I read a lot of posts, and I tried this solution in Mathematica, using NMinimize and not NonLiNearModelFit.

f[a0_, a1_, a2_][t_] := a0*Exp[a1*t + a2*t^2 ];
error[a0_, a1_, a2_][x_, y_] := (f[a0, a1, a2]@x - y);
errorSum[a0_, a1_, a2_] = Total[(error[a0, a1, a2] @@@ data)];
bestFitVals = 
NMinimize[{errorSum[a0, a1, 
a2], {600 < a0 < 1500&& -0.01 < a1 < 0 && -0.0000001 < a2 < 0}}(*,
Thread[0<{a0,a1,a2}]*), {a0, a1, a2}, Method -> "SimulatedAnnealing"]
soln[x_] = f[a0, a1, a2]@x /. bestFitVals[[2]]
limits = Flatten@{x, Through[{Min, Max}@data[[All, 1]]]};
Show@{Plot[soln@x, Evaluate@limits, PlotRange -> All], 
ListPlot[data, PlotStyle -> Red]}

I tried all the methods available to NMinimize, changing starting values or removing them, but no choice works. So finally a collegue suggested that I try another software package: SigmaPlot. I said him that Mathematica is superior, but Sigmaplot gives me an exponential solution, with these coefficients

a0 = 1072.637060511283
a1 = -2.128998307513597e-3
a2 = -4.390602833555603e-6.

and with residuals that are very low.

I tried with another set of data, but the solution is the same: Mathematica is not be able to regress or minimize and Sigmaplot is able to.

I hope my question is not trivial, but I hope more that someone can solve my question: is Sigmaplot superior of Mathematica, or is it the case that I'm not applying Mathematica properly to my problem?

Community
  • 1
Mary
  • 629
  • 3
  • 10
  • 2
    So, it would appear this is an issue of giving reasonable starting values. Obviously, the "preferred" value for a0 is large, while a1 and a2 are close to zero. If you don't specify initial values, Mathematica just chooses 1 by default, hence the difficulty you encountered. Try nlm = NonlinearModelFit[data, func, {{a0, 1000}, {a1, 0}, {a2, 0}}, t] instead and you will get the same result as from SigmaPlot. – Oleksandr R. Feb 26 '13 at 19:16
  • 1
    @s0rce: the data are the same I put inside Sigmaplot and the residuals of exponential fitting is very low compared with linear fit. – Mary Feb 26 '13 at 19:19
  • @Oleksandr R: Very good idea. – Mary Feb 26 '13 at 19:22
  • Maybe nlm = NonlinearModelFit[data, func, {{a0, 10}, {a1, 0.01}, {a2, 0.0001}}, t]; Show[ ListPlot[data, PlotStyle -> Directive[Blue, PointSize[Large]]], Plot[nlm[t], {t, 220, 255}], Frame -> True] and nlm[{"BestFit", "ParameterTable"}] – user1066 Feb 26 '13 at 20:01
  • Are you sure your exponential function is a0*Exp[a1 t + a2 t^2] and not, say, a0*Exp[-(a1 t + a2 t^2)]? – user1066 Feb 26 '13 at 20:17
  • @TomD: I can try! Thanks – Mary Feb 26 '13 at 21:12
  • 2
    Unless you are confident the error variance is constant across different values of $x$, you should instead be fitting the model $\log(y) = \log(a_0) + a_1 t + a_2 t^2$. Use LinearModelFit; it's fast, easy, and can output a huge amount of additional information. – whuber Feb 27 '13 at 00:33
  • @whuber I am a complete novice when it comes to fitting data (at least the theoretical side of it), so, other than speed, why is the constant error variance important in the exponential fit? Doesn't linear fitting assume the same thing? – rcollyer Feb 27 '13 at 04:20
  • @whuber: can you suggest me some link where I can see this important properties – Mary Feb 27 '13 at 08:35
  • 1
    @rc Constant error variance is assumed by (unweighted) least-squares fitting. When that variance changes, intuitively it makes sense to put less weight on the points with larger variance--and in fact using appropriately weighted least squares is one way to fix up the procedure. But there are great advantages to re-expressing the data to produce constant error variances: interpretation is simpler, prediction is simpler and more reliable, and in physical applications often the re-expressions have meaning. These issues are extensively discussed on http://stats.stackexchange.com, Mary. – whuber Feb 27 '13 at 14:07
  • @whuber So, if I'm reading into this correctly, because it is easier to perform a weighted linear fit than a weighted non-linear fit, converting to the linear form is superior if we don't have constant variances, correct? (btw, rc doesn't ping me. (: ) – rcollyer Mar 05 '13 at 13:18
  • @rcollyer Ease of fitting can be a concern, but achieving constancy of variance and/or a linear model are more fundamental than that, because they are related to how one interprets the model and the quality of predictions or estimates made with it. Philosophy (such as Occam's Razor) and experience are part of the basis for such advice, which cannot be fully justified mathematically (nor contradicted mathematically, either). BTW, I am surprised that "_at_rc" doesn't work: pinging is supposed to recognize prefixes. Thanks for the heads-up about that. – whuber Mar 05 '13 at 15:33

1 Answers1

3

As mentioned in the comments, the default initial parameter (1) for NonlinearModelFit is not appropriate for this question, and initial parameters need to be supplied:

nlm = NonlinearModelFit[data, func, {{a0, 1000}, {a1, 0}, {a2, 0}}, t]
nlm["BestFitParameters"]
(* {a0 -> 1072.64, a1 -> -0.002129, a2 -> -4.3906*10^-6} *)

Apparently SigmaPlot has some mechanism for automatically determining initial parameters but those are contained in libraries that are not immediately accessible to those who don't own a copy of SigmaPlot.

I can only speculate on the rules that SigmaPlot uses to determine initial parameters, and since it uses a library of functions (see page 867 of this long PDF), there is likely no one-shoe-fits-all approach.

Instead of initializing all parameters to 1, another trivial approach would be to initialize all parameters to 0. Oddly enough, this approach works in this instance:

nlm = NonlinearModelFit[data, func, {{a0, 0}, {a1, 0}, {a2, 0}}, t]
nlm["BestFitParameters"] 
{a0 -> 1072.64, a1 -> -0.002129, a2 -> -4.3906*10^-6}
bobthechemist
  • 19,693
  • 4
  • 52
  • 138