VariableImportance.cxx

Go to the documentation of this file.

// @(#)root/tmva $Id$
// Author: Omar Zapata and Sergei Gleyzer

#include "TMVA/VariableImportance.h"

#include "TMVA/Config.h"
#include "TMVA/DataSetInfo.h"
#include "TMVA/Envelope.h"
#include "TMVA/Factory.h"
#include "TMVA/OptionMap.h"
#include "TMVA/MethodBase.h"
#include "TMVA/MethodCategory.h"
#include "TMVA/MsgLogger.h"
#include "TMVA/Types.h"
#include "TMVA/VarTransformHandler.h"

#include "TAxis.h"
#include "TGraph.h"
#include "TCanvas.h"
#include "TH1.h"
#include "TRandom3.h"
#include "TStyle.h"
#include "TSystem.h"

#include <bitset>
#include <iostream>
#include <memory>
#include <utility>


//number of bits for bitset
#define NBITS          32

TMVA::VariableImportanceResult::VariableImportanceResult():fImportanceValues("VariableImportance"),
                                                           fImportanceHist(nullptr)
{
    
}

TMVA::VariableImportanceResult::VariableImportanceResult(const VariableImportanceResult &obj)
{
    fImportanceValues = obj.fImportanceValues;
    fImportanceHist   = obj.fImportanceHist;
}


void TMVA::VariableImportanceResult::Print() const
{
    TMVA::MsgLogger::EnableOutput();
    TMVA::gConfig().SetSilent(kFALSE);   
    
    MsgLogger fLogger("VariableImportance");
    if(fType==VIType::kShort)
    {
        fLogger<<kINFO<<"Variable Importance Results (Short)"<<Endl;
    }else if(fType==VIType::kAll)
    {
        fLogger<<kINFO<<"Variable Importance Results (All)"<<Endl;        
    }else{
        fLogger<<kINFO<<"Variable Importance Results (Random)"<<Endl;                
    }
    
    fImportanceValues.Print();
    TMVA::gConfig().SetSilent(kTRUE);   
}


TCanvas* TMVA::VariableImportanceResult::Draw(const TString name) const
{
    TCanvas *c=new TCanvas(name.Data());
    fImportanceHist->Draw("");
    fImportanceHist->GetXaxis()->SetTitle(" Variable Names ");
    fImportanceHist->GetYaxis()->SetTitle(" Importance (%) ");
    c->Draw();
    return c;
}

TMVA::VariableImportance::VariableImportance(TMVA::DataLoader *dataloader):TMVA::Envelope("VariableImportance",dataloader,nullptr),fType(VIType::kShort)
{
    fClassifier=std::unique_ptr<Factory>(new TMVA::Factory("VariableImportance","!V:!ROC:!ModelPersistence:Silent:Color:!DrawProgressBar:AnalysisType=Classification"));
}

TMVA::VariableImportance::~VariableImportance()
{
    fClassifier=nullptr;
}


void TMVA::VariableImportance::Evaluate()
{
    TString methodName    = fMethod.GetValue<TString>("MethodName");
    TString methodTitle   = fMethod.GetValue<TString>("MethodTitle");
    TString methodOptions = fMethod.GetValue<TString>("MethodOptions");
    
    //NOTE: Put the type of VI Algorithm in the results Print
    if(fType==VIType::kShort)
    {
        EvaluateImportanceShort();        
    }else if(fType==VIType::kAll)
    {
        EvaluateImportanceAll();        
    }else{
        UInt_t nbits=fDataLoader->GetDefaultDataSetInfo().GetNVariables();
        if(nbits<10)
            Log()<<kERROR<<"Running variable importance with less that 10 varibales in Random mode "<<
                           "can to produce inconsisten results"<<Endl;
        EvaluateImportanceRandom(pow(2,nbits)); 
    }
    fResults.fType = fType;
    TMVA::MsgLogger::EnableOutput();
    TMVA::gConfig().SetSilent(kFALSE);   
    Log()<<kINFO<<"Evaluation done."<<Endl;
    TMVA::gConfig().SetSilent(kTRUE);   
}

ULong_t TMVA::VariableImportance::Sum(ULong_t i)
{
    ULong_t sum=0;
    for(ULong_t n=0;n<i;n++) sum+=pow(2,n);
    return sum;
}

