# Fitting to the implicit function

I already asked questions about this issue here but there have been some changes that I want to ask a new question.

My data is :

data = {{290, 0.012263719}, {240, 0.007675481}, {220,
0.008038809}, {200, 0.008608707}, {170, 0.010805872}, {150,
0.008832903}, {130, 0.009263129}, {120, 0.011290667}, {100,
0.014344114}, {40, 0.025720622}, {30, 0.028876792}, {20,
0.035088327}, {15, 0.042578946}, {4.2, 0.033039767}};


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

a x ((4 Γ^2 + (y^2-R))/(y^2-R)^3)


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

y == \[Omega]
Sqrt[1 + \[beta] \[Omega]/ y (1/(Exp[y/x] - 1) + 1/2 - y/(2 \[Omega]))]


The parameters that I want to fit are {a, Γ, β,\[Omega],R} and here is my code:

fitfunc[b_?NumericQ, c_?NumericQ, \[Beta]_?NumericQ,
R_?NumericQ, \[Omega]_?NumericQ, \[CapitalGamma]_?NumericQ,
x_?NumericQ] :=
x (10 b + (c/100000) (y^2 - R))/(y^2 - R)^3 /.
FindRoot[
y == \[Omega] Sqrt[1 + \[Beta] \[Omega]/y (1/(Exp[y/x] - 1) + 1/2 - y/2)], {y,
1}];
sol = NonlinearModelFit[
data, {fitfunc[b, c, \[Beta], R, \[Omega], \[CapitalGamma],
x], {b > 0, c > 0, \[Beta] > 0,
R > 0, \[Omega] > 0, \[CapitalGamma] > 0}}, {{b, 8}, {c,
1.3}, {\[Beta], 57}, {\[Omega], 23}, {\[CapitalGamma], 4}}, x];
sol["BestFitParameters"]
(*{b\[Rule]12800,c\[Rule]200000000,\[Beta]\[Rule]57,R \[Rule]3,\
\[Omega] \[Rule]23,\[CapitalGamma]\[Rule]4}*)
(*Original parameters*)
Solve[{a 4 \[CapitalGamma]^2 == 10 b, a == c/100000} /.
sol["BestFitParameters"], {a, \[CapitalGamma]}][]
(*{a\[Rule] 2000,\[CapitalGamma]\[Rule]4}*)

Show[ListPlot[data],
Plot[fitfunc[b, c, \[Beta], R, \[Omega], \[CapitalGamma], x] /.
sol["BestFitParameters"], {x, 1, 300}]]


I know that for this set {a=2000, Γ=4, β=57,\[Omega]=23 ,R=3}} could be one answer ( I found them by changing parameters with hand):

data = {{290, 0.012263719}, {240, 0.007675481}, {220,
0.008038809}, {200, 0.008608707}, {170, 0.010805872}, {150,
0.008832903}, {130, 0.009263129}, {120, 0.011290667}, {100,
0.014344114}, {40, 0.025720622}, {30, 0.028876792}, {20,
0.035088327}, {15, 0.042578946}, {4.2, 0.033039767}};
Show[ListPlot[data],
Plot [{ 2 10^3 x  (4 4^2 + (y^2 - R^2 ))/(y^2 - R^2 )^3 /.
R -> 3} /.
FindRoot[
y == Subscript[\[Omega], 0]
Sqrt[1 +
57 Subscript[\[Omega], 0]/
y (1/(Exp[y/x] - 1) + 1/2 -
y/(2 Subscript[\[Omega], 0]))] /. {Subscript[\[Omega],
0] -> 23}, {y, 1}], {x, 1, 300}], PlotRange -> {0, 0.05}] I don't know why it couldn't fit and find this set at least.

