<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Import
Step1:  Reading initial data
Step2:  Define label
label = signB * signVtx > 0
* same sign of B and vtx -> label = 1
* opposite sign of B and vtx -> label = 0
Step3:  Define B-like events for training and others for prediction
Step4:  Define features
Step5:  Find good vtx to define sign B
trying to guess sign of B based on sign of vtx. If the guess is good, the vtx will be used on next step to train classifier.
2-folding random forest selection for right tagged events
Step7:  Training on all vtx
in this case we don't use preselection with RandomForest
DT full
Step8:  Report for all vtx
Step9:  Calibrating results $p(\text{vrt same sign}|B)$ and combining them
Step10:  Comparison of different models
Step11:  Implementing the best vertex model
and saving its predictions
<ASSISTANT_TASK:>
Python Code:
import pandas
import numpy
from sklearn.ensemble import RandomForestClassifier
from sklearn.metrics import roc_curve, roc_auc_score
from rep.metaml import FoldingClassifier
from rep.data import LabeledDataStorage
from rep.report import ClassificationReport
from rep.report.metrics import RocAuc
from utils import get_N_B_events, get_events_number, get_events_statistics
Explanation: Import
End of explanation
import root_numpy
data = pandas.DataFrame(root_numpy.root2array('datasets/1016_vtx.root'))
Explanation: Reading initial data
End of explanation
event_id_column = 'event_id'
data[event_id_column] = data.runNum.apply(str) + '_' + (data.evtNum.apply(int)).apply(str)
# reconstructing sign of B
data['signB'] = data.tagAnswer * (2 * data.iscorrect - 1)
# assure sign is +1 or -1
data['signVtx'] = (data.signVtx.values > 0) * 2 - 1
data['label'] = (data.signVtx.values * data.signB.values > 0) * 1
data.head()
get_events_statistics(data)
N_pass = get_events_number(data)
tagging_efficiency = 1. * N_pass / get_N_B_events()
tagging_efficiency_delta = sqrt(N_pass) / get_N_B_events()
print tagging_efficiency, tagging_efficiency_delta
Bdata_tracks = pandas.read_csv('models/Bdata_tracks_PID_less.csv')
Bdata_tracks.index = Bdata_tracks.event_id
data['initial_pred'] = Bdata_tracks.ix[data.event_id, 'track_relation_prob'].values
Explanation: Define label
label = signB * signVtx > 0
* same sign of B and vtx -> label = 1
* opposite sign of B and vtx -> label = 0
End of explanation
sweight_threshold = 1.
data_sw_passed = data[data.N_sig_sw > sweight_threshold]
data_sw_not_passed = data[data.N_sig_sw <= sweight_threshold]
get_events_statistics(data_sw_passed)
Explanation: Define B-like events for training and others for prediction
End of explanation
features =  ['mult', 'nnkrec', 'ptB', 'vflag', 'ipsmean', 'ptmean', 'vcharge', 
             'svm', 'svp', 'BDphiDir', 'svtau', 'docamax']
Explanation: Define features
End of explanation
data_sw_passed_lds = LabeledDataStorage(data_sw_passed, data_sw_passed.label, data_sw_passed.N_sig_sw)
Explanation: Find good vtx to define sign B
trying to guess sign of B based on sign of vtx. If the guess is good, the vtx will be used on next step to train classifier.
2-folding random forest selection for right tagged events
End of explanation
from hep_ml.decisiontrain import DecisionTrainClassifier
from hep_ml.losses import LogLossFunction
from hep_ml.losses import HessianLossFunction
from hep_ml.commonutils import check_sample_weight
from scipy.special import expit
class LogLossFunctionTagging(HessianLossFunction):
    Logistic loss function (logloss), aka binomial deviance, aka cross-entropy,
       aka log-likelihood loss.
      
    def fit(self, X, y, sample_weight):
        self.sample_weight = check_sample_weight(y, sample_weight=sample_weight,
                                                 normalize=True, normalize_by_class=True)
        self.initial_pred = numpy.log(X['initial_pred'].values)
        self.y_signed = 2 * y - 1
        self.minus_y_signed = - self.y_signed
        self.y_signed_times_weights = self.y_signed * self.sample_weight
        HessianLossFunction.fit(self, X, y, sample_weight=self.sample_weight)
        return self
    def __call__(self, y_pred):
        y_pred = y_pred + self.initial_pred
        return numpy.sum(self.sample_weight * numpy.logaddexp(0, self.minus_y_signed * y_pred))
    def negative_gradient(self, y_pred):
        y_pred = y_pred + self.initial_pred
        return self.y_signed_times_weights * expit(self.minus_y_signed * y_pred)
    def hessian(self, y_pred):
        y_pred = y_pred + self.initial_pred
        expits = expit(y_pred)
        return self.sample_weight * expits * (1 - expits)
    def prepare_tree_params(self, y_pred):
        y_pred = y_pred + self.initial_pred
        return self.y_signed * expit(self.minus_y_signed * y_pred), self.sample_weight
tt_base = DecisionTrainClassifier(learning_rate=0.02, n_estimators=1500, depth=6, 
                                  max_features=8, loss=LogLossFunctionTagging(regularization=100), train_features=features)
tt_folding = FoldingClassifier(tt_base, n_folds=2, random_state=11,
                               features=features + ['initial_pred'])
%time tt_folding.fit_lds(data_sw_passed_lds)
pass
from scipy.special import expit, logit
data_temp = data_sw_not_passed
print roc_auc_score(data_temp.signB.values, data_temp['initial_pred'].values, sample_weight=data_temp.N_sig_sw.values)
p = tt_folding.predict_proba(data_temp)[:, 1]
print roc_auc_score(data_temp.signB.values,
                    log(data_temp['initial_pred'].values) + logit(p) * data_temp.signB.values,
                    sample_weight=data_temp.N_sig_sw.values)
hist(tt_folding.estimators[0].loss.initial_pred, )
pass
Explanation: Training on all vtx
in this case we don't use preselection with RandomForest
DT full
End of explanation
report = ClassificationReport({'tt': tt_folding}, data_sw_passed_lds)
report.learning_curve(RocAuc())
report.compute_metric(RocAuc())
report.roc()
Explanation: Report for all vtx
End of explanation
models = []
from utils import get_result_with_bootstrap_for_given_part
models.append(get_result_with_bootstrap_for_given_part(tagging_efficiency, tagging_efficiency_delta, tt_folding, 
                                                      [data_sw_passed, data_sw_not_passed], 
                                                      logistic=True, name="tt-log",
                                                      sign_part_column='signVtx', part_name='vertex'))
models.append(get_result_with_bootstrap_for_given_part(tagging_efficiency, tagging_efficiency_delta, tt_folding, 
                                                      [data_sw_passed, data_sw_not_passed], 
                                                      logistic=False, name="tt-iso",
                                                      sign_part_column='signVtx', part_name='vertex'))
Explanation: Calibrating results $p(\text{vrt same sign}|B)$ and combining them
End of explanation
pandas.concat(models)
pandas.concat(models)
Explanation: Comparison of different models
End of explanation
from utils import prepare_B_data_for_given_part
Bdata_prepared = prepare_B_data_for_given_part(tt_folding, [data_sw_passed, data_sw_not_passed], logistic=False,
                                               sign_part_column='signVtx', part_name='vertex')
Bdata_prepared.to_csv('models/Bdata_vertex_new_loss.csv', header=True, index=False)
Explanation: Implementing the best vertex model
and saving its predictions
End of explanation
<END_TASK>