# Fitting to a complicated implicit function

I have this data:

data = {{290, 0.112089478}, {240, 0.072525194}, {220,
0.063863975}, {200, 0.075709626}, {170, 0.076462613}, {150,
0.070010656}, {130, 0.065390206}, {120, 0.064750065}, {100,
0.134944754}, {91, 0.257195474}, {78, 0.470626063}, {60,
0.45755743}, {40, 0.358068721}, {30, 0.307480383}};


and want to fit it using a function of this form:

a x ((4 Γ^2 + y^2)/y^6)


My only variable is x, y is a function of x and the dependency is given implicitly by:

y == Sqrt[1 + β/y (1/(Exp[y/x] - 1) + 1/2 - y/2)]


The parameters that I want to fit are {a, Γ, β}

And this is what I have done so far:

fitfunc[a_?NumericQ, Γ_?NumericQ, β_?NumericQ, x_?NumericQ] :=
y /. FindRoot[
a x ((4 Γ^2 + y^2)/y^6) /.
{y == Sqrt[1 + β/y (1/(Exp[y/x] - 1) + 1/2 - y/2)]}, {y, 1.}
];
FindFit[data, fitfunc[a, Γ, β, x], {a, Γ, β}, x]

• Is it possible that there's a typo in the FindRoot[...] part? Currently, there's an error - it seems like ((4 Γ^2 + y^2)/y^6) /. should go outside the FindRoot expression Sep 14, 2017 at 21:41
• Yeah, it is possible. Actually, I tried this one  {fitfunc[a_?NumericQ, [CapitalGamma]_?NumericQ, x_?NumericQ] := a x ((4 [CapitalGamma]^2 + y^2)/y^6) /. FindRoot[ y == Sqrt[1 + [Beta]/y (1/(Exp[y/x] - 1) + 1/2 - y/2)], {y, 10}];} but it still gives me error Sep 15, 2017 at 0:32
• What is the magnitude order of these parameters? If their values are far from the default value 1.0, the function cannot be easily fitted. Sep 15, 2017 at 3:30

The wildly different orders of magnitude pointed out by @robsondenke is one of the problems. Another is that x and y are highly correlated. But the biggest problem is that the model just doesn't fit the data.

First I reformulate the model so that it is a bit more stable with respect to finding the maximum likelihood estimates of the parameters. And I use NonlinearModelFit as one can extract many more goodies than FindFit with few differences in input.

data = {{290, 0.112089478}, {240, 0.072525194}, {220, 0.063863975},
{200, 0.075709626}, {170, 0.076462613}, {150, 0.070010656},
{130, 0.065390206}, {120, 0.064750065}, {100, 0.134944754},
{91, 0.257195474}, {78, 0.470626063}, {60, 0.45755743},
{40, 0.358068721}, {30, 0.307480383}};

fitfunc[b_?NumericQ, c_?NumericQ, β_?NumericQ, x_?NumericQ] :=
x (10 b + (c/100000) y^2)/y^6 /. FindRoot[y == Sqrt[1 + β/y (1/(Exp[y/x] - 1) + 1/2 - y/2)], {y, 1}];

sol = NonlinearModelFit[data, {fitfunc[b, c, β, x], {b > 0, c > 0, β > 0}}, {{b, 8}, {c, 1.3}, {β, 15}}, x];
sol["BestFitParameters"]
(* {b -> 8.176014208713987,c -> 7.773373668159528,β -> 14.584722480599252} *)

(* Original parameters *)
Solve[{a 4 Γ^2 == 10 b, a == c/100000} /. sol["BestFitParameters"], {a, Γ}][]
(* {a -> 0.00007773373668159529,Γ -> 512.7858713433279} *)

Show[ListPlot[data],
Plot[fitfunc[b, c, β, x] /. sol["BestFitParameters"], {x, 1, 300}]]
` It's a pretty bad fit. Is there some theoretical reason why this particular function?

• "[...] the equation that I need to fit my data into a y Log[1 + [Beta]/(y^3 + ([Beta]/2 - 1) y - [Beta]/2)] ((4 Γ^2 + y^2)/y^6)"
– Kuba
Sep 15, 2017 at 7:46