4
$\begingroup$

Here is my attempt.

data = {{0., 2.61}, {0.1, 2.62}, {0.2, 2.62}, {0.3, 2.62}, {0.4, 
2.63}, {0.5, 2.63}, {0.6, 2.74}, {0.7, 2.98}, {0.8, 3.66}, {0.9, 
5.04}, {1., 7.52}, {1.1, 10.74}, {1.2, 12.62}, {1.3, 10.17}, {1.4,
 5}, {1.5, 2.64}, {1.6, 11.5}, {1.65, 35.4}};
NonlinearModelFit[data,a*Cosh[b*x^c*Sin[d*x^e]], {{a, 3}, {b, 0.2}, {c, 1}, {d, 1}, 
{e, 3}}, x, Method -> "LevenbergMarquardt"]

NonlinearModelFit::nrjnum: The Jacobian is not a matrix of real numbers at {a,b,c,d,e} = {3.,0.2,1.,1.,3.}.

The same issues appear with other initial points. As far as I know it, a good fit with the residuals of quantity $0.01$ exists and can be found with Mathcad (not MATLAB).

Improve this question
$\endgroup$
9
  • 2
    $\begingroup$ The Jacobian contains Log[x] and your data contains x == 0 $\endgroup$ – Coolwater Jan 28 '19 at 17:09
  • $\begingroup$ @Coolwater: There are other methods of fitting. $\endgroup$ – user64494 Jan 28 '19 at 17:43
  • $\begingroup$ It's all about providing good starting values when there are 7 parameters (a, b, c, d, e, and the error variance) to estimate with only 18 (or 17 legitimate) limited precision data points. Following @Coolwater 's approach gets you good starting values for use in NonlinearModelFit. $\endgroup$ – JimB Jan 28 '19 at 18:31
  • $\begingroup$ @JimB: I tried you advice, obtaining "NonlinearModelFit::nrjnum: The Jacobian is not a matrix of real numbers at {a,b,c,d,e} = {2.75699,1.66039,2.92883,4.74005,0.70619}.". $\endgroup$ – user64494 Jan 28 '19 at 18:37
  • $\begingroup$ Drop the observation with x=0 and it will work just fine. $\endgroup$ – JimB Jan 28 '19 at 18:38
7
$\begingroup$

You can minimize the cost function by using NMinimize 

  ClearAll["Global`*"]