TH1F* TMVA::VariableImportance::GetImportance(const UInt_t nbits,std::vector<Float_t> &importances,std::vector<TString> &varNames)
{
    TH1F *vihist  = new TH1F("vihist", "", nbits, 0, nbits);
    
    gStyle->SetOptStat(000000);
    
    Float_t normalization = 0.0;
    for (UInt_t i = 0; i < nbits; i++) normalization += importances[i];
    
    Float_t roc = 0.0;
    
    gStyle->SetTitleXOffset(0.4);
    gStyle->SetTitleXOffset(1.2);
    
    
    for (UInt_t i = 1; i < nbits + 1; i++) {
        roc = 100.0 * importances[i - 1] / normalization;
        vihist->GetXaxis()->SetBinLabel(i, varNames[i - 1].Data());
        vihist->SetBinContent(i, roc);
    }
    
    vihist->LabelsOption("v >", "X");
    vihist->SetBarWidth(0.97);
    vihist->SetFillColor(TColor::GetColor("#006600"));

    vihist->GetXaxis()->SetTitle(" Variable Names ");
    vihist->GetXaxis()->SetTitleSize(0.045);
    vihist->GetXaxis()->CenterTitle();
    vihist->GetXaxis()->SetTitleOffset(1.24);
    
    vihist->GetYaxis()->SetTitle(" Importance (%)");
    vihist->GetYaxis()->SetTitleSize(0.045);
    vihist->GetYaxis()->CenterTitle();
    vihist->GetYaxis()->SetTitleOffset(1.24);
    
    vihist->GetYaxis()->SetRangeUser(-7, 50);
    vihist->SetDirectory(0);
    
    return vihist;
}

void TMVA::VariableImportance::EvaluateImportanceShort()
{    
    TString methodName    = fMethod.GetValue<TString>("MethodName");
    TString methodTitle   = fMethod.GetValue<TString>("MethodTitle");
    TString methodOptions = fMethod.GetValue<TString>("MethodOptions");
    
    uint32_t x = 0;
    uint32_t y = 0;
    //getting number of variables and variable names from loader
    const UInt_t nbits = fDataLoader->GetDefaultDataSetInfo().GetNVariables();
    std::vector<TString> varNames = fDataLoader->GetDefaultDataSetInfo().GetListOfVariables();
    
    ULong_t range = Sum(nbits);
    
    //vector to save importances
    std::vector<Float_t> importances(nbits);
    for (UInt_t i = 0; i < nbits; i++)importances[i] = 0;
    
    Float_t SROC, SSROC; //computed ROC value for every Seed and SubSeed
    
    x = range;
    
    std::bitset<NBITS>  xbitset(x);
    if (x == 0) Log()<<kFATAL<<"Error: need at least one variable."; //dataloader need at least one variable
    
    
    //creating loader for seed
    TMVA::DataLoader *seeddl = new TMVA::DataLoader(xbitset.to_string());
    
    //adding variables from seed
    for (UInt_t index = 0; index < nbits; index++){ 
        if (xbitset[index]) seeddl->AddVariable(varNames[index], 'F');
    }
    
    //Loading Dataset
    DataLoaderCopy(seeddl,fDataLoader.get());
    
    //Booking Seed
    fClassifier->BookMethod(seeddl, methodName, methodTitle, methodOptions);
    
    //Train/Test/Evaluation
    fClassifier->TrainAllMethods();
    fClassifier->TestAllMethods();
    fClassifier->EvaluateAllMethods();
    
    //getting ROC
    SROC = fClassifier->GetROCIntegral(xbitset.to_string(), methodTitle);
    
    delete seeddl;  
    fClassifier->DeleteAllMethods();
    fClassifier->fMethodsMap.clear();
        
    for (uint32_t i = 0; i < NBITS; ++i) {
        if (x & (1 << i)) {
            y = x & ~(1 << i);
            std::bitset<NBITS>  ybitset(y);
            //need at least one variable
            //NOTE: if subssed is zero then is the special case
            //that count in xbitset is 1
            Double_t ny = log(x - y) / 0.693147;
            if (y == 0) {
                importances[ny] = SROC - 0.5;
                continue;
            }
            
            //creating loader for subseed
            TMVA::DataLoader *subseeddl = new TMVA::DataLoader(ybitset.to_string());
            //adding variables from subseed
            for (UInt_t index = 0; index < nbits; index++) {
                if (ybitset[index]) subseeddl->AddVariable(varNames[index], 'F');
            }
            
            //Loading Dataset
            DataLoaderCopy(subseeddl,fDataLoader.get());
            
            //Booking SubSeed
            fClassifier->BookMethod(subseeddl, methodName, methodTitle, methodOptions);
            
            //Train/Test/Evaluation
            fClassifier->TrainAllMethods();
            fClassifier->TestAllMethods();
            fClassifier->EvaluateAllMethods();
            
            //getting ROC
            SSROC = fClassifier->GetROCIntegral(ybitset.to_string(), methodTitle);
            importances[ny] += SROC - SSROC;
            
            delete subseeddl;
            fClassifier->DeleteAllMethods();
       fClassifier->fMethodsMap.clear();
        }
    }
    Float_t normalization = 0.0;
    for (UInt_t i = 0; i < nbits; i++) normalization += importances[i];
    
    for(UInt_t i=0;i<nbits;i++){
        //adding values
        fResults.fImportanceValues[varNames[i]]=(100.0 * importances[i] / normalization);
        //adding sufix
        fResults.fImportanceValues[varNames[i]]=fResults.fImportanceValues.GetValue<TString>(varNames[i])+" % ";
    }
    fResults.fImportanceHist = std::shared_ptr<TH1F>(GetImportance(nbits,importances,varNames));
}

