{"cells":[{"metadata":{"_uuid":"4b77dbff244e1109e1600f5c9b98c4fdf530e855"},"cell_type":"markdown","source":"Full disclosure: I am actually a particle physicist, although I've never worked much with collider experiments.\n\nIn this challenge, we are tasked with predicting whether or not a given reconstructed event in a particle physics experiment is from a \n\n$$ \\tau \\rightarrow \\mu^+\\mu^-\\mu^- $$\n\ndecay.\n\nTau leptons are basically like very heavy electrons, with a mass of 1.78 GeV rather than 511 keV. These are still quite light by the standards of proton colliders. The top quark, which is the heaviest known elementary particle, has a mass about 100 times higher. The LHC was built so that it would be able to create TeV-scale particles such as supersymmetric partners to regular Standard Model particles.\n\nTaus are unstable and decay with a lifetime of $3\\times10^{-13}$ s, or 300 femtoseconds. Due to conservation of lepton number, taus are created in $\\tau^+\\tau^-$ pairs or in conjunction with neutrinos. In a facility like the LHC, this means that the taus will typically have a large amount of kinetic energy in the lab frame. The apparent lifetime in the lab frame is boosted from $\\tau = 300$~fs to $\\gamma\\tau$, where $\\gamma$ is a relativistic factor equal to the ratio between the total energy in the lab frame and the rest mass. So, if $\\gamma=10$ (that is, the $\\tau$ has a total energy 17.8 GeV), it will have a lifetime of 3 picoseconds. Since the tau will be going nearly the speed of light, it will travel on average around 1 mm from the initial vertex before decaying.\n\nIf you can reconstruct the separation between the tau decay and the vertex (at least in the transverse dimension), this will let you do things like compare the estimated $\\gamma$ of the tau (based on muon kinematics) to the distance traveled. This will give you some ability to at least start to select tau-like events."},{"metadata":{"collapsed":true,"trusted":false,"_uuid":"ede39acd6a3b2d53eed982e28e3f6a208f975275"},"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport numpy as np\nimport pandas as pd\n%matplotlib inline","execution_count":null,"outputs":[]},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"98186e3402ca9d4a89e1022b683d07a3a4d3c892"},"cell_type":"code","source":"folder = '../input/'\ntrain = pd.read_csv(folder+'training.csv', index_col=0)\ntest = pd.read_csv(folder+'test.csv',index_col=0)\n#check_correlation = pd.read_csv(folder+'check_correlation.csv', index_col='id')\n#check_agreement = pd.read_csv(folder+'check_agreement.csv', index_col='id')\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"f1bf358a608c590933068a1d46eabf4927499d80"},"cell_type":"code","source":"train.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"39ac7f87529d000cda6c0a3dd1f63f60315abc0a"},"cell_type":"code","source":"train.info()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"fbbcc401247d0303aa6e2402e1070a41c6ceb5d9"},"cell_type":"markdown","source":"We see that there are no null values, which is good. This is all truth or reconstructed information, so we do not expect any null values anywhere.\n\nThere are quite a few different features here. One of the main challenges of this sort of dataset is figuring out how to select or construct features that do not show noticeable differences between simulated data and real data.\n\nThese sorts of tests are very common in particle physics. Most high energy physics experiments rely heavily on detailed simulations of both the physics they are trying to measure and also the interactions of particles in their detectors. The models going into these simulations always have various uncertainties and systematic errors, so you  need to be very careful that you are not training on the idiosyncracies of your simulation but rather on features that are fairly insensitive to data/simulation differences."},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"4c3093db2d8a2f98a369daaf7ec4889d1b5e5ba6"},"cell_type":"code","source":"train.signal.value_counts()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"a7d5df4730b70718a1b1748ea8d8ea679eae14e6"},"cell_type":"markdown","source":"The data aren't perfectly balanced, but this looks like we are not in the regime of finding rare events over a large background. Either there are large weights that need to be applied somewhere to renormalize things or we are seeing a set of events that are preselected to have a decent signal to background ratio."},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"805dfa122c25048c0724d2209f44fa53b9af12d8"},"cell_type":"code","source":"signal = train[train.signal==1]\nnonsignal = train[train.signal==0]","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"01bd92ca2bcf35ecb72f4e75c267e62fb5459003"},"cell_type":"markdown","source":"Now let's plot a few of the variables.\n\nThe transverse momentum ($p_t$) of the 3-muon system would be the magnitude of the vector sum of the transverse momenta of the muons. This shows how much momentum is going away from axis defining the beam direction. This is important in many collider analyses because when protons collide at high energies, we are really seeing subatomic particles withing the protons (i.e. quarks and gluons). These particles, called partons, carry only a fraction of the total proton momentum. As a result, there can be a significant overall boost in the longitudinal direction. There is going to be little transverse momentum in the collision, so high $p_t$ particles somehow got their transverse momentum from the collision. If you could measure all the particles from the collision, the total $p_t$ would be very close to 0. In practice, there can be missing $p_t$ from things like neutrinos, fragments of the protons that are too close to the beam to measure, and even exotic theoretical particles. "},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"bd541a705137f47d11275b07558970f4fc4c9a12"},"cell_type":"code","source":"plt.hist(signal.pt/1000, range=(0,25), bins=100, label='Signal', alpha=0.7,normed=True)\nplt.hist(nonsignal.pt/1000, range=(0,25), bins=100, label='Background',alpha=0.7,normed=True)\nplt.xlabel('Transverse Momentum [GeV?]')\nplt.ylabel('Normalized Fraction')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"7fd8417b820b983274f70a735ab57a2e5ff09e3a"},"cell_type":"markdown","source":"I already explained what the flight distance represents. This is the distance between the initial vertex and the location where the proposed $\\tau$ lepton decayed."},{"metadata":{"scrolled":true,"trusted":false,"collapsed":true,"_uuid":"70b5269560225f1af0234db8fd28542a93f13b8b"},"cell_type":"code","source":"plt.hist(signal.FlightDistance, range=(0,100), bins=100, label='Signal', alpha=0.7,normed=True)\nplt.hist(nonsignal.FlightDistance, range=(0,100), bins=100, label='Background',alpha=0.7,normed=True)\nplt.xlabel('Flight Distance')\nplt.ylabel('Normalized Fraction')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"c55f4d67a012fd309b577ba77bdda15f74b924be"},"cell_type":"markdown","source":"Dira gives the angle between the collision vertex/decay vertex vector and the direction of the 3-muon system. In $\\tau\\rightarrow 3\\mu$ decays, all of the $\\tau$'s momentum goes into the muons, so the momentum vector does not change. Too large an angle and we know that there is something not quite right with our proposed event type."},{"metadata":{"scrolled":true,"trusted":false,"collapsed":true,"_uuid":"ad90dbf22a5926aac8ca5e86aacea74ed9c221c0"},"cell_type":"code","source":"plt.hist(signal.dira, range=(0.998,1), bins=100, label='Signal', alpha=0.7,normed=True)\nplt.hist(nonsignal.dira,range=(0.999,1),bins=100, label='Background',alpha=0.7,normed=True)\nplt.xlabel('Dira')\nplt.ylabel('Normalized Fraction')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"53fcba9bef2ae68bef08be2efd1b9a9336bdc16f"},"cell_type":"markdown","source":"The vertex $\\chi^2$ is a measure of the goodness of fit of the $\\tau$ vertex. I believe this means the goodness of fit of the three muons converging on the same location. This will typically be low for true $\\tau\\rightarrow3\\mu$ events and high for events where the 3 muons are not really from the same decay.\n\nHere we see a sudden cutoff in the histograms around $\\chi^2=15$. I did not put that there. It's in the data, so this looks like one of the variables used in the preselection process used to collect the data.\n\nIt's quite likely that this is not going to be a good feature to use in this challenge."},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"d8407825560a846e3b04818429e887639003fbed"},"cell_type":"code","source":"plt.hist(signal.VertexChi2, range=(0,16), bins=100, label='Signal', alpha=0.7,normed=True)\nplt.hist(nonsignal.VertexChi2,range=(0,16),bins=100, label='Background',alpha=0.7,normed=True)\nplt.xlabel(r'Vertex $\\chi^2$')\nplt.ylabel('Normalized Fraction')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d6750ba701b801bd865a023bb846d2a765400995"},"cell_type":"markdown","source":"I'm not 100% sure of what the impact parameter is but my guess is that this is the estimated closest approach of the tau candidate to the vertex position. That is, if we take the momentum vector of the three muons as the tau momentum and then propagate the tau backward in time, this is the closest it gets to the vertex."},{"metadata":{"scrolled":true,"trusted":false,"collapsed":true,"_uuid":"eb0087ac0c062a4d5a263057eba1a9b052858ec2"},"cell_type":"code","source":"plt.hist(signal.IP, range=(0,0.5), bins=100, label='Signal', alpha=0.7,normed=True)\nplt.hist(nonsignal.IP,range=(0,0.5),bins=100, label='Background',alpha=0.7,normed=True)\nplt.xlabel('Impact Paraneter')\nplt.ylabel('Normalized Fraction')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"13c5f8f3811acade0663a8efdb5c1bc7f2a3e650"},"cell_type":"markdown","source":"There are a lot of track isolation variables. This likely show how far the muon tracks are from other tracks. Basically, these tell you how clean the muon tracks are. This is important because if we want precision measurements of where things are near the vertex, we don't want our tracks to be contaminated from stray hits from other particles. If the tracks aren't sufficiently isolated, you might assign a hit to the wrong track and throw off your whole reconstruction model."},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"78c5494a401d13cb25d4a9ab3d02f447717259db"},"cell_type":"code","source":"plt.hist(signal.iso, range=(0,20), bins=20, label='Signal', alpha=0.7,normed=True)\nplt.hist(nonsignal.iso,range=(0,20),bins=20, label='Background',alpha=0.7,normed=True)\nplt.xlabel('Track Isolation Variable')\nplt.ylabel('Normalized Fraction')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"2dccd1e546552e15152d40a0dd5e3ea7eb868232"},"cell_type":"markdown","source":"Finally, the lifetime should show the estimated lifetime of each particle. We expect this took exponential (not including things like efficiency and resolution). This will be much too short to directly measure, so it is likely extrapolated by taking the flight distance, the 3-muon momentum, and the known tau mass to get a lifetime. This is calculated as:\n\n$$\\tau = \\frac{\\ell}{\\sqrt{1+p^2/m_\\tau^2}}$$\n\nwhere $\\ell$ is the flight distance, $p$ is the magnitude of the total momentum and $m_\\tau$  is the tau mass."},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"baaffd8fc694aaff72335a2deeb3e1e849189aa3"},"cell_type":"code","source":"plt.hist(signal.LifeTime, range=(0,1e-2), bins=50, label='Signal', alpha=0.7,normed=True)\nplt.hist(nonsignal.LifeTime,range=(0,1e-2),bins=50, label='Background',alpha=0.7,normed=True)\nplt.xlabel('Lifetime')\nplt.ylabel('Normalized Fraction')\nplt.legend()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"909f1b30629e9347d0fc42a5788d80e74043ca64"},"cell_type":"markdown","source":"To start out, we know from the examples that hte flight distance, lifetime, and transverse momentum should be safe variables to use that should typically pass the various tests. I've also checked that the impact parameter works as well. I'll train a simple untuned LightGBM model model."},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"8586dab08f36a8b11f93da0471584a10cbe3dbe9"},"cell_type":"code","source":"import lightgbm as lgb\nfeature_names = ['FlightDistance', 'LifeTime', 'pt', 'IP']\nfeatures = train[feature_names]\ntarget = train['signal']\ntrain_set = lgb.Dataset(features,train.signal)\nparams = {\n    'task': 'train',\n    'boosting_type': 'gbdt',\n    'objective': 'binary',\n    'num_leaves': 31,\n    'metric': {'auc'},\n    'learning_rate': 0.01,\n#    'feature_fraction': 0.8,\n#    'bagging_fraction': 0.8     \n    \n}\ncv_output = lgb.cv(\n    params,\n    train_set,\n    num_boost_round=400,\n    nfold=10,\n)","execution_count":null,"outputs":[]},{"metadata":{"scrolled":false,"trusted":false,"collapsed":true,"_uuid":"6c5c4fe07d68051b6a436a42241c48c01733e2e1"},"cell_type":"code","source":"best_niter = np.argmax(cv_output['auc-mean'])\nbest_score = cv_output['auc-mean'][best_niter]\nprint('Best number of iterations: {}'.format(best_niter))\nprint('Best CV score: {}'.format(best_score))\nmodel = lgb.train(params, train_set, num_boost_round=best_niter)","execution_count":null,"outputs":[]},{"metadata":{"collapsed":true,"trusted":false,"_uuid":"5681b5e26c92ea5ca4017bb31e3c1349c5a1fac3"},"cell_type":"code","source":"test_features = test[feature_names]","execution_count":null,"outputs":[]},{"metadata":{"trusted":false,"collapsed":true,"_uuid":"fd9530e52a402e239eaabdcfc3eb347e9882fbef"},"cell_type":"code","source":"predictions = model.predict(test_features)\ntest['prediction'] = predictions\ntest[['prediction']].to_csv('lightgbm_starter.csv')","execution_count":null,"outputs":[]},{"metadata":{"collapsed":true,"trusted":false,"_uuid":"b3d3a7390060bf997b24db74abb4b767cfc61b3c"},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.5","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}