0
$\begingroup$

I have not tried programming this myself yet, so I do not have any example code. Mostly, I am looking for advice on whether I am going to be on the right track if I try something.

Imagine a sophisticated Dynamic interface where equations, parameter values, and so on can be entered. Data is plotted and a model is overlayed, along with residuals. I then specify which parameters I want to allow to vary, and start an optimization by pressing a button (either via NonlinearModelFit, or by manually calculating the sum of the squares of the residuals and running NMinimize on that -- yes, there are some reasons why that latter approach may be preferable). I then use StepMonitor to update the plots and the parameter values as they change.

Occasionally, however, one will see the fit wandering away from what is reasonable. Or you realize that you should have fixed a value instead of letting it vary. So I would like a button to abort the optimization process, where the values from the last-calculated step would be returned, but I do not want the rest of the interface aborted. I want to smoothly catch the results of the optimization abort in order to update values and then let the user keep going with the interactive editing and calculations.

The first thing that came to mind was to use something like a Break, Interrupt, Abort, or Throw statement in an EvaluationMonitor, along with a CheckAbort, or Catch or something. But I have no idea if I am on the right track at all. Any suggestions?

$\endgroup$
1
$\begingroup$

It looks like Throw/Catch may be the easiest way to do this. Here is some sample code. First, setting up a model and some sample data to fit.

sample[t_] = (0.002 + 101 t - 461000 t^2 + 2.218 10^9 t^3 - 3.64 10^12 t^4 + 
              3.17 10^15 t^5) Exp[-8653 t];
data = SetPrecision[Table[{t, sample[t] + RandomVariate[NormalDistribution[0, 0.00001]]}, 
             {t, 0, 0.002, 0.000004}],30];
rateeqs = {a'[t] == k1b b[t] + ksqb b[t] a[t] + kttb b[t]^2 + kbd b[t] c[t] - kdb a[t] d[t], 
           b'[t] == -k1b b[t] - ksqb b[t] a[t] - kttb b[t]^2 - kbd b[t] c[t] + kdb a[t] d[t], 
           c'[t] == k1d d[t] + ksqd d[t] c[t] + kttd d[t]^2 + kdb a[t] d[t] - kbd b[t] c[t], 
           d'[t] == -k1d d[t] - ksqd d[t] c[t] - kttd d[t]^2 - kdb a[t] d[t] + kbd b[t] c[t]};
initconc = {a[0] == a0, b[0] == b0, c[0] == c0, d[0] == d0};
additionaltdeps = {abs60[t] == 5 eps60 b[t], abs70[t] == 5 eps70 d[t], 
                   abs[t] == abs60[t] + abs70[t]};
additionalinitcond = {abs60[0] == 5 eps60 b[0], abs70[0] == 5 eps70 d[0], 
                      abs[0] == abs60[0] + abs70[0]};
tdepvars = {a, b, c, d, abs60, abs70, abs};
fixedparams = {k1b -> 6000, k1d -> 100, ksqb -> 10^6, ksqd -> 10^6, 
               kttb -> 10^9, kttd -> 10^9, a0 -> 4 10^-5, c0 -> 2 10^-5, 
               eps60 -> 3500, eps70 -> 12000};
varparams = {kbd, kdb, b0, d0};
initguesses = {kbd -> 5 10^8, kdb -> 10^8, b0 -> 10^-7, d0 -> 10^-8};
solution = ParametricNDSolve[(Join[rateeqs, initconc, additionaltdeps, 
           additionalinitcond] /. fixedparams), tdepvars, {t, 0, 0.002}, 
           varparams, WorkingPrecision -> 30];

Now assigning the initial guesses to a dynamic monitoring variable.

tmp = varparams /. initguesses;
Column[{Show[ListPlot[data, ImageSize -> Automatic -> 400, ImagePadding -> {{50, 1}, {1, 1}}, 
                      Frame -> True], 
             Plot[((abs /. solution) @@ tmp)[t], {t, 0, 0.002}, PlotStyle -> Red,
                  PlotRange -> Full]], 
        ListPlot[{#1, #2 - ((abs /. solution) @@ tmp)[#1]} & @@@ data, 
                 ImageSize -> Automatic -> 400, ImagePadding -> {{50, 1}, {50, 1}},
                 Frame -> True, AspectRatio -> 0.2]}]

enter image description here

And now a DynamicModule to display the fit as it happens as well as a button to start the fit and another to abort it.

DynamicModule[{fitInProgress, abort, nlm}, 
  fitInProgress = False; 
  abort = False; 
  nlm = Null; 
  Column[{Dynamic@Show[ListPlot[data, ImageSize -> Automatic -> 400, 
                           ImagePadding -> {{50, 1}, {1, 1}}, Frame -> True], 
                       Plot[((abs /. solution) @@ tmp)[t], {t, 0, 0.002}, 
                           PlotStyle -> Red, PlotRange -> Full]], 
          Dynamic@ListPlot[{#1, #2 - ((abs /. solution) @@ tmp)[#1]} & @@@ data, 
                           ImageSize -> Automatic -> 400, 
                           ImagePadding -> {{50, 1}, {50, 1}}, Frame -> True, 
                           AspectRatio -> 0.2], 
  Button["Perform Fit", 
          fitInProgress = True; 
          tmp = Catch[nlm = NonlinearModelFit[data, ((abs /. solution) @@ varparams)[t], 
                                  Evaluate[{#, # /. initguesses} & /@ varparams], t, 
                                  Method -> "LevenbergMarquardt", StepMonitor :> (tmp = varparams), 
                                  EvaluationMonitor :> If[abort, Throw[varparams]], 
                                  Gradient -> "FiniteDifference", WorkingPrecision -> 30]; 
                      Throw[varparams /. nlm["BestFitParameters"]]]; 
          abort = False; 
          fitInProgress = False, ImageSize -> Automatic, Method -> "Queued"], 
 Button["Abort Fit", abort = True, ImageSize -> Automatic, Method -> "Preemptive"], 
 Row[{"nlm = ", Dynamic@nlm}], 
 Row[{"tmp = ", Dynamic@SetPrecision[tmp, 5]}]}]]

enter image description here

$\endgroup$

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.