void TMVA::VariableImportance::EvaluateImportanceRandom(UInt_t seeds)
{
    TString methodName    = fMethod.GetValue<TString>("MethodName");
    TString methodTitle   = fMethod.GetValue<TString>("MethodTitle");
    TString methodOptions = fMethod.GetValue<TString>("MethodOptions");

    TRandom3 *rangen = new TRandom3(0);  //Random Gen.
    
    uint32_t x = 0;
    uint32_t y = 0;
    
    //getting number of variables and variable names from loader
    const UInt_t nbits = fDataLoader->GetDefaultDataSetInfo().GetNVariables();
    std::vector<TString> varNames = fDataLoader->GetDefaultDataSetInfo().GetListOfVariables();
    
    ULong_t range = pow(2, nbits);
    
    //vector to save importances
    std::vector<Float_t> importances(nbits);
    Float_t importances_norm = 0;
    
    for (UInt_t i = 0; i < nbits; i++)importances[i] = 0;
    
    Float_t SROC, SSROC; //computed ROC value for every Seed and SubSeed
    
    x = range;
        
    for (UInt_t n = 0; n < seeds; n++) {
        x = rangen -> Integer(range);
        
        std::bitset<NBITS>  xbitset(x);
        if (x == 0) continue; //dataloader need at least one variable
        
        
        //creating loader for seed
        TMVA::DataLoader *seeddl = new TMVA::DataLoader(xbitset.to_string());
        
        //adding variables from seed
        for (UInt_t index = 0; index < nbits; index++) {
            if (xbitset[index]) seeddl->AddVariable(varNames[index], 'F');
        }
        
        //Loading Dataset
        DataLoaderCopy(seeddl,fDataLoader.get());
        
        //Booking Seed
        fClassifier->BookMethod(seeddl, methodName, methodTitle, methodOptions);
        
        //Train/Test/Evaluation
        fClassifier->TrainAllMethods();
        fClassifier->TestAllMethods();
        fClassifier->EvaluateAllMethods();
        
        //getting ROC
        SROC = fClassifier->GetROCIntegral(xbitset.to_string(), methodTitle);
        
        delete seeddl;
        fClassifier->DeleteAllMethods();
   fClassifier->fMethodsMap.clear();
        
        for (uint32_t i = 0; i < 32; ++i) {
            if (x & (1 << i)) {
                y = x & ~(1 << i);
                std::bitset<NBITS>  ybitset(y);
                //need at least one variable
                //NOTE: if subssed is zero then is the special case
                //that count in xbitset is 1
                Double_t ny = log(x - y) / 0.693147;
                if (y == 0) {
                    importances[ny] = SROC - 0.5;
                    importances_norm += importances[ny];
                    continue;
                }
                
                //creating loader for subseed
                TMVA::DataLoader *subseeddl = new TMVA::DataLoader(ybitset.to_string());
                //adding variables from subseed
                for (UInt_t index = 0; index < nbits; index++) {
                    if (ybitset[index]) subseeddl->AddVariable(varNames[index], 'F');
                }
                
                //Loading Dataset
                DataLoaderCopy(subseeddl,fDataLoader.get());
                
                //Booking SubSeed
                fClassifier->BookMethod(subseeddl, methodName, methodTitle, methodOptions);
                
                //Train/Test/Evaluation
                fClassifier->TrainAllMethods();
                fClassifier->TestAllMethods();
                fClassifier->EvaluateAllMethods();
                
                //getting ROC
                SSROC = fClassifier->GetROCIntegral(ybitset.to_string(), methodTitle);
                importances[ny] += SROC - SSROC;
                
                delete subseeddl;
                fClassifier->DeleteAllMethods();
      fClassifier->fMethodsMap.clear();
            }
        }
    }
    
    Float_t normalization = 0.0;
    for (UInt_t i = 0; i < nbits; i++) normalization += importances[i];
    
    for(UInt_t i=0;i<nbits;i++){
        //adding values
        fResults.fImportanceValues[varNames[i]]=(100.0 * importances[i] / normalization);
        //adding sufix
        fResults.fImportanceValues[varNames[i]]=fResults.fImportanceValues.GetValue<TString>(varNames[i])+" % ";
    }
    fResults.fImportanceHist = std::shared_ptr<TH1F>(GetImportance(nbits,importances,varNames));
    delete rangen;

}


