Documentation / Manuel utilisateur en C++ :
The example script below uses the THSIC
class to compute and display the indices. Note that HSIC indexes do not require a specific design of experiments, they can be computed on a given sample.
void sensitivityHSICFunctionFlowrate(){
gROOT->LoadMacro("UserFunctions.C");
// Define the DataServer
TDataServer *tds = new TDataServer("tdsflowreate", "DataBase flowreate");
tds->addAttribute( new TUniformDistribution("rw", 0.05, 0.15));
tds->addAttribute( new TUniformDistribution("r", 100.0, 50000.0));
tds->addAttribute( new TUniformDistribution("tu", 63070.0, 115600.0));
tds->addAttribute( new TUniformDistribution("tl", 63.1, 116.0));
tds->addAttribute( new TUniformDistribution("hu", 990.0, 1110.0));
tds->addAttribute( new TUniformDistribution("hl", 700.0, 820.0));
tds->addAttribute( new TUniformDistribution("l", 1120.0, 1680.0));
tds->addAttribute( new TUniformDistribution("kw", 9855.0, 12045.0));
// Generation of the sample (it can be a given sample).
Int_t nS = 500;
TSampling *sampling = new TSampling(tds, "lhs", nS);
sampling->generateSample();
TLauncherFunction * tlf = new TLauncherFunction(tds, "flowrateModel");
tlf->setDrawProgressBar(kFALSE);
tlf->run();
// Create a THSIC object, compute indexes and print results
THSIC * thsic = new THSIC(tds, "rw:r:tu:tl:hu:hl:l:kw","flowrateModel");
thsic->computeIndexes("quiet");
thsic->getResultTuple()->SetScanField(60);
thsic->getResultTuple()->Scan("Out:Inp:Method:Order:Value:CILower:CIUpper");
// Print HSIC indexes
TCanvas *can = new TCanvas("c1", "Graph for the Macro sensitivityHSICFunctionFlowrate",5,64,1270,667);
thsic->drawIndexes("Flowrate", "", "hist,first,nonewcanv");
}
There is one kind of constructor to build a THSIC
object based on provided data. The HSIC compute indices on given data, the TDataServer
must be filled. The constructor does not require a model to run (this has to be done).
// Create a THSIC object data
THSIC(TDataServer *tds, const char *inp, const char *out , const char *Option="")
This constructor takes up to four arguments:
- a pointer to a
TDataServer
object, - a string to specify the names of the input factors separated by ':' (ex. "rw:r:tu:tl:hu:hl:l:kw"), no default is accepted,
- a string to specify the names of the output variables of the model, no default is accepted.
- a string to specify the options, by default empty "". Possibility are "empiri" or "median" to specify the method to estimate the variance (see Section VI.8.1.4).
Here is an example of how to use the constructor with provided data:
THSIC * thsic = new THSIC(tds, "rw:r:tu:tl:hu:hl:l:kw", "flowrateModel");
To compute the indices, run the method computeIndexes
:
thsic->computeIndexes();
Note this method computes the 'V-stat' indices. If the "unbiased" option is provided, it computes the 'U-stat' indices. It computes three values for each variable:
- the HSIC measures,
- the R2HSIC indices, which are considered as the sensitivity indices
- the p-value of the independence test
setEstimatedVariance
: Set the method used to estimate the variance of the gaussian kernel ( kUnknown
| kStdEmpirical
| kMedianeDelta
)
kStdEmpirical
: the parameters is estimated by the emperical standatd deviation of the attributeskMedianeDelta
: use the median distance between points
thsic->setEstimatedVariance(URANIE::Sensitivity::THSIC::kStdEmpirical);
Two technics are available to compute the p-values to test the independence between inputs and outputs:
- By default, the gamma approximation is used.
- The following command allows using the bootstrap
thsic->computeIndexes("nperm=100"); // specify the number of permutation
computeIndexes("nperm=200")
: Specify the number of bootstrap permutation using the option "nperm=" incomputeIndexes
The measures, once computed, are stored in a TTree
. To get this
TTree
, use the method TSensitivity::getResultTuple()
:
TTree * resultsHSIC = thsic->getResultTuple();
Several methods exist in ROOT to extract data from a TTree
, it is advised to look for them
into the ROOT documentation. We propose two ways of extracting the value of each coefficient from the
TTree
.
The method use the method getValue
of the THSIC
object specifying the
order of the extract value ("R2HSic", "HSic" or "pValues"), the related input and possibly more selected options.
double hl_r2hsic_index = thsic->getValue("R2HSic","hl");
double hl_hsic_index = thsic->getValue("HSic","hl");
double hl_pvalues_index = thsic->getValue("pValues","hl");