{"cells":[{"metadata":{},"cell_type":"markdown","source":"## Introduction\nThe aim of this competition is to predict the **scalar coupling coefficient** of atoms in a molecule using certain chemical features of the molecule provided in the dataset. The scalar coupling coefficient is simply a measure of the magnetic interaction between two atoms in a molecule.\n\nIn this kernel, I will explore different features provided in the data and give a brief chemical background of each feature. I will also visualize the distributions of these features in the dataset."},{"metadata":{},"cell_type":"markdown","source":""},{"metadata":{},"cell_type":"markdown","source":"## Acknowledgements\n\nThe courtesy for the chemistry images in this kernel goes to the book [Chemistry for the IB Diploma by Steve Owen](https://www.ibdocuments.com/IB%20BOOKS/Group%204%20-%20Sciences/Chemistry/CAMBRIDGE/Chemistry%20-%20Steve%20Owen%20-%20Second%20Edition%20-%20Cambridge%202014.pdf). I would like to thank [funkyboy](https://www.kaggle.com/super13579) for [this kernel](https://www.kaggle.com/super13579/simple-eda-and-lightgbm) from which I borrowed many ideas for visulization. I would also like to thank [Chemistry LibreTexts](https://chem.libretexts.org) from which I picked up many Chemistry concepts explained in this kernel."},{"metadata":{},"cell_type":"markdown","source":"### Import necessary libraries"},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import os\nimport gc\nimport numpy as np\nimport pandas as pd\n\nimport seaborn as sns\nimport matplotlib.pyplot as plt\nimport warnings\nwarnings.filterwarnings('ignore')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Check the files provided in the dataset"},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":true},"cell_type":"code","source":"os.listdir('../input')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Load training data and get list of molecule types"},{"metadata":{"trusted":true},"cell_type":"code","source":"train = pd.read_csv('../input/train.csv')\ntrain.head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"typelist = list(train['type'].value_counts().index)\ntypelist","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Visualize distribution of scalar coupling coefficient"},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.distplot(train['scalar_coupling_constant'], color='orangered')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"From the above graph, we can see that the distribution of *scalar_coupling_coefficient* is skewed to the left, but the distribution is not perfectly unimodal. There is a significantly smaller peak close to 100 making it bimodal. The mode and mean is approximately 0."},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.figure(figsize=(26, 24))\nfor i, col in enumerate(typelist):\n plt.subplot(4,2, i + 1)\n sns.distplot(train[train['type']==col]['scalar_coupling_constant'],color ='indigo')\n plt.title(col)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Above, there are the distributions of the target for different molecule categories. Some distributions are unimodal, some heavily skewed to the left and others perfectly bimodal."},{"metadata":{},"cell_type":"markdown","source":"## Dipole Moment"},{"metadata":{},"cell_type":"markdown","source":"When two electrical charges, of opposite sign and equal magnitude, are separated by a distance, an electric dipole is established. In the case of a molecule, there is a dipole established between two atoms when there is large difference in [electronegativity](https://en.wikipedia.org/wiki/Electronegativity). This causes the electron density between the atoms to be imbalanced, creating a dipole (more electrons are attracted to the higher electronegativity molecule).\n\nThe size of a dipole is measured by its **dipole moment** (μ). Dipole moment is measured in Debye units, which is equal to the distance between the charges multiplied by the charge. The dipole moment is calculated using the formula provided below where :\n* μ is the dipole moment vector\n* qi is the magnitude of the ith charge and\n* ri is the position vector of the ith charge\n\n"},{"metadata":{},"cell_type":"markdown","source":"The dipole moment acts in the direction of the vector quantity. An example of a polar molecule is H2O. Because of the lone pair on oxygen, the structure of H2O is bent (via [VSEPR theory](https://en.wikipedia.org/wiki/VSEPR_theory)), which that the vectors representing the dipole moment of each bond do not cancel each other out. Hence, water is a [polar molecule](https://en.wikipedia.org/wiki/Chemical_polarity). The case of water is shown below. As one can see, the two dipole moment vectors do not cancel out and instead, these vectors add up to form a net dipole moment in the upward direction (shown by the arrow in the figure)."},{"metadata":{},"cell_type":"markdown","source":"
Water molecule
"},{"metadata":{},"cell_type":"markdown","source":""},{"metadata":{},"cell_type":"markdown","source":"But, when the dipole moment vectors in a molecule cancel out, the molecule becomes non-polar. For example, in CO2, which has a linear structure (180 degrees angle between the two carbon-oxygen bonds) due to the lack of a lone pair on carbon, the dipole moment vectors of each bond cancel out as they are in exactly opposite directions. Therefore, carbon dioxide is a non-polar molecule. The case of carbon dioxide is shown below."},{"metadata":{},"cell_type":"markdown","source":"
Carbon dioxide molecule
"},{"metadata":{},"cell_type":"markdown","source":""},{"metadata":{},"cell_type":"markdown","source":"### Load dipole moment data"},{"metadata":{"trusted":true},"cell_type":"code","source":"dipole_moments = pd.read_csv('../input/dipole_moments.csv')\ndipole_moments.head(10)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Visualize the distribution of dipole moments in X, Y and Z directions"},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.distplot(dipole_moments.X, color='mediumseagreen')\nplt.title('Dipole moment along X-axis')\nplt.show()\nsns.distplot(dipole_moments.Y, color='seagreen')\nplt.title('Dipole moment along Y-axis')\nplt.show()\nsns.distplot(dipole_moments.Z, color='green')\nplt.title('Dipole moment along Z-axis')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The distributions of dipole moment along the X and Y axes are approximately normal with a mean of 0, with the X-axis distribution having a greater standard deviation and range. On the other hand, the dipole moment along the Z-axis has a slightly skewed distribution (skewed to the right), with a secondary peak around 1 in addition to the primary peak (mode) above 0."},{"metadata":{},"cell_type":"markdown","source":"### Visualize the distribution of dipole moments in all directions for each molecule type"},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.figure(figsize=(26, 24))\nfor i, col in enumerate(typelist):\n plt.subplot(4,2, i + 1)\n sns.distplot(dipole_moments[train['type']==col]['X'],color = 'orange', kde=False)\n sns.distplot(dipole_moments[train['type']==col]['Y'],color = 'red', kde=False)\n sns.distplot(dipole_moments[train['type']==col]['Z'],color = 'blue', kde=False)\n plt.title(col)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"In the above figures, the orange, red and blue distributions represent the dipole moment distributions along the X, Y and Z axes respectively. They all are normal distributions with a mean of 0, but the standard deviation (\"spread\") increase form Z to Y to X."},{"metadata":{},"cell_type":"markdown","source":"## Potential Energy"},{"metadata":{},"cell_type":"markdown","source":"The [potential energy](https://en.wikipedia.org/wiki/Potential_energy) of a molecule (or any object for that matter) is the energy contained in the molecule by virtue of its position, composition or arrangement in a force field. Potential energy in the case of molecules is proportional to the sum of several forces between atoms and molecules in the substance. These forces include [covalent bonds](https://en.wikipedia.org/wiki/Covalent_bond) between atoms, [electrostatic forces](https://en.wikipedia.org/wiki/Electrostatic_forces) between oppositely charged ions and [nuclear forces](https://en.wikipedia.org/wiki/Nuclear_force) that hold the protons and neutrons in atomic nuclei together. It is these forces that constitute the potential energy of a molecule and as one can see, all these forces depend on the relative positions of all particles in the substance (atoms, protons, neutrons etc), as well as the charges and masses of these particles.\n\nFor example, the potential energy of liquid bromine is higher than that of bromine vapor. This is because the molecules of bromine vapor are further apart from each other as compared to liquid bromine molecules. This results in lower forces of attraction between the molecules of bromine vapor (**force is inversely proportional to distance2**), and thus bromine vapor molecules have a lower potential energy as compared to liquid bromine molecules. Here, **the relative positions of the molecules affect the potential energy**.\n\n
Bromine liquid
\n"},{"metadata":{},"cell_type":"markdown","source":"
Bromine vapor
\n\n\nAlso, water has a higher potential energy than carbon dioxide because there is [hydrogen bonding](https://en.wikipedia.org/wiki/Hydrogen_bond) between the hydrogen and oxygen atoms of different H2O molecules (due to the polar nature of H2O). This hydrogen bond is a strong intermolecular force that brings the molecules of water closer together. But, on the other hand, carbon dioxide only has [London forces](https://www.chem.purdue.edu/gchelp/liquids/disperse.html) between atoms of different CO2 molecules. These London forces are a lot weaker than the hydrogen bonds between water molecules, and thus water has a greater potential energy than carbon dioxide. Here, **the intermolecular forces affect the potential energy**.\n\n\n\nPotential energy need not always refer to that in molecules. It may also refer to [gravitational](https://en.wikipedia.org/wiki/Potential_energy#Gravitational_potential_energy), [electrical](https://en.wikipedia.org/wiki/Potential_energy#Electric_potential_energy) or [magnetic](https://en.wikipedia.org/wiki/Potential_energy#Magnetic_potential_energy) potential energy."},{"metadata":{},"cell_type":"markdown","source":"### Load potential energy data"},{"metadata":{"trusted":true},"cell_type":"code","source":"potential_energy = pd.read_csv('../input/potential_energy.csv')\npotential_energy.head(10)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Visualize the distribution of potential energy"},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.distplot(potential_energy.potential_energy, color='darkblue', kde=False)\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The distribution of potential energy of the molecules is approximately normal with a mean of around -400."},{"metadata":{},"cell_type":"markdown","source":"### Visualize the distribution of potential energy for each molecule type"},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.figure(figsize=(26, 24))\nfor i, col in enumerate(typelist):\n plt.subplot(4,2, i + 1)\n sns.distplot(potential_energy[train['type']==col]['potential_energy'], color = 'orangered')\n plt.title(col)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Above are the distributions of potential energy for each molecule type. One can see the distributions are very different for each molecule type."},{"metadata":{},"cell_type":"markdown","source":"## Magnetic Shielding\n"},{"metadata":{},"cell_type":"markdown","source":"Magnetic shielding is a complicated quantum mechanical concept which is difficult to explain in plain terms. But, [here](https://www.sciencedirect.com/topics/physics-and-astronomy/magnetic-shielding) is a good article explaining what magnetic shielding tensors are.\n\nEssentially, we are given 2-D magnetic shielding tensors for each molecule. Each 2-D tensor is a 3-by-3 matrix with quantities representing magnetic shielding in all possible direction combinations (XX, XY, XZ, YX, YY, YZ, ZX, ZY and ZZ)."},{"metadata":{},"cell_type":"markdown","source":"### Load magnetic shielding tensor data"},{"metadata":{"trusted":true},"cell_type":"code","source":"magnetic_shielding_tensors = pd.read_csv('../input/magnetic_shielding_tensors.csv')\nmagnetic_shielding_tensors.head(10)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Define helper function to remove outliers"},{"metadata":{"trusted":true},"cell_type":"code","source":"def is_outlier(points, thresh=3.5):\n \"\"\"\n Returns a boolean array with True if points are outliers and False \n otherwise.\n\n Parameters:\n -----------\n points : An numobservations by numdimensions array of observations\n thresh : The modified z-score to use as a threshold. Observations with\n a modified z-score (based on the median absolute deviation) greater\n than this value will be classified as outliers.\n\n Returns:\n --------\n mask : A numobservations-length boolean array.\n\n References:\n ----------\n Boris Iglewicz and David Hoaglin (1993), \"Volume 16: How to Detect and\n Handle Outliers\", The ASQC Basic References in Quality Control:\n Statistical Techniques, Edward F. Mykytka, Ph.D., Editor. \n \"\"\"\n if len(points.shape) == 1:\n points = points[:,None]\n median = np.median(points, axis=0)\n diff = np.sum((points - median)**2, axis=-1)\n diff = np.sqrt(diff)\n med_abs_deviation = np.median(diff)\n\n modified_z_score = 0.6745 * diff / med_abs_deviation\n\n return modified_z_score > thresh","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Visualize the magnetic shielding in each direction combination"},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.distplot(magnetic_shielding_tensors.XX[~is_outlier(magnetic_shielding_tensors.XX)], color='red')\nplt.title('Magnetic Shielding (XX)')\nplt.show()\nsns.distplot(magnetic_shielding_tensors.XY[~is_outlier(magnetic_shielding_tensors.XY)], color='orangered')\nplt.title('Magnetic Shielding (XY)')\nplt.show()\nsns.distplot(magnetic_shielding_tensors.XZ[~is_outlier(magnetic_shielding_tensors.XZ)], color='orange')\nplt.title('Magnetic Shielding (XZ)')\nplt.show()\nsns.distplot(magnetic_shielding_tensors.YX[~is_outlier(magnetic_shielding_tensors.YX)], color='yellow')\nplt.title('Magnetic Shielding (YX)')\nplt.show()\nsns.distplot(magnetic_shielding_tensors.YY[~is_outlier(magnetic_shielding_tensors.YY)], color='green')\nplt.title('Magnetic Shielding (YY)')\nplt.show()\nsns.distplot(magnetic_shielding_tensors.YZ[~is_outlier(magnetic_shielding_tensors.YZ)], color='blue')\nplt.title('Magnetic Shielding (YZ)')\nplt.show()\nsns.distplot(magnetic_shielding_tensors.ZX[~is_outlier(magnetic_shielding_tensors.ZX)], color='darkblue')\nplt.title('Magnetic Shielding (ZX)')\nplt.show()\nsns.distplot(magnetic_shielding_tensors.ZY[~is_outlier(magnetic_shielding_tensors.ZY)], color='indigo')\nplt.title('Magnetic Shielding (ZY)')\nplt.show()\nsns.distplot(magnetic_shielding_tensors.ZZ[~is_outlier(magnetic_shielding_tensors.ZZ)], color='darkviolet')\nplt.title('Magnetic Shielding (ZZ)')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The distributions of magnetic shielding in each direction seems to be roughly linear with a mean of 0. But, some of them have steep, sharp slopes from the peak to the tails (YX for example), while others have smooth, bulgy tails (ZY for example). All of them have unque shapes, but most of them fall under these two categories."},{"metadata":{},"cell_type":"markdown","source":"## Mulliken Charges\nMulliken charges arise from [Mulliken population analysis](http://iqc.udg.es/articles/pdf/iqc413.pdf), and are calculated using the methods of computational chemistry. They provide a means to calculate the partial charge in an atom.\n\nA partial charge is a non-integer charge value when measured in elementary charge units. Partial charge is more commonly called net atomic charge. It is represented by the Greek lowercase letter δ, namely δ− or δ+.\n\nPartial charges are created due to the asymmetric distribution of electrons in chemical bonds. For example, in a polar covalent bond like HCl, the shared electrons are drawn more towards the chlorine atom due to its higher electronegativity (greater tendency to attract electrons) as compared to the hydrogen atom. This creates a higher electron density around the chlorine atom, and thus chlorine gains a partial negative charge and hydrogen gains a partial positive charge due to the lower electron density surrounding it. These charges are relative, which means that the chlorine atom has a negative charge only because it is more negative than hydrogen. The case of HCl is shown below.\n\n
Hydrogen chloride
\n"},{"metadata":{},"cell_type":"markdown","source":"### Load mulliken charge data"},{"metadata":{"trusted":true},"cell_type":"code","source":"mulliken_charges = pd.read_csv('../input/mulliken_charges.csv')\nmulliken_charges.head(10)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Visualize distribution of mulliken charge"},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.distplot(mulliken_charges.mulliken_charge, color = 'seagreen')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The distribution of mulliken charges of the molecules peaks at around 0.175 and it is clearly unimodal. But, the distribution has a very uneven (small peaks and valleys) appearance on the tails. Additionally, there is a clear rightward skew."},{"metadata":{},"cell_type":"markdown","source":"### Visualize distribution of mulliken charge for each atom index"},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.distplot(mulliken_charges.loc[mulliken_charges.atom_index == 0].mulliken_charge, color = 'blue')\nplt.title('Atom index 0')\nplt.show()\nsns.distplot(mulliken_charges.loc[mulliken_charges.atom_index == 1].mulliken_charge, color = 'darkblue')\nplt.title('Atom index 1')\nplt.show()\nsns.distplot(mulliken_charges.loc[mulliken_charges.atom_index == 2].mulliken_charge, color = 'blueviolet')\nplt.title('Atom index 2')\nplt.show()\nsns.distplot(mulliken_charges.loc[mulliken_charges.atom_index == 3].mulliken_charge, color = 'purple')\nplt.title('Atom index 3')\nplt.show()\nsns.distplot(mulliken_charges.loc[mulliken_charges.atom_index == 3].mulliken_charge, color = 'indigo')\nplt.title('Atom index 4')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"As one can see from the distributions above, there are four peaks (quadmodal ?). The means of the distributions tend to decrease as the atom index increases (the mulliken charges are higher for lower atom indices). I am not sure about the reason for this trend."},{"metadata":{},"cell_type":"markdown","source":"That's it ! Thanks for reading this kernel. Hope you found it useful :)"}],"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}