• The code you've written doesn't line up with the problem you've stated. Where did b and c come from? Also a x ((4 Γ^2 + (y^2-R)/(y^2-R)^3) has mismatched parentheses and contains a potentially degenerate expression (y^2-R)/(y^2-R)^3. Oct 18, 2017 at 2:59
• Well considering the above plot it should have an answer around this set '{a=2000, Γ=4, β=57,[Omega]=23 ,R=3}}' and I guess it shouldn't be degenerate at least around this set. Oct 18, 2017 at 4:34
• You've left off R in the NonlinearModelFit statement and there's no $\Gamma$ in the definition of fitfunc.
– JimB
Oct 18, 2017 at 4:48
• Clarification: You've left off R in the parameter and starting value list in the NonlinearModelFit statement and there's no $\Gamma$ in the definition of fitfunc.
– JimB
Oct 18, 2017 at 5:26

I get a dissimilar set of values that fit adequately well by developing an objective function and using NMinimize on it:

fitfunc[a_?NumericQ, \[Beta]_?NumericQ, \[CapitalGamma]_?NumericQ, R_?NumericQ, \[Omega]_?NumericQ, x_?NumericQ] := a x ((4 \[CapitalGamma]^2 + (y^2 - R))/(y^2 - R)^3)
/. FindRoot[y == \[Omega] Sqrt[1 + \[Beta] \[Omega]/y (1/(Exp[y/x] - 1)) + 1/2 - y/(2 \[Omega])], {y, 1}];

of = Total[(Norm@fitfunc[a,\[Beta],\[CapitalGamma],R,\[Omega],#]&/@(Transpose[data][])-Transpose[data][])^2];

mres = NMinimize[{of,a>0,\[Beta]>0,\[CapitalGamma]>0,R>0,\[Omega]>0},{a,\[Beta],\[CapitalGamma],R,\[Omega]}];

Show[ListPlot[data],Plot[fitfunc[a,\[Beta],\[CapitalGamma],R,\[Omega],x]/.mres[],{x,1,300}, PlotRange->{0,0.05}]]


NMinimize has trouble with it too, but at least on my machine it has a much more reasonable fit than NonlinearModelFit gets. The Norm in the objective function is to avoid issues caused by complex numbers appearing during NMinimize's operation.

Edit: Since adding the extra division pointed out in the comments, the objective function above may fail to fit the first few points as well. If that's an issue, consider replacing the entire line with something like:

of = Total[(1000/Transpose[data][])(Norm@fitfunc[a,\[Beta],\[CapitalGamma],R,\[Omega],#]&/@(Transpose[data][])-Transpose[data][])^2];


The (1000/Transpose[data][]) term will weight the accuracy of the terms near the origin higher, and so this process is more likely to achieve a good fit there.

• Well there is a mismatch in your code with mine, Oct 18, 2017 at 5:53
• Two points: there is a mismatch in the equation in your code with mine, 'FindRoot[y == [Omega] Sqrt[1 + [Beta] [Omega] y (1/(Exp[y/x] - 1)) + 1/2 - y/(2 [Omega])]' should be 'FindRoot[y == [Omega] Sqrt[1 + [Beta] [Omega]/ y ((1/(Exp[y/x] - 1) + 1/2 - y/(2 [Omega]))]'. Second, after performing this revision, the plot becomes disconnected. It doesn't say what are the values of these parameters. Oct 18, 2017 at 6:05
• The former is an easy fix. The parameters it fits are the second argument of NMinimize, so {a,\[Beta],\[CapitalGamma],R,\[Omega]}. The parameters are stored as mres[] above. Oct 18, 2017 at 6:08
• @Kassik: Regarding the plot becoming disconnected, adding PlotRange->{0,0.05} to the Plot options helped for me. Plot doesn't normally plot points that would be outside of its graphic, but Show doesn't inform Plot of what the regions of the other graphics are, so Plot still does its own thing. Oct 18, 2017 at 6:20
• 'PlotRange->{0,0.05}' is not fixing it because for the set of data stored in 'mres[]', the plot is not convergent with 'NMinimize' Oct 18, 2017 at 6:47