0
$\begingroup$

I had asked a question regarding simultaneous plotting here.

My Program will fit real part and imaginary part of a data simultaneously and gives out the value something like this

enter image description here

I column has frequency, II column has Amplitud, III column has Phase.

Here is my modified program from the responses I got from previous question.

Qfit[data1_, fmin_, fmax_] := Module[{da1 = {}, da2 = {}, L, A, B, kint, kport, f, f0, index, mini,
res, i, r, myF, nlm, data = {}, allData}, Clear[A, B, kint, kport, f, f0, index, mini, res, L]; data = Select[data1, fmin <= #[[1]] <= fmax &];da1 = Transpose[{data[[All, 1]], 
 data[[All, 2]]*Cos[data[[All, 3]]]}];da2 = Transpose[{data[[All, 
  1]], (data[[All, 2]]*Sin[data[[All, 3]]])}]; mini = Flatten[Position[da1[[All, 2]], Min[da1[[All, 2]]]]][[1]];res = da1[[mini, 1]];allData =Join[{1, Sequence @@ #} & /@ da1, {2, Sequence @@ #} & /@ da2]; i[f_] := A*B*(4 (f - f0)^2 + (kint - kport) (kint + kport))/(4 (f - f0)^2 + (kint + kport)^2);r[f_] := A*(1 - B (4 kport (f - f0))/(4 (f - f0)^2 + (kint + kport)^2));myF[index_, f_] := KroneckerDelta[index - 1] i[f] + KroneckerDelta[index - 2] r[f]; nlm = NonlinearModelFit[allData, myF[index, 
 f], {{kint, 0.001}, {kport, 0.001}, {f0, res}, {A, 0.1}, {B, 
  0.05}}, {index, f}, MaxIterations -> Infinity];fitparams = nlm["BestFitParameters"]; Print[fitparams];Set @@@ fitparams;Print["Internal Quality factor (Loss inside cavity)  =  ", f0/kint, "    External Quality factor (Input port 1 loss)  =  ", f0/kport, "    Total loss or Quality factor   =  ", f0/(kint + kport)]; Print["The graphs are ", Show[ListPlot[{da1, da2}], 
Plot[{nlm[1, f], nlm[2, f]}, {f, fmin, fmax}], 
ImageSize -> Large],Show[ListPlot[Transpose@{Last /@ da1, Last /@ da2}, Joined -> True,
  AxesLabel -> {"Re", "Im"}], 
ParametricPlot[{nlm[1, f], nlm[2, f]}, {f, Min@(First /@ da1), 
  Max@(First /@ da1)}, PlotStyle -> Red], PlotRange -> All, 
ImageSize -> Large]]]

This usually work for all most all the values, but for some data like this,here is a data from experiment.

For this we have to compute as follows Qfit[data,10.45,10.475].

It is not converging at all. I do not why it is happening I appreciate any quick responses, and sorry for my ignorance.

|improve this question
$\endgroup$
  • $\begingroup$ Lack of convergence can happen when the model is overparameterized or poorly parameterized or there just isn't a very good fit between the data and the model. Because you have so much data, you can fit the two pieces separately. Then if they really have parameters in common, the estimates of those parameters should be near each other. But in this case they are not and for the r[f] piece the resulting standard errors are near zero for all parameters and for the i[f] piece the resulting standard errors are also all near zero with a correlation matrix of all 1's or -1's. $\endgroup$ – JimB Jan 14 '17 at 7:31
  • $\begingroup$ Thanks Jim Baldwin for your comment, here is the same data with less amount of points pastebin.com/hAvuJY37. But you see the points is , inorder to extract the values correctly, I have to plot both at the same time. But I am really stuck now, I had used this same for more 1000 fits never had any problems so far. it would be great if you could help me out with this $\endgroup$ – TM90 Jan 14 '17 at 7:39
  • $\begingroup$ Can you give a dataset where it does work correctly? $\endgroup$ – JimB Jan 14 '17 at 7:40
  • $\begingroup$ Here is a data set, pastebin.com/q4bx6d4m , and i used the following to plot and got the correct results-Qfit[data, 10.45, 10.6]. Would you like to see the result as well ?I do not how to copy the image file here $\endgroup$ – TM90 Jan 14 '17 at 7:42
  • $\begingroup$ here is the result, {kint$68514->-0.00287178,kport$68514->-0.00147122,f0$68514->10.509,A$68514->0.00374764,B$68514->158.606} Internal Quality factor (Loss inside cavity) = -3659.4 External Quality factor (Input port 1 loss) = -7143.06 Total loss or Quality factor = -2419.76 $\endgroup$ – TM90 Jan 14 '17 at 7:46
