Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TMultiLayerPerceptron.cxx
Go to the documentation of this file.
1// @(#)root/mlp:$Id$
2// Author: Christophe.Delaere@cern.ch 20/07/03
3
4/*************************************************************************
5 * Copyright (C) 1995-2003, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12/** \class TMultiLayerPerceptron
13
14
15This class describes a neural network.
16There are facilities to train the network and use the output.
17
18The input layer is made of inactive neurons (returning the
19optionally normalized input) and output neurons are linear.
20The type of hidden neurons is free, the default being sigmoids.
21(One should still try to pass normalized inputs, e.g. between [0.,1])
22
23The basic input is a TTree and two (training and test) TEventLists.
24Input and output neurons are assigned a value computed for each event
25with the same possibilities as for TTree::Draw().
26Events may be weighted individually or via TTree::SetWeight().
276 learning methods are available: kStochastic, kBatch,
28kSteepestDescent, kRibierePolak, kFletcherReeves and kBFGS.
29
30This implementation, written by C. Delaere, is *inspired* from
31the mlpfit package from J.Schwindling et al. with some extensions:
32
33 - the algorithms are globally the same
34 - in TMultilayerPerceptron, there is no limitation on the number of
35 layers/neurons, while MLPFIT was limited to 2 hidden layers
36 - TMultilayerPerceptron allows you to save the network in a root file, and
37 provides more export functionalities
38 - TMultilayerPerceptron gives more flexibility regarding the normalization of
39 inputs/outputs
40 - TMultilayerPerceptron provides, thanks to Andrea Bocci, the possibility to
41 use cross-entropy errors, which allows to train a network for pattern
42 classification based on Bayesian posterior probability.
43
44### Introduction
45
46Neural Networks are more and more used in various fields for data
47analysis and classification, both for research and commercial
48institutions. Some randomly chosen examples are:
49
50 - image analysis
51 - financial movements predictions and analysis
52 - sales forecast and product shipping optimisation
53 - in particles physics: mainly for classification tasks (signal
54 over background discrimination)
55
56More than 50% of neural networks are multilayer perceptrons. This
57implementation of multilayer perceptrons is inspired from the
58<A HREF="http://schwind.home.cern.ch/schwind/MLPfit.html">MLPfit
59package</A> originally written by Jerome Schwindling. MLPfit remains
60one of the fastest tool for neural networks studies, and this ROOT
61add-on will not try to compete on that. A clear and flexible Object
62Oriented implementation has been chosen over a faster but more
63difficult to maintain code. Nevertheless, the time penalty does not
64exceed a factor 2.
65
66### The MLP
67
68The multilayer perceptron is a simple feed-forward network with
69the following structure:
70
71\image html mlp.png
72
73It is made of neurons characterized by a bias and weighted links
74between them (let's call those links synapses). The input neurons
75receive the inputs, normalize them and forward them to the first
76hidden layer.
77
78Each neuron in any subsequent layer first computes a linear
79combination of the outputs of the previous layer. The output of the
80neuron is then function of that combination with <I>f</I> being
81linear for output neurons or a sigmoid for hidden layers. This is
82useful because of two theorems:
83
84 1. A linear combination of sigmoids can approximate any
85 continuous function.
86 2. Trained with output = 1 for the signal and 0 for the
87 background, the approximated function of inputs X is the probability
88 of signal, knowing X.
89
90### Learning methods
91
92The aim of all learning methods is to minimize the total error on
93a set of weighted examples. The error is defined as the sum in
94quadrature, divided by two, of the error on each individual output
95neuron.
96In all methods implemented, one needs to compute
97the first derivative of that error with respect to the weights.
98Exploiting the well-known properties of the derivative, especially the
99derivative of compound functions, one can write:
100
101 - for a neuron: product of the local derivative with the
102 weighted sum on the outputs of the derivatives.
103 - for a synapse: product of the input with the local derivative
104 of the output neuron.
105
106This computation is called back-propagation of the errors. A
107loop over all examples is called an epoch.
108Six learning methods are implemented.
109
110#### Stochastic minimization:
111
112is the most trivial learning method. This is the Robbins-Monro
113stochastic approximation applied to multilayer perceptrons. The
114weights are updated after each example according to the formula:
115\f$w_{ij}(t+1) = w_{ij}(t) + \Delta w_{ij}(t)\f$
116
117with
118
119\f$\Delta w_{ij}(t) = - \eta(d e_p / d w_{ij} + \delta) + \epsilon \Delta w_{ij}(t-1)\f$
120
121The parameters for this method are Eta, EtaDecay, Delta and
122Epsilon.
123
124#### Steepest descent with fixed step size (batch learning):
125
126It is the same as the stochastic
127minimization, but the weights are updated after considering all the
128examples, with the total derivative dEdw. The parameters for this
129method are Eta, EtaDecay, Delta and Epsilon.
130
131#### Steepest descent algorithm:
132
133Weights are set to the minimum along the line defined by the gradient. The
134only parameter for this method is Tau. Lower tau = higher precision =
135slower search. A value Tau = 3 seems reasonable.
136
137#### Conjugate gradients with the Polak-Ribiere updating formula:
138
139Weights are set to the minimum along the line defined by the conjugate gradient.
140Parameters are Tau and Reset, which defines the epochs where the direction is
141reset to the steepest descent.
142
143#### Conjugate gradients with the Fletcher-Reeves updating formula:
144
145Weights are set to the minimum along the line defined by the conjugate gradient. Parameters
146are Tau and Reset, which defines the epochs where the direction is
147reset to the steepest descent.
148
149#### Broyden, Fletcher, Goldfarb, Shanno (BFGS) method:
150
151 Implies the computation of a NxN matrix
152computation, but seems more powerful at least for less than 300
153weights. Parameters are Tau and Reset, which defines the epochs where
154the direction is reset to the steepest descent.
155
156### How to use it...
157
158TMLP is build from 3 classes: TNeuron, TSynapse and
159TMultiLayerPerceptron. Only TMultiLayerPerceptron should be used
160explicitly by the user.
161
162TMultiLayerPerceptron will take examples from a TTree
163given in the constructor. The network is described by a simple
164string: The input/output layers are defined by giving the expression for
165each neuron, separated by comas. Hidden layers are just described
166by the number of neurons. The layers are separated by colons.
167In addition, input/output layer formulas can be preceded by '@' (e.g "@out")
168if one wants to also normalize the data from the TTree.
169Input and outputs are taken from the TTree given as second argument.
170Expressions are evaluated as for TTree::Draw(), arrays are expended in
171distinct neurons, one for each index.
172This can only be done for fixed-size arrays.
173If the formula ends with "!", softmax functions are used for the output layer.
174One defines the training and test datasets by TEventLists.
175
176Example:
177~~~ {.cpp}
178TMultiLayerPerceptron("x,y:10:5:f",inputTree);
179~~~
180
181Both the TTree and the TEventLists can be defined in
182the constructor, or later with the suited setter method. The lists
183used for training and test can be defined either explicitly, or via
184a string containing the formula to be used to define them, exactly as
185for a TCut.
186
187The learning method is defined using the TMultiLayerPerceptron::SetLearningMethod() .
188Learning methods are :
189
190 - TMultiLayerPerceptron::kStochastic,
191 - TMultiLayerPerceptron::kBatch,
192 - TMultiLayerPerceptron::kSteepestDescent,
193 - TMultiLayerPerceptron::kRibierePolak,
194 - TMultiLayerPerceptron::kFletcherReeves,
195 - TMultiLayerPerceptron::kBFGS
196
197A weight can be assigned to events, either in the constructor, either
198with TMultiLayerPerceptron::SetEventWeight(). In addition, the TTree weight
199is taken into account.
200
201Finally, one starts the training with
202TMultiLayerPerceptron::Train(Int_t nepoch, Option_t* options). The
203first argument is the number of epochs while option is a string that
204can contain: "text" (simple text output) , "graph"
205(evoluting graphical training curves), "update=X" (step for
206the text/graph output update) or "+" (will skip the
207randomisation and start from the previous values). All combinations
208are available.
209
210Example:
211~~~ {.cpp}
212net.Train(100,"text, graph, update=10");
213~~~
214
215When the neural net is trained, it can be used
216directly ( TMultiLayerPerceptron::Evaluate() ) or exported to a
217standalone C++ code ( TMultiLayerPerceptron::Export() ).
218
219Finally, note that even if this implementation is inspired from the mlpfit code,
220the feature lists are not exactly matching:
221
222 - mlpfit hybrid learning method is not implemented
223 - output neurons can be normalized, this is not the case for mlpfit
224 - the neural net is exported in C++, FORTRAN or PYTHON
225 - the drawResult() method allows a fast check of the learning procedure
226
227In addition, the paw version of mlpfit had additional limitations on the number of
228neurons, hidden layers and inputs/outputs that does not apply to TMultiLayerPerceptron.
229*/
230
231
233#include "TSynapse.h"
234#include "TNeuron.h"
235#include "TClass.h"
236#include "TTree.h"
237#include "TEventList.h"
238#include "TRandom3.h"
239#include "TTimeStamp.h"
240#include "TRegexp.h"
241#include "TCanvas.h"
242#include "TH2.h"
243#include "TGraph.h"
244#include "TLegend.h"
245#include "TMultiGraph.h"
246#include "TDirectory.h"
247#include "TSystem.h"
248#include <iostream>
249#include <fstream>
250#include "TMath.h"
251#include "TTreeFormula.h"
252#include "TTreeFormulaManager.h"
253#include "TMarker.h"
254#include "TLine.h"
255#include "TText.h"
256#include "TObjString.h"
257#include <cstdlib>
258
259
260////////////////////////////////////////////////////////////////////////////////
261/// Default constructor
262
264{
265 if(!TClass::GetClass("TTreePlayer")) gSystem->Load("libTreePlayer");
266 fNetwork.SetOwner(true);
267 fFirstLayer.SetOwner(false);
268 fLastLayer.SetOwner(false);
269 fSynapses.SetOwner(true);
270 fData = nullptr;
271 fCurrentTree = -1;
273 fStructure = "";
274 fWeight = "1";
275 fTraining = nullptr;
276 fTrainingOwner = false;
277 fTest = nullptr;
278 fTestOwner = false;
279 fEventWeight = nullptr;
280 fManager = nullptr;
282 fEta = .1;
283 fEtaDecay = 1;
284 fDelta = 0;
285 fEpsilon = 0;
286 fTau = 3;
287 fLastAlpha = 0;
288 fReset = 50;
291 fextF = "";
292 fextD = "";
293}
294
295////////////////////////////////////////////////////////////////////////////////
296/// The network is described by a simple string:
297/// The input/output layers are defined by giving
298/// the branch names separated by comas.
299/// Hidden layers are just described by the number of neurons.
300/// The layers are separated by colons.
301///
302/// Ex: "x,y:10:5:f"
303///
304/// The output can be prepended by '@' if the variable has to be
305/// normalized.
306/// The output can be followed by '!' to use Softmax neurons for the
307/// output layer only.
308///
309/// Ex: "x,y:10:5:c1,c2,c3!"
310///
311/// Input and outputs are taken from the TTree given as second argument.
312/// training and test are the two TEventLists defining events
313/// to be used during the neural net training.
314/// Both the TTree and the TEventLists can be defined in the constructor,
315/// or later with the suited setter method.
316
321 const char* extF, const char* extD)
322{
323 if(!TClass::GetClass("TTreePlayer")) gSystem->Load("libTreePlayer");
324 fNetwork.SetOwner(true);
325 fFirstLayer.SetOwner(false);
326 fLastLayer.SetOwner(false);
327 fSynapses.SetOwner(true);
329 fData = data;
330 fCurrentTree = -1;
333 fTrainingOwner = false;
334 fTest = test;
335 fTestOwner = false;
336 fWeight = "1";
337 fType = type;
339 fextF = extF;
340 fextD = extD;
341 fEventWeight = nullptr;
342 fManager = nullptr;
343 if (data) {
344 BuildNetwork();
345 AttachData();
346 }
348 fEta = .1;
349 fEpsilon = 0;
350 fDelta = 0;
351 fEtaDecay = 1;
352 fTau = 3;
353 fLastAlpha = 0;
354 fReset = 50;
355}
356
357////////////////////////////////////////////////////////////////////////////////
358/// The network is described by a simple string:
359/// The input/output layers are defined by giving
360/// the branch names separated by comas.
361/// Hidden layers are just described by the number of neurons.
362/// The layers are separated by colons.
363///
364/// Ex: "x,y:10:5:f"
365///
366/// The output can be prepended by '@' if the variable has to be
367/// normalized.
368/// The output can be followed by '!' to use Softmax neurons for the
369/// output layer only.
370///
371/// Ex: "x,y:10:5:c1,c2,c3!"
372///
373/// Input and outputs are taken from the TTree given as second argument.
374/// training and test are the two TEventLists defining events
375/// to be used during the neural net training.
376/// Both the TTree and the TEventLists can be defined in the constructor,
377/// or later with the suited setter method.
378
380 const char * weight, TTree * data,
384 const char* extF, const char* extD)
385{
386 if(!TClass::GetClass("TTreePlayer")) gSystem->Load("libTreePlayer");
387 fNetwork.SetOwner(true);
388 fFirstLayer.SetOwner(false);
389 fLastLayer.SetOwner(false);
390 fSynapses.SetOwner(true);
392 fData = data;
393 fCurrentTree = -1;
396 fTrainingOwner = false;
397 fTest = test;
398 fTestOwner = false;
399 fWeight = weight;
400 fType = type;
402 fextF = extF;
403 fextD = extD;
404 fEventWeight = nullptr;
405 fManager = nullptr;
406 if (data) {
407 BuildNetwork();
408 AttachData();
409 }
411 fEta = .1;
412 fEtaDecay = 1;
413 fDelta = 0;
414 fEpsilon = 0;
415 fTau = 3;
416 fLastAlpha = 0;
417 fReset = 50;
418}
419
420////////////////////////////////////////////////////////////////////////////////
421/// The network is described by a simple string:
422/// The input/output layers are defined by giving
423/// the branch names separated by comas.
424/// Hidden layers are just described by the number of neurons.
425/// The layers are separated by colons.
426///
427/// Ex: "x,y:10:5:f"
428///
429/// The output can be prepended by '@' if the variable has to be
430/// normalized.
431/// The output can be followed by '!' to use Softmax neurons for the
432/// output layer only.
433///
434/// Ex: "x,y:10:5:c1,c2,c3!"
435///
436/// Input and outputs are taken from the TTree given as second argument.
437/// training and test are two cuts (see TTreeFormula) defining events
438/// to be used during the neural net training and testing.
439///
440/// Example: "Entry$%2", "(Entry$+1)%2".
441///
442/// Both the TTree and the cut can be defined in the constructor,
443/// or later with the suited setter method.
444
446 const char * training,
447 const char * test,
449 const char* extF, const char* extD)
450{
451 if(!TClass::GetClass("TTreePlayer")) gSystem->Load("libTreePlayer");
452 fNetwork.SetOwner(true);
453 fFirstLayer.SetOwner(false);
454 fLastLayer.SetOwner(false);
455 fSynapses.SetOwner(true);
457 fData = data;
458 fCurrentTree = -1;
460 {
462 fTraining = new TEventList(Form("fTrainingList_%zu",(size_t)this));
463 }
464 fTrainingOwner = true;
465 {
467 fTest = new TEventList(Form("fTestList_%zu",(size_t)this));
468 }
469 fTestOwner = true;
470 fWeight = "1";
472 if(testcut=="") testcut = Form("!(%s)",training);
473 fType = type;
475 fextF = extF;
476 fextD = extD;
477 fEventWeight = nullptr;
478 fManager = nullptr;
479 if (data) {
480 BuildNetwork();
481 data->Draw(Form(">>fTrainingList_%zu",(size_t)this),training,"goff");
482 data->Draw(Form(">>fTestList_%zu",(size_t)this),(const char *)testcut,"goff");
483 AttachData();
484 }
485 else {
486 Warning("TMultiLayerPerceptron::TMultiLayerPerceptron","Data not set. Cannot define datasets");
487 }
489 fEta = .1;
490 fEtaDecay = 1;
491 fDelta = 0;
492 fEpsilon = 0;
493 fTau = 3;
494 fLastAlpha = 0;
495 fReset = 50;
496}
497
498////////////////////////////////////////////////////////////////////////////////
499/// The network is described by a simple string:
500/// The input/output layers are defined by giving
501/// the branch names separated by comas.
502/// Hidden layers are just described by the number of neurons.
503/// The layers are separated by colons.
504///
505/// Ex: "x,y:10:5:f"
506///
507/// The output can be prepended by '@' if the variable has to be
508/// normalized.
509/// The output can be followed by '!' to use Softmax neurons for the
510/// output layer only.
511///
512/// Ex: "x,y:10:5:c1,c2,c3!"
513///
514/// Input and outputs are taken from the TTree given as second argument.
515/// training and test are two cuts (see TTreeFormula) defining events
516/// to be used during the neural net training and testing.
517///
518/// Example: "Entry$%2", "(Entry$+1)%2".
519///
520/// Both the TTree and the cut can be defined in the constructor,
521/// or later with the suited setter method.
522
524 const char * weight, TTree * data,
525 const char * training,
526 const char * test,
528 const char* extF, const char* extD)
529{
530 if(!TClass::GetClass("TTreePlayer")) gSystem->Load("libTreePlayer");
531 fNetwork.SetOwner(true);
532 fFirstLayer.SetOwner(false);
533 fLastLayer.SetOwner(false);
534 fSynapses.SetOwner(true);
536 fData = data;
537 fCurrentTree = -1;
539 {
541 fTraining = new TEventList(Form("fTrainingList_%zu",(size_t)this));
542 }
543 fTrainingOwner = true;
544 {
546 fTest = new TEventList(Form("fTestList_%zu",(size_t)this));
547 }
548 fTestOwner = true;
549 fWeight = weight;
551 if(testcut=="") testcut = Form("!(%s)",training);
552 fType = type;
554 fextF = extF;
555 fextD = extD;
556 fEventWeight = nullptr;
557 fManager = nullptr;
558 if (data) {
559 BuildNetwork();
560 data->Draw(Form(">>fTrainingList_%zu",(size_t)this),training,"goff");
561 data->Draw(Form(">>fTestList_%zu",(size_t)this),(const char *)testcut,"goff");
562 AttachData();
563 }
564 else {
565 Warning("TMultiLayerPerceptron::TMultiLayerPerceptron","Data not set. Cannot define datasets");
566 }
568 fEta = .1;
569 fEtaDecay = 1;
570 fDelta = 0;
571 fEpsilon = 0;
572 fTau = 3;
573 fLastAlpha = 0;
574 fReset = 50;
575}
576
577////////////////////////////////////////////////////////////////////////////////
578/// Destructor
579
585
586////////////////////////////////////////////////////////////////////////////////
587/// Set the data source
588
590{
591 if (fData) {
592 std::cerr << "Error: data already defined." << std::endl;
593 return;
594 }
595 fData = data;
596 if (data) {
597 BuildNetwork();
598 AttachData();
599 }
600}
601
602////////////////////////////////////////////////////////////////////////////////
603/// Set the event weight
604
606{
608 if (fData) {
609 if (fEventWeight) {
611 delete fEventWeight;
612 }
613 fManager->Add((fEventWeight = new TTreeFormula("NNweight",fWeight.Data(),fData)));
614 }
615}
616
617////////////////////////////////////////////////////////////////////////////////
618/// Sets the Training dataset.
619/// Those events will be used for the minimization
620
622{
623 if(fTraining && fTrainingOwner) delete fTraining;
624 fTraining = train;
625 fTrainingOwner = false;
626}
627
628////////////////////////////////////////////////////////////////////////////////
629/// Sets the Test dataset.
630/// Those events will not be used for the minimization but for control
631
633{
634 if(fTest && fTestOwner) delete fTest;
635 fTest = test;
636 fTestOwner = false;
637}
638
639////////////////////////////////////////////////////////////////////////////////
640/// Sets the Training dataset.
641/// Those events will be used for the minimization.
642/// Note that the tree must be already defined.
643
645{
646 if(fTraining && fTrainingOwner) delete fTraining;
647 {
649 fTraining = new TEventList(Form("fTrainingList_%zu",(size_t)this));
650 }
651 fTrainingOwner = true;
652 if (fData) {
653 fData->Draw(Form(">>fTrainingList_%zu",(size_t)this),train,"goff");
654 }
655 else {
656 Warning("TMultiLayerPerceptron::TMultiLayerPerceptron","Data not set. Cannot define datasets");
657 }
658}
659
660////////////////////////////////////////////////////////////////////////////////
661/// Sets the Test dataset.
662/// Those events will not be used for the minimization but for control.
663/// Note that the tree must be already defined.
664
666{
667 if(fTest && fTestOwner) {delete fTest; fTest=nullptr;}
668 if(fTest) if(strncmp(fTest->GetName(),Form("fTestList_%zu",(size_t)this),10)) delete fTest;
669 {
671 fTest = new TEventList(Form("fTestList_%zu",(size_t)this));
672 }
673 fTestOwner = true;
674 if (fData) {
675 fData->Draw(Form(">>fTestList_%zu",(size_t)this),test,"goff");
676 }
677 else {
678 Warning("TMultiLayerPerceptron::TMultiLayerPerceptron","Data not set. Cannot define datasets");
679 }
680}
681
682////////////////////////////////////////////////////////////////////////////////
683/// Sets the learning method.
684/// Available methods are: kStochastic, kBatch,
685/// kSteepestDescent, kRibierePolak, kFletcherReeves and kBFGS.
686/// (look at the constructor for the complete description
687/// of learning methods and parameters)
688
693
694////////////////////////////////////////////////////////////////////////////////
695/// Sets Eta - used in stochastic minimisation
696/// (look at the constructor for the complete description
697/// of learning methods and parameters)
698
700{
701 fEta = eta;
702}
703
704////////////////////////////////////////////////////////////////////////////////
705/// Sets Epsilon - used in stochastic minimisation
706/// (look at the constructor for the complete description
707/// of learning methods and parameters)
708
713
714////////////////////////////////////////////////////////////////////////////////
715/// Sets Delta - used in stochastic minimisation
716/// (look at the constructor for the complete description
717/// of learning methods and parameters)
718
720{
721 fDelta = delta;
722}
723
724////////////////////////////////////////////////////////////////////////////////
725/// Sets EtaDecay - Eta *= EtaDecay at each epoch
726/// (look at the constructor for the complete description
727/// of learning methods and parameters)
728
733
734////////////////////////////////////////////////////////////////////////////////
735/// Sets Tau - used in line search
736/// (look at the constructor for the complete description
737/// of learning methods and parameters)
738
740{
741 fTau = tau;
742}
743
744////////////////////////////////////////////////////////////////////////////////
745/// Sets number of epochs between two resets of the
746/// search direction to the steepest descent.
747/// (look at the constructor for the complete description
748/// of learning methods and parameters)
749
751{
752 fReset = reset;
753}
754
755////////////////////////////////////////////////////////////////////////////////
756/// Load an entry into the network
757
759{
760 if (!fData) return;
762 if (fData->GetTreeNumber() != fCurrentTree) {
763 ((TMultiLayerPerceptron*)this)->fCurrentTree = fData->GetTreeNumber();
764 fManager->Notify();
765 ((TMultiLayerPerceptron*)this)->fCurrentTreeWeight = fData->GetWeight();
766 }
768 for (Int_t i=0;i<nentries;i++) {
769 TNeuron *neuron = (TNeuron *)fNetwork.UncheckedAt(i);
770 neuron->SetNewEvent();
771 }
772}
773
774////////////////////////////////////////////////////////////////////////////////
775/// Train the network.
776/// nEpoch is the number of iterations.
777/// option can contain:
778/// - "text" (simple text output)
779/// - "graph" (evoluting graphical training curves)
780/// - "update=X" (step for the text/graph output update)
781/// - "+" will skip the randomisation and start from the previous values.
782/// - "current" (draw in the current canvas)
783/// - "minErrorTrain" (stop when NN error on the training sample gets below minE
784/// - "minErrorTest" (stop when NN error on the test sample gets below minE
785/// All combinations are available.
786
788{
789 Int_t i;
790 TString opt = option;
791 opt.ToLower();
792 // Decode options and prepare training.
793 Int_t verbosity = 0;
794 Bool_t newCanvas = true;
795 Bool_t minE_Train = false;
796 Bool_t minE_Test = false;
797 if (opt.Contains("text"))
798 verbosity += 1;
799 if (opt.Contains("graph"))
800 verbosity += 2;
802 if (opt.Contains("update=")) {
803 TRegexp reg("update=[0-9]*");
804 TString out = opt(reg);
805 displayStepping = atoi(out.Data() + 7);
806 }
807 if (opt.Contains("current"))
808 newCanvas = false;
809 if (opt.Contains("minerrortrain"))
810 minE_Train = true;
811 if (opt.Contains("minerrortest"))
812 minE_Test = true;
813 TVirtualPad *canvas = nullptr;
814 TMultiGraph *residual_plot = nullptr;
815 TGraph *train_residual_plot = nullptr;
816 TGraph *test_residual_plot = nullptr;
817 if ((!fData) || (!fTraining) || (!fTest)) {
818 Error("Train","Training/Test samples still not defined. Cannot train the neural network");
819 return;
820 }
821 Info("Train","Using %d train and %d test entries.",
822 fTraining->GetN(), fTest->GetN());
823 // Text and Graph outputs
824 if (verbosity % 2)
825 std::cout << "Training the Neural Network" << std::endl;
826 if (verbosity / 2) {
828 if(newCanvas)
829 canvas = new TCanvas("NNtraining", "Neural Net training");
830 else {
831 canvas = gPad;
832 if(!canvas) canvas = new TCanvas("NNtraining", "Neural Net training");
833 }
836 canvas->SetLeftMargin(0.14);
837 train_residual_plot->SetLineColor(4);
838 test_residual_plot->SetLineColor(2);
841 residual_plot->Draw("LA");
842 if (residual_plot->GetXaxis()) residual_plot->GetXaxis()->SetTitle("Epoch");
843 if (residual_plot->GetYaxis()) residual_plot->GetYaxis()->SetTitle("Error");
844 }
845 // If the option "+" is not set, one has to randomize the weights first
846 if (!opt.Contains("+"))
847 Randomize();
848 // Initialisation
849 fLastAlpha = 0;
851 Double_t *buffer = new Double_t[els];
852 Double_t *dir = new Double_t[els];
853 for (i = 0; i < els; i++)
854 buffer[i] = 0;
857 TMatrixD gamma(matrix_size, 1);
858 TMatrixD delta(matrix_size, 1);
859 // Epoch loop. Here is the training itself.
862 for (Int_t iepoch = 0; (iepoch < nEpoch) && (!minE_Train || training_E>minE) && (!minE_Test || test_E>minE) ; iepoch++) {
863 switch (fLearningMethod) {
865 {
866 MLP_Stochastic(buffer);
867 break;
868 }
870 {
871 ComputeDEDw();
872 MLP_Batch(buffer);
873 break;
874 }
876 {
877 ComputeDEDw();
878 SteepestDir(dir);
879 if (LineSearch(dir, buffer))
880 MLP_Batch(buffer);
881 break;
882 }
884 {
885 ComputeDEDw();
886 if (!(iepoch % fReset)) {
887 SteepestDir(dir);
888 } else {
889 Double_t norm = 0;
890 Double_t onorm = 0;
891 for (i = 0; i < els; i++)
892 onorm += dir[i] * dir[i];
893 Double_t prod = 0;
894 Int_t idx = 0;
895 TNeuron *neuron = nullptr;
896 TSynapse *synapse = nullptr;
898 for (i=0;i<nentries;i++) {
899 neuron = (TNeuron *) fNetwork.UncheckedAt(i);
900 prod -= dir[idx++] * neuron->GetDEDw();
901 norm += neuron->GetDEDw() * neuron->GetDEDw();
902 }
904 for (i=0;i<nentries;i++) {
906 prod -= dir[idx++] * synapse->GetDEDw();
907 norm += synapse->GetDEDw() * synapse->GetDEDw();
908 }
909 ConjugateGradientsDir(dir, (norm - prod) / onorm);
910 }
911 if (LineSearch(dir, buffer))
912 MLP_Batch(buffer);
913 break;
914 }
916 {
917 ComputeDEDw();
918 if (!(iepoch % fReset)) {
919 SteepestDir(dir);
920 } else {
921 Double_t norm = 0;
922 Double_t onorm = 0;
923 for (i = 0; i < els; i++)
924 onorm += dir[i] * dir[i];
925 TNeuron *neuron = nullptr;
926 TSynapse *synapse = nullptr;
928 for (i=0;i<nentries;i++) {
929 neuron = (TNeuron *) fNetwork.UncheckedAt(i);
930 norm += neuron->GetDEDw() * neuron->GetDEDw();
931 }
933 for (i=0;i<nentries;i++) {
935 norm += synapse->GetDEDw() * synapse->GetDEDw();
936 }
938 }
939 if (LineSearch(dir, buffer))
940 MLP_Batch(buffer);
941 break;
942 }
944 {
945 SetGammaDelta(gamma, delta, buffer);
946 if (!(iepoch % fReset)) {
947 SteepestDir(dir);
948 bfgsh.UnitMatrix();
949 } else {
950 if (GetBFGSH(bfgsh, gamma, delta)) {
951 SteepestDir(dir);
952 bfgsh.UnitMatrix();
953 } else {
954 BFGSDir(bfgsh, dir);
955 }
956 }
957 if (DerivDir(dir) > 0) {
958 SteepestDir(dir);
959 bfgsh.UnitMatrix();
960 }
961 if (LineSearch(dir, buffer)) {
962 bfgsh.UnitMatrix();
963 SteepestDir(dir);
964 if (LineSearch(dir, buffer)) {
965 Error("TMultiLayerPerceptron::Train()","Line search fail");
966 iepoch = nEpoch;
967 }
968 }
969 break;
970 }
971 }
972 // Security: would the learning lead to non real numbers,
973 // the learning should stop now.
975 Error("TMultiLayerPerceptron::Train()","Stop.");
976 iepoch = nEpoch;
977 }
978 // Process other ROOT events. Time penalty is less than
979 // 1/1000 sec/evt on a mobile AMD Athlon(tm) XP 1500+
983 // Intermediate graph and text output
984 if ((verbosity % 2) && ((!(iepoch % displayStepping)) || (iepoch == nEpoch - 1))) {
985 std::cout << "Epoch: " << iepoch
986 << " learn=" << training_E
987 << " test=" << test_E
988 << std::endl;
989 }
990 if (verbosity / 2) {
993 if (!iepoch) {
996 for (i = 1; i < nEpoch; i++) {
997 train_residual_plot->SetPoint(i, i, trp);
998 test_residual_plot->SetPoint(i, i, tep);
999 }
1000 }
1001 if ((!(iepoch % displayStepping)) || (iepoch == nEpoch - 1)) {
1002 if (residual_plot->GetYaxis()) {
1003 residual_plot->GetYaxis()->UnZoom();
1004 residual_plot->GetYaxis()->SetTitleOffset(1.4);
1005 residual_plot->GetYaxis()->SetDecimals();
1006 }
1007 canvas->Modified();
1008 canvas->Update();
1009 }
1010 }
1011 }
1012 // Cleaning
1013 delete [] buffer;
1014 delete [] dir;
1015 // Final Text and Graph outputs
1016 if (verbosity % 2)
1017 std::cout << "Training done." << std::endl;
1018 if (verbosity / 2) {
1019 TLegend *legend = new TLegend(.75, .80, .95, .95);
1020 legend->AddEntry(residual_plot->GetListOfGraphs()->At(0),
1021 "Training sample", "L");
1022 legend->AddEntry(residual_plot->GetListOfGraphs()->At(1),
1023 "Test sample", "L");
1024 legend->Draw();
1025 }
1026}
1027
1028////////////////////////////////////////////////////////////////////////////////
1029/// Computes the output for a given event.
1030/// Look at the output neuron designed by index.
1031
1033{
1034 GetEntry(event);
1035 TNeuron *out = (TNeuron *) (fLastLayer.At(index));
1036 if (out)
1037 return out->GetValue();
1038 else
1039 return 0;
1040}
1041
1042////////////////////////////////////////////////////////////////////////////////
1043/// Error on the output for a given event
1044
1046{
1047 GetEntry(event);
1048 Double_t error = 0;
1049 // look at 1st output neuron to determine type and error function
1051 if (nEntries == 0) return 0.0;
1052 switch (fOutType) {
1053 case (TNeuron::kSigmoid):
1054 error = GetCrossEntropyBinary();
1055 break;
1056 case (TNeuron::kSoftmax):
1057 error = GetCrossEntropy();
1058 break;
1059 case (TNeuron::kLinear):
1060 error = GetSumSquareError();
1061 break;
1062 default:
1063 // default to sum-of-squares error
1064 error = GetSumSquareError();
1065 }
1066 error *= fEventWeight->EvalInstance();
1067 error *= fCurrentTreeWeight;
1068 return error;
1069}
1070
1071////////////////////////////////////////////////////////////////////////////////
1072/// Error on the whole dataset
1073
1075{
1076 TEventList *list =
1078 Double_t error = 0;
1079 Int_t i;
1080 if (list) {
1081 Int_t nEvents = list->GetN();
1082 for (i = 0; i < nEvents; i++) {
1083 error += GetError(list->GetEntry(i));
1084 }
1085 } else if (fData) {
1086 Int_t nEvents = (Int_t) fData->GetEntries();
1087 for (i = 0; i < nEvents; i++) {
1088 error += GetError(i);
1089 }
1090 }
1091 return error;
1092}
1093
1094////////////////////////////////////////////////////////////////////////////////
1095/// Error on the output for a given event
1096
1098{
1099 Double_t error = 0;
1100 for (Int_t i = 0; i < fLastLayer.GetEntriesFast(); i++) {
1101 TNeuron *neuron = (TNeuron *) fLastLayer[i];
1102 error += neuron->GetError() * neuron->GetError();
1103 }
1104 return (error / 2.);
1105}
1106
1107////////////////////////////////////////////////////////////////////////////////
1108/// Cross entropy error for sigmoid output neurons, for a given event
1109
1111{
1112 Double_t error = 0;
1113 for (Int_t i = 0; i < fLastLayer.GetEntriesFast(); i++) {
1114 TNeuron *neuron = (TNeuron *) fLastLayer[i];
1115 Double_t output = neuron->GetValue(); // sigmoid output and target
1116 Double_t target = neuron->GetTarget(); // values lie in [0,1]
1117 if (target < DBL_EPSILON) {
1118 if (output == 1.0)
1119 error = DBL_MAX;
1120 else
1121 error -= TMath::Log(1 - output);
1122 } else
1123 if ((1 - target) < DBL_EPSILON) {
1124 if (output == 0.0)
1125 error = DBL_MAX;
1126 else
1127 error -= TMath::Log(output);
1128 } else {
1129 if (output == 0.0 || output == 1.0)
1130 error = DBL_MAX;
1131 else
1132 error -= target * TMath::Log(output / target) + (1-target) * TMath::Log((1 - output)/(1 - target));
1133 }
1134 }
1135 return error;
1136}
1137
1138////////////////////////////////////////////////////////////////////////////////
1139/// Cross entropy error for a softmax output neuron, for a given event
1140
1142{
1143 Double_t error = 0;
1144 for (Int_t i = 0; i < fLastLayer.GetEntriesFast(); i++) {
1145 TNeuron *neuron = (TNeuron *) fLastLayer[i];
1146 Double_t output = neuron->GetValue(); // softmax output and target
1147 Double_t target = neuron->GetTarget(); // values lie in [0,1]
1148 if (target > DBL_EPSILON) { // (target == 0) => dE = 0
1149 if (output == 0.0)
1150 error = DBL_MAX;
1151 else
1152 error -= target * TMath::Log(output / target);
1153 }
1154 }
1155 return error;
1156}
1157
1158////////////////////////////////////////////////////////////////////////////////
1159/// Compute the DEDw = sum on all training events of dedw for each weight
1160/// normalized by the number of events.
1161
1163{
1164 Int_t i,j;
1167 for (i=0;i<nentries;i++) {
1169 synapse->SetDEDw(0.);
1170 }
1171 TNeuron *neuron;
1173 for (i=0;i<nentries;i++) {
1174 neuron = (TNeuron *) fNetwork.UncheckedAt(i);
1175 neuron->SetDEDw(0.);
1176 }
1177 Double_t eventWeight = 1.;
1178 if (fTraining) {
1179 Int_t nEvents = fTraining->GetN();
1180 for (i = 0; i < nEvents; i++) {
1182 eventWeight = fEventWeight->EvalInstance();
1183 eventWeight *= fCurrentTreeWeight;
1185 for (j=0;j<nentries;j++) {
1187 synapse->SetDEDw(synapse->GetDEDw() + (synapse->GetDeDw()*eventWeight));
1188 }
1190 for (j=0;j<nentries;j++) {
1191 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
1192 neuron->SetDEDw(neuron->GetDEDw() + (neuron->GetDeDw()*eventWeight));
1193 }
1194 }
1196 for (j=0;j<nentries;j++) {
1198 synapse->SetDEDw(synapse->GetDEDw() / (Double_t) nEvents);
1199 }
1201 for (j=0;j<nentries;j++) {
1202 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
1203 neuron->SetDEDw(neuron->GetDEDw() / (Double_t) nEvents);
1204 }
1205 } else if (fData) {
1206 Int_t nEvents = (Int_t) fData->GetEntries();
1207 for (i = 0; i < nEvents; i++) {
1208 GetEntry(i);
1209 eventWeight = fEventWeight->EvalInstance();
1210 eventWeight *= fCurrentTreeWeight;
1212 for (j=0;j<nentries;j++) {
1214 synapse->SetDEDw(synapse->GetDEDw() + (synapse->GetDeDw()*eventWeight));
1215 }
1217 for (j=0;j<nentries;j++) {
1218 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
1219 neuron->SetDEDw(neuron->GetDEDw() + (neuron->GetDeDw()*eventWeight));
1220 }
1221 }
1223 for (j=0;j<nentries;j++) {
1225 synapse->SetDEDw(synapse->GetDEDw() / (Double_t) nEvents);
1226 }
1228 for (j=0;j<nentries;j++) {
1229 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
1230 neuron->SetDEDw(neuron->GetDEDw() / (Double_t) nEvents);
1231 }
1232 }
1233}
1234
1235////////////////////////////////////////////////////////////////////////////////
1236/// Randomize the weights
1237
1239{
1241 Int_t j;
1243 TNeuron *neuron;
1244 TTimeStamp ts;
1245 TRandom3 gen(ts.GetSec());
1246 for (j=0;j<nentries;j++) {
1248 synapse->SetWeight(gen.Rndm() - 0.5);
1249 }
1251 for (j=0;j<nentries;j++) {
1252 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
1253 neuron->SetWeight(gen.Rndm() - 0.5);
1254 }
1255}
1256
1257////////////////////////////////////////////////////////////////////////////////
1258/// Connects the TTree to Neurons in input and output
1259/// layers. The formulas associated to each neuron are created
1260/// and reported to the network formula manager.
1261/// By default, the branch is not normalised since this would degrade
1262/// performance for classification jobs.
1263/// Normalisation can be requested by putting '@' in front of the formula.
1264
1266{
1267 Int_t j = 0;
1268 TNeuron *neuron = nullptr;
1269 Bool_t normalize = false;
1271
1272 // Set the size of the internal array of parameters of the formula
1276
1277 //first layer
1278 const TString input = TString(fStructure(0, fStructure.First(':')));
1279 const TObjArray *inpL = input.Tokenize(", ");
1281 // make sure nentries == entries in inpL
1282 R__ASSERT(nentries == inpL->GetLast()+1);
1283 for (j=0;j<nentries;j++) {
1284 normalize = false;
1285 const TString brName = ((TObjString *)inpL->At(j))->GetString();
1286 neuron = (TNeuron *) fFirstLayer.UncheckedAt(j);
1287 if (brName[0]=='@')
1288 normalize = true;
1289 fManager->Add(neuron->UseBranch(fData,brName.Data() + (normalize?1:0)));
1290 if(!normalize) neuron->SetNormalisation(0., 1.);
1291 }
1292 delete inpL;
1293
1294 // last layer
1296 fStructure(fStructure.Last(':') + 1,
1297 fStructure.Length() - fStructure.Last(':')));
1298 const TObjArray *outL = output.Tokenize(", ");
1300 // make sure nentries == entries in outL
1301 R__ASSERT(nentries == outL->GetLast()+1);
1302 for (j=0;j<nentries;j++) {
1303 normalize = false;
1304 const TString brName = ((TObjString *)outL->At(j))->GetString();
1305 neuron = (TNeuron *) fLastLayer.UncheckedAt(j);
1306 if (brName[0]=='@')
1307 normalize = true;
1308 fManager->Add(neuron->UseBranch(fData,brName.Data() + (normalize?1:0)));
1309 if(!normalize) neuron->SetNormalisation(0., 1.);
1310 }
1311 delete outL;
1312
1313 fManager->Add((fEventWeight = new TTreeFormula("NNweight",fWeight.Data(),fData)));
1314 //fManager->Sync();
1315
1316 // Set the old values
1318}
1319
1320////////////////////////////////////////////////////////////////////////////////
1321/// Expand the structure of the first layer
1322
1324{
1326 const TObjArray *inpL = input.Tokenize(", ");
1327 Int_t nneurons = inpL->GetLast()+1;
1328
1330 fStructure(fStructure.First(':') + 1,
1331 fStructure.Length() - fStructure.First(':')));
1333 Int_t i = 0;
1334 // loop on input neurons
1335 for (i = 0; i<nneurons; i++) {
1336 const TString name = ((TObjString *)inpL->At(i))->GetString();
1337 TTreeFormula f("sizeTestFormula",name,fData);
1338 // Variable size arrays are unreliable
1339 if(f.GetMultiplicity()==1 && f.GetNdata()>1) {
1340 Warning("TMultiLayerPerceptron::ExpandStructure()","Variable size arrays cannot be used to build implicitly an input layer. The index 0 will be assumed.");
1341 }
1342 // Check if we are coping with an array... then expand
1343 // The array operator used is {}. It is detected in TNeuron, and
1344 // passed directly as instance index of the TTreeFormula,
1345 // so that complex compounds made of arrays can be used without
1346 // parsing the details.
1347 else if(f.GetNdata()>1) {
1348 for(Int_t j=0; j<f.GetNdata(); j++) {
1349 if(i||j) newInput += ",";
1350 newInput += name;
1351 newInput += "{";
1352 newInput += j;
1353 newInput += "}";
1354 }
1355 continue;
1356 }
1357 if(i) newInput += ",";
1358 newInput += name;
1359 }
1360 delete inpL;
1361
1362 // Save the result
1364}
1365
1366////////////////////////////////////////////////////////////////////////////////
1367/// Instantiates the network from the description
1368
1370{
1373 TString hidden = TString(
1374 fStructure(fStructure.First(':') + 1,
1375 fStructure.Last(':') - fStructure.First(':') - 1));
1377 fStructure(fStructure.Last(':') + 1,
1378 fStructure.Length() - fStructure.Last(':')));
1379 Int_t bll = atoi(TString(
1380 hidden(hidden.Last(':') + 1,
1381 hidden.Length() - (hidden.Last(':') + 1))).Data());
1382 if (input.Length() == 0) {
1383 Error("BuildNetwork()","malformed structure. No input layer.");
1384 return;
1385 }
1386 if (output.Length() == 0) {
1387 Error("BuildNetwork()","malformed structure. No output layer.");
1388 return;
1389 }
1391 BuildHiddenLayers(hidden);
1393}
1394
1395////////////////////////////////////////////////////////////////////////////////
1396/// Instantiates the neurons in input
1397/// Inputs are normalised and the type is set to kOff
1398/// (simple forward of the formula value)
1399
1401{
1402 const TObjArray *inpL = input.Tokenize(", ");
1403 const Int_t nneurons =inpL->GetLast()+1;
1404 TNeuron *neuron = nullptr;
1405 Int_t i = 0;
1406 for (i = 0; i<nneurons; i++) {
1407 const TString name = ((TObjString *)inpL->At(i))->GetString();
1408 neuron = new TNeuron(TNeuron::kOff, name);
1409 fFirstLayer.AddLast(neuron);
1410 fNetwork.AddLast(neuron);
1411 }
1412 delete inpL;
1413}
1414
1415////////////////////////////////////////////////////////////////////////////////
1416/// Builds hidden layers.
1417
1419{
1420 Int_t beg = 0;
1421 Int_t end = hidden.Index(":", beg + 1);
1422 Int_t prevStart = 0;
1424 Int_t layer = 1;
1425 while (end != -1) {
1426 BuildOneHiddenLayer(hidden(beg, end - beg), layer, prevStart, prevStop, false);
1427 beg = end + 1;
1428 end = hidden.Index(":", beg + 1);
1429 }
1430
1431 BuildOneHiddenLayer(hidden(beg, hidden.Length() - beg), layer, prevStart, prevStop, true);
1432}
1433
1434////////////////////////////////////////////////////////////////////////////////
1435/// Builds a hidden layer, updates the number of layers.
1436
1440{
1441 TNeuron *neuron = nullptr;
1442 TSynapse *synapse = nullptr;
1443 TString name;
1444 if (!sNumNodes.IsAlnum() || sNumNodes.IsAlpha()) {
1445 Error("BuildOneHiddenLayer",
1446 "The specification '%s' for hidden layer %d must contain only numbers!",
1447 sNumNodes.Data(), layer - 1);
1448 } else {
1449 Int_t num = atoi(sNumNodes.Data());
1450 for (Int_t i = 0; i < num; i++) {
1451 name.Form("HiddenL%d:N%d",layer,i);
1452 neuron = new TNeuron(fType, name, "", (const char*)fextF, (const char*)fextD);
1453 fNetwork.AddLast(neuron);
1454 for (Int_t j = prevStart; j < prevStop; j++) {
1455 synapse = new TSynapse((TNeuron *) fNetwork[j], neuron);
1457 }
1458 }
1459
1460 if (!lastLayer) {
1461 // tell each neuron which ones are in its own layer (for Softmax)
1463 for (Int_t i = prevStop; i < nEntries; i++) {
1464 neuron = (TNeuron *) fNetwork[i];
1465 for (Int_t j = prevStop; j < nEntries; j++)
1466 neuron->AddInLayer((TNeuron *) fNetwork[j]);
1467 }
1468 }
1469
1472 layer++;
1473 }
1474}
1475
1476////////////////////////////////////////////////////////////////////////////////
1477/// Builds the output layer
1478/// Neurons are linear combinations of input, by default.
1479/// If the structure ends with "!", neurons are set up for classification,
1480/// ie. with a sigmoid (1 neuron) or softmax (more neurons) activation function.
1481
1483{
1484 Int_t nneurons = output.CountChar(',')+1;
1485 if (fStructure.EndsWith("!")) {
1486 fStructure = TString(fStructure(0, fStructure.Length() - 1)); // remove "!"
1487 if (nneurons == 1)
1489 else
1491 }
1493 Int_t prevStart = prevStop - prev;
1494 Ssiz_t pos = 0;
1495 TNeuron *neuron;
1497 TString name;
1498 Int_t i,j;
1499 for (i = 0; i<nneurons; i++) {
1500 Ssiz_t nextpos=output.Index(",",pos);
1501 if (nextpos!=kNPOS)
1502 name=output(pos,nextpos-pos);
1503 else name=output(pos,output.Length());
1504 pos+=nextpos+1;
1505 neuron = new TNeuron(fOutType, name);
1506 for (j = prevStart; j < prevStop; j++) {
1507 synapse = new TSynapse((TNeuron *) fNetwork[j], neuron);
1509 }
1510 fLastLayer.AddLast(neuron);
1511 fNetwork.AddLast(neuron);
1512 }
1513 // tell each neuron which ones are in its own layer (for Softmax)
1515 for (i = prevStop; i < nEntries; i++) {
1516 neuron = (TNeuron *) fNetwork[i];
1517 for (j = prevStop; j < nEntries; j++)
1518 neuron->AddInLayer((TNeuron *) fNetwork[j]);
1519 }
1520
1521}
1522
1523////////////////////////////////////////////////////////////////////////////////
1524/// Draws the neural net output
1525/// It produces an histogram with the output for the two datasets.
1526/// Index is the number of the desired output neuron.
1527/// "option" can contain:
1528/// - test or train to select a dataset
1529/// - comp to produce a X-Y comparison plot
1530/// - nocanv to not create a new TCanvas for the plot
1531
1533{
1534 TString opt = option;
1535 opt.ToLower();
1536 TNeuron *out = (TNeuron *) (fLastLayer.At(index));
1537 if (!out) {
1538 Error("DrawResult()","no such output.");
1539 return;
1540 }
1541 //TCanvas *canvas = new TCanvas("NNresult", "Neural Net output");
1542 if (!opt.Contains("nocanv"))
1543 new TCanvas("NNresult", "Neural Net output");
1544 const Double_t *norm = out->GetNormalisation();
1545 TEventList *events = nullptr;
1547 Int_t i;
1548 if (opt.Contains("train")) {
1549 events = fTraining;
1550 setname = Form("train%d",index);
1551 } else if (opt.Contains("test")) {
1552 events = fTest;
1553 setname = Form("test%d",index);
1554 }
1555 if ((!fData) || (!events)) {
1556 Error("DrawResult()","no dataset.");
1557 return;
1558 }
1559 if (opt.Contains("comp")) {
1560 //comparison plot
1561 TString title = "Neural Net Output control. ";
1562 title += setname;
1563 setname = "MLP_" + setname + "_comp";
1564 TH2D *hist = ((TH2D *) gDirectory->Get(setname.Data()));
1565 if (!hist)
1566 hist = new TH2D(setname.Data(), title.Data(), 50, -1, 1, 50, -1, 1);
1567 hist->Reset();
1568 Int_t nEvents = events->GetN();
1569 for (i = 0; i < nEvents; i++) {
1570 GetEntry(events->GetEntry(i));
1571 hist->Fill(out->GetValue(), (out->GetBranch() - norm[1]) / norm[0]);
1572 }
1573 hist->Draw();
1574 } else {
1575 //output plot
1576 TString title = "Neural Net Output. ";
1577 title += setname;
1578 setname = "MLP_" + setname;
1579 TH1D *hist = ((TH1D *) gDirectory->Get(setname.Data()));
1580 if (!hist)
1581 hist = new TH1D(setname, title, 50, 1, -1);
1582 hist->Reset();
1583 Int_t nEvents = events->GetN();
1584 for (i = 0; i < nEvents; i++)
1585 hist->Fill(Result(events->GetEntry(i), index));
1586 hist->Draw();
1587 if (opt.Contains("train") && opt.Contains("test")) {
1588 events = fTraining;
1589 setname = "train";
1590 hist = ((TH1D *) gDirectory->Get("MLP_test"));
1591 if (!hist)
1592 hist = new TH1D(setname, title, 50, 1, -1);
1593 hist->Reset();
1594 nEvents = events->GetN();
1595 for (i = 0; i < nEvents; i++)
1596 hist->Fill(Result(events->GetEntry(i), index));
1597 hist->Draw("same");
1598 }
1599 }
1600}
1601
1602////////////////////////////////////////////////////////////////////////////////
1603/// Dumps the weights to a text file.
1604/// Set filename to "-" (default) to dump to the standard output
1605
1607{
1609 std::ostream * output;
1610 if (filen == "") {
1611 Error("TMultiLayerPerceptron::DumpWeights()","Invalid file name");
1612 return kFALSE;
1613 }
1614 if (filen == "-")
1615 output = &std::cout;
1616 else
1617 output = new std::ofstream(filen.Data());
1618 TNeuron *neuron = nullptr;
1619 *output << "#input normalization" << std::endl;
1621 Int_t j=0;
1622 for (j=0;j<nentries;j++) {
1623 neuron = (TNeuron *) fFirstLayer.UncheckedAt(j);
1624 *output << neuron->GetNormalisation()[0] << " "
1625 << neuron->GetNormalisation()[1] << std::endl;
1626 }
1627 *output << "#output normalization" << std::endl;
1629 for (j=0;j<nentries;j++) {
1630 neuron = (TNeuron *) fLastLayer.UncheckedAt(j);
1631 *output << neuron->GetNormalisation()[0] << " "
1632 << neuron->GetNormalisation()[1] << std::endl;
1633 }
1634 *output << "#neurons weights" << std::endl;
1636 while ((neuron = (TNeuron *) it->Next()))
1637 *output << neuron->GetWeight() << std::endl;
1638 delete it;
1640 TSynapse *synapse = nullptr;
1641 *output << "#synapses weights" << std::endl;
1642 while ((synapse = (TSynapse *) it->Next()))
1643 *output << synapse->GetWeight() << std::endl;
1644 delete it;
1645 if (filen != "-") {
1646 ((std::ofstream *) output)->close();
1647 delete output;
1648 }
1649 return kTRUE;
1650}
1651
1652////////////////////////////////////////////////////////////////////////////////
1653/// Loads the weights from a text file conforming to the format
1654/// defined by DumpWeights.
1655
1657{
1659 Double_t w;
1660 if (filen == "") {
1661 Error("TMultiLayerPerceptron::LoadWeights()","Invalid file name");
1662 return kFALSE;
1663 }
1664 char *buff = new char[100];
1665 std::ifstream input(filen.Data());
1666 // input normalzation
1667 input.getline(buff, 100);
1669 Float_t n1,n2;
1670 TNeuron *neuron = nullptr;
1671 while ((neuron = (TNeuron *) it->Next())) {
1672 input >> n1 >> n2;
1673 neuron->SetNormalisation(n2,n1);
1674 }
1675 input.getline(buff, 100);
1676 // output normalization
1677 input.getline(buff, 100);
1678 delete it;
1680 while ((neuron = (TNeuron *) it->Next())) {
1681 input >> n1 >> n2;
1682 neuron->SetNormalisation(n2,n1);
1683 }
1684 input.getline(buff, 100);
1685 // neuron weights
1686 input.getline(buff, 100);
1687 delete it;
1689 while ((neuron = (TNeuron *) it->Next())) {
1690 input >> w;
1691 neuron->SetWeight(w);
1692 }
1693 delete it;
1694 input.getline(buff, 100);
1695 // synapse weights
1696 input.getline(buff, 100);
1698 TSynapse *synapse = nullptr;
1699 while ((synapse = (TSynapse *) it->Next())) {
1700 input >> w;
1701 synapse->SetWeight(w);
1702 }
1703 delete it;
1704 delete[] buff;
1705 return kTRUE;
1706}
1707
1708////////////////////////////////////////////////////////////////////////////////
1709/// Returns the Neural Net for a given set of input parameters
1710/// #%parameters must equal #%input neurons
1711
1713{
1715 TNeuron *neuron;
1716 while ((neuron = (TNeuron *) it->Next()))
1717 neuron->SetNewEvent();
1718 delete it;
1720 Int_t i=0;
1721 while ((neuron = (TNeuron *) it->Next()))
1722 neuron->ForceExternalValue(params[i++]);
1723 delete it;
1724 TNeuron *out = (TNeuron *) (fLastLayer.At(index));
1725 if (out)
1726 return out->GetValue();
1727 else
1728 return 0;
1729}
1730
1731////////////////////////////////////////////////////////////////////////////////
1732/// Exports the NN as a function for any non-ROOT-dependant code
1733/// Supported languages are: only C++ , FORTRAN and Python (yet)
1734/// This feature is also useful if you want to plot the NN as
1735/// a function (TF1 or TF2).
1736
1738{
1740 lg.ToUpper();
1741 Int_t i;
1743 Warning("TMultiLayerPerceptron::Export","Request to export a network using an external function");
1744 }
1745 if (lg == "C++") {
1747 Int_t slash = basefilename.Last('/')+1;
1749
1750 TString classname = basefilename;
1751 TString header = filename;
1752 header += ".h";
1754 source += ".cxx";
1755 std::ofstream headerfile(header);
1756 std::ofstream sourcefile(source);
1757 headerfile << "#ifndef " << basefilename << "_h" << std::endl;
1758 headerfile << "#define " << basefilename << "_h" << std::endl << std::endl;
1759 headerfile << "class " << classname << " { " << std::endl;
1760 headerfile << "public:" << std::endl;
1761 headerfile << " " << classname << "() {}" << std::endl;
1762 headerfile << " ~" << classname << "() {}" << std::endl;
1763 sourcefile << "#include \"" << header << "\"" << std::endl;
1764 sourcefile << "#include <cmath>" << std::endl << std::endl;
1765 headerfile << " double Value(int index";
1766 sourcefile << "double " << classname << "::Value(int index";
1767 for (i = 0; i < fFirstLayer.GetEntriesFast(); i++) {
1768 headerfile << ",double in" << i;
1769 sourcefile << ",double in" << i;
1770 }
1771 headerfile << ");" << std::endl;
1772 sourcefile << ") {" << std::endl;
1773 for (i = 0; i < fFirstLayer.GetEntriesFast(); i++)
1774 sourcefile << " input" << i << " = (in" << i << " - "
1775 << ((TNeuron *) fFirstLayer[i])->GetNormalisation()[1] << ")/"
1776 << ((TNeuron *) fFirstLayer[i])->GetNormalisation()[0] << ";"
1777 << std::endl;
1778 sourcefile << " switch(index) {" << std::endl;
1779 TNeuron *neuron;
1781 Int_t idx = 0;
1782 while ((neuron = (TNeuron *) it->Next()))
1783 sourcefile << " case " << idx++ << ":" << std::endl
1784 << " return neuron" << neuron << "();" << std::endl;
1785 sourcefile << " default:" << std::endl
1786 << " return 0.;" << std::endl << " }"
1787 << std::endl;
1788 sourcefile << "}" << std::endl << std::endl;
1789 headerfile << " double Value(int index, double* input);" << std::endl;
1790 sourcefile << "double " << classname << "::Value(int index, double* input) {" << std::endl;
1791 for (i = 0; i < fFirstLayer.GetEntriesFast(); i++)
1792 sourcefile << " input" << i << " = (input[" << i << "] - "
1793 << ((TNeuron *) fFirstLayer[i])->GetNormalisation()[1] << ")/"
1794 << ((TNeuron *) fFirstLayer[i])->GetNormalisation()[0] << ";"
1795 << std::endl;
1796 sourcefile << " switch(index) {" << std::endl;
1797 delete it;
1799 idx = 0;
1800 while ((neuron = (TNeuron *) it->Next()))
1801 sourcefile << " case " << idx++ << ":" << std::endl
1802 << " return neuron" << neuron << "();" << std::endl;
1803 sourcefile << " default:" << std::endl
1804 << " return 0.;" << std::endl << " }"
1805 << std::endl;
1806 sourcefile << "}" << std::endl << std::endl;
1807 headerfile << "private:" << std::endl;
1808 for (i = 0; i < fFirstLayer.GetEntriesFast(); i++)
1809 headerfile << " double input" << i << ";" << std::endl;
1810 delete it;
1812 idx = 0;
1813 while ((neuron = (TNeuron *) it->Next())) {
1814 if (!neuron->GetPre(0)) {
1815 headerfile << " double neuron" << neuron << "();" << std::endl;
1816 sourcefile << "double " << classname << "::neuron" << neuron
1817 << "() {" << std::endl;
1818 sourcefile << " return input" << idx++ << ";" << std::endl;
1819 sourcefile << "}" << std::endl << std::endl;
1820 } else {
1821 headerfile << " double input" << neuron << "();" << std::endl;
1822 sourcefile << "double " << classname << "::input" << neuron
1823 << "() {" << std::endl;
1824 sourcefile << " double input = " << neuron->GetWeight()
1825 << ";" << std::endl;
1826 TSynapse *syn = nullptr;
1827 Int_t n = 0;
1828 while ((syn = neuron->GetPre(n++))) {
1829 sourcefile << " input += synapse" << syn << "();" << std::endl;
1830 }
1831 sourcefile << " return input;" << std::endl;
1832 sourcefile << "}" << std::endl << std::endl;
1833
1834 headerfile << " double neuron" << neuron << "();" << std::endl;
1835 sourcefile << "double " << classname << "::neuron" << neuron << "() {" << std::endl;
1836 sourcefile << " double input = input" << neuron << "();" << std::endl;
1837 switch(neuron->GetType()) {
1838 case (TNeuron::kSigmoid):
1839 {
1840 sourcefile << " return ((input < -709. ? 0. : (1/(1+exp(-input)))) * ";
1841 break;
1842 }
1843 case (TNeuron::kLinear):
1844 {
1845 sourcefile << " return (input * ";
1846 break;
1847 }
1848 case (TNeuron::kTanh):
1849 {
1850 sourcefile << " return (tanh(input) * ";
1851 break;
1852 }
1853 case (TNeuron::kGauss):
1854 {
1855 sourcefile << " return (exp(-input*input) * ";
1856 break;
1857 }
1858 case (TNeuron::kSoftmax):
1859 {
1860 sourcefile << " return (exp(input) / (";
1861 Int_t nn = 0;
1862 TNeuron* side = neuron->GetInLayer(nn++);
1863 sourcefile << "exp(input" << side << "())";
1864 while ((side = neuron->GetInLayer(nn++)))
1865 sourcefile << " + exp(input" << side << "())";
1866 sourcefile << ") * ";
1867 break;
1868 }
1869 default:
1870 {
1871 sourcefile << " return (0.0 * ";
1872 }
1873 }
1874 sourcefile << neuron->GetNormalisation()[0] << ")+" ;
1875 sourcefile << neuron->GetNormalisation()[1] << ";" << std::endl;
1876 sourcefile << "}" << std::endl << std::endl;
1877 }
1878 }
1879 delete it;
1880 TSynapse *synapse = nullptr;
1882 while ((synapse = (TSynapse *) it->Next())) {
1883 headerfile << " double synapse" << synapse << "();" << std::endl;
1884 sourcefile << "double " << classname << "::synapse"
1885 << synapse << "() {" << std::endl;
1886 sourcefile << " return (neuron" << synapse->GetPre()
1887 << "()*" << synapse->GetWeight() << ");" << std::endl;
1888 sourcefile << "}" << std::endl << std::endl;
1889 }
1890 delete it;
1891 headerfile << "};" << std::endl << std::endl;
1892 headerfile << "#endif // " << basefilename << "_h" << std::endl << std::endl;
1893 headerfile.close();
1894 sourcefile.close();
1895 std::cout << header << " and " << source << " created." << std::endl;
1896 }
1897 else if(lg == "FORTRAN") {
1898 TString implicit = " implicit double precision (a-h,n-z)\n";
1899 std::ofstream sigmoid("sigmoid.f");
1900 sigmoid << " double precision FUNCTION SIGMOID(X)" << std::endl
1901 << implicit
1902 << " IF(X.GT.37.) THEN" << std::endl
1903 << " SIGMOID = 1." << std::endl
1904 << " ELSE IF(X.LT.-709.) THEN" << std::endl
1905 << " SIGMOID = 0." << std::endl
1906 << " ELSE" << std::endl
1907 << " SIGMOID = 1./(1.+EXP(-X))" << std::endl
1908 << " ENDIF" << std::endl
1909 << " END" << std::endl;
1910 sigmoid.close();
1912 source += ".f";
1913 std::ofstream sourcefile(source);
1914
1915 // Header
1916 sourcefile << " double precision function " << filename
1917 << "(x, index)" << std::endl;
1919 sourcefile << " double precision x(" <<
1920 fFirstLayer.GetEntriesFast() << ")" << std::endl << std::endl;
1921
1922 // Last layer
1923 sourcefile << "C --- Last Layer" << std::endl;
1924 TNeuron *neuron;
1926 Int_t idx = 0;
1927 TString ifelseif = " if (index.eq.";
1928 while ((neuron = (TNeuron *) it->Next())) {
1929 sourcefile << ifelseif.Data() << idx++ << ") then" << std::endl
1930 << " " << filename
1931 << "=neuron" << neuron << "(x);" << std::endl;
1932 ifelseif = " else if (index.eq.";
1933 }
1934 sourcefile << " else" << std::endl
1935 << " " << filename << "=0.d0" << std::endl
1936 << " endif" << std::endl;
1937 sourcefile << " end" << std::endl;
1938
1939 // Network
1940 sourcefile << "C --- First and Hidden layers" << std::endl;
1941 delete it;
1943 idx = 0;
1944 while ((neuron = (TNeuron *) it->Next())) {
1945 sourcefile << " double precision function neuron"
1946 << neuron << "(x)" << std::endl
1947 << implicit;
1948 sourcefile << " double precision x("
1949 << fFirstLayer.GetEntriesFast() << ")" << std::endl << std::endl;
1950 if (!neuron->GetPre(0)) {
1951 sourcefile << " neuron" << neuron
1952 << " = (x(" << idx+1 << ") - "
1953 << ((TNeuron *) fFirstLayer[idx])->GetNormalisation()[1]
1954 << "d0)/"
1955 << ((TNeuron *) fFirstLayer[idx])->GetNormalisation()[0]
1956 << "d0" << std::endl;
1957 idx++;
1958 } else {
1959 sourcefile << " neuron" << neuron
1960 << " = " << neuron->GetWeight() << "d0" << std::endl;
1961 TSynapse *syn;
1962 Int_t n = 0;
1963 while ((syn = neuron->GetPre(n++)))
1964 sourcefile << " neuron" << neuron
1965 << " = neuron" << neuron
1966 << " + synapse" << syn << "(x)" << std::endl;
1967 switch(neuron->GetType()) {
1968 case (TNeuron::kSigmoid):
1969 {
1970 sourcefile << " neuron" << neuron
1971 << "= (sigmoid(neuron" << neuron << ")*";
1972 break;
1973 }
1974 case (TNeuron::kLinear):
1975 {
1976 break;
1977 }
1978 case (TNeuron::kTanh):
1979 {
1980 sourcefile << " neuron" << neuron
1981 << "= (tanh(neuron" << neuron << ")*";
1982 break;
1983 }
1984 case (TNeuron::kGauss):
1985 {
1986 sourcefile << " neuron" << neuron
1987 << "= (exp(-neuron" << neuron << "*neuron"
1988 << neuron << "))*";
1989 break;
1990 }
1991 case (TNeuron::kSoftmax):
1992 {
1993 Int_t nn = 0;
1994 TNeuron* side = neuron->GetInLayer(nn++);
1995 sourcefile << " div = exp(neuron" << side << "())" << std::endl;
1996 while ((side = neuron->GetInLayer(nn++)))
1997 sourcefile << " div = div + exp(neuron" << side << "())" << std::endl;
1998 sourcefile << " neuron" << neuron ;
1999 sourcefile << "= (exp(neuron" << neuron << ") / div * ";
2000 break;
2001 }
2002 default:
2003 {
2004 sourcefile << " neuron " << neuron << "= 0.";
2005 }
2006 }
2007 sourcefile << neuron->GetNormalisation()[0] << "d0)+" ;
2008 sourcefile << neuron->GetNormalisation()[1] << "d0" << std::endl;
2009 }
2010 sourcefile << " end" << std::endl;
2011 }
2012 delete it;
2013
2014 // Synapses
2015 sourcefile << "C --- Synapses" << std::endl;
2016 TSynapse *synapse = nullptr;
2018 while ((synapse = (TSynapse *) it->Next())) {
2019 sourcefile << " double precision function " << "synapse"
2020 << synapse << "(x)\n" << implicit;
2021 sourcefile << " double precision x("
2022 << fFirstLayer.GetEntriesFast() << ")" << std::endl << std::endl;
2023 sourcefile << " synapse" << synapse
2024 << "=neuron" << synapse->GetPre()
2025 << "(x)*" << synapse->GetWeight() << "d0" << std::endl;
2026 sourcefile << " end" << std::endl << std::endl;
2027 }
2028 delete it;
2029 sourcefile.close();
2030 std::cout << source << " created." << std::endl;
2031 }
2032 else if(lg == "PYTHON") {
2033 TString classname = filename;
2035 pyfile += ".py";
2036 std::ofstream pythonfile(pyfile);
2037 pythonfile << "from math import exp" << std::endl << std::endl;
2038 pythonfile << "from math import tanh" << std::endl << std::endl;
2039 pythonfile << "class " << classname << ":" << std::endl;
2040 pythonfile << "\tdef value(self,index";
2041 for (i = 0; i < fFirstLayer.GetEntriesFast(); i++) {
2042 pythonfile << ",in" << i;
2043 }
2044 pythonfile << "):" << std::endl;
2045 for (i = 0; i < fFirstLayer.GetEntriesFast(); i++)
2046 pythonfile << "\t\tself.input" << i << " = (in" << i << " - "
2047 << ((TNeuron *) fFirstLayer[i])->GetNormalisation()[1] << ")/"
2048 << ((TNeuron *) fFirstLayer[i])->GetNormalisation()[0] << std::endl;
2049 TNeuron *neuron;
2051 Int_t idx = 0;
2052 while ((neuron = (TNeuron *) it->Next()))
2053 pythonfile << "\t\tif index==" << idx++
2054 << ": return self.neuron" << neuron << "();" << std::endl;
2055 pythonfile << "\t\treturn 0." << std::endl;
2056 delete it;
2058 idx = 0;
2059 while ((neuron = (TNeuron *) it->Next())) {
2060 pythonfile << "\tdef neuron" << neuron << "(self):" << std::endl;
2061 if (!neuron->GetPre(0))
2062 pythonfile << "\t\treturn self.input" << idx++ << std::endl;
2063 else {
2064 pythonfile << "\t\tinput = " << neuron->GetWeight() << std::endl;
2065 TSynapse *syn;
2066 Int_t n = 0;
2067 while ((syn = neuron->GetPre(n++)))
2068 pythonfile << "\t\tinput = input + self.synapse"
2069 << syn << "()" << std::endl;
2070 switch(neuron->GetType()) {
2071 case (TNeuron::kSigmoid):
2072 {
2073 pythonfile << "\t\tif input<-709. : return " << neuron->GetNormalisation()[1] << std::endl;
2074 pythonfile << "\t\treturn ((1/(1+exp(-input)))*";
2075 break;
2076 }
2077 case (TNeuron::kLinear):
2078 {
2079 pythonfile << "\t\treturn (input*";
2080 break;
2081 }
2082 case (TNeuron::kTanh):
2083 {
2084 pythonfile << "\t\treturn (tanh(input)*";
2085 break;
2086 }
2087 case (TNeuron::kGauss):
2088 {
2089 pythonfile << "\t\treturn (exp(-input*input)*";
2090 break;
2091 }
2092 case (TNeuron::kSoftmax):
2093 {
2094 pythonfile << "\t\treturn (exp(input) / (";
2095 Int_t nn = 0;
2096 TNeuron* side = neuron->GetInLayer(nn++);
2097 pythonfile << "exp(self.neuron" << side << "())";
2098 while ((side = neuron->GetInLayer(nn++)))
2099 pythonfile << " + exp(self.neuron" << side << "())";
2100 pythonfile << ") * ";
2101 break;
2102 }
2103 default:
2104 {
2105 pythonfile << "\t\treturn 0.";
2106 }
2107 }
2108 pythonfile << neuron->GetNormalisation()[0] << ")+" ;
2109 pythonfile << neuron->GetNormalisation()[1] << std::endl;
2110 }
2111 }
2112 delete it;
2113 TSynapse *synapse = nullptr;
2115 while ((synapse = (TSynapse *) it->Next())) {
2116 pythonfile << "\tdef synapse" << synapse << "(self):" << std::endl;
2117 pythonfile << "\t\treturn (self.neuron" << synapse->GetPre()
2118 << "()*" << synapse->GetWeight() << ")" << std::endl;
2119 }
2120 delete it;
2121 pythonfile.close();
2122 std::cout << pyfile << " created." << std::endl;
2123 }
2124}
2125
2126////////////////////////////////////////////////////////////////////////////////
2127/// Shuffle the Int_t index[n] in input.
2128///
2129/// Input:
2130/// - index: the array to shuffle
2131/// - n: the size of the array
2132///
2133/// Output:
2134/// - index: the shuffled indexes
2135///
2136/// This method is used for stochastic training
2137
2139{
2140 TTimeStamp ts;
2141 TRandom3 rnd(ts.GetSec());
2142 Int_t j, k;
2143 Int_t a = n - 1;
2144 for (Int_t i = 0; i < n; i++) {
2145 j = (Int_t) (rnd.Rndm() * a);
2146 k = index[j];
2147 index[j] = index[i];
2148 index[i] = k;
2149 }
2150 return;
2151}
2152
2153////////////////////////////////////////////////////////////////////////////////
2154/// One step for the stochastic method
2155/// buffer should contain the previous dw vector and will be updated
2156
2158{
2159 Int_t nEvents = fTraining->GetN();
2160 Int_t *index = new Int_t[nEvents];
2161 Int_t i,j,nentries;
2162 for (i = 0; i < nEvents; i++)
2163 index[i] = i;
2164 fEta *= fEtaDecay;
2165 Shuffle(index, nEvents);
2166 TNeuron *neuron;
2168 for (i = 0; i < nEvents; i++) {
2170 // First compute DeDw for all neurons: force calculation before
2171 // modifying the weights.
2173 for (j=0;j<nentries;j++) {
2174 neuron = (TNeuron *) fFirstLayer.UncheckedAt(j);
2175 neuron->GetDeDw();
2176 }
2177 Int_t cnt = 0;
2178 // Step for all neurons
2180 for (j=0;j<nentries;j++) {
2181 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
2182 buffer[cnt] = (-fEta) * (neuron->GetDeDw() + fDelta)
2183 + fEpsilon * buffer[cnt];
2184 neuron->SetWeight(neuron->GetWeight() + buffer[cnt++]);
2185 }
2186 // Step for all synapses
2188 for (j=0;j<nentries;j++) {
2190 buffer[cnt] = (-fEta) * (synapse->GetDeDw() + fDelta)
2191 + fEpsilon * buffer[cnt];
2192 synapse->SetWeight(synapse->GetWeight() + buffer[cnt++]);
2193 }
2194 }
2195 delete[]index;
2196}
2197
2198////////////////////////////////////////////////////////////////////////////////
2199/// One step for the batch (stochastic) method.
2200/// DEDw should have been updated before calling this.
2201
2203{
2204 fEta *= fEtaDecay;
2205 Int_t cnt = 0;
2207 TNeuron *neuron = nullptr;
2208 // Step for all neurons
2209 while ((neuron = (TNeuron *) it->Next())) {
2210 buffer[cnt] = (-fEta) * (neuron->GetDEDw() + fDelta)
2211 + fEpsilon * buffer[cnt];
2212 neuron->SetWeight(neuron->GetWeight() + buffer[cnt++]);
2213 }
2214 delete it;
2216 TSynapse *synapse = nullptr;
2217 // Step for all synapses
2218 while ((synapse = (TSynapse *) it->Next())) {
2219 buffer[cnt] = (-fEta) * (synapse->GetDEDw() + fDelta)
2220 + fEpsilon * buffer[cnt];
2221 synapse->SetWeight(synapse->GetWeight() + buffer[cnt++]);
2222 }
2223 delete it;
2224}
2225
2226////////////////////////////////////////////////////////////////////////////////
2227/// Sets the weights to a point along a line
2228/// Weights are set to [origin + (dist * dir)].
2229
2231{
2232 Int_t idx = 0;
2233 TNeuron *neuron = nullptr;
2234 TSynapse *synapse = nullptr;
2236 while ((neuron = (TNeuron *) it->Next())) {
2237 neuron->SetWeight(origin[idx] + (dir[idx] * dist));
2238 idx++;
2239 }
2240 delete it;
2242 while ((synapse = (TSynapse *) it->Next())) {
2243 synapse->SetWeight(origin[idx] + (dir[idx] * dist));
2244 idx++;
2245 }
2246 delete it;
2247}
2248
2249////////////////////////////////////////////////////////////////////////////////
2250/// Sets the search direction to steepest descent.
2251
2253{
2254 Int_t idx = 0;
2255 TNeuron *neuron = nullptr;
2256 TSynapse *synapse = nullptr;
2258 while ((neuron = (TNeuron *) it->Next()))
2259 dir[idx++] = -neuron->GetDEDw();
2260 delete it;
2262 while ((synapse = (TSynapse *) it->Next()))
2263 dir[idx++] = -synapse->GetDEDw();
2264 delete it;
2265}
2266
2267////////////////////////////////////////////////////////////////////////////////
2268/// Search along the line defined by direction.
2269/// buffer is not used but is updated with the new dw
2270/// so that it can be used by a later stochastic step.
2271/// It returns true if the line search fails.
2272
2274{
2275 Int_t idx = 0;
2277 TNeuron *neuron = nullptr;
2278 TSynapse *synapse = nullptr;
2279 // store weights before line search
2283 for (j=0;j<nentries;j++) {
2284 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
2285 origin[idx++] = neuron->GetWeight();
2286 }
2288 for (j=0;j<nentries;j++) {
2290 origin[idx++] = synapse->GetWeight();
2291 }
2292 // try to find a triplet (alpha1, alpha2, alpha3) such that
2293 // Error(alpha1)>Error(alpha2)<Error(alpha3)
2295 Double_t alpha1 = 0.;
2297 if (alpha2 < 0.01)
2298 alpha2 = 0.01;
2299 if (alpha2 > 2.0)
2300 alpha2 = 2.0;
2304 Double_t err3 = err2;
2305 Bool_t bingo = false;
2306 Int_t icount;
2307 if (err1 > err2) {
2308 for (icount = 0; icount < 100; icount++) {
2309 alpha3 *= fTau;
2312 if (err3 > err2) {
2313 bingo = true;
2314 break;
2315 }
2316 alpha1 = alpha2;
2317 err1 = err2;
2318 alpha2 = alpha3;
2319 err2 = err3;
2320 }
2321 if (!bingo) {
2323 delete[]origin;
2324 return true;
2325 }
2326 } else {
2327 for (icount = 0; icount < 100; icount++) {
2328 alpha2 /= fTau;
2331 if (err1 > err2) {
2332 bingo = true;
2333 break;
2334 }
2335 alpha3 = alpha2;
2336 err3 = err2;
2337 }
2338 if (!bingo) {
2340 delete[]origin;
2341 fLastAlpha = 0.05;
2342 return true;
2343 }
2344 }
2345 // Sets the weights to the bottom of parabola
2346 fLastAlpha = 0.5 * (alpha1 + alpha3 -
2347 (err3 - err1) / ((err3 - err2) / (alpha3 - alpha2)
2348 - (err2 - err1) / (alpha2 - alpha1)));
2349 fLastAlpha = fLastAlpha < 10000 ? fLastAlpha : 10000;
2352 // Stores weight changes (can be used by a later stochastic step)
2353 idx = 0;
2355 for (j=0;j<nentries;j++) {
2356 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
2357 buffer[idx] = neuron->GetWeight() - origin[idx];
2358 idx++;
2359 }
2361 for (j=0;j<nentries;j++) {
2363 buffer[idx] = synapse->GetWeight() - origin[idx];
2364 idx++;
2365 }
2366 delete[]origin;
2367 return false;
2368}
2369
2370////////////////////////////////////////////////////////////////////////////////
2371/// Sets the search direction to conjugate gradient direction
2372/// beta should be:
2373///
2374/// \f$||g_{(t+1)}||^2 / ||g_{(t)}||^2\f$ (Fletcher-Reeves)
2375///
2376/// \f$g_{(t+1)} (g_{(t+1)}-g_{(t)}) / ||g_{(t)}||^2\f$ (Ribiere-Polak)
2377
2379{
2380 Int_t idx = 0;
2382 TNeuron *neuron = nullptr;
2383 TSynapse *synapse = nullptr;
2385 for (j=0;j<nentries;j++) {
2386 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
2387 dir[idx] = -neuron->GetDEDw() + beta * dir[idx];
2388 idx++;
2389 }
2391 for (j=0;j<nentries;j++) {
2393 dir[idx] = -synapse->GetDEDw() + beta * dir[idx];
2394 idx++;
2395 }
2396}
2397
2398////////////////////////////////////////////////////////////////////////////////
2399/// Computes the hessian matrix using the BFGS update algorithm.
2400/// from gamma (g_{(t+1)}-g_{(t)}) and delta (w_{(t+1)}-w_{(t)}).
2401/// It returns true if such a direction could not be found
2402/// (if gamma and delta are orthogonal).
2403
2405{
2406 TMatrixD gd(gamma, TMatrixD::kTransposeMult, delta);
2407 if ((Double_t) gd[0][0] == 0.)
2408 return true;
2412 Double_t a = 1 / (Double_t) gd[0][0];
2413 Double_t f = 1 + ((Double_t) gHg[0][0] * a);
2414 TMatrixD res( TMatrixD(delta, TMatrixD::kMult,
2416 res *= f;
2417 res -= (TMatrixD(delta, TMatrixD::kMult, tmp) +
2420 res *= a;
2421 bfgsh += res;
2422 return false;
2423}
2424
2425////////////////////////////////////////////////////////////////////////////////
2426/// Sets the gamma \f$(g_{(t+1)}-g_{(t)})\f$ and delta \f$(w_{(t+1)}-w_{(t)})\f$ vectors
2427/// Gamma is computed here, so ComputeDEDw cannot have been called before,
2428/// and delta is a direct translation of buffer into a TMatrixD.
2429
2431 Double_t * buffer)
2432{
2434 Int_t idx = 0;
2436 TNeuron *neuron = nullptr;
2437 TSynapse *synapse = nullptr;
2439 for (j=0;j<nentries;j++) {
2440 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
2441 gamma[idx++][0] = -neuron->GetDEDw();
2442 }
2444 for (j=0;j<nentries;j++) {
2446 gamma[idx++][0] = -synapse->GetDEDw();
2447 }
2448 for (Int_t i = 0; i < els; i++)
2449 delta[i].Assign(buffer[i]);
2450 //delta.SetElements(buffer,"F");
2451 ComputeDEDw();
2452 idx = 0;
2454 for (j=0;j<nentries;j++) {
2455 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
2456 gamma[idx++][0] += neuron->GetDEDw();
2457 }
2459 for (j=0;j<nentries;j++) {
2461 gamma[idx++][0] += synapse->GetDEDw();
2462 }
2463}
2464
2465////////////////////////////////////////////////////////////////////////////////
2466/// scalar product between gradient and direction
2467/// = derivative along direction
2468
2470{
2471 Int_t idx = 0;
2473 Double_t output = 0;
2474 TNeuron *neuron = nullptr;
2475 TSynapse *synapse = nullptr;
2477 for (j=0;j<nentries;j++) {
2478 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
2479 output += neuron->GetDEDw() * dir[idx++];
2480 }
2482 for (j=0;j<nentries;j++) {
2484 output += synapse->GetDEDw() * dir[idx++];
2485 }
2486 return output;
2487}
2488
2489////////////////////////////////////////////////////////////////////////////////
2490/// Computes the direction for the BFGS algorithm as the product
2491/// between the Hessian estimate (bfgsh) and the dir.
2492
2494{
2496 TMatrixD dedw(els, 1);
2497 Int_t idx = 0;
2499 TNeuron *neuron = nullptr;
2500 TSynapse *synapse = nullptr;
2502 for (j=0;j<nentries;j++) {
2503 neuron = (TNeuron *) fNetwork.UncheckedAt(j);
2504 dedw[idx++][0] = neuron->GetDEDw();
2505 }
2507 for (j=0;j<nentries;j++) {
2509 dedw[idx++][0] = synapse->GetDEDw();
2510 }
2512 for (Int_t i = 0; i < els; i++)
2513 dir[i] = -direction[i][0];
2514 //direction.GetElements(dir,"F");
2515}
2516
2517////////////////////////////////////////////////////////////////////////////////
2518/// Draws the network structure.
2519/// Neurons are depicted by a blue disk, and synapses by
2520/// lines connecting neurons.
2521/// The line width is proportional to the weight.
2522
2524{
2525#define NeuronSize 2.5
2526
2528 Float_t xStep = 1./(nLayers+1.);
2529 Int_t layer;
2530 for(layer=0; layer< nLayers-1; layer++) {
2532 if(layer==0) {
2534 nNeurons_this = input.CountChar(',')+1;
2535 }
2536 else {
2537 Int_t cnt=0;
2538 TString hidden = TString(fStructure(fStructure.First(':') + 1,fStructure.Last(':') - fStructure.First(':') - 1));
2539 Int_t beg = 0;
2540 Int_t end = hidden.Index(":", beg + 1);
2541 while (end != -1) {
2542 Int_t num = atoi(TString(hidden(beg, end - beg)).Data());
2543 cnt++;
2544 beg = end + 1;
2545 end = hidden.Index(":", beg + 1);
2546 if(layer==cnt) nNeurons_this = num;
2547 }
2548 Int_t num = atoi(TString(hidden(beg, hidden.Length() - beg)).Data());
2549 cnt++;
2550 if(layer==cnt) nNeurons_this = num;
2551 }
2553 if(layer==nLayers-2) {
2555 nNeurons_next = output.CountChar(',')+1;
2556 }
2557 else {
2558 Int_t cnt=0;
2559 TString hidden = TString(fStructure(fStructure.First(':') + 1,fStructure.Last(':') - fStructure.First(':') - 1));
2560 Int_t beg = 0;
2561 Int_t end = hidden.Index(":", beg + 1);
2562 while (end != -1) {
2563 Int_t num = atoi(TString(hidden(beg, end - beg)).Data());
2564 cnt++;
2565 beg = end + 1;
2566 end = hidden.Index(":", beg + 1);
2567 if(layer+1==cnt) nNeurons_next = num;
2568 }
2569 Int_t num = atoi(TString(hidden(beg, hidden.Length() - beg)).Data());
2570 cnt++;
2571 if(layer+1==cnt) nNeurons_next = num;
2572 }
2576 TSynapse *theSynapse = nullptr;
2577 Float_t maxWeight = 0;
2578 while ((theSynapse = (TSynapse *) it->Next()))
2579 maxWeight = maxWeight < theSynapse->GetWeight() ? theSynapse->GetWeight() : maxWeight;
2580 delete it;
2585 synapse->Draw();
2586 theSynapse = (TSynapse *) it->Next();
2587 if (!theSynapse) continue;
2588 synapse->SetLineWidth(Int_t((theSynapse->GetWeight()/maxWeight)*10.));
2589 synapse->SetLineStyle(1);
2590 if(((TMath::Abs(theSynapse->GetWeight())/maxWeight)*10.)<0.5) synapse->SetLineStyle(2);
2591 if(((TMath::Abs(theSynapse->GetWeight())/maxWeight)*10.)<0.25) synapse->SetLineStyle(3);
2592 }
2593 }
2594 delete it;
2595 }
2596 for(layer=0; layer< nLayers; layer++) {
2597 Float_t nNeurons = 0;
2598 if(layer==0) {
2600 nNeurons = input.CountChar(',')+1;
2601 }
2602 else if(layer==nLayers-1) {
2604 nNeurons = output.CountChar(',')+1;
2605 }
2606 else {
2607 Int_t cnt=0;
2608 TString hidden = TString(fStructure(fStructure.First(':') + 1,fStructure.Last(':') - fStructure.First(':') - 1));
2609 Int_t beg = 0;
2610 Int_t end = hidden.Index(":", beg + 1);
2611 while (end != -1) {
2612 Int_t num = atoi(TString(hidden(beg, end - beg)).Data());
2613 cnt++;
2614 beg = end + 1;
2615 end = hidden.Index(":", beg + 1);
2616 if(layer==cnt) nNeurons = num;
2617 }
2618 Int_t num = atoi(TString(hidden(beg, hidden.Length() - beg)).Data());
2619 cnt++;
2620 if(layer==cnt) nNeurons = num;
2621 }
2622 Float_t yStep = 1./(nNeurons+1.);
2623 for(Int_t neuron=0; neuron<nNeurons; neuron++) {
2624 TMarker* m = new TMarker(xStep*(layer+1),yStep*(neuron+1),20);
2625 m->SetMarkerColor(4);
2627 m->Draw();
2628 }
2629 }
2630 const TString input = TString(fStructure(0, fStructure.First(':')));
2631 const TObjArray *inpL = input.Tokenize(" ,");
2632 const Int_t nrItems = inpL->GetLast()+1;
2633 Float_t yStep = 1./(nrItems+1);
2634 for (Int_t item = 0; item < nrItems; item++) {
2635 const TString brName = ((TObjString *)inpL->At(item))->GetString();
2636 TText* label = new TText(0.5*xStep,yStep*(item+1),brName.Data());
2637 label->Draw();
2638 }
2639 delete inpL;
2640
2642 yStep=1./(numOutNodes+1);
2644 TNeuron* neuron=(TNeuron*)fLastLayer[outnode];
2645 if (neuron && neuron->GetName()) {
2646 TText* label = new TText(xStep*nLayers,
2647 yStep*(outnode+1),
2648 neuron->GetName());
2649 label->Draw();
2650 }
2651 }
2652}
#define f(i)
Definition RSha256.hxx:104
#define a(i)
Definition RSha256.hxx:99
int Int_t
Signed integer 4 bytes (int)
Definition RtypesCore.h:59
float Float_t
Float 4 bytes (float)
Definition RtypesCore.h:71
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
double Double_t
Double 8 bytes.
Definition RtypesCore.h:73
constexpr Ssiz_t kNPOS
The equivalent of std::string::npos for the ROOT class TString.
Definition RtypesCore.h:131
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
#define gDirectory
Definition TDirectory.h:385
#define R__ASSERT(e)
Checks condition e and reports a fatal error if it's false.
Definition TError.h:125
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
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 filename
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 target
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
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
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void reg
char name[80]
Definition TGX11.cxx:110
int nentries
TMatrixT< Double_t > TMatrixD
Definition TMatrixDfwd.h:23
#define NeuronSize
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2495
R__EXTERN TSystem * gSystem
Definition TSystem.h:572
#define gPad
static void SetMaxima(Int_t maxop=1000, Int_t maxpar=1000, Int_t maxconst=1000)
static function to set the maximum value of 3 parameters
static void GetMaxima(Int_t &maxop, Int_t &maxpar, Int_t &maxconst)
static function to get the maximum value of 3 parameters -maxop : maximum number of operations -maxpa...
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition TAttMarker.h:39
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition TAttMarker.h:46
virtual void SetLeftMargin(Float_t leftmargin)
Set Pad left margin in fraction of the pad width.
Definition TAttPad.cxx:108
The Canvas class.
Definition TCanvas.h:23
static TClass * GetClass(const char *name, Bool_t load=kTRUE, Bool_t silent=kFALSE)
Static method returning pointer to TClass of the specified class name.
Definition TClass.cxx:2973
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
TDirectory::TContext keeps track and restore the current directory.
Definition TDirectory.h:89
<div class="legacybox"><h2>Legacy Code</h2> TEventList is a legacy interface: there will be no bug fi...
Definition TEventList.h:31
virtual Long64_t GetEntry(Int_t index) const
Return value of entry at index in the list.
virtual Int_t GetN() const
Definition TEventList.h:56
A TGraph is an object made of two arrays X and Y with npoints each.
Definition TGraph.h:41
1-D histogram with a double per channel (see TH1 documentation)
Definition TH1.h:927
void Reset(Option_t *option="") override
Reset.
Definition TH1.cxx:10504
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3326
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3048
2-D histogram with a double per channel (see TH1 documentation)
Definition TH2.h:356
void Reset(Option_t *option="") override
Reset this histogram: contents, errors, etc.
Definition TH2.cxx:4202
Int_t Fill(Double_t) override
Invalid Fill method.
Definition TH2.cxx:356
This class displays a legend box (TPaveText) containing several legend entries.
Definition TLegend.h:23
Use the TLine constructor to create a simple line.
Definition TLine.h:22
Manages Markers.
Definition TMarker.h:22
void Draw(Option_t *option="") override
Draw this marker with its current attributes.
Definition TMarker.cxx:198
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition TMultiGraph.h:34
This class describes a neural network.
TTreeFormula * fEventWeight
! formula representing the event weight
void BuildOneHiddenLayer(const TString &sNumNodes, Int_t &layer, Int_t &prevStart, Int_t &prevStop, Bool_t lastLayer)
Builds a hidden layer, updates the number of layers.
void SteepestDir(Double_t *)
Sets the search direction to steepest descent.
void BuildNetwork()
Instantiates the network from the description.
TObjArray fNetwork
Collection of all the neurons in the 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.
TEventList * fTest
! EventList defining the events in the test dataset
bool GetBFGSH(TMatrixD &, TMatrixD &, TMatrixD &)
Computes the hessian matrix using the BFGS update algorithm.
void BuildHiddenLayers(TString &)
Builds hidden layers.
void BuildFirstLayer(TString &)
Instantiates the neurons in input Inputs are normalised and the type is set to kOff (simple forward o...
void SetTau(Double_t tau)
Sets Tau - used in line search (look at the constructor for the complete description of learning meth...
TMultiLayerPerceptron()
Default constructor.
Double_t GetSumSquareError() const
Error on the output for a given event.
void ConjugateGradientsDir(Double_t *, Double_t)
Sets the search direction to conjugate gradient direction beta should be:
Double_t fTau
! Tau - used in line search - Default=3.
TTree * fData
! pointer to the tree used as datasource
Double_t Result(Int_t event, Int_t index=0) const
Computes the output for a given event.
void SetGammaDelta(TMatrixD &, TMatrixD &, Double_t *)
Sets the gamma and delta vectors Gamma is computed here, so ComputeDEDw cannot have been called bef...
TEventList * fTraining
! EventList defining the events in the training dataset
TString fStructure
String containing the network structure.
Int_t fReset
! number of epochs between two resets of the search direction to the steepest descent - Default=50
Bool_t LoadWeights(Option_t *filename="")
Loads the weights from a text file conforming to the format defined by DumpWeights.
void MLP_Batch(Double_t *)
One step for the batch (stochastic) method.
TNeuron::ENeuronType fOutType
Type of output neurons.
Double_t fCurrentTreeWeight
! weight of the current tree in a chain
ELearningMethod fLearningMethod
! The Learning Method
Double_t fLastAlpha
! internal parameter used in line search
Int_t fCurrentTree
! index of the current tree in a chain
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++ ,...
Double_t fEpsilon
! Epsilon - used in stochastic minimisation - Default=0.
void Train(Int_t nEpoch, Option_t *option="text", Double_t minE=0)
Train the network.
TNeuron::ENeuronType GetType() const
void BFGSDir(TMatrixD &, Double_t *)
Computes the direction for the BFGS algorithm as the product between the Hessian estimate (bfgsh) and...
void SetTestDataSet(TEventList *test)
Sets the Test dataset.
Bool_t fTrainingOwner
! internal flag whether one has to delete fTraining or not
void SetLearningMethod(TMultiLayerPerceptron::ELearningMethod method)
Sets the learning method.
void SetTrainingDataSet(TEventList *train)
Sets the Training dataset.
void BuildLastLayer(TString &, Int_t)
Builds the output layer Neurons are linear combinations of input, by default.
Double_t fDelta
! Delta - used in stochastic minimisation - Default=0.
TTreeFormulaManager * fManager
! TTreeFormulaManager for the weight and neurons
void Randomize() const
Randomize the weights.
Bool_t LineSearch(Double_t *, Double_t *)
Search along the line defined by direction.
void ExpandStructure()
Expand the structure of the first layer.
Double_t fEta
! Eta - used in stochastic minimisation - Default=0.1
Double_t GetError(Int_t event) const
Error on the output for a given event.
Double_t fEtaDecay
! EtaDecay - Eta *= EtaDecay at each epoch - Default=1.
void SetEtaDecay(Double_t ed)
Sets EtaDecay - Eta *= EtaDecay at each epoch (look at the constructor for the complete description o...
void AttachData()
Connects the TTree to Neurons in input and output layers.
void SetData(TTree *)
Set the data source.
void SetEventWeight(const char *)
Set the event weight.
Bool_t DumpWeights(Option_t *filename="-") const
Dumps the weights to a text file.
TString fWeight
String containing the event weight.
void SetDelta(Double_t delta)
Sets Delta - used in stochastic minimisation (look at the constructor for the complete description of...
~TMultiLayerPerceptron() override
Destructor.
Double_t GetCrossEntropy() const
Cross entropy error for a softmax output neuron, for a given event.
void SetReset(Int_t reset)
Sets number of epochs between two resets of the search direction to the steepest descent.
Bool_t fTestOwner
! internal flag whether one has to delete fTest or not
void Shuffle(Int_t *, Int_t) const
Shuffle the Int_t index[n] in input.
Double_t DerivDir(Double_t *)
scalar product between gradient and direction = derivative along direction
void MLP_Stochastic(Double_t *)
One step for the stochastic method buffer should contain the previous dw vector and will be updated.
void Draw(Option_t *option="") override
Draws the network structure.
TObjArray fSynapses
Collection of all the synapses in the network.
void MLP_Line(Double_t *, Double_t *, Double_t)
Sets the weights to a point along a line Weights are set to [origin + (dist * dir)].
TNeuron::ENeuronType fType
Type of hidden neurons.
TObjArray fLastLayer
Collection of the output neurons; subset of fNetwork.
TString fextD
String containing the derivative name.
void ComputeDEDw() const
Compute the DEDw = sum on all training events of dedw for each weight normalized by the number of eve...
Double_t GetCrossEntropyBinary() const
Cross entropy error for sigmoid output neurons, for a given event.
void DrawResult(Int_t index=0, Option_t *option="test") const
Draws the neural net output It produces an histogram with the output for the two datasets.
void SetEta(Double_t eta)
Sets Eta - used in stochastic minimisation (look at the constructor for the complete description of l...
TObjArray fFirstLayer
Collection of the input neurons; subset of fNetwork.
void GetEntry(Int_t) const
Load an entry into the network.
void SetEpsilon(Double_t eps)
Sets Epsilon - used in stochastic minimisation (look at the constructor for the complete description ...
TString fextF
String containing the function name.
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
This class describes an elementary neuron, which is the basic element for a Neural Network.
Definition TNeuron.h:25
Double_t GetWeight() const
Definition TNeuron.h:48
void SetWeight(Double_t w)
Sets the neuron weight to w.
Definition TNeuron.cxx:1143
Double_t GetDEDw() const
Definition TNeuron.h:53
Double_t GetValue() const
Computes the output using the appropriate function and all the weighted inputs, or uses the branch as...
Definition TNeuron.cxx:943
void SetDEDw(Double_t in)
Sets the derivative of the total error wrt the neuron weight.
Definition TNeuron.cxx:1163
Double_t GetDeDw() const
Computes the derivative of the error wrt the neuron weight.
Definition TNeuron.cxx:1079
TNeuron * GetInLayer(Int_t n) const
Definition TNeuron.h:37
Double_t GetError() const
Computes the error for output neurons.
Definition TNeuron.cxx:1058
TTreeFormula * UseBranch(TTree *, const char *)
Sets a formula that can be used to make the neuron an input.
Definition TNeuron.cxx:873
TSynapse * GetPre(Int_t n) const
Definition TNeuron.h:35
void ForceExternalValue(Double_t value)
Uses the branch type to force an external value.
Definition TNeuron.cxx:1120
Double_t GetTarget() const
Computes the normalized target pattern for output neurons.
Definition TNeuron.cxx:1069
const Double_t * GetNormalisation() const
Definition TNeuron.h:50
ENeuronType GetType() const
Returns the neuron type.
Definition TNeuron.cxx:862
void SetNewEvent() const
Inform the neuron that inputs of the network have changed, so that the buffered values have to be rec...
Definition TNeuron.cxx:1152
void SetNormalisation(Double_t mean, Double_t RMS)
Sets the normalization variables.
Definition TNeuron.cxx:1132
void AddInLayer(TNeuron *)
Tells a neuron which neurons form its layer (including itself).
Definition TNeuron.cxx:852
ENeuronType
Definition TNeuron.h:29
@ kLinear
Definition TNeuron.h:29
@ kSigmoid
Definition TNeuron.h:29
@ kExternal
Definition TNeuron.h:29
@ kSoftmax
Definition TNeuron.h:29
@ kGauss
Definition TNeuron.h:29
@ kOff
Definition TNeuron.h:29
@ kTanh
Definition TNeuron.h:29
Iterator of object array.
Definition TObjArray.h:117
TObject * Next() override
Return next object in array. Returns 0 when no more objects in array.
An array of TObjects.
Definition TObjArray.h:31
Int_t GetEntriesFast() const
Definition TObjArray.h:58
TIterator * MakeIterator(Bool_t dir=kIterForward) const override
Returns an array iterator.
TObject * At(Int_t idx) const override
Definition TObjArray.h:164
void AddLast(TObject *obj) override
Add object in the next empty slot in the array.
TObject * UncheckedAt(Int_t i) const
Definition TObjArray.h:84
Collectable string class.
Definition TObjString.h:28
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition TObject.cxx:293
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition TObject.cxx:1045
Random number generator class based on M.
Definition TRandom3.h:27
Regular expression class.
Definition TRegexp.h:31
Basic string class.
Definition TString.h:138
Ssiz_t Length() const
Definition TString.h:425
void ToLower()
Change string to lower-case.
Definition TString.cxx:1189
Bool_t EndsWith(const char *pat, ECaseCompare cmp=kExact) const
Return true if string ends with the specified string.
Definition TString.cxx:2250
Ssiz_t First(char c) const
Find first occurrence of a character c.
Definition TString.cxx:545
const char * Data() const
Definition TString.h:384
Ssiz_t Last(char c) const
Find last occurrence of a character c.
Definition TString.cxx:938
Int_t CountChar(Int_t c) const
Return number of times character c occurs in the string.
Definition TString.cxx:522
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:640
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition TString.h:659
This is a simple weighted bidirectional connection between two neurons.
Definition TSynapse.h:20
virtual int Load(const char *module, const char *entry="", Bool_t system=kFALSE)
Load a shared library.
Definition TSystem.cxx:1868
virtual Bool_t ProcessEvents()
Process pending events (GUI, timers, sockets).
Definition TSystem.cxx:414
Base class for several text objects.
Definition TText.h:22
The TTimeStamp encapsulates seconds and ns since EPOCH.
Definition TTimeStamp.h:45
Used to coordinate one or more TTreeFormula objects.
bool Notify() override
This method must be overridden to handle object notification (the base implementation is no-op).
virtual void Add(TTreeFormula *)
Add a new formula to the list of formulas managed The manager of the formula will be changed and the ...
virtual void Remove(TTreeFormula *)
Remove a formula from this manager.
Used to pass a selection expression to the Tree drawing routine.
T EvalInstance(Int_t i=0, const char *stringStack[]=nullptr)
Evaluate this treeformula.
A TTree represents a columnar dataset.
Definition TTree.h:89
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:5718
virtual Double_t GetWeight() const
Definition TTree.h:623
void Draw(Option_t *opt) override
Default Draw method for all objects.
Definition TTree.h:470
virtual Long64_t GetEntries() const
Definition TTree.h:502
virtual Int_t GetTreeNumber() const
Definition TTree.h:598
TVirtualPad is an abstract base class for the Pad and Canvas classes.
Definition TVirtualPad.h:51
virtual void Modified(Bool_t flag=1)=0
virtual void Update()=0
const Int_t n
Definition legend1.C:16
Bool_t IsNaN(Double_t x)
Definition TMath.h:903
Double_t Log(Double_t x)
Returns the natural logarithm of x.
Definition TMath.h:767
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition TMath.h:673
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:124
TCanvas * slash()
Definition slash.C:1
TMarker m
Definition textangle.C:8
static void output()