{"cells":[{"metadata":{"_uuid":"a47869a8604420767b21511e0969381cb54f8f27"},"cell_type":"markdown","source":"# 1 . Read file"},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"collapsed":true},"cell_type":"code","source":"import pandas as pd, xgboost\n\n# Read the training dataset:\ntrain = pd.read_csv('../input/training.csv', index_col='id')\ntrain = train[train['min_ANNmuon'] > 0.4]","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","collapsed":true,"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":false},"cell_type":"markdown","source":"# 2. Feature Engineering"},{"metadata":{"trusted":true,"_uuid":"94b6d479e2011d74d3827d899ebe338f103e8927","collapsed":true},"cell_type":"code","source":"import numpy as np\n\n# Physical constants:\nc = 299.792458     # Speed of light\nm_mu = 105.6583715 # Muon mass (in MeV)\nm_tau = 1776.82    # Tau mass (in MeV)\n\n# List of the features for the first booster:\nlist1 = [\n# Original features:\n         'FlightDistance',\n         'FlightDistanceError',\n         'LifeTime',\n         'IP',\n         'IPSig',\n         'VertexChi2',\n         'dira',\n         'pt',\n         'DOCAone',\n         'DOCAtwo',\n         'DOCAthree',\n         'IP_p0p2',\n         'IP_p1p2',\n         'isolationa',\n         'isolationb',\n         'isolationc',\n         'isolationd',\n         'isolatione',\n         'isolationf',\n         'iso',\n         'CDF1',\n         'CDF2',\n         'CDF3',\n         'ISO_SumBDT',\n         'p0_IsoBDT',\n         'p1_IsoBDT',\n         'p2_IsoBDT',\n         'p0_track_Chi2Dof',\n         'p1_track_Chi2Dof',\n         'p2_track_Chi2Dof',\n         'p0_IP',\n         'p0_IPSig',\n         'p1_IP',\n         'p1_IPSig',\n         'p2_IP',\n         'p2_IPSig',\n# Extra features:\n         'E',\n         'FlightDistanceSig',\n         'DOCA_sum',\n         'isolation_sum',\n         'IsoBDT_sum',\n         'track_Chi2Dof',\n         'IP_sum',\n         'IPSig_sum',\n         'CDF_sum'\n        ]\n\n# List of the features for the second booster:\nlist2 = [\n# Original features:\n         'dira',\n         'pt',\n         'p0_pt',\n         'p0_p',\n         'p0_eta',\n         'p1_pt',\n         'p1_p',\n         'p1_eta',\n         'p2_pt',\n         'p2_p',\n         'p2_eta',\n# Extra features:\n         'E',\n         'pz',\n         'beta',\n         'gamma',\n         'beta_gamma',\n         'Delta_E',\n         'Delta_M',\n         'flag_M',\n         'E0',\n         'E1',\n         'E2',\n         'E0_ratio',\n         'E1_ratio',\n         'E2_ratio',\n         'p0_pt_ratio',\n         'p1_pt_ratio',\n         'p2_pt_ratio',\n         'eta_01',\n         'eta_02',\n         'eta_12',\n         't_coll'\n         ]\n\n# Function to add extra features:\ndef add_features(df):\n    # Number of events:\n    N = len(df)\n\n    # Internal arrays:\n    p012_p = np.zeros(3)\n    p012_pt = np.zeros(3)\n    p012_z = np.zeros(3)\n    p012_eta = np.zeros(3)\n    p012_IsoBDT = np.zeros(3)\n    p012_track_Chi2Dof = np.zeros(3)\n    p012_IP = np.zeros(3)\n    p012_IPSig = np.zeros(3)\n    CDF123 = np.zeros(3)\n    isolation = np.zeros(6)\n\n    # Kinematic features related to the mother particle:\n    E = np.zeros(N)\n    pz = np.zeros(N)\n    beta = np.zeros(N)\n    gamma = np.zeros(N)\n    beta_gamma = np.zeros(N)\n    M_lt = np.zeros(N)\n    M_inv = np.zeros(N)\n    Delta_E = np.zeros(N)\n    Delta_M = np.zeros(N)\n    flag_M = np.zeros(N)\n\n    # Kinematic features related to the final-state particles p0, p1, and p2:\n    E012 = np.zeros((N,3))\n    E012_ratio = np.zeros((N,3))\n    p012_pt_ratio = np.zeros((N,3))\n    eta_01 = np.zeros(N)\n    eta_02 = np.zeros(N)\n    eta_12 = np.zeros(N)\n    t_coll = np.zeros(N)\n\n    # Other extra features:\n    FlightDistanceSig = np.zeros(N)\n    DOCA_sum = np.zeros(N)\n    isolation_sum = np.zeros(N)\n    IsoBDT_sum = np.zeros(N)\n    track_Chi2Dof = np.zeros(N)\n    IP_sum = np.zeros(N)\n    IPSig_sum = np.zeros(N)\n    CDF_sum = np.zeros(N)\n\n    for i in range(N):\n        # Read some of the original features:  \n        pt = df['pt'].values[i]\n        dira = df['dira'].values[i]\n        LifeTime = df['LifeTime'].values[i]\n        FlightDistance = df['FlightDistance'].values[i]\n        FlightDistanceError = df['FlightDistanceError'].values[i]\n        DOCAone = df['DOCAone'].values[i]\n        DOCAtwo = df['DOCAtwo'].values[i]\n        DOCAthree = df['DOCAthree'].values[i]\n        isolation[0] = df['isolationa'].values[i]\n        isolation[1] = df['isolationb'].values[i]\n        isolation[2] = df['isolationc'].values[i]\n        isolation[3] = df['isolationd'].values[i]\n        isolation[4] = df['isolatione'].values[i]\n        isolation[5] = df['isolationf'].values[i]\n\n    for j in range(3):\n        p012_p[j] = df['p'+str(j)+'_p'].values[i]\n        p012_pt[j] = df['p'+str(j)+'_pt'].values[i]\n        p012_eta[j] = df['p'+str(j)+'_eta'].values[i]\n        p012_IsoBDT[j] = df['p'+str(j)+'_IsoBDT'].values[i]\n        p012_track_Chi2Dof[j] = df['p'+str(j)+'_track_Chi2Dof'].values[i]\n        p012_IP[j] = df['p'+str(j)+'_IP'].values[i]\n        p012_IPSig[j] = df['p'+str(j)+'_IPSig'].values[i]\n        CDF123[j] = df['CDF'+str(j+1)].values[i]\n\n    # Differences between pseudorapidities of the final-state particles:\n    eta_01[i] = p012_eta[0] - p012_eta[1]\n    eta_02[i] = p012_eta[0] - p012_eta[2]\n    eta_12[i] = p012_eta[1] - p012_eta[2]\n\n    # Transverse collinearity of the final-state particles (equals to 1 if they are collinear):\n    t_coll[i] = sum(p012_pt[:])/pt\n\n    # Longitudinal momenta of the final-state particles:\n    p012_z[:] = p012_pt[:]*np.sinh(p012_eta[:])\n\n    # Energies of the final-state particles:\n    E012[i,:] = np.sqrt(np.square(m_mu) + np.square(p012_p[:]))\n\n    # Energy and momenta of the mother particle:\n    E[i] = sum(E012[i,:])\n    pz[i] = sum(p012_z[:])\n    p = np.sqrt(np.square(pt) + np.square(pz[i]))\n\n    # Energies and momenta of the final-state particles relative to those of the mother particle:\n    E012_ratio[i,:] = E012[i,:]/E[i]\n    p012_pt_ratio[i,:] = p012_pt[:]/pt\n\n    # Mass of the mother particle calculated from FlightDistance and LifeTime:\n    beta_gamma[i] = FlightDistance/(LifeTime*c)\n    M_lt[i] = p/beta_gamma[i]\n\n    # If M_lt is around the tau mass then flag_M = 1 (otherwise 0):\n    if np.fabs(M_lt[i] - m_tau - 1.44) < 17: flag_M[i] = 1\n\n    # Invariant mass of the mother particle calculated from its energy and momentum:\n    M_inv[i] = np.sqrt(np.square(E[i]) - np.square(p))\n\n    # Relativistic gamma and beta of the mother particle:\n    gamma[i] = E[i]/M_inv[i]\n    beta[i] = np.sqrt(np.square(gamma[i]) - 1.)/gamma[i]\n\n    # Difference between M_lt and M_inv:\n    Delta_M[i] = M_lt[i] - M_inv[i]\n\n    # Difference between energies of the mother particle calculated in two different ways:\n    Delta_E[i] = np.sqrt(np.square(M_lt[i]) + np.square(p)) - E[i]\n\n    # Other extra features:\n    FlightDistanceSig[i] = FlightDistance/FlightDistanceError\n    DOCA_sum[i] = DOCAone + DOCAtwo + DOCAthree\n    isolation_sum[i] = sum(isolation[:])\n    IsoBDT_sum[i] = sum(p012_IsoBDT[:])\n    track_Chi2Dof[i] = np.sqrt(sum(np.square(p012_track_Chi2Dof[:] - 1.)))\n    IP_sum[i] = sum(p012_IP[:])\n    IPSig_sum[i] = sum(p012_IPSig[:])\n    CDF_sum[i] = sum(CDF123[:])\n\n    # Kinematic features related to the mother particle:\n    df['E'] = E\n    df['pz'] = pz\n    df['beta'] = beta\n    df['gamma'] = gamma\n    df['beta_gamma'] = beta_gamma\n    df['M_lt'] = M_lt\n    df['M_inv'] = M_inv\n    df['Delta_E'] = Delta_E\n    df['Delta_M'] = Delta_M\n    df['flag_M'] = flag_M\n\n    # Kinematic features related to the final-state particles:\n    df['E0'] = E012[:,0]\n    df['E1'] = E012[:,1]\n    df['E2'] = E012[:,2]\n    df['E0_ratio'] = E012_ratio[:,0]\n    df['E1_ratio'] = E012_ratio[:,1]\n    df['E2_ratio'] = E012_ratio[:,2]\n    df['p0_pt_ratio'] = p012_pt_ratio[:,0]\n    df['p1_pt_ratio'] = p012_pt_ratio[:,1]\n    df['p2_pt_ratio'] = p012_pt_ratio[:,2]\n    df['eta_01'] = eta_01\n    df['eta_02'] = eta_02\n    df['eta_12'] = eta_12\n    df['t_coll'] = t_coll\n\n    # Other extra features:\n    df['FlightDistanceSig'] = FlightDistanceSig\n    df['DOCA_sum'] = DOCA_sum\n    df['isolation_sum'] = isolation_sum\n    df['IsoBDT_sum'] = IsoBDT_sum\n    df['track_Chi2Dof'] = track_Chi2Dof\n    df['IP_sum'] = IP_sum\n    df['IPSig_sum'] = IPSig_sum\n    df['CDF_sum'] = CDF_sum\n\n    return df","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"53573a4ca7890241fa3459a58ce3be6bc8bdd6b5"},"cell_type":"markdown","source":"# 3.Model Training"},{"metadata":{"trusted":true,"_uuid":"45eeae978244819d0f3689989c09634ca29f1c65","collapsed":true},"cell_type":"code","source":"# Add extra features:\ntrain = add_features(train)\n\n# Random seed:\nrandom_state = 42\n\n# Weight for the first classifier:\nw1 = 0.82\n\n# Numbers of trees:\nnum_trees1 = 200 # Classifier 1\nnum_trees2 = 100 # Classifier 2\n\n# Parameters of the classifiers:\nparams = {'objective': 'binary:logistic',\n          'eta': 0.05,\n          'max_depth': 4,\n          'scale_pos_weight': 5.,\n          'silent': 1,\n          'seed': random_state}\n\n# Train the first (geometric) XGBoost classifier:\nbst1 = xgboost.train(params,\n                     xgboost.DMatrix(train[list1],\n                     train['signal']), num_trees1)\nbst1.save_model('bst1.model')\n\n# Train the second (kinematic) XGBoost classifier:\nbst2 = xgboost.train(params,\n                     xgboost.DMatrix(train[list2],\n                     train['signal']), num_trees2)\nbst2.save_model('bst2.model')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"06564f8fc594d3ff75469cb5c8811cb079af29a4"},"cell_type":"markdown","source":"# 4. Predict the Data"},{"metadata":{"trusted":true,"_uuid":"3ab041062065287b1435d3845ce9bdcf503be04a","collapsed":true},"cell_type":"code","source":"import csv\n# Load the first booster:\nbst1 = xgboost.Booster()\nbst1.load_model(\"bst1.model\")\n\n# Load the second booster:\nbst2 = xgboost.Booster()\nbst2.load_model(\"bst2.model\")\n\n# Create a submission file:\nwith open('submission.csv', 'w') as csvfile:\n    csv.writer(csvfile, delimiter=',').writerow(['id', 'prediction'])\n\n# Prediction and output:\nfor chunk in pd.read_csv(\"../input/test.csv\", index_col='id', chunksize=100000):\n    # Add extra features:\n    chunk = add_features(chunk)\n\n    # Predict probabilities:\n    probs1 = bst1.predict(xgboost.DMatrix(chunk[list1])) # Booster 1\n    probs2 = bst2.predict(xgboost.DMatrix(chunk[list2])) # Booster 2\n\n    # Weighted average of the predictions:\n    result = pd.DataFrame({'id': chunk.index})\n    result['prediction'] = 0.5*(w1*probs1 + (1 - w1)*probs2)\n\n    # Write to the submission file:\n    result.to_csv('submission.csv', index=False, header=False, sep=',', mode='a')","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.4","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}