2
$\begingroup$

This is an extended comment rather than an answer.

I've downloaded the problem data and I believe the trouble comes from the common parameters not having anywhere near the same values between the two subsets of the data. (However, there is a confounding issue described below.)

I took your code and modified it so that it estimates the two regressions separately:

Qfit[data1_, fmin_, fmax_] := 
 Module[{da1 = {}, da2 = {}, L, A, B, kint, kport, f, f0, index, mini,
    res, i, r, myF, nlm, data = {}, allData}, 
  Clear[A, B, kint, kport, f, f0, index, mini, res, L]; 
  data = Select[data1, fmin <= #[[1]] <= fmax &]; 
  da1 = Transpose[{data[[All, 1]], 
     data[[All, 2]]*Cos[data[[All, 3]]]}]; 
  da2 = Transpose[{data[[All, 1]], (data[[All, 2]]*
       Sin[data[[All, 3]]])}]; 
  mini = Flatten[Position[da1[[All, 2]], Min[da1[[All, 2]]]]][[1]]; 
  res = da1[[mini, 1]];
  allData = 
   Join[{1, Sequence @@ #} & /@ da1, {2, Sequence @@ #} & /@ da2]; 
  ListPointPlot3D[allData]; 
  i[f_] := A*
    B*(4 (f - f0)^2 + (kint - kport) (kint + 
          kport))/(4 (f - f0)^2 + (kint + kport)^2); 
  r[f_] := A*(1 - 
      B (4 kport (f - f0))/(4 (f - f0)^2 + (kint + kport)^2)); 
  myF[index_, f_] := 
   KroneckerDelta[index - 1] i[f] + KroneckerDelta[index - 2] r[f]; 
  allData1 = Select[allData, #[[1]] == 1 &]; 
  allData2 = Select[allData, #[[1]] == 2 &];
  allData1 = allData1[[All, {2, 3}]];
  allData2 = allData2[[All, {2, 3}]];
  nlm1 = NonlinearModelFit[allData1, 
    myF[1, f], {{kint, 0.001}, {kport, 0.001}, {f0, res}, {A, 
      0.1}, {B, 0.05}}, f, MaxIterations -> Infinity];
  nlm2 = NonlinearModelFit[allData2, 
    myF[2, f], {{kint, 0.001}, {kport, 0.001}, {f0, res}, {A, 
      0.1}, {B, 0.05}}, f, MaxIterations -> Infinity];

  fitparams1 = nlm1["BestFitParameters"];
  fitparams2 = nlm2["BestFitParameters"]; 
  Print[{fitparams1, fitparams2}];
  Print[nlm1["CorrelationMatrix"]];
  Print[nlm2["CorrelationMatrix"]];
  Print["The graphs are ", 
   Show[ListPlot[{da1, da2}], 
    Plot[{nlm1[f], nlm2[f]}, {f, fmin, fmax}], ImageSize -> Large], 
   Show[ListPlot[Transpose@{Last /@ da1, Last /@ da2}, Joined -> True,
      AxesLabel -> {"Re", "Im"}], 
    ParametricPlot[{nlm1[f], nlm2[f]}, {f, Min@(First /@ da1), 
      Max@(First /@ da1)}, PlotStyle -> Red], PlotRange -> All, 
    ImageSize -> Large]]]
Qfit[data, 10.45, 10.6]

While there was a warning about "The step size in the search has become less than the tolerance prescribed by the PrecisionGoal option, but the gradient is larger than the tolerance specified by the AccuracyGoal option. There is a possibility that the method has stalled at a point that is not a local minimum.", the fits look fine:

Separate fits

But the (cleaned-up) solutions look very different:

sol1 = {kint -> -0.0027306463823735722`, 
  kport -> -0.0013994544277047404`,
  f0 -> 10.508996991971562`, 
  A -> 0.9737627424533896`, 
  B -> 0.6095062568309071`};
sol2 = {kint -> 0.032106337286351976`, 
  kport -> -0.02754995227500653`,
  f0 -> 10.50900015588828`, 
  A -> 0.0036249688207138334`,
  B -> 9.185735216137683`};

Now fitting i[f] we can't estimate both A and B but only their product. When the product A B is examined for both fits, the products are very different.

The confounding factor is that the estimated correlations among the parameter estimators are many times close to +1 or -1 which can cause trouble with the estimation procedure. It's possible that a reformulation of the model with maybe ksum = kint + kport and kdiff = kint - kport might make the estimation process more stable.

If the data sets that do work fine have high correlations among the parameter estimators, then I would highly recommend some sort of reparameterization of the model.

|improve this answer
$\endgroup$
2
$\begingroup$

Maybe this is the answer as to why there is instability and why the extreme correlation between some of the parameter estimators.

Look at the limits of the two functions being fit (i[f] and r[f]) as $f \rightarrow \infty$ and $f \rightarrow -\infty$:

Limit[A*B*(4 (f - f0)^2 + (kint - kport) (kint + kport))/(4 (f - f0)^2 + (kint + kport)^2),
  f -> ∞]
(* A B *)

Limit[A*B*(4 (f - f0)^2 + (kint - kport) (kint + kport))/(4 (f - f0)^2 + (kint + kport)^2), 
  f -> -∞]
(* A B *)

Limit[A*(1 - B (4 kport (f - f0))/(4 (f - f0)^2 + (kint + kport)^2)), f -> ∞]
(* A *)

Limit[A*(1 - B (4 kport (f - f0))/(4 (f - f0)^2 + (kint + kport)^2)), f -> -∞]
(* A *)

If we look at the data, the value for A in the r[f] curve would appear to be zero or very close to zero. That forces B to be very large.

Is there a fix? A re-parameterization might do it. Rather than having A and B we use a and ab:

i[f_] := ab*(4 (f - f0)^2 + (kint - kport) (kint + kport))/(4 (f - f0)^2 + (kint + kport)^2); 
r[f_] := a - ab (4 kport (f - f0))/(4 (f - f0)^2 + (kint + kport)^2);

Here's the full updated code:

Qfit[data1_, fmin_, fmax_] := 
 Module[{da1 = {}, da2 = {}, L, A, B, kint, kport, f, f0, index, mini,
    res, i, r, myF, nlm, data = {}, allData}, 
  Clear[A, B, kint, kport, f, f0, index, mini, res, L]; 
  data = Select[data1, fmin <= #[[1]] <= fmax &]; 
  da1 = Transpose[{data[[All, 1]], 
     data[[All, 2]]*Cos[data[[All, 3]]]}]; 
  da2 = Transpose[{data[[All, 1]], (data[[All, 2]]*
       Sin[data[[All, 3]]])}]; 
  mini = Flatten[Position[da1[[All, 2]], Min[da1[[All, 2]]]]][[1]]; 
  res = da1[[mini, 1]]; 
  allData = 
   Join[{1, Sequence @@ #} & /@ da1, {2, Sequence @@ #} & /@ da2]; 
  i[f_] := ab*(4 (f - f0)^2 + (kint - kport) (kint + 
          kport))/(4 (f - f0)^2 + (kint + kport)^2); 
  r[f_] := a - 
    ab (4 kport (f - f0))/(4 (f - f0)^2 + (kint + kport)^2); 
  myF[index_, f_] := 
   KroneckerDelta[index - 1] i[f] + KroneckerDelta[index - 2] r[f]; 
  nlm = NonlinearModelFit[allData, 
    myF[index, 
     f], {{kint, 0.001}, {kport, 0.001}, {f0, res}, {a, 0.1}, {ab, 
      0.05}}, {index, f}, MaxIterations -> Infinity]; 
  fitparams = nlm["BestFitParameters"]; Print[fitparams]; 
  Set @@@ fitparams; 
  Print["Internal Quality factor (Loss inside cavity)  =  ", f0/kint, 
   "    External Quality factor (Input port 1 loss)  =  ", f0/kport, 
   "    Total loss or Quality factor   =  ", f0/(kint + kport)]; 
  Print["The graphs are ", 
   Show[ListPlot[{da1, da2}], 
    Plot[{nlm[1, f], nlm[2, f]}, {f, fmin, fmax}], 
    ImageSize -> Large], 
   Show[ListPlot[Transpose@{Last /@ da1, Last /@ da2}, Joined -> True,
      AxesLabel -> {"Re", "Im"}], 
    ParametricPlot[{nlm[1, f], nlm[2, f]}, {f, Min@(First /@ da1), 
      Max@(First /@ da1)}, PlotStyle -> Red], PlotRange -> All, 
    ImageSize -> Large]]]

Now when we run the analysis we get the following:

Qfit[data, 10.45, 10.475]

Data and fit

|improve this answer
$\endgroup$
  • $\begingroup$ Thank you very much Jim Baldwin. I appreciate your help. You were right about the parameterization issues $\endgroup$ – TM90 Jan 14 '17 at 23:12

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.