Minimizing a functor - hand holding needed

VIVE L’AMOUR!
I’ve got a class: class chi_2_functor { public: int NDF; chi_2_functor(int ndf) { NDF = ndf; } double operator() (double *x, double *p) { double v = 0.0; if (p) for (int i = 0; i < NDF; i++) v += (x[i] - p[i]) * (x[i] - p[i]); else for (int i = 0; i < NDF; i++) v += (x[i] - i) * (x[i] - i); return v; } };
I’ve successfully been using it with TF1, TF2 and TF3: [code]{
double x, y, z, v;
double p[11] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10};

chi_2_functor *cf1 = new chi_2_functor(1);

TF1 *f1 = new TF1(“f1”, cf1, -1, 1, 11, “chi_2_functor”);
f1->SetParameters§;
x = f1->GetMinimumX(); v = f1->GetMinimum();
std::cout << "f1( " << x << " ) = " << v << std::endl;

TF1 *f10 = new TF1(“f10”, cf1, -1, 1, 0, “chi_2_functor”);
x = f10->GetMinimumX(); v = f10->GetMinimum();
std::cout << "f10( " << x << " ) = " << v << std::endl;

chi_2_functor *cf2 = new chi_2_functor(2);

TF2 *f2 = new TF2(“f2”, cf2, -1, 1, 0, 2, 11, “chi_2_functor”);
f2->SetParameters§;
f2->GetMinimumXY(x, y); v = f2->Eval(x, y);
std::cout << "f2( " << x << " , " << y << " ) = " << v << std::endl;

TF2 *f20 = new TF2(“f20”, cf2, -1, 1, 0, 2, 0, “chi_2_functor”);
f20->GetMinimumXY(x, y); v = f20->Eval(x, y);
std::cout << "f20( " << x << " , " << y << " ) = " << v << std::endl;

chi_2_functor *cf3 = new chi_2_functor(3);

TF3 *f3 = new TF3(“f3”, cf3, -1, 1, 0, 2, 1, 3, 11, “chi_2_functor”);
f3->SetParameters§;
f3->GetMinimumXYZ(x, y, z); v = f3->Eval(x, y, z);
std::cout << "f3( " << x << " , " << y << " , " << z << " ) = " << v << std::endl;

TF3 *f30 = new TF3(“f30”, cf3, -1, 1, 0, 2, 1, 3, 0, “chi_2_functor”);
f30->GetMinimumXYZ(x, y, z); v = f30->Eval(x, y, z);
std::cout << “f30( " << x << " , " << y << " , " << z << " ) = " << v << std::endl;
}[/code]
Well, now I need to find minima for 4 to 11 parameters … and there’s no TF4, …, TF11.
I think I’d like to use Minuit2 (in the first approximation, later maybe GSL) … but I cannot find any suitable example.
In general, I would need to “set limits” on (some) parameters and I would need to “fix” (some) others.
The real problem is - how can I re-use my “chi_2_functor” class and its ''operator()” when defining / creating the “FCN” to be minimized?
Note: I don’t really want to change the layout of the “chi_2_functor” class (a lot of my existing source code depends on it). I could, however, say “Double_t operator() (Double_t *x, Double_t *p = 0);”, if that makes things easier (i.e. for the current purposes, I could live fine without “p”).
If possible, I would prefer to stay away from RooFit - I’d like to do everything “manually”, following “conventional” procedures.
A pitiful case, am I not?
Pepe Le Pew.

Hi,

For minimizing directly a function (using Minuit2 or whatever minimizer you like) see the tutorial
tutorials/fit/NumericalMinimization.C

or also
root.cern.ch/drupal/content/nume … idim_minim

Concerning your chi_2_functor class, you need just to have the signature as

Double_t operator() (Double_t *x )

where x is the vector of your parameter you want to minimize.

Cheers
Lorenzo

Well, I have two ordinary c++ functions now:
double chi2(double a, double b, double c, …);
double chi2(const double *x);
if I try:
ROOT::Math::Functor f(&chi2, 11);
I get:

error: no matching function for call to ‘ROOT::Math::Functor::Functor(<unresolved overloaded function type>, int)’ /.../Functor.h:389: note: candidates are: ROOT::Math::Functor::Functor(const ROOT::Math::Functor&) /.../Functor.h:376: note: ROOT::Math::Functor::Functor(void*, unsigned int, const char*, const char*) /.../Functor.h:351: note: ROOT::Math::Functor::Functor()
If I do not want to rename the second “chi2” function, how can I enforce its usage in the “ROOT::Math::Functor” constructor?

