# Evaluation of PLR test statistic for discovery case

Dear ROOTers,

As far as I understand, the ROOT macro StandardHypoTestInvDemo.C is specially designed to calculate a confidence interval for the particular case in which the NULL hypothesis is S+B (i.e. you are testing at which level of signal your model is incompatible with the data you observed).

I am trying to use the same code to test the B-only hypothesis, so that if if the null hypothesis is ruled out at a given confidence level, one could claim a discovery. In some tutorials online I’ve seen that in principle that should be as easy as swaping the order of the hypotheses in the FrequentistCalculator:
FrequentistCalculator(*data, *sbModel, *bModel); //ALT=S+B, NULL=B-only

After having looked carefully at what the ProfileLikelihoodTestStat class is doing, I don’t think that only doing the thing above would be entirely correct. For discovery studies we want to evaluate the so-called q_0 test statistic, see for instance (12) in G.Cowan et.al. Hence, no matter if the pseudo-data to construct the test statistic distributions has been generated using \mu=0 (B-only) or \mu=\mu_test (S+B), the evaluation of the test statistic should always be done considering \mu=0. Correct me if I’m wrong, but I think that right now this method always evaluates the test statistic assuming \mu=\mu_test (for the conditional likelihood).

I think that this functionality could be simply added to the ProfileLikelihoodTestStat class by adding something like the following lines of code close to the beginning of ProfileLikelihoodTestStat::EvaluateProfileLikelihood:

if(fLimitType==oneSidedDiscovery){
firstPOI->setVal(0.);
initial_mu_value = 0.;
}


In this way, the value for the parameter of interest used in the evaluation of the conditional likelihood will always be 0.

Am I missing anything obvious in here? This is the only way I have found to evaluate the test statistic for discovery q_0, but it requires modifying one of the original RooStats classes.

Best,
Ibles

Hi,

For discovery significance you should use the similar macro, StandardHypoTestDemo.C, see https://root.cern.ch/doc/master/StandardHypoTestDemo_8C.html

The tested value of \mu to compute the profiled likelihood ratio is provided by the used by defining a snapshot in the ModelConfig with the desired value (e.g. \mu = 0).

Best Regards

Lorenzo

I have already provided a ModelConfig with a snapshop in which the POI is 0, and I think that it has been correctly picked up by the code. My concern is that I am not entirely sure that the ProfileLikelihoodTestStat class is actually using \mu=0 for the evaluation of the likelihood (as opposed to pseudo-data generation).

From lines 104 and 164 of ProfileLikelihoodTestStat.cxx I get the impression that the EvaluateProfileLikelihood function is taking the initial value of the POI to evaluate the conditional likelihood. Printing the variable initial_mu_value on screen tells me that in the evaluation of the test statistic under both the NULL (pseudo-data(\mu=0)) and ALT (pseudo-data(\mu=\mu_test)) models, the tested \mu is always equal to “\mu_test”, rather than 0.

I have also tried to use the B-only model in the definition of the test statictic object, like in the example you sent me:
ProfileLikelihoodTestStatMod profll(*bModel->GetPdf())
However, this did not change the value of initial_mu_value.

Any thoughts of what might be happening? Also, could you confirm that for both the NULL and ALT distributions you are really testing \mu=0, and not \mu=\mu_test, for each step in the inverted hypothesis testing?

Thanks and kind regards,
Ibles

Hi,

It is correct that the ProfileLikelihoodTestStat https://root.cern.ch/doc/v606/ProfileLikelihoodTestStat_8cxx_source.html class will use the provided test value of mu to evaluate the test statistic. But when computing a discovery significance, this test value is the one of the B model (i.e. mu = 0).

In case of a discovery significance there is no scan to different value of mu, only one test (Null hypothesis mu = 0 vs Alt hypothesis mu = standard value).

I am sure the test statistic , shown for example in the plot in

https://root.cern.ch/doc/master/StandardHypoTestDemo_8C.html

is evaluated at mu = 0 for the B toys otherwise you cannot get a half central chi2 dsitribution for the red histogram. If it is evaluated at a different test value than the one used to generate toys you will get a non-central chi-square, as it is the case for the blue histogram (S+B).

Lorenzo

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