# Output of NonlinearModelFit differs from the correct result

I'm having a bad time dealing with the NonlinearModelFit in Mathematica 8, since the result given is a bit imprecise. An example is given on potential regression, as follows:

data = {{0.9, 1.1}, {2, 2}, {6, 4.1}, {6.9, 4.5}, {9.9, 5.5},
{10.7, 5.9}, {14, 6.7}, {15.9, 7.3}}

NonlinearModelFit[data, b*t^a, {a, b}, t]


The result I get is: 1.34232*t^0.615024, however it should be 1.23...*t^0.6548.

What am I doing wrong? I experience the same thing when I try to do a exponential regression (like this: NonlinearModelFit[data, b*a^t, {a, b}, t]).

-

## migrated from stackoverflow.comApr 30 '12 at 19:39

This question came from our site for professional and enthusiast programmers.

Why do you think "[..] the result should be 1.23...*t^0.6548 [..]" in the first place? The results you get for a non-linear fit depend on many factors, for instance the optimization/minimization algorithm used, starting values, convergence criteria etc. – Frank Apr 30 '12 at 17:05
@FrankNiemeyer I agree, i just plotted the data and fit; seem "pretty good" to me. There's a lot of details in there to worry about... – tkott Apr 30 '12 at 19:56
Well, it was because of the test we made, where I put the first result as the answer, and I got an error. However, I've spoken to my teacher about the fact that Mathematic presents a more valid answer (than the calculator TI-84) and I'll be able to put that as an answer in the future, without getting an error. Thanks for your help! – Frederik Brinck Jensen May 1 '12 at 7:17

This is the fit returned by Mathematica :

nlm = NonlinearModelFit[data, b*t^a, {a, b}, t]


The model is

nlm[t]

(* 1.34232 t^0.615024 *)


The residuals (model(x)-y) are :

nlm["FitResiduals"]

(* {-0.158097,-0.0558756,0.059473,0.0967975,0.00211152,0.132971,-0.103818,-0.0577413} *)


The sum of the squares of the residuals is :

res1 = Total[nlm["FitResiduals"]^2]

(* 0.0728215 *)


An explicit check :

Total[(nlm[#[[1]]] - #[[2]])^2 & /@ data]

(* 0.0728215 *)


The sum of the squares of the residual relative to the other model is

y[t_] = 1.23 t^0.6548

res2 = Total[(y[#[[1]]] - #[[2]])^2 & /@ data]

(* 0.152818 *)


Since res1 < res2 I would say that Mathematica returns a better fit. Why do you say the best fit is the second function ?

-
Would you please elaborate a little on your answer, since I'm not yet that tough going on Mathematica? – Frederik Brinck Jensen Apr 30 '12 at 15:26
@FrederikBrinckJensen Tried to be more clear. – b.gatessucks Apr 30 '12 at 15:35
Thanks b.gatessucks. The thing is, that I'm using Mathematica for my exams coming up at high school, however, the results I get – though they're closer – are not valid in terms of correction. So I usually miss some important points. Is it only due to the fact that Mathematica is more precise, that I get a different result. It has nothing to do with the function I'm regression as exponential (ba^t) or potential (bt^a)? Thanks for your patience. – Frederik Brinck Jensen Apr 30 '12 at 16:05
@FrederikBrinckJensen what do you mean "are not valid in terms of correction"? – tkott Apr 30 '12 at 20:00
I meant that the result didn't correspond with the right answer, however as I said above, I presented the results to my teacher and she will let me use this result in the future. – Frederik Brinck Jensen May 1 '12 at 7:19

You will want to run:

FindFit[data, b*t^a, {a, b}, t, NormFunction -> (Norm[#, 1] &)]

{a -> 0.66955, b -> 1.20678}


Normally in statistics you try to reduce the 2-Norm. In some cases however you want to reduce the 1 norm of the residuals instead. Reducing the 1-Norm is more robust when there are large outliers that would have undue influence on the values of the parameters.

Long story short, Mathematica is giving you the correct answer. The answer above is just another correct answer and the one you wanted to get.

Some software packages will choose this norm instead of the 2-norm for you. That's irresponsible silly

-
can you specify the norm function for the nonlinearfit? BTW, I love that this is a 3rd "Correct" answer :) mmmm fitting – tkott Apr 30 '12 at 20:55
Nope. NonlinearModelFit doesn't have this option. You'll want to use FindFit for more general problems such as this. – Searke Apr 30 '12 at 21:28
Hmm, ok. Thanks! I always considered NonlinearModelFit to be more general than FindFit – tkott Apr 30 '12 at 21:51
In the end, they are all just fancy ways of minimizing a function. – Searke Apr 30 '12 at 22:03

Ah, I know where the "correct" answer came from:

Exp[a] u^b /. FindFit[Log[data], a + b u, {a, b}, u]
1.230270887036395*u^0.6547870743058625


As has already been shown to you, the results of a true nonlinear fit (which is what Mathematica does through the Levenberg-Marquardt algorithm) are almost always better than the results obtained through an initial linearization of the model, which in this case is $\log\,y=\log\,a+b\log\,x$. This is because the application of the logarithm distorts the errors in the data as well, so in truth you are not doing a true "least-squares" fit. There is a way to modify the linearization approach to make things slightly better, which I discussed in this math.SE answer, but since FindFit[] and NonlinearModelFit[] are already there, it is preferable to use them.

-