Hello.
I do some minimisations using ROOT. Everything was fine until I tried to do the same thing with IDL using their CURVEFIT() procedure. The results are the same, but IDL finishes in about the same time as the C++ code. That was strange at first because C++ should easily beat any IDL procedure. The problem is that IDL needs only 3 iterations (they don’t really explain what they mean by iteration) to find the minimum. All Minuit2 algorithms on the other hand need about 30 NCalls() to get the same result. For GSL algorithms NCalls() returns 0, but they seem to be slower than Minuit2.
Could this be explained by wrong setting of the minimizer?
Hello.
I wrote a dirty version that I hope demonstrates the problem very well. The archive can be downloaded from http://pmal.cust.aspone.cz/minuit2_vs_idl.zip (there are two directories in the archive. Test is the C++ version and idl_test contains IDL routines that do the same thing. data.txt in Test/Debug is the data file). I am using Minuit2 because I have been working with threads and Minuit2 is thread safe. The program was built in Eclipse CDT with the following flags (not all of them are necessary for this little demonstration).
Hi,
I don’t know how to run the IDL version, but what I get in Minuit2 is 4 iterations and 45 function calls. Here is the log obtained with the verbose option:
How many function calls do you get in IDL ? What really counts is the the total number function calls, because is there where the majority of the time is normally spent.
In the comparison you should also be careful, on what is used as tolerance for stopping the iterations. In the given example, you see from the log, that after 2 iterations the edm is ~ 1E-6. However, we require an effective tolerance of 10^-8 (there is a factor 10^-4) to take into account which is applied inside Minuit2.
Also in the last iteration, which requires 12 function calls, the second derivative matrix Hesse is computed. I am not sure that this is done in your IDL minimizer.
Hello and thanks.
One of the reasons why IDL seems to be so fast is this:
The model function is a*f(x) + b
In my C++ code f(x) was evaluated every time DoEval() called chiSquared(). IDL evaluates f(x) once and stores the results in an array. This way it does not have to evaluate f(x) for every point and repetition. I updated the code to do the same thing and got from 0.6s to 0.04s per minimize(). I tried to play with the step size and got to 0.03s with the C++ code. Strategy/Edmval/SetValidError don’t seem to make much difference… Unfortunately, IDL is still at 0.02s. I also tried to get the number of calls of the minimized function made by IDL using a simple print statement and it looks like IDL needs ~12 calls.
The difference in time can then be due to the difference in function calls. Looking at the documentation of the IDL routine (astro.virginia.edu/class/oco … /idl50.htm) I noticed that
the derivatives are computed using forward difference. Minuit uses a more complex method which requires more
function calls (at least 2 for each derivative).
Probably the comparison with a linear function, like in your case, does not make much sense. It is also known that Minuit is not the most efficient method for fitting. Methods like the Levenberg–Marquardt (implemented in ROOT as “GSLMultiFit”) or similar ones (“Fumili” or “Fumili2” in ROOT) are more efficient, but often less robust.
Fumili2 is not making any difference as far as I can tell. Not sure how to implement the function for GSLMultiFit. GSL’s BFGS2 does not want to converge… The updated C++ version is at: http://pmal.cust.aspone.cz/minuit2_vs_idl_2.zip in case you are interested.
If your problem is a simple linear regression, I suggest to visit the class TLinearFitter at root.cern.ch/root/html/TLinearFitter.html
and example in $ROOTSYS/tutorials/fit/fitLinear.C, etc
Thanks.
I tried the TLinearFitter, but NR3 simple linear regression (Fitab) is slightly faster. It is probably also more appropriate for the given problem.
The parameter values and chi2 I am getting from NR3 and TLinearFitter are very much the same, but errors returned by TLinearFitter seem to be ~3 times smaller than that from the NR3 fitter.
Both methods are called without weights.
The way it is done in the NR3 method is that they compute the parameters (a, b) and errors (siga, sigb) using sigma_i = sigma = 1.0 and then recompute sigma^2 = (X^2)/DOF to get a “perfect fit”. Finally, rescaling the parameter errors: siga = sigasigma, sigb = sigbsigma.
Is this expected behaviour?
Thanks.
PS
I also tried with NR3 Fitlin (linear fit using normal equations) that should be the same method as TLinearFitter (with sigma_i = 1.0) and the errors are ~1000 times smaller that the errors returned by TLinearFitter.