Then, the page http://root.cern.ch/root/html/ROOT__Math__Functor.html says … “a free C function of type double ()(double *)”. I have found, however, that it MUST be “double ()(const double *)”, otherwise one gets:

/.../Functor.h: In member function ‘double ROOT::Math::FunctorHandler<ParentFunctor, Func>::DoEval(const double*) const [with ParentFunctor = ROOT::Math::Functor, Func = double (*)(double*)]’: /.../chi2_fit_cxx_ACLiC_dict.cxx:1497: instantiated from here /.../Functor.h:86: error: invalid conversion from ‘const double*’ to ‘double*’
Maybe it would be good if one corrected the description on the page mentioned above.
(In your first post here you say “Double_t operator() (Double_t *x)”, but I believe this should also be “Double_t operator() (const Double_t *x)” and that’s how it is described on the page mentioned above.)

Could you, please, have a look at the three fit results below and tell me what you think about the messages that appear there.
I am trying to use “Minuit2 / Combined” with “SetTolerance(0.000001);” - the source code is based on the “NumericalMinimization.C” tutorial.

FIT @ p0 = 1 (“Valid minimum - status = 1” means “QUITE OK”? -> “CovMatrixStatus()” returns 2): Minuit2Minimizer: Minimize with max iterations 1000000 edmval 1e-06 strategy 1 Info in <Minuit2>: matrix forced pos-def by adding to diagonal : padd = 0.00631114 Info in <Minuit2>: MnHesse: matrix was forced pos. def. Info in <Minuit2>: Minuit2Minimizer::Minimize : Covar was made pos def Minuit2Minimizer : Valid minimum - status = 1 FVAL = 12.0074262361659585 Edm = 3.8272731535285797e-11 Nfcn = 6169 p0 = 1 (fixed) p1 = 0.045557455 +/- 0.040940145 (limited) p2 = -0.47076808 +/- 0.016720148 (limited) p3 = 0 (fixed) p4 = 0.31532817 +/- 0.14644255 (limited) p5 = 0 (fixed) p6 = 1.2280403 +/- 0.052268778 (limited) p7 = 0.10554927 +/- 0.089167039 (limited) p8 = 0.058886162 +/- 0.23076764 (limited) p9 = -0.22283159 +/- 0.017421336 (limited) p10 = 1 (fixed)

FIT @ p0 = 4: Minuit2Minimizer: Minimize with max iterations 1000000 edmval 1e-06 strategy 1 Info in <Minuit2>: matrix forced pos-def by adding to diagonal : padd = 0.00638057 Info in <Minuit2>: MnHesse: matrix was forced pos. def. Info in <Minuit2>: VariableMetricBuilder: Tolerance is not sufficient, continue the minimization Info in <Minuit2>: edm = 1.32528e-08 Info in <Minuit2>: required : edmval = 2e-10 Info in <Minuit2>: Minuit2Minimizer::Minimize : Covar was made pos def Minuit2Minimizer : Valid minimum - status = 1 FVAL = 12.0074659910234587 Edm = 1.32527662495515725e-08 Nfcn = 1904 p0 = 4 (fixed) p1 = 0.085535234 +/- 0.058582369 (limited) p2 = -0.87955545 +/- 0.029964089 (limited) p3 = 0 (fixed) p4 = 0.58402282 +/- 0.24903038 (limited) p5 = 0 (fixed) p6 = 2.2943462 +/- 0.088552805 (limited) p7 = 0.1642453 +/- 0.058619207 (limited) p8 = 0.1422711 +/- 0.065355502 (limited) p9 = 0.43458384 +/- 0.014148478 (limited) p10 = 1 (fixed)

