Ratio between TKDE

Dear ROOT experts,
I’d like to use the TKDE class in one of my analysis to model a ratio between two dataset. I would like to know if there is some way to make a ratio between two KD functions.
For example I’d like to do something like this:

TKDE * kdeS = new TKDE(n, &dataS[0], xmin,xmax, "", rho);
TKDE * kdeB = new TKDE(n, &dataB[0], xmin,xmax, "", rho);
TKDE * kdeW;
kdeW = kdeS/(kdeS+kdeB);
TF1 * const hk = kdeW->GetFunction(1000); \\ <----- :)

Is it possible?


ROOT Version: 6.12/07
Platform: linuxx8664gcc


Hi,

It is not possible directly, but it can be done easily using a lambda function

auto hk = new TF1("knew","kkdeW",[&](double *x, double *) { double yS = (*kdeS)(x[0]);  double yB = (*kdeB)(x[0]); return yS/(yS+yB); }, xmin, xmax, 0);

And ini case it hit is too slow, you can do the same using the corresponding TF1 that you get from TKDE::GetFunction

Lorenzo

Hi Lorenzo,
thank you for the suggestion. Can you please expand a bit the second part of your post?

Did you mean retrieve the functions of the two KDE with TKDE::GetFunction and then make the ratio?
I would like to use the TF1 only at the end (to model the ratio) and make the ratio between the KDE, in order to avoid as much bias as possible

Alberto

Hello Alberto

Yes I was meaning retrieving the TF1 objects from TKDE and then use those ones. This is equivalent to binning the KDE values and then use some interpolation.

Yes if you want to avoid this then the code I have put before should be fine.

Cheers

Lorenzo

I’m trying unsuccessfully to make it work. I’ll copy paste my actual code, so we can avoid naming confusion. I’m running the macro on the fly with CLING. The TKDE class are created correctly since I can plot them.

---------- As first I tried:

    // create TKDE class
    float xMin = 0.0;
    float xMax = 1.0;
    float rhoRT = 0.2;
    float rhoWT = 0.2;

    TKDE *kdeRT = new TKDE(nRT, &vKDERT[0],xMin,xMax,"", rhoRT);
    TKDE *kdeWT = new TKDE(nWT, &vKDEWT[0],xMin,xMax,"", rhoWT);

    TF1 *pdfRT = new TF1("pdfRT",kdeRT,xMin,xMax,0);
    TF1 *pdfWT = new TF1("pdfWT",kdeWT,xMin,xMax,0);

    TCanvas *c10 = new TCanvas();
    mva_RT->Scale(1./mva_RT->Integral(),"width" );
    mva_RT->Draw("");
    pdfRT->Draw("SAME"); // THIS WORKS

    TCanvas *c11 = new TCanvas();
    mva_WT->Scale(1./mva_WT->Integral(),"width" );
    mva_WT->Draw("");
    pdfWT->Draw("SAME"); // THIS WORKS

    auto pdfW = new TF1("knew", "kkdeW",
    [&](double *x, double *p) { double yRT = (*kdeRT)(x[0]); double yWT = (*kdeWT)(x[0]); return yWT/(yWT+yRT); }, 
    xMin,xMax,0);

    TCanvas *c12 = new TCanvas();
    pdfW->Draw();

---------- The code does not even starts and gives the following error:

IIn file included from /build/cmsbld/auto-builds/CMSSW_10_3_0_pre1-slc6_amd64_gcc700/build/CMSSW_10_3_0_pre1-build/BUILD/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/build/input_line_12:24:
In file included from /cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/etc/../etc/dictpch/allHeaders.h:522:
In file included from /cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/etc/../include/TAdvancedGraphicsDialog.h:35:
In file included from /cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/etc/../include/TF1.h:35:
/cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/etc/../include/Math/ParamFunctor.h:182:9: error: array initializer must be an initializer list or string literal
      : fObj(pObj), fMemFn(pMemFn)
        ^
/cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/etc/../include/Math/ParamFunctor.h:294:19: note: in instantiation of member function 'ROOT::Math::ParamMemFunHandler<ROOT::Math::ParamFunctorTempl<double>, char [6], (lambda at
      /lustre/cmswork/abragagn/BPH/BTag/OSMuon/src/PDAnalysis/Ntu/bin/ntuples/MVAtraining/fitMVA.C:176:5)>::ParamMemFunHandler' requested here
      : fImpl(new ParamMemFunHandler<ParamFunctorTempl<T>, PtrObj, MemFn>(p, memFn))
                  ^
/cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/etc/../include/TF1.h:383:137: note: in instantiation of function template specialization 'ROOT::Math::ParamFunctorTempl<double>::ParamFunctorTempl<char [6], (lambda at
      /lustre/cmswork/abragagn/BPH/BTag/OSMuon/src/PDAnalysis/Ntu/bin/ntuples/MVAtraining/fitMVA.C:176:5)>' requested here
      TF1(EFType::kTemplScalar, name, xmin, xmax, npar, ndim, addToGlobList, new TF1Parameters(npar), new TF1FunctorPointerImpl<double>(ROOT::Math::ParamFunctor(p...
                                                                                                                                        ^
/lustre/cmswork/abragagn/BPH/BTag/OSMuon/src/PDAnalysis/Ntu/bin/ntuples/MVAtraining/fitMVA.C:175:21: note: in instantiation of function template specialization 'TF1::TF1<char [6], (lambda at
      /lustre/cmswork/abragagn/BPH/BTag/OSMuon/src/PDAnalysis/Ntu/bin/ntuples/MVAtraining/fitMVA.C:176:5)>' requested here
    auto pdfW = new TF1("knew", "kkdeW",

---------- Then looking at the TF1 documentation I tried to reproduce the same syntax so I changed pdfW declaration to:

    auto pdfW = new TF1("pdfW", //I removed the second " "
        [&](double *x, double *p) { double yRT = (*kdeRT)(x[0]); double yWT = (*kdeWT)(x[0]); return yWT/(yWT+yRT); }, 
        xMin,xMax,0);

---------- Now the code runs up until the KDE block, then crashes with the following output:

#0  0x000000335c6ac7be in waitpid () from /lib64/libc.so.6
#1  0x000000335c63e5c9 in do_system () from /lib64/libc.so.6
#2  0x00007f23a744b931 in TUnixSystem::StackTrace() () from /cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/lib/libCore.so
#3  0x00007f23a41f2715 in cling::MultiplexInterpreterCallbacks::PrintStackTrace() () from /cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/lib/libCling.so
#4  0x00007f23a41f219b in cling_runtime_internal_throwIfInvalidPointer () from /cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/lib/libCling.so
#5  0x00007f23a6f6418f in ?? ()
#6  0x00007f23a72a7bb0 in ?? () from /cvmfs/cms.cern.ch/slc6_amd64_gcc700/lcg/root/6.12.07-ogkkac/lib/libCore.so
#7  0x0000000007301860 in ?? ()
#8  0x000000000000001a in ?? ()
#9  0x0000000007301400 in ?? ()
#10 0x0000000003b3f348 in ?? ()
#11 0x0000000001554910 in ?? ()
#12 0x0000000007ecac48 in ?? ()
#13 0x0000000000000000 in ?? ()
terminate called after throwing an instance of 'cling::InvalidDerefException'
  what():  Trying to dereference null pointer or trying to call routine taking non-null arguments

---------- I even tried with

    auto pdfW = new TF1("pdfW",
        "[&](double *x, double *p) { return pdfRT; }", 
        xMin,xMax,0);

---------- to no success:

input_line_40:2:41: error: no viable conversion from '(lambda at input_line_40:2:73)' to 'std::function<double (double *, double *)>'
 std::function<double(double*,double*)> lambda__id4567358677425919240 = [&](double *x, double *p) { return pdfRT; } ;
                                        ^                               ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
/build/cmsbld/auto-builds/CMSSW_10_3_0_pre1-slc6_amd64_gcc700/build/CMSSW_10_3_0_pre1-build/slc6_amd64_gcc700/external/gcc/7.0.0-omkpbe2/bin/../lib/gcc/x86_64-unknown-linux-gnu/7.3.1/../../../../include/c++/7.3.1/bits/std_function.h:421:7: note: candidate constructor not viable: no known conversion from '(lambda at input_line_40:2:73)' to 'std::nullptr_t' (aka 'nullptr_t') for 1st argument
      function(nullptr_t) noexcept
      ^
/build/cmsbld/auto-builds/CMSSW_10_3_0_pre1-slc6_amd64_gcc700/build/CMSSW_10_3_0_pre1-build/slc6_amd64_gcc700/external/gcc/7.0.0-omkpbe2/bin/../lib/gcc/x86_64-unknown-linux-gnu/7.3.1/../../../../include/c++/7.3.1/bits/std_function.h:432:7: note: candidate constructor not viable: no known conversion from '(lambda at input_line_40:2:73)' to 'const std::function<double (double *, double *)> &' for 1st
      argument
      function(const function& __x);
      ^
/build/cmsbld/auto-builds/CMSSW_10_3_0_pre1-slc6_amd64_gcc700/build/CMSSW_10_3_0_pre1-build/slc6_amd64_gcc700/external/gcc/7.0.0-omkpbe2/bin/../lib/gcc/x86_64-unknown-linux-gnu/7.3.1/../../../../include/c++/7.3.1/bits/std_function.h:441:7: note: candidate constructor not viable: no known conversion from '(lambda at input_line_40:2:73)' to 'std::function<double (double *, double *)> &&' for 1st argument
      function(function&& __x) noexcept : _Function_base()
      ^
/build/cmsbld/auto-builds/CMSSW_10_3_0_pre1-slc6_amd64_gcc700/build/CMSSW_10_3_0_pre1-build/slc6_amd64_gcc700/external/gcc/7.0.0-omkpbe2/bin/../lib/gcc/x86_64-unknown-linux-gnu/7.3.1/../../../../include/c++/7.3.1/bits/std_function.h:465:2: note: candidate template ignored: requirement '_Callable<(lambda at input_line_40:2:73), TF1 *>::value' was not satisfied [with _Functor =
      (lambda at input_line_40:2:73), $1 = void]
        function(_Functor);
        ^
input_line_42:2:3: error: use of undeclared identifier 'lambda__id4567358677425919240'
 (lambda__id4567358677425919240)
  ^
Error in <HandleInterpreterException>: Error evaluating expression (lambda__id4567358677425919240).
Execution of your code was aborted.
Error in <TFormula::TFormula>: Syntax error in building the lambda expression [&](double *x, double *p) { return pdfRT; }
Error in <TFormula::Eval>: Formula is invalid and not ready to execute 

Warning in <TCanvas::ResizePad>: Inf/NaN propagated to the pad. Check drawn objects.
Warning in <TCanvas::ResizePad>: c1_n3 height changed from 0 to 10

root [1] Warning in <TCanvas::ResizePad>: Inf/NaN propagated to the pad. Check drawn objects.
Warning in <TCanvas::ResizePad>: c1_n3 height changed from 0 to 10

Warning in <TCanvas::ResizePad>: Inf/NaN propagated to the pad. Check drawn objects.
Warning in <TCanvas::ResizePad>: c1_n3 height changed from 0 to 10

Alberto

Hello,

Your second case, the one that crashes, has the correct syntax. When you pass a lambda or a function pointer to a TF1, you should not pass a second string for the title.

The crashes occurs when drawing because the lambda has been deleted at drawing time. This is a known problem and the solution is to use clone when drawing the TF1 object, by using DrawClone().
So just do at the end:

pdfW->DrawClone()

Lorenzo

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.