void TMVA::VariableImportance::EvaluateImportanceAll()
{

    TString methodName    = fMethod.GetValue<TString>("MethodName");
    TString methodTitle   = fMethod.GetValue<TString>("MethodTitle");
    TString methodOptions = fMethod.GetValue<TString>("MethodOptions");
    
    uint32_t x = 0;
    uint32_t y = 0;
    
    //getting number of variables and variable names from loader
    const UInt_t nbits = fDataLoader->GetDefaultDataSetInfo().GetNVariables();
    std::vector<TString> varNames = fDataLoader->GetDefaultDataSetInfo().GetListOfVariables();
    
    ULong_t range = pow(2, nbits);
    
    //vector to save importances
    std::vector<Float_t> importances(nbits);
    
    //vector to save ROC-Integral values
    std::vector<Float_t> ROC(range);
    ROC[0]=0.5;
    for (UInt_t i = 0; i < nbits; i++) importances[i] = 0;
    
    Float_t SROC, SSROC; //computed ROC value
    for ( x = 1; x <range ; x++) {
        
        std::bitset<NBITS>  xbitset(x);
        if (x == 0) continue; //dataloader need at least one variable
        
        //creating loader for seed
        TMVA::DataLoader *seeddl = new TMVA::DataLoader(xbitset.to_string());
        
        //adding variables from seed
        for (UInt_t index = 0; index < nbits; index++) {
            if (xbitset[index]) seeddl->AddVariable(varNames[index], 'F');
        }
        
        DataLoaderCopy(seeddl,fDataLoader.get());
        
        seeddl->PrepareTrainingAndTestTree(fDataLoader->GetDefaultDataSetInfo().GetCut("Signal"), fDataLoader->GetDefaultDataSetInfo().GetCut("Background"), fDataLoader->GetDefaultDataSetInfo().GetSplitOptions());
        
        //Booking Seed
        fClassifier->BookMethod(seeddl, methodName, methodTitle, methodOptions);
        
        //Train/Test/Evaluation
        fClassifier->TrainAllMethods();
        fClassifier->TestAllMethods();
        fClassifier->EvaluateAllMethods();
        
        //getting ROC
        ROC[x] = fClassifier->GetROCIntegral(xbitset.to_string(), methodTitle);
        
        delete seeddl;    
        fClassifier->DeleteAllMethods();
   fClassifier->fMethodsMap.clear();
    }
    
    
    for ( x = 0; x <range ; x++) 
    {
        SROC=ROC[x];
        for (uint32_t i = 0; i < NBITS; ++i) {
            if (x & (1 << i)) {
                y = x & ~(1 << i);
                std::bitset<NBITS>  ybitset(y);
                
                Float_t ny = log(x - y) / 0.693147;
                if (y == 0) {
                    importances[ny] = SROC - 0.5;
                    continue;
                }
                
                //getting ROC
                SSROC = ROC[y];
                importances[ny] += SROC - SSROC;
            }
            
        }
    }
    Float_t normalization = 0.0;
    for (UInt_t i = 0; i < nbits; i++) normalization += importances[i];
    
    for(UInt_t i=0;i<nbits;i++){
        //adding values
        fResults.fImportanceValues[varNames[i]]=(100.0 * importances[i] / normalization);
        //adding sufix
        fResults.fImportanceValues[varNames[i]]=fResults.fImportanceValues.GetValue<TString>(varNames[i])+" % ";
    }
    fResults.fImportanceHist = std::shared_ptr<TH1F>(GetImportance(nbits,importances,varNames));
}