# (Non)LinearModelFit hangs on "trivial" linear fit (Fit succeeds)

A simplified version of my function and a simplified subset of my data (real data also has weights, if that matters, but I tested both with and without those):

myfunc[x_, a_] = 8*a*Exp[-8/Sqrt[x]]/(Sqrt[x]*Abs[WhittakerW[E*I/Sqrt[x],
1/2, 8/10*I*Sqrt[x]*a]]^2)
mydata = {{4.5, 195}, {2.9, 175}, {2.1, 95}}


myfunc also needs a global rescaling. If I want to fit both parameter a and the scale, NonlinearModelFit works great out of the box and quickly (a Timing on a fresh kernel returned 0.6 seconds):

NonlinearModelFit[mydata, scale*myfunc[x, a], {scale, a}, {x}]


If, however, I only want to find the best scaling for an a value of my choice, this:

NonlinearModelFit[mydata, scale*myfunc[x, 4], {scale}, {x}]


takes forever in my computer (hangs until I abort). I am surprised: why the "complicated" fit works so easily, and this "trivial" linear one not?

LinearModelFit is supposedly the right tool here, but:

LinearModelFit[mydata, {myfunc[x, 4]}, {x}, IncludeConstantBasis -> False]


hangs in the same way (I let it run for about 1.5 hours, then aborted). There must be something wrong I am doing. I tried using 4.0 instead of 4 for a, enclosing with an N[], adding WorkingPrecision -> MachinePrecision, suggesting a good starting value for scale, but I get the same result.

NonlinearModelFit works adding Method -> "NMinimize", although it is slow (Timing gave 19 seconds), and raises a warning

NonlinearModelFit::lmnl: The model (32 E^(-(8/Sqrt[x])) scale)/(Sqrt[x] Abs[WhittakerW[(I E)/Sqrt[x],1/2,(16 I Sqrt[x])/5]]^2) is linear in the parameters {scale}, but a nonlinear method or non-Euclidean norm was specified, so nonlinear methods will be used.

I do not really know what NMinimize is doing here, I guess it applies the "right" simplification to myfunc, so the computation is easier?

Finally, after wasting one day, Fit just works:

Fit[mydata, {myfunc[x, 4]}, {x}]


This reinforces my impression that I just need to apply some magical simplification, that Fit always employs, that NonlinearModelFit uses only for the two-parameter fit, and that Method -> "NMinimize" will cause. Is LinearModelFit evaluating myfunc much more times than needed? I am clueless.

LinearModelFit has several niceties that I would rather not give up for what appears to be such a stupid issue. Furthermore, I am probably going to learn a new important Mathematica quirk from this. Thus:

Why the above code, using (Non)LinearModelFit for a simple scale*f[x] model, hangs, and how should I modify it to make it work?

For completeness, I am adding some tests on my system (12.2.0.0) after @xzczd suggestions in the comments. This one:

Clear[undefinedname]
FindFit[{undefinedname}, b WhittakerW[I, 1, I], b, {x}]


hangs, regardless of whether some data is passed. LinearModelFit[{undefinedname},...] instead won't get tricked and stop immediately. However, passing real data, those do not hang:

FindFit[mydata, b WhittakerW[I, 1, I * 1.0], ... ]
LinearModelFit[mydata, Abs[WhittakerW[I, 1, I * 1.0]], ... ]
LinearModelFit[mydata, Abs[WhittakerW[I, 1, I]], ...] (*5x slower than the previous one*)

NonlinearModelFit[mydata, scale*myfunc[x, 4], {scale}, {x},
Method -> "ConjugateGradient"] (*Hangs with InteriorPoint*)


It sounds like both the linear fitter and InteriorPoint get confused when they cannot explicitly evaluate WhittakerW (as it always is when the function depends on x).

• Looks like a bug to me. LinearModelFit doesn't hang in v9. A simpler example that hangs: Clear[aha, myfunc]; myfunc[x_?NumericQ, a_?NumericQ] := WhittakerW[I, 1, I a]; FindFit[{aha}, scale myfunc[x, 1], scale, {x}]. BTW, this sampe doesn't hang in v9.0.1, but if I change myfunc[x, 1] to myfunc[x, 3.47303], FindFit hangs. (Yeah myfunc[x, 3.47302 doesn't hang in v9. ) Have you reported this to WRI? Jan 6 at 15:48
• This seems to be related to the method InteriorPoint. All other available methods give an output within a minute: {"Gradient", "ConjugateGradient"(*, "InteriorPoint"*), "QuasiNewton", "Newton", "NMinimize", "LevenbergMarquardt"}. Jan 6 at 16:09
• @xzczd Very interesting. It took me a while before realizing you really meant to pass an undefined {aha}, but it does hang! I updated the question with a couple of tests. Good catch about simply changing method, although I have to give up linear fitting that way. No I haven't reported it yet, I was assuming there was some simple obscure option to pass to solve the problem. I will now. Jan 6 at 17:29

I hope that your real data consists of more than 3 data points when fitting a function with 3 parameters (a, b, and an error variance).

While I don't know why the functions hang, for your function there is a simple fix. Because myfunc[x, 4] does not change with the scaling parameter, just replace the predictor with the values of the function. Then both NonlinearModelFit and 'LinearModelFit work fine.

myfunc[x_, a_] = 8*a*Exp[-8/Sqrt[x]]/(Sqrt[x]*
Abs[WhittakerW[E*I/Sqrt[x], 1/2, 8/10*I*Sqrt[x]*a]]^2)

(* Replace the predictor with the result of the function *)
mydata = {{4.5, 195}, {2.9, 175}, {2.1, 95}};
mydata[[All, 1]] = myfunc[#, 4] & /@ mydata[[All, 1]];

(* Fit the model to the data *)
NonlinearModelFit[mydata, b x, {{b, 16}}, x] // Normal
(* 18.4684 x *)
LinearModelFit[mydata, x, x, IncludeConstantBasis -> False] // Normal
(* 18.4684 x *)
`