data = {{0., 2.61}, {0.1, 2.62}, {0.2, 2.62}, {0.3, 2.62}, {0.4, 
2.63}, {0.5, 2.63}, {0.6, 2.74}, {0.7, 2.98}, {0.8, 3.66}, {0.9, 
5.04}, {1., 7.52}, {1.1, 10.74}, {1.2, 12.62}, {1.3, 10.17}, {1.4,
 5}, {1.5, 2.64}, {1.6, 11.5}, {1.65, 35.4}};

        model[x_] := a Cosh[b x^c Sin[d x^e]]
        cost = Total[(#2 - model[#1])^2 & @@@ data];
        fit = NMinimize[{cost, 1 < a < 3, 1 < b < 3, 1 < c < 3, 1 < d < 3, 1 < e < 3},
     {a, b, c, d, e}, Method -> "DifferentialEvolution"]

{0.00112965, {a -> 2.61748, b -> 1.71949, c -> 2.30924, d -> 1.50333, e -> 1.84597}}

Thread[{a, b, c, d, e} = {a, b, c, d, e} /. Last@fit];

Show[ListPlot[data], Plot[model[k], {k, 0, 1.7}]]

enter image description here

You can also use NonlinearModelFit with appropriate constrain on parameters and method choice. 

ClearAll["Global`*"]

data = {{0., 2.61}, {0.1, 2.62}, {0.2, 2.62}, {0.3, 2.62}, {0.4, 
2.63}, {0.5, 2.63}, {0.6, 2.74}, {0.7, 2.98}, {0.8, 3.66}, {0.9, 
5.04}, {1., 7.52}, {1.1, 10.74}, {1.2, 12.62}, {1.3, 10.17}, {1.4,
 5}, {1.5, 2.64}, {1.6, 11.5}, {1.65, 35.4}};

model[x_] := a Cosh[b x^c Sin[d x^e]]
nlm = NonlinearModelFit[
  data, {model[x], 1 < a < 3, 1 < b < 3, 1 < c < 3, 1 < d < 3, 
   1 < e < 3}, {a, b, c, d, e}, x, 
  Method -> {"NMinimize", {Method -> "DifferentialEvolution"}}, MaxIterations -> 1000]

nlm // Normal

2.61749 Cosh[1.71949 x^2.30926 Sin[1.50334 x^1.84596]]

nlm["BestFitParameters"]

{a -> 2.61749, b -> 1.71949, c -> 2.30926, d -> 1.50334, e -> 1.84596}

Same picture

Improve this answer
$\endgroup$
2
  • $\begingroup$ Why was my comment " ClearAll kills the data so the data should be inputted after ClearAll" deleted by an unknown to me person? This is unfair. $\endgroup$ – user64494 Jan 29 '19 at 12:29
  • $\begingroup$ Maybe somebody deleted because I modified the code and now it included the data in the code after ClearAll. $\endgroup$ – OkkesDulgerci Jan 29 '19 at 12:34
8
$\begingroup$

Below err is defined as the residual sum of squares. Using ?NumericQ on one of the arguments prevents FindMinimum from exact differentiation. In that way no problem arises even though the point with x == 0 is included.

err[a_, b_, c_, d_, e_?NumericQ] = Total[(a*Cosh[b*#^c*Sin[d*#^e]] - #2)^2 & @@@ data];

res[x_] = a*Cosh[b*x^c*Sin[d*x^e]] /. Last[FindMinimum[err[a, b, c, d, e], {
                            {a, 3, 1/100, 25}, {b, 2, 1/100, 25}, {c, 1, 1/100, 25},
                            {d, 23/10, 1/100, 25}, {e, 2, 1/100, 25}}, Method -> "InteriorPoint"]]

Show[Plot[res[x], {x, 0, 1.65}, PlotRange -> All], ListPlot[data]]

Alternative start values result in a smaller sum of squares:

FindMinimum[err[a, b, c, d, e],
  {{a, 2.6174, 1/100, 25}, {b, 1.7195, 1/100, 25}, {c, 2.3092, 1/100, 25},
   {d, 1.5033, 1/100, 25}, {e, 1.845, 1/100, 25}}, Method -> "InteriorPoint"]

{0.0011296543, {a -> 2.6174825, b -> 1.7194932, c -> 2.3092448, d -> 1.5033314, e -> 1.845972}}

These are highly likely optimal: When x == 0 the regression formula equals just a and because of the first data point a should presumably be close to 2.61.

If the last factor and Cosh are cancelled from the data values (using the right branch of ArcCosh) we get something we 2 peaks at high x-values:

data2 = Thread[{data[[All, 1]], PadLeft[{-1, -1}, 18, 1] ArcCosh[data[[All, 2]]/2.61]}]
ListPlot[data2]

We can solve for b and c such that the peaks are interpolated leaving only d and e for estimation:

sol = First[Solve[b #1^c Sin[d #1^e] == #2 & @@@ data2[[{-1, -6}]], {b, c}] /. C[1] -> 0 // Chop]

For different values of d and e let's look at the residual sum of squares excluding x == 0, because sol can't be evaluated in that instance

err0[a_, b_, c_, d_, e_?NumericQ] = Total[(a*Cosh[b*#^c*Sin[d*#^e]] - #2)^2 & @@@ Rest[data]];
search = Table[{d, e, Log[If[Im[#] == Im[#2] == 0, 
            Quiet[Min[err0[2.61, #, #2, d, e], 10^5.]], 10^5.]]} & @@
              ({b, c} /. sol), {d, 1/50, 7, 1/50}, {e, 1/50, 7, 1/50}] // Catenate;
ListPointPlot3D[search, PlotRange -> All, AxesLabel -> {x, y, z}]

the 2 smallest local minima of which coincide with the previous fits.

Improve this answer
$\endgroup$
7
  • $\begingroup$ Thank you. Likely a somewhat better fit than your $2.75699 \cosh \left(1.66039 x^{2.92883} \sin \left(4.74005 x^{0.70619}\right)\right) $ is $[a= 2.61748239217892,b= 1.71949328471199,c= 2.30924398627606,d= 1.50333104215966,e= 1.84597270613482]$. Can you kindly compare the results? $\endgroup$ – user64494 Jan 28 '19 at 17:42
  • $\begingroup$ Colleagues, let us call things by their proper names. The answer of @Coolwater is not bad, but is not optimal. $\endgroup$ – user64494 Jan 28 '19 at 18:40
  • $\begingroup$ That's interesting, but for $[a= 2.61748239217892,b= 1.71949328471199,c= 2.30924398627606,d= 1.50333104215966,e= 1.84597270613482]$ the sum of squared residuals equals $0.0000627585700895793 $. Hope you feel the difference. $\endgroup$ – user64494 Jan 28 '19 at 19:01
  • $\begingroup$ You can see it at dropbox.com/s/o08c8x0cm0dutcl/deep%20fit.pdf?dl=0 . $\endgroup$ – user64494 Jan 28 '19 at 19:50
  • $\begingroup$ From where "Alternative start values {a,}, {b, 1.7195,}, {c, 2.3092,}, {d, 1.5033}, {e, 1.845}" which lead to {a -> 2.6174825, b -> 1.7194932, c -> 2.3092448, d -> 1.5033314, e -> 1.845972}} are taken? $\endgroup$ – user64494 Jan 28 '19 at 20:17
3
$\begingroup$

Just an extended comment that for this particular data/model combination it is difficult to find the values of parameters that minimize the residual sum of squares. I think there are two main issues:

  1. The surface is extremely bumpy which can get any iterative search algorithms to get lost. (This is not necessarily a Mathematic, Maple, R, SAS, or *MathCAD" issue.)

  2. Some of the parameter estimators are highly correlated with each other. That also causes problems for finding the optimal parameter values.

Below is a contour plot of the sum of squares evaluated at @user64494 's optimal values of $a$, $b$, and $c$ in the neighborhoods of $d$ and $e$. One can see the extreme bumpiness and the highly correlated nature of the estimators of $d$ and $e$. The red dot indicates the location of the optimal values of $d$ and $e$.

Bumpy surface near optimal values of d and e

Improve this answer
$\endgroup$
1
  • $\begingroup$ +1. Thank you for your valuable info. I think the increase of the data size to 30 changes nothing here. $\endgroup$ – user64494 Jan 28 '19 at 21:00
0
$\begingroup$

NonlinearModelFit works without additional assumptions(Thanks @JimB 's comment) with Method->"NMinimize"

mod = NonlinearModelFit[Rest[data],a*Cosh[b*x^c*Sin[d*x^e]], { a , b , c, d,  e }, x,Method -> "NMinimize" ]
Show[{ListPlot[data], Plot[Normal[mod], {x, 0, 1.65}]}]

enter image description here

Improve this answer
$\endgroup$
3
  • $\begingroup$ 1. You changed the data by Rest. 2. Your fit $ 2.76855 \cosh \left(1.65609 x^{2.93246} \sin \left(4.74049 x^{0.705923}\right)\right)$ is far away from optimal ones (see dropbox.com/s/o08c8x0cm0dutcl/deep%20fit.pdf?dl=0 ). $\endgroup$ – user64494 Jan 28 '19 at 20:10
  • $\begingroup$ Yes I omitted the first data point(as @JimB recommended). My plot seems to fit the data quite well , not more. $\endgroup$ – Ulrich Neumann Jan 28 '19 at 20:16
  • $\begingroup$ One of the aims of my question is to pay attention of Mathematica developers to the fact that NonlinearModelFit faces problems with some model functions. $\endgroup$ – user64494 Jan 28 '19 at 20:21

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.