FIT @ p0 = 0.4 (“Valid minimum - status = 0” means “PERFECT”? -> “CovMatrixStatus()” returns 3 -> but the errors on fitted parameters are much larger than in the two examples above): Minuit2Minimizer: Minimize with max iterations 1000000 edmval 1e-06 strategy 1 Info in <Minuit2>: VariableMetricBuilder: Tolerance is not sufficient, continue the minimization Info in <Minuit2>: edm = 9.21189e-10 Info in <Minuit2>: required : edmval = 2e-10 Minuit2Minimizer : Valid minimum - status = 0 FVAL = 12.0074262384092538 Edm = 9.21188502830194783e-10 Nfcn = 6624 p0 = 0.4 (fixed) p1 = 0.023773689 +/- 0.18619995 (limited) p2 = -0.24556061 +/- 0.16167679 (limited) p3 = 0 (fixed) p4 = 0.16447292 +/- 0.25682105 (limited) p5 = 0 (fixed) p6 = 0.64057292 +/- 0.40286584 (limited) p7 = 0.0733344 +/- 0.25451507 (limited) p8 = 0.012423947 +/- 0.15619706 (limited) p9 = -0.59461343 +/- 0.13577978 (limited) p10 = 1 (fixed)

Hi,

In the first two fits the covariance matrix was forced to be positive defined, by increasing the diagonal elements. This makes the error artificially smaller. Probably the minimum is fine, but the errors are not reliable. The last fit looks perfect !

Concerning your other posts, unfortunately you cannot have the two function with same name and different signature. If you want to keep the same name you could bring the one you want to use for the Functor in a namespace or in a class.
You are right also about the signature, it must be double operator() (const double *). I will fix it in the reference documentation.

Thanks for the comments,

Lorenzo

1 Like

Hi Pepe,

For the case with 2 C++ function with the same name but different prototype, you could try:double chi2(double a, double b, double c, ...); double chi2(const double *x); ... typedef double (*fitting_function)(const double *x); .... fitting_function fitter = chi2; ROOT::Math::Functor f(fitter, 11);

Cheers,
Philippe.

VIVE L’AMOUR!
well, in almost all cases, I get “Valid minimum - status = 1” and the “CovMatrixStatus()” returns 2: Info in <Minuit2>: VariableMetricBuilder: no improvement in line search Info in <Minuit2>: VariableMetricBuilder: machine accuracy limits further improvement. Info in <Minuit2>: matrix forced pos-def by adding to diagonal : padd = 0.00634671 Info in <Minuit2>: MnHesse: matrix was forced pos. def. Info in <Minuit2>: VariableMetricBuilder: Tolerance is not sufficient, continue the minimization Info in <Minuit2>: edm = 4.7934e-11 Info in <Minuit2>: required : edmval = 1.19209e-11 Info in <Minuit2>: Minuit2Minimizer::Minimize : Covar was made pos def Minuit2Minimizer : Valid minimum - status = 1
As you write in your previous post, the “minimum” is fine but the “errors” are not reliable in this case.
In another thread, you once advised me … "But the best test for the errors is to run Minos. If it succeeds normally the errors returned are reliable."
I haven’t found any “Minos()” method in the “ROOT::Math::Minimizer” interface (there exists “Hesse()”, though).
If I simply try “minimizer->GetMinosError(8, eL, eH);”, I get: Info in <Minuit2>: MnMinos Parameter is at Lower limit. : par = 8 Minos: Parameter is at Lower limit. (bool)1
Well, I think that after I find the minimum with all variables “limited” (i.e. they were created with “minimizer->SetLimitedVariable(…)”), I could remove all limits … and … re-run the error analysis or even simply minimize it again.
The question is … what is the easiest way to do this?
How can I change the “type” of the “variable” (from “limited” to “free”) leaving it’s current value and step untouched?
Or maybe, how can I set new limits for the “variable” leaving it’s current value and step untouched (well, I could set huge “artificial” limits, which would make the “variable” “virtually-free”)?
How can I re-run the error analysis (without limits)?
I am stupid. No?
Pepe Le Pew.

Hi Pepe,

At the moment you need to call the Minimizer::Clear() method to clear the parameter states and re-set them as free ones. You can save their values and step sizes (using the errors) from the previous minimization) when you re-set them.
The opposite, setting first a parameter free and then limited is possible, just by calling Minimizer::SetLimitedVariable the second time.
I need to fix this to make possible to set a parameter bounded and then free.

Best Regards

Lorenzo

[quote]… values and step sizes (using the errors) …[/quote] Well, I assume I get the values using the “X()” method and the errors using the “Errors()” method, but how do I recalculate errors into steps sizes?

Hi,

Just use as step sizes the error values. The step sizes is an indication on the scale of the parameter and it is required for numerical calculation of the derivatives. Minuit internally computes the optimal step size for the particular derivative computation

Lorenzo