Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches

Detailed Description

View in nbviewer Open in SWAN
Example of a Multi Layer Perceptron For a LEP search for invisible Higgs boson, a neural network was used to separate the signal from the background passing some selection cuts.

Here is a simplified version of this network, taking into account only WW events.

Training the Neural Network
Epoch: 0 learn=0.12676 test=0.125856
Epoch: 10 learn=0.101039 test=0.0959826
Epoch: 20 learn=0.0991537 test=0.0944351
Epoch: 30 learn=0.0954485 test=0.0897768
Epoch: 40 learn=0.0919049 test=0.0889158
Epoch: 50 learn=0.0912535 test=0.0883491
Epoch: 60 learn=0.0906147 test=0.0883948
Epoch: 70 learn=0.0901361 test=0.0873496
Epoch: 80 learn=0.0895824 test=0.0870785
Epoch: 90 learn=0.0891055 test=0.0855251
Epoch: 99 learn=0.0887071 test=0.08477
Training done.
test.py created.
Network with structure: @msumf,@ptsumf,@acolin:5:3:type
inputs with low values in the differences plot may not be needed
@msumf -> 0.0177375 +/- 0.0163272
@ptsumf -> 0.0247588 +/- 0.0383736
@acolin -> 0.0312125 +/- 0.0364969
void mlpHiggs(Int_t ntrain=100) {
const char *fname = "mlpHiggs.root";
TFile *input = 0;
if (!gSystem->AccessPathName(fname)) {
input = TFile::Open(fname);
} else if (!gSystem->AccessPathName(Form("%s/legacy/mlp/%s", TROOT::GetTutorialDir().Data(), fname))) {
input = TFile::Open(Form("%s/legacy/mlp/%s", TROOT::GetTutorialDir().Data(), fname));
} else {
printf("accessing %s file from http://root.cern/files\n",fname);
input = TFile::Open(Form("http://root.cern/files/%s",fname));
}
if (!input) return;
TTree *sig_filtered = (TTree *) input->Get("sig_filtered");
TTree *bg_filtered = (TTree *) input->Get("bg_filtered");
TTree *simu = new TTree("MonteCarlo", "Filtered Monte Carlo Events");
Float_t ptsumf, qelep, nch, msumf, minvis, acopl, acolin;
sig_filtered->SetBranchAddress("ptsumf", &ptsumf);
sig_filtered->SetBranchAddress("qelep", &qelep);
sig_filtered->SetBranchAddress("nch", &nch);
sig_filtered->SetBranchAddress("msumf", &msumf);
sig_filtered->SetBranchAddress("minvis", &minvis);
sig_filtered->SetBranchAddress("acopl", &acopl);
sig_filtered->SetBranchAddress("acolin", &acolin);
bg_filtered->SetBranchAddress("ptsumf", &ptsumf);
bg_filtered->SetBranchAddress("qelep", &qelep);
bg_filtered->SetBranchAddress("nch", &nch);
bg_filtered->SetBranchAddress("msumf", &msumf);
bg_filtered->SetBranchAddress("minvis", &minvis);
bg_filtered->SetBranchAddress("acopl", &acopl);
bg_filtered->SetBranchAddress("acolin", &acolin);
simu->Branch("ptsumf", &ptsumf, "ptsumf/F");
simu->Branch("qelep", &qelep, "qelep/F");
simu->Branch("nch", &nch, "nch/F");
simu->Branch("msumf", &msumf, "msumf/F");
simu->Branch("minvis", &minvis, "minvis/F");
simu->Branch("acopl", &acopl, "acopl/F");
simu->Branch("acolin", &acolin, "acolin/F");
simu->Branch("type", &type, "type/I");
type = 1;
Int_t i;
for (i = 0; i < sig_filtered->GetEntries(); i++) {
sig_filtered->GetEntry(i);
simu->Fill();
}
type = 0;
for (i = 0; i < bg_filtered->GetEntries(); i++) {
bg_filtered->GetEntry(i);
simu->Fill();
}
// Build and train the NN ptsumf is used as a weight since we are primarily
// interested by high pt events.
// The datasets used here are the same as the default ones.
new TMultiLayerPerceptron("@msumf,@ptsumf,@acolin:5:3:type",
"ptsumf",simu,"Entry$%2","(Entry$+1)%2");
mlp->Train(ntrain, "text,graph,update=10");
mlp->Export("test","python");
// Use TMLPAnalyzer to see what it looks for
TCanvas* mlpa_canvas = new TCanvas("mlpa_canvas","Network analysis");
mlpa_canvas->Divide(2,2);
TMLPAnalyzer ana(mlp);
// Initialisation
ana.GatherInformations();
// output to the console
ana.CheckNetwork();
mlpa_canvas->cd(1);
// shows how each variable influences the network
ana.DrawDInputs();
mlpa_canvas->cd(2);
// shows the network structure
mlp->Draw();
mlpa_canvas->cd(3);
// draws the resulting network
ana.DrawNetwork(0,"type==1","type==0");
mlpa_canvas->cd(4);
// Use the NN to plot the results for each sample
// This will give approx. the same result as DrawNetwork.
// All entries are used, while DrawNetwork focuses on
// the test sample. Also the xaxis range is manually set.
TH1F *bg = new TH1F("bgh", "NN output", 50, -.5, 1.5);
TH1F *sig = new TH1F("sigh", "NN output", 50, -.5, 1.5);
bg->SetDirectory(0);
sig->SetDirectory(0);
Double_t params[3];
for (i = 0; i < bg_filtered->GetEntries(); i++) {
bg_filtered->GetEntry(i);
params[0] = msumf;
params[1] = ptsumf;
params[2] = acolin;
bg->Fill(mlp->Evaluate(0, params));
}
for (i = 0; i < sig_filtered->GetEntries(); i++) {
sig_filtered->GetEntry(i);
params[0] = msumf;
params[1] = ptsumf;
params[2] = acolin;
sig->Fill(mlp->Evaluate(0,params));
}
bg->SetFillStyle(3008); bg->SetFillColor(kBlue);
sig->SetFillStyle(3003); sig->SetFillColor(kRed);
bg->SetStats(0);
sig->SetStats(0);
bg->Draw();
sig->Draw("same");
TLegend *legend = new TLegend(.75, .80, .95, .95);
legend->AddEntry(bg, "Background (WW)");
legend->AddEntry(sig, "Signal (Higgs)");
legend->Draw();
mlpa_canvas->cd(0);
delete input;
}
int Int_t
Definition RtypesCore.h:45
float Float_t
Definition RtypesCore.h:57
double Double_t
Definition RtypesCore.h:59
@ kRed
Definition Rtypes.h:66
@ kBlue
Definition Rtypes.h:66
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void input
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2489
R__EXTERN TSystem * gSystem
Definition TSystem.h:555
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition TAttFill.h:37
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition TAttFill.h:39
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition TAttLine.h:40
The Canvas class.
Definition TCanvas.h:23
TVirtualPad * cd(Int_t subpadnumber=0) override
Set current canvas & pad.
Definition TCanvas.cxx:716
A ROOT file is composed of a header, followed by consecutive data records (TKey instances) with a wel...
Definition TFile.h:53
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=ROOT::RCompressionSetting::EDefaults::kUseCompiledDefault, Int_t netopt=0)
Create / open a file.
Definition TFile.cxx:4067
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:621
virtual void SetDirectory(TDirectory *dir)
By default, when a histogram is created, it is added to the list of histogram objects in the current ...
Definition TH1.cxx:8901
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3340
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3062
virtual void SetStats(Bool_t stats=kTRUE)
Set statistics option on/off.
Definition TH1.cxx:8954
This class displays a legend box (TPaveText) containing several legend entries.
Definition TLegend.h:23
TLegendEntry * AddEntry(const TObject *obj, const char *label="", Option_t *option="lpf")
Add a new entry to this legend.
Definition TLegend.cxx:317
void Draw(Option_t *option="") override
Draw this legend with its current attributes.
Definition TLegend.cxx:422
This utility class contains a set of tests useful when developing a neural network.
This class describes a neural network.
Double_t Evaluate(Int_t index, Double_t *params) const
Returns the Neural Net for a given set of input parameters #parameters must equal #input neurons.
void Export(Option_t *filename="NNfunction", Option_t *language="C++") const
Exports the NN as a function for any non-ROOT-dependant code Supported languages are: only C++ ,...
void Train(Int_t nEpoch, Option_t *option="text", Double_t minE=0)
Train the network.
void Draw(Option_t *option="") override
Draws the network structure.
void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0) override
Automatic pad generation by division.
Definition TPad.cxx:1153
static const TString & GetTutorialDir()
Get the tutorials directory in the installation. Static utility function.
Definition TROOT.cxx:3068
virtual Bool_t AccessPathName(const char *path, EAccessMode mode=kFileExists)
Returns FALSE if one can access a file using the specified access mode.
Definition TSystem.cxx:1296
A TTree represents a columnar dataset.
Definition TTree.h:79
virtual Int_t Fill()
Fill all branches.
Definition TTree.cxx:4603
virtual Int_t GetEntry(Long64_t entry, Int_t getall=0)
Read all branches of entry and return total number of bytes read.
Definition TTree.cxx:5638
virtual Int_t SetBranchAddress(const char *bname, void *add, TBranch **ptr=nullptr)
Change branch address, dealing with clone trees properly.
Definition TTree.cxx:8380
virtual Long64_t GetEntries() const
Definition TTree.h:463
TBranch * Branch(const char *name, T *obj, Int_t bufsize=32000, Int_t splitlevel=99)
Add a new branch, and infer the data type from the type of obj being passed.
Definition TTree.h:353
Author
Christophe Delaere

Definition in file mlpHiggs.C.