Home | News | Documentation | Download

Changing constraint on nuisance parameters in HistFactory


#1

Dear HistFactory users,

I am performing a fit using a model constructed to HistFactory. In the documentation it is written that the standard constraint term for a nuisance parameter (in particular I am adding a HistoSys variation) is a unit gaussian Pdf. I wanted to change this constraint term with a Uniform Constraint. So as a first attempt I have taken the model generated with the Measurement::MakeModelAndMeasurementFast() function, and substituted by hand the RooGaussian constraint with a RooUniform Pdf, using a RooCustomizer.

I subsequently realized that there is a method of the HistFactory::Measurement class meant to change the constraint term of a systematics called Measurement::AddUniformSyst(). Using this metohod, though, does not give me back the same result on the nuisance parameters during the fit.

Is this something understood? Am I missing some part of the substitutions when using the RooCustomizer?
Thanks a lot in advance,
Simone


#2

Hello @simeloni,

try
pdf->Print("t")
first to print a tree of the PDF. Or use the graphViz interface to see what’s different between the two PDFs:
https://root.cern.ch/doc/master/classRooAbsArg.html#a77b2077d0a0fa6a260f272ccf4a19a8f

When we have more information, we can investigate further.


#3

Hello @Hageboeck,
thanks for offering some help.

I have used your suggestion and had a look at the model tree. Given that the model is very complicated, I have tried to reduce a little bit the tree by including only the sample on to which the systematic has been added.

Below I am attaching quasi-code snippet that explains what I am doing

#define CUSTOMIZER_STRATEGY kTrue

//Define the templates file
TFile* fin = new TFile("Templates_File.root");

//----Define the measurement 
RooStats::HistFactory::Measurement meas("R","R");
meas.SetLumi(1.0);

//---Define the channel
RooStats::HistFactory::Channel chan("chan");
chan.SetStatErrorConfig(1.e-5,"Poisson");

//----Define the channel
RooStats::HistFactory::Sample* Bd2DD = new RooStats::HistFactory::Sample("hNBd2DD", "hNBd2DD", "Templates_file.root")
Bd2DD->SetNormalizeByTheory(kFalse);

Bd2DD.AddNormFactor("NDD");
Bd2DD.AddNormFactor("MCBd2DD");
Bd2DD.AddNormFactor("fBd2DD");

Bd2DD.AddHistoSys("variation",      "hNBd2DD_m1sigma",      "Templates_file.root", "hNBd2DD_p1sigma",      "Templates_file.root");
Bd2DD.AddHistoSys("quad_variation", "hNBd2DD_mqsigma_quad", "Templates_file.root", "hNBd2DD_p1sigma_quad", "Templates_file.root");

chan.AddSample(*Bd2DD);
chan.SetData("hNData", "Templates_file.root");

meas.AddChannel(chan);
int num_channels = 1;

//---Set some normalization parameters constant
meas.AddConstantParam("MCBd2DD");
meas.AddConstantParam("fBd2DD");

if(!CUSTOMIZER_STRATEGY){
	//---Changing the Gaussian constraint to a uniform constraint
	meas.AddUniformSyst("variation");
	meas.AddUniformSyst("quad_variation");
}

//---Define the workspace of the measurement, and the ModelConfig attached
RooWorkSpace* w;
w = RooStats::HistFactory::MakeModelAndMeasurementFast(meas);
RooStats::ModelConfig *mc  (RooStats::ModelConfig*) w->obj("ModelConfig");

//----Get the pdf and define a new pdf that will contain the modified model if CUSTOMIZER_STRATEGY is true
RooArgSet* obs = (RooArgSet*) mc->GetObservables();
Roocategory* idx = (RooCategory*) obs->find("channelCat");

RooSimultaneous* model  = (RooSimultaneous*)mc->GetPdf();
RooSimultanoeus* model_cust = new RooSimultanoeus("simPdf_modified", "simPdf_modified", *idx); 

if(CUSTOMIZER_STRATEGY){
	//----Take the old gaussian constraint
	RooRealVar* alpha_variation      = (RooRealVar*) mc->GetNuisanceParameters()->find("alpha_variation");
	RooRealVar* alpha_quad_variation = (RooRealVar*) mc->GetNuisanceParameters()->fing("alpha_quad_variation");
	RooGaussian* alpha_variationConstraint        = (RooGaussian*) w->pdf("alpha_variationConstraint");
	RooGaussian* alpha__quad_variationConstraint  = (RooGaussian*) w->pdf("alpha_quad_variationConstraint");

	//---Define the new constraint
	RooUniform* alpha_variationConstraint_Uniform      = new RooUniform("alpha_variationConstraint_Uniform",      "alpha_variationConstraint_Uniform",      RooArgSet(*alpha_variation));
	RooUniform* alpha_quad_variationConstraint_Unifofm = new RooUniform("alpha_quad_variationConstraint_Uniform", "alpha_quad_variationConstraint_Uniform", RooArgSet(*alpha_quad_variation)); 

	//---Customize the model
	RooCustomizer* cust_;
	RooAbsPdf* pdf_;
	for (unsigned int j=0; j<num_channels; j++){
		idx->setIndex(j);
		pdf_ = model->getPdf(idx->getLabel());
		cust_ = new RooCustomizer(*pdf_, "cust_");
		cust_->replaceArg(*alpha_variatonConstraint,       *alpha_variationConstraint_Uniform);
		cust_->replaceArg(*alpha_quad_variationConstraint, *alpha_quad_variationConstraint_Uniform);
		model_cust->addPadf((RooAbsPdf&)*cust_->build(),idx->getLabel());
	}
}

//---Print the model in a .dot file to be processed with graphViz
if(!CUSTOMIZER_STRATEGY){
	model->graphVizTree("TreeUniform.dot");
}
else{
	model_cust->graphVizTree("TreeUniform_RooCustomizer.dot");
}

The only channel that has been added has two systematics attached to it, called variaton and quad_variation. I then change the gaussian constraints with a uniform constraint, choosing between the two methods with the variable CUSTOMIZER_STRATEGY.

I have also attached the two graph trees i obtain with the two methods. I am sorry if the names do not fully coincide with the pseudo-code from above, but I have tried to change the names from my original code to keep the tree image compact.
TreeUniform_RooCustomizer.pdf (55.9 KB)
TreeUniform.pdf (57.2 KB)

The main difference I can see by eye is that when using the function AddUniformSyst() from the HistFactory::Measurement class, the PieceWiseInterpolation class that performs the template morphing does not take the alpha values directly from the alpha RooRealVar’s, but from a RooPolyVar with constant parameters +1 and -1 which x values are the alpha variables . Can you please explain why is this happening and if this can be the source of the different results?

Thanks a lot,
Simone


#4

Dear all,

We would really need any feedback from the HistFactory experts about this.
Any help is very much appreciated!

Thanks a lot and sorry for upping the topic.
Simone


#5

Hello Simone,

it looks like there are no HistFactory experts around.

Looking into the code of HistoToWorkspaceFactoryFast::AddConstraintTerms, I found that when adding a uniform constraint, it actually adds a Gaussian with sigma=100. That should be quite close to uniform, but it’s obviously not. I assume that the level of disagreement you find is more than that, correct?

When making a model first and replacing the constraints afterwards, it’s indeed using a RooPolyVar that transforms the value beta [0, 10] starting at 1 with a uniform constraint into alpha [-1, 1]. If I see this correctly, that would allow the nuisance parameter to go from -1 to 9 with equal probability. Could this explain the results you see?
If the replacement function is used, histfactory should print something like
edit for <nuisance parameter> with rel uncert = <uncert>.
Try to search for this in the output.